Qualitative assessment of optimizing the well spacings based on the economic analysis

The design of well spacing significantly influences the sustainability and economic benefit of geothermal energy extraction. However, most studies have predominantly employed heat production-related parameters as indicators of well spacing, and a comprehensive analysis of well spacing design based on an economic model is necessary for practical implementation. In this study, an economic indicator considering the benefits derived from heat production and operating costs is proposed and applied in the Caofeidian, a typical abandoned oilfield in the Bohai Bay Basin. It offers a refined portrayal of directional wells, moving beyond rudimentary representations, to capture their appropriate degree of complexity and behavior in drilling configurations. First, by integrating thermophysical information and site investigation data from previous oil investigations, a heterogeneous 3D model is constructed to forecast the 30-year temperature and pressure evolution. Then, a modified levelized cost of heat (LCOH-HT) is proposed to perform economic analysis in optimizing the well spacing, revealing an optimal range of 300–600 m for the different selected wells. In comparison with results derived solely from heat production considerations, drilling and pumping costs contribute to a 300 m reduction in the optimal well spacing based on the proposed approach, as a larger well spacing leads to increased hydraulic losses and drilling cost, necessitating greater pumping efforts and costs. This finding under-scores the need to balance economic and thermal considerations. In addition, we found the difference in the optimal well spacing in space is also caused by the porosity variations. Porosity affects fluid temperature and pressure, leading to changes in the benefits and costs associated with pressure fluctuations. Notably, this novel economic analysis method is not limited to spacing optimization; it can also be used to optimize operating parameters, such as the flow rate, which could provide practical strategies for geothermal energy extraction.


