Geothermal reservoir simulation of hot sedimentary aquifer system using FEFLOW®

The study presents the simulation of hot sedimentary aquifer for geothermal utilization. Hot sedimentary aquifer (HSA) is a conduction-dominated hydrothermal play type utilizing deep aquifer, which is heated by near normal heat flow. One of the examples of HSA is Bavarian Molasse Basin in South Germany. This system typically uses doublet wells: an injection and production well. The simulation was run for 3650 days of simulation time. The technical feasibility and performance are analysed in regards to the extracted energy from this concept. Several parameters are compared to determine the model performance. Parameters such as reservoir characteristics, temperature information and well information are defined. Several assumptions are also defined to simplify the simulation process. The main results of the simulation are heat period budget or total extracted heat energy, and heat rate budget or heat production rate. Qualitative approaches for sensitivity analysis are conducted by using five parameters in which assigned lower and higher value scenarios.


Background
As geothermal industry, especially the EGS technologies', growing in the last couple of decades, low enthalpy geothermal system is gained more attention. One of the plays in the low enthalpy geothermal system is hot sedimentary aquifer (HSA). Moeck [1] classified the geothermal plays by its geological setting. HSA is a conduction-dominated hydrothermal play type utilizing deep aquifer, which is heated by near normal heat flow.
HSA play type typically is found in passive plate tectonic where tectonic activity is relatively rare and no near active igneous body acted as a heat source. HSA play type is mainly located in intracratonic basin, orogenic belts and associated foreland basins [1]. The target reservoirs are typically high porosity/high permeability or high porosity/low permeability formations that need to be enhanced or stimulated using EGS to gain high productivity. Figure 1 shows the typical lithologies or formations which is targeted in HSA. With temperature gradient 32˚C / km as an assumption, category A is less than 3 km depth geothermal play and mostly used for district heating; while category B geothermal plays that located at more than 3 km deep and is suitable for heating and electricity, and category C is potential HDR system. Reservoir fluid of HSA is rarely coming from meteoric water, instead the water is the result of deep circulation system within the basin or major crustal faults.
Typically, doublet well is used in order to utilize HSA system. Doublet well is an open loop system, which uses an injection and production well to circulate water to and from the reservoir. One of the successful HSA geothermal plays is Bavarian Molasse Basin of the Alps, which used as heating and electricity generation in Southern Germany. Bavarian Molasse Basin is a foreland basin associated with orogeny belt. As the sedimentary sequences of foreland basin are subsiding towards the orogen due to the weight of the formations, the basin formed wedge shape and normal faults regime, thus caused the positive geothermal gradient within the basin.

Objective
The purpose of this assignment are to generate a transient HSA model using FEFLOW software package and analyse technical feasibility and the performance regarding the extracted energy from this concept. The procedures of the 3D modelling is briefly described in the next chapter. The input data for the base and modified model are also represented in the next chapter. Furthermore, the heat recovery from different scenarios (base and modified model) is compared in order to obtain the significant influences of given parameters.

FeFlow
The HSA modeling in this project assignment is performed by FEFLOW, a numerical and computation software. FEFLOW is a commercial simulation software for groundwater flow and heat transfer which is fully associated by density dependent and finite element [2]. The related processes in this software are operated in one software environment, one user interface, and one simulation model. As an intuitive program, FEFLOW is a relatively user-friendly simulation which is concurrently powerful and comprehensive for complicated subsurface modelling. The stages of FEFLOWmodelling consist of preprocessing, simulation, and post-processing [3]. Complex fluid flow and geometry are also possible to be handled as an excellence in this software.
The example of utilizations related to groundwater application are mine water management, tunneling, agriculture, and geothermal energy [3]. FEFLOW provides several built-in modules those quite applicable for Geothermal modelling such as Borehole Exchanger and Multilayer well. The approaching method for this HSA modelling is using a built-in module or a special program interface that included multilayer well modelling [4].

