Modeling thermal front dynamics in geothermal reservoirs using an open-source MRST–MATLAB simulation framework

1. Introduction

The depletion of oil and gas fields becomes evident when extraction is no longer economically viable, often resulting in the abandonment of wells and a transition in industry practices. As fields mature and production declines, the cost of extracting the remaining resources increases, further reducing economic feasibility and leading to well abandonment. In response, the repurposing of oil and gas wells for geothermal energy has emerged as a modern practice, utilizing existing infrastructure and reservoirs for low-temperature direct applications as well as electricity generation. This conversion enables the continued use of established wells, drilling equipment, and support facilities, thereby lowering upfront investment and minimizing environmental impact. Such an approach offers dual advantages: it mitigates the environmental risks associated with abandoned wells while simultaneously advancing a sustainable energy alternative to fossil fuels. In early 1990s the exploration in Cambrian sandstone reservoirs in Baltic Basin including Lithuania resulted in oil production from small fields. Now these reservoirs are depleting and had an opportunity for Enhanced Oil Recovery (EOR). Moreover, the plans are to repurpose these reservoirs for Geothermal Energy (GE) [1]. Although not a major production source, these sandstones remain important for exploration, hydrocarbon recovery, and geothermal development. Geothermal energy depends on the hot water production rate and in turn on the borehole temperature. Hence, knowledge of accurate underground temperature is very important for planning geothermal plants because it improves performance and make energy production more cost-efficient. Thus, measuring subsurface temperatures is important for both science and practice which helps guide the design and optimization of geothermal energy extraction. This increases efficiency and reduces environmental impacts. Therefore, numerical simulation is a key tool for evaluating geothermal reservoirs, but it often simplifies complex Multiphysics processes such as fluid flow, heat transfer, and geochemical interactions. While effective, these simplifications limit predictive accuracy, especially in heterogeneous systems. Recent advances in machine learning (ML) and artificial intelligence (AI) offer new possibilities for upscaling reservoir properties, modeling subsurface flow, and forecasting geothermal performance. This research combines physics-based simulations with ML methods to develop hybrid workflows for improved geothermal modeling. The MATLAB Reservoir Simulation Toolbox (MRST), with its geothermal extension, provides an open-source, flexible, and transparent platform to achieve these objectives. MRST allows not only to construct homogenous and heterogenous geothermal models but also helps in generation of synthetic datasets for ML training and validation. Moreover, it also helps in benchmarking of surrogate models against established numerical solvers.

2. Purpose and methodology

The primary objective of this study is to evaluate the potential for repurposing depleted oil and gas reservoirs for geothermal energy production. To address this objective, data from 19 reservoirs were analyzed, and 5 of the most promising sites were identified based on their capacity to deliver significant heat output and maximize power generation [1-4]. The selected sites are described in the Table 1. Our previous studies demonstrated that geothermal power output can be improved by optimizing well design and operating conditions [1-2]. Using a commercial reservoir simulator, vertical wells were shown to provide baseline production, while horizontal wells significantly increased water yield and power generation at the same temperature difference [3-4]. Further enhancement was achieved by lowering the reinjection temperature, which increased the temperature difference and resulted in a substantial gain in power output compared to previous models. The use of a commercial simulator ensures that these findings reflect realistic reservoir behavior (True behavior), supporting their applicability to field-scale geothermal development.

Similarly, the use of MRST can model both vertical and horizontal wells by defining well trajectories in the grid. Moreover, it also allows setting of different completions, well indices, and injection or production rates to compare well performance. The MRST can calculate production rates from reservoir pressure and permeability and can combine flow and temperature to estimate enthalpy changes and power output. The thermal module of MRST makes it possible to simulate heat transport and test different reinjection temperatures. By lowering reinjection temperature, MRST can show how the temperature difference increases and how this affects energy production. Moreover, it also has the capability to predict the temperature differences by varying the permeability both in homogenous and heterogenous scenario. Furthermore, it also supports scenario studies, making it possible to compare vertical and horizontal wells or different reinjection strategies and quantify improvements in flow and power [4-5]. While MRST is research-oriented (not a commercial simulator), it can capture the same physics of Multiphase flow, heat transfer, well control, and reinjection. It is also flexible and transparent, so we can validate our findings before moving to a commercial simulator.