Introduction
Geothermal energy, as a clean energy resource, is experiencing rapid growth in China (Wang et al. 2015;Zhang and Hu 2018) and facilitating the nation's journey towards carbon neutrality.In this context, determining how to effectively and sustainably use geothermal energy is a common concern for scientists (Anderson and Rezaie 2019;Kaygusuz and Kaygusuz 2004;Shortall et al. 2015;Soltani et al. 2021).The utilization of geothermal resources typically involves extraction and reinjection processes.Reinjection is required not only for extracting geothermal energy but also for maintaining reservoir pressure and plays a pivotal role in determining the sustainability of a geothermal system (Kaya et al. 2011;Liu et al. 2018;Quinao and Zarrouk 2018).Therefore, considerable efforts have been devoted to exploring optimal reinjection strategies (Aliyu and Chen 2017; Chen and Jiang 2015; Kong et al. 2017;Liang et al. 2018;Liu et al. 2020), including the well spacing design and the optimization of the reinjection flow rate and temperature.The design of well spacing takes precedence and demands careful attention; the significance of this factor lies in its substantial impact on the operational efficiency and longevity of geothermal systems, as well as its direct connection to the initial investment, because it is unchangeable once the drilling process is completed.
Many optimization objectives have been proposed to determine the optimal well spacing.And most related studies have focused on heat performance-related indicators (Ren et al. 2023), such as thermal breakthrough and reservoir temperature changes, which influence the heat performance of a geothermal system.While heat performance is an important factor to consider, an economic analysis that accounts for costs is also necessary for obtaining integrated information related to the production of geothermal energy.Kong et al. (2017) performed cost estimation and accounted for the increased costs resulting from temperature changes and water table drawdown during the extraction process, whereas the drilling costs were ignored due to their use of a 2D model.Even in the 3D model, geothermal wells are sometimes treated as simple straight wells (Wang et al. 2019).However, in practice the surface area, access, environment, license, and permit are limited, so that many wells are designed directionally to reach the reservoir (Marbun et al. 2021;Wang et al. 2020).The application of directional wells significantly affects reservoir parameters and therefore should not be treated simply as point or straight wells in the model (Blöcher et al. 2010).In addition, the costs of geothermal projects are mainly focused on drilling (Allahvirdizadeh 2020;Kipsang 2015), which cannot be neglected in the economic indicators.In comparison, the levelized cost of heat (LCOH) employed in technoeconomic simulations provides a comprehensive indicator for cost calculation.The LCOH thoroughly accounts for three primary aspects: the total upfront capital investment, the operation and maintenance (O&M) costs and the quantity of heat generated (Armstead and Tester 1987;OECD/IEA, 2015).This model incorporates various parameters, including the cost of wells, the expenses for the field data collection system, and the O&M costs for make-up water and pumping, etc. (Beckers and McCabe 2019).However, it is important to note that these parameters are predominantly tailored for enhanced geothermal systems.For instance, the O&M costs for make-up water are specific to these systems.Therefore, when applying the LCOH to hydrothermal systems, adjustments are crucial, especially in the context of altered pumping costs resulting from pressure variations.
Since models are generally established to predict changes in key parameters (temperature, pressure, etc.), the accuracy of these models is pivotal for determining the optimal well spacing.Many researchers have predicted temperature changes after long-term operation based on the assumption of reservoir homogeneity.However, it has been demonstrated that numerous reservoirs, even those conventionally considered homogeneous, such as sandstones, exhibit non homogeneity, significantly influencing the heat transport process (Babaei and Nick 2019;García-Valladares et al. 2006).For example, Willems et al. (2017) showed that reservoir heterogeneity has a significant effect on connectivity and heat production in geothermal systems.The limited availability of reservoir characterization stagnates accurate modelling.Furthermore, a considerable number of models are focused primarily on temperature variations, neglecting the impact of well spacing variations on pressure.Contrary to this simplified perspective, the actual geothermal mining process involves complex, coupled thermohydro-mechanical phenomena (Guo et al. 2021;Huang et al. 2021;Pandey et al. 2018;Xu et al. 2015).Pressure variations, a critical component of this process, can directly impact the cost of subsequent operation and maintenance.
Therefore, in this work, the Caofeidian geothermal oilfield, an abandoned oilfield, was chosen as the study area.Sufficient well logging information and monitoring data were obtained to establish a complex three-dimensional heterogeneous numerical model.Then, a novel optimization indicator based on the LCOH, encompassing the benefits of heat production and operating costs, was proposed to optimize the well spacing.Finally, the economic objective proposed by previous model was used for comparative analysis, and conclusions were obtained regarding the diverse factors affecting the optimal well spacing.

Methods
Our procedure comprised two steps, as illustrated in Fig. 1.Initially, 3D numerical models were established to simulate heat extraction and reinjection processes and predict hydraulic and temperature variations.Subsequently, a novel economic analysis method was introduced and applied to the study area to optimize the well spacing.This section provides an overview of the theoretical model, including the governing equations and the economic analysis method.

Numerical model
The OpenGeoSys (OGS) simulator was used to establish 3D models to predict hydraulic and temperature changes during fluid extraction and reinjection.OGS can couple the fluid flow and heat transport processes in a reservoir based on a flexible finite element method.In the OGS simulator (Huang et al. 2023), the wells themselves are treated as line elements, and the water and heat flow processes in wells, including fluid circulation and the associated heat transport, are simulated based on governing equations for these 1D vertical line elements.For the reservoir part, 3D prism elements are used to discretize the different sediment layers.
The water flow in porous media was described by the Darcy's continuum equations, expressed as where φ is the porosity, ρ l is the density of the fluid (kg/m3), Q is the fluid source/sink term (kg/m 3 /s), v l is the fluid velocity (m/s) can be expressed by Darcy's law: where K is the intrinsic permeability tensor (m2), µ is the fluid dynamic viscosity (Pas), p is the liquid pressure (Pa).
The heat transport equation takes into account both advective and diffusive fluxes: where h l is the specific enthalpy (J/kg), is the heat conductivity of the porous medium (W/m•K), Q heat is the heat source/sink term (W/m 3 ), H t is the total enthalpy (J/m 3 ), including both fluids and rocks: (1) where ρ R is the density of the rock (kg/m 3 ), c pR is the specific heat capacity of the porous medium (J/kg•K).