Procedure
In the FEFLOW software, defining reservoir area by plotting the coordinates is the first stage of the modelling. Afterwards, the well locations are defined by plotting points in the reservoir area. After generating and smoothing the mesh, defining the hydraulic head and conductivity, 3D layer configuration is created in order to define the properties of each layer. Subsequently, the transient state is chosen for this modelling instead of steady-state condition. Moreover, the injection and production are created in certain coordinates those have been determined at the previous stage. Before running the simulation, the temperature difference and simulation time have to be defined by using open loop plug-in and simulation time control respectively.

Input Data
Four basic information of reservoir, aquifer, temperature, and well are involved as significant parameters in developing this HSA model. The basic reservoir information offers four components such as area, hydraulic head, hydraulic conductivity, and porosity as represented in table 1. The values of aquifer are also represented in this table; thereby, the properties of aquifer in fourth layer can be described with these parameters. The temperature characteristic is organized into two sections for whole layers. The temperature of reservoir section is determined by temperature gradient and surface temperature as represented in table 1. Otherwise, the aquifer section is fixed by 95°C. The different temperature of injection-production well is furthermore specified by the open loop. The procedure of multilayer well implementation is described in previous chapter.
The values of longitudinal dispersivity, volumetric heat capacity of fluid, thermal conductivity of fluid are inputted identically as the FEFLOW's default values. Moreover, the values of radius, well depth, well distance, and flow rate have to be defined respectively as shown in the table 1. These data are further simulated and used as basic model, which is compared to other models with modified parameter. In order to simplify the modelling process, following assumptions were made such as: 1) reservoir is considered to be homogenous and isotropic; 2) the inflow temperature at a certain depth is assumed to be equal to the corresponding initial temperature; 3) only the heat transport is simulated and the mass transport is considered in steady state; 4) the flow rate of well is pumped constantly for 10 simulation years.

Heat Period Budget
Heat period budget represents produced thermal energy of HSA system after 10 simulation years. Negative values of heat period budget are obtained after running the simulation, since FEFLOW defines the extracted energy as negative value and the storage energy as positive value. Otherwise, these values were expressed in Joule, which were inappropriate due to the huge numbers. Thus, the authors decided to convert the Joule to Gigawatthour unit. The heat period budget for base model showed approximately 15.5 GWh after 10 years as represented in figure 2.

Figure 2. Heat Period Budget after 10 Simulation Years
Flow rate was the most significant parameter, that influenced the produced heat period budget. The heat period budget reached almost 20 GWh, when the high scenario was applied in this HSA system. Inversely, low scenario of flow rate for this system caused a less heat period budget. Therefore, the flow rate was direct proportional to the heat period budget. The depth played also a significant role for this HSA system. Based on lower and higher scenario in figure 2, it could be concluded that the deeper the aquifer depth could generate more thermal energy of HSA system. The hydraulic conductivity and thermal conductivity of solid were also direct proportional to the heat period budget. The relationship of those two parameter is represented also in figure 3. In contrast, the heat period budget was inverse proportional to the porosity. By comparing between the low and high scenario to the base model, the highest heat period budget was generated by the low scenario.

Heat Rate Budget
One of the results obtained from the simulation is heat rate budget. Heat rate budget is the amount of thermal power extracted by the production wells. FEFLOW plotted the result in a graph of thermal power against the simulation time. As mentioned in previous chapter, this study simulated five different parameters by three different scenarios, the heat rate budget of each parameter are described as follows. Graphs of the comparisons are attached in the appendix of the paper.

Comparison by The Depth
The first parameter that was investigated is the depth of aquifer. Appendix 1 shows the comparison between three depths of aquifer scenarios. The base model scenario was aquifer with thickness of 500 m and bottom depth at 2000 m, the second scenario was 1500 m bottom depth with 500 m thick aquifer layer, and the last scenario was model with bottom depth of aquifer at 2500 m. At 10 years of simulation time, the production well at base model scenario performed at around 170 kW, while the low scenario extracted 154 kW thermal power and the high scenario obtained higher value at 196 kW. The high scenario with 0.0034 kW / day was the smallest decline rate of thermal power from the production well than other scenarios.