In our study the numerical model is defined over a two-dimensional (2D) vertical domain of LENGTH × BREADTH (m²). Cold water at 40 °C is injected along the left boundary (INJECTOR), creating a thermal front that propagates toward the right-hand production side (PRODUCER). Permeability and porosity are first prescribed as homogeneous values, then extended to heterogeneous cases using external text files (.txt). Thermal properties are added to the rock object, which may be defined as scalars, arrays, or function handles. The reservoir is initialized with uniform pressure and temperature, extended by an explicit temperature field (T). Flow is driven through forcing terms such as wells, boundary conditions, and sources. Cold-water inflow is imposed via fixed pressure and temperature at the left boundary. Outflow is prescribed at the right boundary with pressure and temperature equal to initial reservoir values. Simulation is advanced through timesteps increasing up to 30 days, for a total simulated period of 50 years as depicted in Fig. 1. The properties used in the model is described in Table 1.

3. Results

The preceding section outlines the model setup implemented in both T-Navigator and MRST. The predicted thermal front as an output obtained from T-Navigator was compared with that generated by MRST for the Vilkyciai reservoir depicted in Fig. 1. As shown in Fig. 1, both toolboxes predicted that the thermal front reached the same distance (500 m) within the 50-year simulation period. The results confirm that the MRST model accurately reproduces the thermal front obtained from the commercial simulator. Thus, this analysis benchmarks the MRST workflow as a reliable tool for thermal front generation for other geothermal sites. Following this workflow, additional scenarios were developed for the remaining sites as shown in Fig. 2. Fig. 2 presents the evolution of the thermal front under both homogeneous and heterogeneous conditions for all screened sites, while Fig. 3 illustrates a comparison of the true temperature profiles for both scenarios over a 50-year period. From Fig. 3 it is demonstrated that the MRST workflow produces similar temperature profiles in homogeneous and heterogeneous systems, and the differences observed are attributes due to reservoir heterogeneity.

Fig. 1Model setup in T-Navigator (3D) and MRST (2D)

4. Discussion

This study focuses only on Lithuanian aquifer complexes, so the results are based on local geological and thermal conditions. However, the developed MRST workflow is universal and based on basic principle properties. The same workflow can be used for any geological site by changing the porosity and permeability distributions. The thermal front behavior mainly depends on fluid flow and heat transport processes that occur in most subsurface systems. Therefore, the results represent a general modeling approach rather than conclusions limited to one site. In future work, this workflow can be applied to other available geological sites in Lithuania with different reservoir properties. This will help evaluate the robustness of the workflow under diverse subsurface conditions.

Fig. 2Thermal front determination (2D) for homogenous and heterogenous scenarios for screened sites

a) Diegliai

b) Genciai

c) Nausodis

d) Siuparai

e) Vilkyciai

Fig. 3Comparison of temperature prediction from homogenous and heterogenous scenarios for screened sites

a) Diegliai

b) Genciai

c) Nausodis

d) Siuparai

e) Vilkyciai

5. Conclusions

From the above study, it can be concluded that the 2D geothermal model developed in MRST effectively captures the evolution of the thermal front, providing insights into the reservoir’s thermal response. The results demonstrate that the thermal front propagates toward the production side in both homogeneous and heterogeneous reservoir property scenarios represented within the model. Moreover, slight variation in true temperature prediction is observed from both homogenous and heterogenous scenarios for all the screened sites, when comparing with the commercial simulator T-Navigator. Furthermore, the workflow is able to demonstrate the ability of MRST to generate reliable synthetic geothermal datasets, and such simulations provide valuable insights into geothermal production and storage efficiency through modeling framework which is well-suited for benchmarking and coupling with machine learning approaches.