Economic analysis
A new economic evaluation indicator was proposed to analyze the economics of geothermal extraction based on the standard LCOH (OECD/IEA, 2015, Beckers and McCabe, 2019).LCOH refers to the cost of heat production per unit of heat (Olasolo et al. 2016) and is expressed as where CCAP is the total capital cost, COM t is the annual O&M cost,E Pro,t is the annual direct energy generated, i is the discount rate, and LT is the geothermal field lifetime.CCAP is calculated as follows: where C d , C f , C s and C h are the well drilling cost, field data collection system cost, sur- face plant cost and heat exchange station cost, respectively.E pro is the produced thermal energy over time and is expressed as where T p,t and T i,t are the temperature of the production fluid and the reinjection fluid, respectively, at time t (℃), Q pro is the geothermal fluid production rate (m 3 /yr), c l is the specific heat capacity of the fluid (J/kg•K), t l is the operation period of the well (yr).
Here, we extended and modified the model so that it can be applied to estimate the production cost of hydrothermal resources.Given the conditions of hydrothermal resources and key factors affecting the sustainable development of geothermal exploitation, the following factors were included in the cost estimation process.
This study offers a refined portrayal of directional wells, moving beyond rudimentary representations, to capture their appropriate degree of complexity and behavior in drilling configurations.
(1) Consequently, the drilling costs escalate with both the increase in drilling depth and the expansion of well distance, as expressed in Eq. ( 8): where n is the drilling cost per meter (CNY/m), L p is the length of production well (m), L s is the spacing between production and reinjection points (m).
(2) COM t is the sum of annual labor expenses, the annual depreciation and mainte- nance of equipment and plants, and pumping O&M costs and is written as follows: (4) where C l is the annual labor expenses cost, Cp is the annual pumping O&M costs and can be quantitatively calculated, and Ce is the annual equipment and plant depreciation and maintenance cost.The annual equipment and plant depreciation and maintenance cost is assumed to be the 2.5% of heat exchange station cost and 1% of drilling well cost.
(3) The pumping costs C p in Eq. ( 9) are strongly dependent on the pressure changes in the wells, where the effects of friction are taken into account and are applicable to inclined shafts the relationship can be expressed as: where D is the depth of the production well (m), θ is the angle of inclination of well to the vertical, P is the pressure at the bottom of the production well (Pa), P 0 is standard atmospheric pressure (Pa), f D is the Darcy friction factor, L is the length of the well (m), d is the hydraulic diameter of the well (m), g is gravity, pp is the price of electrical power and the η is the efficiency of pump use and set 0.7.
Finally, the modified LCOH, defined as LCOH-HT, is expressed as follows: Notably, the actual proportionality constant relating the LCOT-HT with the hydraulic pressure and temperature changes, underscoring the necessity for accurate predictive models.

Study case
In our application case, the Caofeidian geothermal oilfield located north of Tangshan, Hebei Province, China, was selected.In terms of geological structure, the Caofeidian geothermal oilfield is located north of the Nanpu Sag in the northcentral part of the Bohai Bay Basin (Fig. 2) (Dong et al. 2021;Li et al. 2015;Wang et al. 2022).According to the project report, there were 20 wells in the study area for heating during the heating period.Each production well was equipped with a reinjection well for recharge during the non-heating period.The depth of the target reservoir ranged from 1500 to 3300 m, and the lithology in the area is sandstone from the Guantao Formation.
A previous site investigation revealed the porosity distributions (Fig. 3a), which were derived from the petrophysical interpretation of well log data (Huang et al. 2022).The permeability model was also built based on permeability log data (Fig. 3b).Heterogeneous porosity and permeability distributions were applied in the model to closely reflect real conditions and produce accurate results.
Groundwater monitoring was used to determine the well temperatures and pressures, and which the temperature and pressure distributions in the study area were obtained via the interpolation method (please refer to Fig. 3c, d).The groundwater flowed from the northwest to southeast at an average flow rate of 0.2 -0.4 m/yr.In this case, the impact of the natural groundwater flow was negligible.Geological investigation revealed that the heat flow in the case study varied from 50 to 74 mW/m 2 , while the thermal (9) gradient ranged between 2.6 and 3.6 ℃/100 m (Huang et al. 2022).The above situation is summarized in Table 1.
The six wells, namely, GR2, GR3, GR4, GR5, GR6 and GR8, were chosen to simulate and optimize the spacing between the production and injection wells.The locations of the wells are shown in Fig. 2d, and the wells were far from each other and temperature and pressure changes around wells did not affect each other during the simulated extraction and reinjection processes.