Comparison by Porosity
Porosity is the second parameter that was analysed. Appendix 2 shows the comparison between three scenarios of porosity of the aquifer. The porosity in the aquifer was applied throughout the layer since the aquifer layer considered as homogeneous and isotropic. The base model scenario was aquifer with porosity of 0.3, the value of porosity in the low scenario was 0.2 and the last scenario was the high scenario with porosity at 0.4. There were only small differences, around 2 kW, between one scenario and the other. After 10 years, the thermal power produced at low, base and high scenarios were 173, 171 and 169 kW, respectively. The production well performed relatively similar with decline rate difference of 0.03 x 10-5 kW / day to each scenario.

Comparison by Flow Rate
The next parameter that was examined in this study is flow rate of the well. As mentioned in modelling section, the flow rates on the injection well and production well were applied at constant and same value for 10 years of simulation time but at different direction. Appendix 3 shows the comparison between flow rates of the production well at three different scenarios. The base model scenario was 60 m3 / day, the low scenario was 40 m3 / day, and the high scenario was 80 m3 / day. After simulated for 10 years of simulation time, the production well at base model scenario extracted 171 kW, while for the low scenario was 123 kW and the high scenario obtained thermal power at 216 kW. However, the high scenario performed relatively high decline rate at 2.563 x 10 -5 kW / day, meanwhile the low scenario had decline rate of 1.934 x 10 -5 kW / day and the base scenario continuously declined at 2.303 x 10 -5 kW / day.

Comparison by Hydraulic Conductivity
The fourth parameter which was focused on in this study is the hydraulic conductivity of the aquifer. Since the model was simplified and the layers considered as homogeneous and isotropic, the hydraulic conductivity of the aquifer had the same value for all axis. Appendix 4 shows the comparison between hydraulic conductivity of the aquifer at three different scenarios. The base model scenario was 2 m / day, the low scenario was 1 m / day, and the high scenario was 3 m / day. After simulated for 10 years of simulation time, the production well at base model scenario obtained 171 kW thermal power, while for the low scenario and the high scenario were 161 kW and 178 kW, respectively. Even though the decline rate of all the scenarios had only small differences, around 0.03 x 10 -5 kW / day, the high scenario performed relatively small decline rate at 2.067 x 10 -5 kW / day, meanwhile the low scenario and the base scenario had decline rate of 2.623 x 10 -5 kW / day and 2.303 x 10 -5 kW / day, subsequently.

Comparison by Thermal Conductivity of Solid
The last parameter that was examined in this study is the thermal conductivity of solid of the aquifer. This parameter was applied homogeneously on the aquifer layer to simplify the model. Appendix 5 shows the comparison between thermal conductivity of solid of the aquifer at three different scenarios. The base model scenario was 3 J / (m.s.K), the low scenario was 2 J / (m.s.K), and the high scenario was 4 J / (m.s.K). After simulated for 3650 days of simulation time, the production well at base model scenario obtained 171 kW thermal power, the low scenario gained 164 kW, and the high scenario extracted 177 kW thermal power. The production well thermal power declined at a small difference rates with base scenario at 2.303 x 10 -5 kW / day, low scenario at 2.446 x 10 -5 kW / day and high scenario at 2.164 x 10 -5 kW / day.