Initial and boundary conditions
As mentioned above, the initial pressure and temperature distributions were derived from well logs.For the upper model boundary, we applied a constant temperature and pressure distributions as Dirichlet boundary conditions.For the bottom boundary, a Neumann boundary condition was applied to account for the heat flow in the study area.The heat flow boundary value was 0.062 W/m 2 , which was the average heat flow in the study area.All the lateral boundaries were assumed to be thermally isolated.In addition, the flow rate and temperature changes associated with reinjection were considered when setting the boundary conditions.

Model setup
The reservoir domain was designed to be 7000 m long, 5000 m wide and 396 m deep (Z: −1929 -−2325 m asl).Several preliminary simulations were conducted with different thermal loads to ensure that the domain size was sufficiently large and that the thermal plume produced by the production and reinjection activities would not interfere with  The size of each typical grid element was 100 m*100 m*10 m; however, the mesh surrounding each well area was further refined, and these elements were designed to be 50 m in length, 50 m wide, and 3 m deep.Several simulations were conducted with even more refined meshes, and the convergence of the simulated groundwater level and temperature distributions was confirmed.The authors are confident that mesh independence and time-step independence were achieved for the simulation results presented in this work.The thermophysical parameters of the subsurface are shown in Table 2. Due to the high temperature and pressure in the reservoir, variations in the fluid density and viscosity resulting from the temperature and pressure changes were accounted for in the simulations with an equation of state according to the IAPWS (Wagner and Kretzschmar, 2008), aiming to better simulate groundwater flow conditions and improve the accuracy of the model.

Model validation
Prior to applying the model for the prediction of long-term thermal impacts, a model validation procedure was implemented to verify the correctness of the mathematical model.

Temperature and pressure changes of different well spacings
To analyze the economics of different well spacings, the temperature and pressure changes for different distances between the production and reinjection wells were investigated first.The temperature and pressure changes for distances from 100 to 900 m are shown in Figs. 5 and 6.
The temperature evolution of all the wells exhibited periodic fluctuations.During the heat extraction period, the water temperature decreased until heat extraction ceased, followed by recovery upon receiving recharge from the surrounding reservoir and a subsurface heat flow.The temperature drop decreased as the well spacing increased for all wells.However, variations occurred among the different wells, with GR6 experiencing a maximum decrease of 37.8 ℃ at a spacing of 100 m, while GR3 exhibited a minimum decrease, reaching 12 ℃.As the well spacing increased to 500 m, the temperature drop of GR6 decreased to ~ 7 ℃, and for well GR2, it slightly decreased to 0.6 ℃.As the distance increased from 600 to 900 m, the temperature drop was basically unchanged, and at 900 m, the temperature drop was basically negligible.A large well spacing allows recharged water more time to absorb heat in the thermal reservoir during the reinjection process, facilitating a more complete exchange of heat with geothermal fluid (Kaya et al. 2011;Kong et al. 2017).The different temperature trends between different wells were attributed to the positional differences between the wells.Due to the heterogeneity of the reservoir, the groundwater flow conditions vary for each well, leading to different outcomes.
The pressure changing with time at different distances from the production wells also demonstrated a periodic evolution.Similar to the temperature variation, the pressure experienced a decline during extraction phases followed by a recuperative trend during reinjection cycles.The average annual pressure in production wells continued to decrease over long-term operation (Huang et al. 2023).With a well spacing of 100 m, the pressure drop in each production well ranged between 0.03 and 0.14 Mpa (equivalent to a water table drop of 3.1 -14.3 m) after 30 years of operation.As the spacing was increased to 900 m, the pressure drop in each extraction well ranged of 0.04 -0.26 MPa (equivalent to a water table drop of 4.1 -26.5 m).With increasing well spacing, the hydraulic loss continued to increase, and the injection water overcame the increasing hydraulic resistance to reach the production wells.
In brief, the model results revealed that well spacing significantly influences both the thermal and hydraulic performance of geothermal systems.Expanding the well spacing elevates the outlet temperature, which gradually approaches the reservoir temperature.As a result, the heat production and thermal efficiency of the geothermal system are enhanced.However, enlarging the well spacing decreases the outlet pressure and some studies have shown that this approach is not conducive to maintaining reservoir pressure and can cause damage to reservoirs (Kaya et al. 2011;Liu et al. 2018).