Temperature Distribution
One of the results of the simulation is the temperature distribution within the model and the wells. As discussed on modelling section, the temperature boundary condition and initial temperature inside the model was defined before the model was simulated. After 3650 days of simulation, the temperature distribution within the model was altered, especially near the wellbore region. Only temperature distribution of base model is presented in this study since there are only small differences of temperature distribution between all of the scenarios that were simulated. Figure 4 is 2D view of the model at 2000 m depth slice. This figure shows temperature distribution of the model at the deepest part of the model. After 10 years of simulation time, there were only small area, around maximum radius of 75 m from the center of the wells, was affected by the injection wells. The lowest temperature was at the very center of the well at 8.8˚C. Meanwhile, the lowest temperature on the production well was around 80 ˚C and did not affected the near well bore region of the well. Figure 5 shows the 3D view of the temperature distribution of the base model. This figure represents that after 3650 days of simulation time, the injection wells caused more disturbances in temperature in the model than the production well. However, the disturbance developed only on the near region of the well while the rest of the model was at steady state condition. Figure 6 shows the base model's average temperature of production and injection well. After 10 years of simulation time, the average temperature of production well was 69 ˚C and the injection well average temperature was 9˚C. The production well was able to maintain high enough temperature to be extracted for longer time.

Discussion
As represented in previous chapter, the greater depth of the aquifer showed higher heat production rate with lower decline rate than the base scenario, which had a direct correlation with thermal gradient. Higher flow rate also generated a higher heat production rate similar to the depth correlation. On that note, the decline rate had to be considered, since the higher flow rate showed a higher decline rate too. In contrast, the porosity expressed a different proportionality to the heat production rate. The lower porosity showed higher heat production rate with lower decline rate. This behaviour could be occurred due to low thermal conductivity of the fluid compared to thermal conductivity of solid in the aquifer. Thus, conduction had higher effect on heat transfer from the solid in the aquifer than the convection of the fluid.
Thermal conductivity of solid and hydraulic conductivity was indeed dependent on the lithology of the aquifer. The higher value of both parameters could conduct more thermal energy to the surface, thus higher conductivity could generate higher heat production rate. The sensitivity analysis was performed only to identify the qualitative relationship between the heat production rate and those parameters.
The temperature distribution in the reservoir is only disturbed for a short distance, around 75 meters, after 10 simulation years. The production well could still generate for longer time until this well had to be shut off, in order to replenish the heat in the aquifer.
The generated thermal energy from this HSA system was not able to be converted into electrical energy, since the average temperature was only circa 70°C. The binary power plant cycle requires the minimum temperature of 85°C at the surface in order to generate electrical energy. If the average temperature exceeded 85°C [5] , the extracted energy was also not economically beneficial, since the thermal energy and the efficiency of binary cycle were to low (ca. 10-13%) [3]. Therefore, this HSA system is more suitable to utilize for heating system, such as district heating. However, the surface temperature from production well cannot be defined clearly. Hence, the extracted heat power from this system cannot be estimated accurately.
The mass balance of this fluid is not considered in this simulation, due to produced geothermal water is injected back to the reservoir, therefore there is an uncertainty regarding the fluid movement in the aquifer which correlates to the heat transfer. Otherwise, lack of published paper related to HSA simulation may cause difficulty in finding reliable data for this study's references.

Conclusion
The base model of transient HSA simulation was built by using FEFLOW and generated almost 15,5 GWh with heat production rate at 171 kW, which declined exponentially at 2.303 x 10 -5 kW / day. The generation of electrical energy from the thermal energy was not feasible, since the heat and temperature could not fulfil the minimum condition. Therefore, only heat utilization might be applied from this HSA system. This study has examined the influence of several significant parameters to the HSA system. The production well in this case could produce higher heat production rate by having better aquifer, which has the following criteria: deeper depth, higher thermal conductivity of solid, higher hydraulic conductivity, and lower porosity. These parameters are empirical properties which are possessed by the aquifer formation or lithology layer. When the one of these parameters could not be improved, increasing the flow rate could be an alternative to generate higher heat production rate. However, the increasing flow rate could cause faster decline rate.
The simulation of HSA simulation can be investigated further in term of using more reliable data or performed sensitivity analysis by using quantitative approaches.