LCOH-HTs of different well spacings
The LCOG-HT method was employed to identify the optimal well spacing, and Fig. 8 illustrates the changes in the LCOH-HT, total heat production, and cost at different well distances.The total heat production and cost were selected for comparative analysis.With increasing well spacing, the total heat production over 30 years exhibited an initial upwards trend followed by a nearly constant phase.This is because larger well spacing allows for more efficient absorption of heat by the reinjection water.Regarding costs, there was an increase associated with increased well spacing, attributed to the decrease in pressure and increase in drilling length, leading to elevated pumping and drilling costs.In addition to heat production, the LCOH-HT emphasizes the importance of economic considerations, aiming to strike a balance between benefits and costs.As depicted in Fig. 7, the LCOH-HT value reached a minimum at certain well spacings, corresponding to the optimal well spacings.Specifically, the optimal well spacings were 400 m for GR3 and GR5 and 600 m and 500 m for GR6 and GR8, respectively.If the total heat production alone was used as an indicator, the optimal well spacing aligned with the minimum distance at which the thermal output essentially remained constant.It was evident that the optimal well spacing determined based on heat production for each well was greater than or equal to the spacing determined based on the LCOH, given the omission of drilling and pump costs.
Furthermore, based on the LCOH-HT, notable disparities existed in the optimal well spacing among the different wells, extending up to 300 m.These differences could be attributed to variances in the deployment locations of wells, primarily stemming from variations in porosity.An analysis of the porosity distribution revealed a high-permeability zone between the production and injection wells of GR2, with an average porosity of 0.20 and a minimum optimal spacing of 300 m.In contrast, GR6 displayed the lowest average porosity of 0.13, coupled with the highest optimal well spacing.In formations with high porosity, groundwater flow velocity accelerates, leading to decreased heat absorption during the transfer from reinjection wells to production wells, consequently affecting heat production.However, faster flow velocity aids in reservoir pressure restoration, reducing the increase in pumping cost due to pressure drop.Therefore, despite the study area being composed of relatively homogeneous sandstone, the porosity varies from 0.13 to 0.20, resulting in a substantial discrepancy in results for an optimal well spacing of 300 m.

Comparison of LCOH-HT with previous model
To verify whether our model will have a positive impact on production, the economic objective function proposed by Kong et al. was used for comparison.The previous model considers the additional cost of electricity due to hydraulic head drawdown in the production well and the extra cost caused by thermal breakthrough and the optimum well spacing yielded the lowest cost.Figure 8 shows a comparison of the results.The optimal well spacings determined with the previous model were 500 m, 900 m, 900 m, 600 m, 800 m, and 900 m for each well.The optimal well spacings determined based on the LCOH-HT model were smaller than those reported by previous model for all wells.Compared to the previous model, the LCOH-HT model accounts for the drilling cost, which increases with well spacing.This result implied that, when utilizing horizontal and inclined wells, it is imperative to account for drilling cost.Because in the case of using inclined wells, increasing the well spacing is achieved by extending the drilling length of the inclined wells, which inevitably leads to increased cost.In conclusion, when optimizing well spacing, it is crucial to fully consider the heat production and associated costs at different well spacings, particularly those related to well spacing.

Optimal flow rate and temperature
Using the LCOH-HT enables the optimization of not only the well spacing but also various other production parameters, such as the reinjection temperature and flow rate.In this section, the reinjection temperature and flow rate of GR2 were individually optimized via LCOH-HT analysis.The well spacing was the optimal well spacing obtained from the analysis, and the other parameters remained constant.Based on actual production needs, reinjection temperatures were set at 20, 25, 30, 35, and 40 ℃. Figure 9a indicates that as the reinjection temperature decreases, the total heat production increases steadily, while the costs remain relatively stable.And the LCOH-HT value continued to increase from 28.2 to 41.2.Notably, the lowest reinjection temperature of 20 ℃ yielded the highest economic benefits.Figure 9c demonstrates that despite a 20 °C decrease in the reinjection temperature, the resultant change in the production temperature was remarkably small.Over 30 years, the observed variation in the production temperature was a mere 0.5 ℃.This phenomenon was attributed to the favorable thermal background in the study area, facilitating rapid heat recovery via the reinjection of water.It also indicated that lower reinjection temperatures are unlikely to cause cooling of the reservoir.
Figure 9b, d displays the changes in the LCOH-HT and temperature when the flow rate was set at 90, 100, 110, 120, and 130 m 3 /h.The total heat production increased with increasing flow rate, accompanied by a rise in costs.This is attributed to the increased pumping costs due to the larger pressure drop resulting from the increased flow rate.However, the increase in costs is relatively minor, leading to a continuous decrease in the LCOH-HT (see Fig. 9b).However, increasing the flow rate leads to significant drops in temperature and pressure.Specifically, as the flow rate increases from 90 to 130 m 3 /h, the temperature decrease escalates from 0.6 to 1.7 ℃ after 30 years (refer to Fig. 9b), and the pressure variation became more pronounced.Although total heat production benefited from the increased flow rate, a high flow rate is not advisable due to the potential adverse effects on the subsurface environment and the reduced lifespan of the geothermal field, caused by the induced temperature and pressure drops.These results indicated the versatility of the LCOH-HT method for broader parameter optimization.

Fig. 1
Fig. 1 Procedure to optimize the well spacings

Fig. 2 a
Fig. 2 a Location of Bohai Bay Basin in China; b location of study area in the Bohai Bay Basin; c location and geological map of study area; d location of the wells

Fig. 9
Fig. 9 LCOH-HT a and b at different reinjection temperatures (a), flow rate (b), associated temperature changes (c) and flow rate changes (d)

Table 1
Information of study area the model boundaries.The domain was discretized into 420,000 tetrahedral elements.
The water table data for two wells (M1 and M2) for the 2023 off-heating season, i.e., March 23 -October 30, were used for model calibration and validation.The simulated pressure at the bottom of the well was converted to the water table level and compared with the measured water table elevations.The model was run numerous times during calibration to provide the best match between the simulation and measured static pressure and temperature data.A comparison of the hydrological results of the model and the measured water table elevations is shown in Fig.4.The red scattered points represent the simulated results, and the blue scattered points represent the observed water

Table 2
Parameters of Caofeidian geothermal reservoir Fig. 4 Comparison of water table result of model and measurement water table for two wells: M1 (a) and M2 (b)

table level .
A comparison of the simulated and measured water table data sets consistently revealed relative water table elevation errors of less than 3%.Our analysis indicated that this difference is acceptable and that the developed model is reasonable and can be applied in further work.