**Numerical Investigation of Techno-Economic Multiobjective Optimization of Geothermal Water Reservoir Development: A Case Study of China**

#### **Luyi Zhang 1, Ruifei Wang 2, Hongqing Song 1,\*, Hui Xie 1,\*, Huifang Fan 1, Pengguang Sun <sup>3</sup> and Li Du <sup>3</sup>**


Received: 20 September 2019; Accepted: 4 November 2019; Published: 6 November 2019

**Abstract:** Renewable geothermal utilization is a significant approach for residential heating with two principal modes, direct geothermal district heating systems (DGDHSs) and indirect geothermal district heating systems (IGDHSs). The key principle of geothermal development design is to prevent premature thermal breakthrough, which could result in low efficiency of geothermal heating systems. In this paper, a new approach considering building heating demand, geothermal water resource protection, and optimal economic benefits is presented systematically. The results simulated by OGS software show that well spacing, reinjection temperature, and production rate are the most significant parameters affecting thermal breakthrough in geothermal reservoirs. In addition, production rate and reinjection temperature have a huge effect on the payback period of investment. Comparing IGDHS to DGDHS, the investment in construction of geothermal wells and the annual water consumption decrease by up to 10% and 50%, respectively. Additionally, electricity costs increase by 5% to 30%. The indirect geothermal district heating system with a well spacing of 300 m, a production rate of 100 m3/h, and a reinjection temperature of 301.15 K is much better for this case, both technically and economically. The systematic calculation approach can be reasonably applied to other regions with geothermal energy utilization.

**Keywords:** geothermal water reservoir; well spacing; direct geothermal district heating system; indirect geothermal district heating system; multiobjective optimization; technical and economic evaluation

#### **1. Introduction**

With the development of the economy and improvement of people's living standards, primary energy consumption has maintained a rapid growth. In northern China, energy consumption of building heating and space cooling account for 20% and 14%, respectively [1]. One measure for energy conservation and emission reduction is to improve building insulation levels and the air tightness of building envelopes. Another measure is to use geothermal resources. Geothermal energy is a stable resource from the earth's interior [2,3]. Theoretically, it can be utilized continuously [4]. As long as the upper limit of geothermal utilization is properly controlled, the service life of a geothermal reservoir can be extended [5,6]. Besides, geothermal energy is generally unaffected by geography, climate, and season, unlike solar, water, wind, biomass, and ocean energy [7,8].

Geothermal heating technology is an effective way to fully develop and utilize geothermal energy [9,10]. It has been applied in some countries because of its high operational efficiency and

significant environmental benefits [11]. However, the high initial investment and long payback period have hindered the popularization of this technology [12,13]. Geothermal well infrastructure investment can exceed 50% of system construction costs for all types of geothermal utilization systems [14,15]. Statistically, the expense of drilling a geothermal well is about \$300 per meter [16]. The cost increases with deepening of the wells [17]. In addition, all extracted groundwater needs to be reinjected to maintain the pressure of the geothermal reservoir. The reinjection rate is far lower than the production rate, resulting in a waste of ground water [18]. Cold reinjection water also leads premature thermal breakthrough, which affects the system's service behavior and operating efficiency.

Researchers have studied the energy production of geothermal temperature fields. Zhang et al. developed a numerical model to evaluate geothermal reservoirs and optimize the production performance of geothermal systems [19]. Held et al. derived a reservoir mathematical model that considered fluid–rock interaction in the utilization of geothermal resources [20]. Nam carried out 3D numerical simulations and verified the reliability by comparing results with experimental results [21]. Salimzadeh et al. set up a fully coupled thermal–hydraulic–mechanical (THM) finite element model [22]. Pandey et al. found that thermochemical and thermomechanical effects were the most important consequences of heterogeneity [23]. Nikhil et al. simplified the crack model by a statistical linear congruent generator method, and carried out sensitivity analysis of fracture aperture size and fluid flow [24].

The previous efforts offered many insights into modeling heat extraction from geothermal reservoirs. Rock deformation should be emphasized, which may lead to unwanted sand components in the water flow [25]. However, previous studies ignored the economic benefits of geothermal projects. The disadvantages of geothermal projects without economic benefit analysis were huge. It not only increased the cost, but also caused wasted geothermal resources. High investment and operating cost hindered the development of geothermal heating technology.

In this study, based on a coupled mathematical and multiobjective optimization model, numerical investigation is carried out to evaluate the heat extraction and economic efficiency of a geothermal project in Xinji County. First, coupling of the thermal–hydraulic–mechanical model and multiobjective function of the geothermal heating system is established. Second, the main factors causing thermal breakthrough are analyzed based on data of a geothermal reservoir in Xinji City. Third, production parameters of direct and indirect geothermal district heating systems are optimized based on the ground heating load and system operation strategies. Finally, the initial investment and operation cost of direct and indirect geothermal district heating systems are analyzed, and the optimal method and production parameters are thus obtained. This optimization and economic evaluation method can be applied to geothermal heating systems in other regions, which is beneficial to the further development of geothermal heating technology.

#### **2. Target Project**

#### *2.1. Background*

The area studied in this paper is located in Xinji City, Hebei, China. The heat-providing area is 250,100 m2, which includes residential and commercial areas. The designed heating load is 9040 kW. After this project is put into operation, the annual heating consumption will be 65,600 GJ.

The sedimentary Miocene Guantao Formation is the most favorable for geothermal heating systems [26,27]. According to a geological investigation conducted by China Petroleum and Chemical Corporation, the geothermal resources of Xinji is regarded as consisting of homogeneous porous media with high permeability. The reservoir thickness is 512 m and the average porosity of this geothermal reservoir is 30%. The depth of the bottom boundary varies from 1600 to 2100 m. The single-well water rate is between 78 and 120 m3/h, and the production temperature is in the range of 333.15–3337.15 K. Laboratory tests show that the reservoir rock permeability is around 160 mD. The average reservoir

temperature at target layers around 1300 m deep is about 333.15 K. The ground water of this geothermal reservoir can be directly used for building heating due to the low relative temperature [26].

To ensure the stable operation of a geothermal system, the geological conditions and well test data need to be considered in the design of the heating system [27]. Reinjection operation data of Xinji oil production well no. 5 are shown in Figure 1. Except during the filter element replacement period, the reinjection rate is mostly maintained at 120 m3/h, which ensures production water being totally reinjected underground.

**Figure 1.** Reinjection rate of Xinji oil production well no. 5.

#### *2.2. Geothermal Heating Strategy*

There are two kinds of geothermal heating modes designed for Xinji. One is a direct geothermal district heating system (DGDHS). As shown in Figure 2a, ground water extracted from the geothermal production well is transported to the heat exchanger. The heating circulation water absorbs heat from the geothermal water. Then, the heat is piped to residential areas. The other one is an indirect geothermal district heating system (IGDHS), shown in Figure 2b. The heat of geothermal water is partially transferred to the heating circulation water at the primary heat exchanger with higher returning water temperature, and the heat is piped to some residential areas. After the initial heat transfer, the geothermal water is transported to the sequential secondary heat exchanger, where the heat is transferred from the geothermal water to the intermediate circulation water. Then, the heat pump extracts the heat from the intermediate circulation water and discharges it into the heating circulation water. Finally, the heat is transported to the remaining residential areas.

**Figure 2.** Operation strategies of geothermal heating systems: (**a**) direct geothermal district heating system (DGDHS) and (**b**) indirect geothermal district heating system (IGDHS).

#### **3. Mathematical Model and Optimization Methodology**

#### *3.1. Thermal–Hydraulic–Mechanical Model*

In this paper, the problems of water flow, heat transport, and rock deformation in the Xinji geothermal reservoir are studied jointly. The process of groundwater extraction is considered to be a flowing process of liquid in porous media consisting of a liquid phase and a solid phase [28]. The process of modeling and theoretical analysis consists of the following assumptions [28–33]:


The mass balance equation, energy balance equation, and rock deformation in porous media are shown in Equations (1)–(3) [31,32]:

$$\frac{\partial(\phi\rho\_f)}{\partial t} + \nabla \cdot (\phi\rho\_f u\_{fs}) = Q\_f \tag{1}$$

$$
\rho c\_p \frac{\partial (T)}{\partial t} + \nabla \cdot (-\lambda \nabla T) + \rho\_f c\_f u\_{fs} \nabla T = 0 \tag{2}
$$

$$
\nabla \cdot \sigma = 0 \tag{3}
$$

where φ is rock porosity; ρ*<sup>f</sup>* is water density (kg/m3); *t* is time (s); *Qf* is the mass flow of sink/source item (kg/(m3·s)); <sup>ρ</sup>*cp* is the heat capacity of porous media (J/(m3·K)); *<sup>T</sup>* is temperature (K); <sup>λ</sup> is the thermal conductivity of porous media (W/(m·K)); σ is the typical Cauchy stress tensor; and *uf s* is fluid flow velocity between rock and water flow (m/s). Both the heat capacity and thermal conductivity of porous media can be estimated as an arithmetic mean of each phase property weighted by its volume fraction.

Then, the flow velocity *u* is divided into the rock velocity *us* and the relative velocity of fluid *uf s*. Huo et al. and Salimzadeh et al. proved that *us* is associated with the expansion coefficient of rock and fluid ∂ε <sup>∂</sup>*<sup>t</sup>* , and ε is calculated from the displacement solved by Equation (3) [32,33]. The term *us* is used to couple the water flow with rock deformation. Equation (1) is rewritten into Equation (4):

$$
\omega\_t \phi \rho \frac{\partial p}{\partial t} + \nabla \cdot (\phi \rho\_f (\mathbf{u}\_{fs} + \mathbf{u}\_s)) = Q\_f \tag{4}
$$

where *p* is fluid pressure (Pa); *ct* is total compressibility (Pa<sup>−</sup>1), *ct* = *c <sup>f</sup>* + *cs*; *c <sup>f</sup>* is water compressibility (Pa<sup>−</sup>1); and *cs* is rock compressibility (Pa<sup>−</sup>1).

In Equation (4), the effect of rock deformation caused by geothermal water extraction is considered. Thus, the flow velocity *uf s* + *us* depicts the movement of water in porous media. This is a very important step. Although *us* is usually a very small value, it cannot be neglected, as it accounts for the rock deformation mechanism.

As shown in Equation (5), *ufs* further expressed by Darcy's law that shows the heat–fluid coupling process:

$$
\Delta u\_{fs} = -\frac{k}{\phi \mu} (\nabla p + \rho\_f g \nabla z) \tag{5}
$$

where *g* is gravitational acceleration (m/s2); *z* is the depth of the geothermal reservoir (m); and *k* is the intrinsic permeability tensor (m2).

In order to specify the solution for the above equations, one needs to prescribe appropriate boundary conditions. The details are shown in [29,30].

(1) Prescribed pressure (Dirichlet condition)

$$p = \overline{p} \text{ on } \Gamma\_p \tag{6}$$

for imposing hydrostatic conditions at domain boundaries or pressure at wellbores.

(2) Prescribed fluid flux (Neumann condition)

$$q\_n = q \cdot n \text{ on } \Gamma\_q \tag{7}$$

with the normal vector *n* for imposing inflow or outflow rate at wellbores.

(3) Prescribed temperature (Dirichlet condition)

$$T = T \text{ on } \Gamma\_T \tag{8}$$

for imposing natural thermal gradients at model boundaries and the temperature of the injected fluid.

The finite element method is used in Open-GeoSys (OGS) to provide numerical solutions to the aforementioned coupled formulation. More information can obtained from [29–33].

#### *3.2. Calculation of Well Number based on Heating Load*

Based on the heating load of residential buildings, the number of geothermal production wells is calculated, which is used to calculate construction costs. The building heat load is estimated based on local weather conditions, the function of the building, and characteristics of the protected exterior construction [34]. The calculation formula is shown in Equation (9):

$$Q\_l = \mathcal{W}\_b \cdot \mathcal{S}\_b + \mathcal{W}\_c \cdot \mathcal{S}\_c \tag{9}$$

where *Wc* is heat load per unit area of residential building (W/m2); *Wb* is heat load per unit area of commercial building (W/m2); *Sc* is residential heat-providing area (m2); *Sb* is commercial heat-providing area (m2); and *Ql* is total heat load (W).

The number of geothermal production wells is determined by the production rate and heating mode, shown in Equations (10) and (11):

$$m\_1 = \frac{(1+\alpha)Q\_{\rm l}}{c\_{\rm lv} \cdot q\_p \cdot \rho\_f \cdot (t\_\mathcal{g} - t\_\mathcal{h})} \tag{10}$$

where *n*<sup>1</sup> is the number of production wells for DGDHS (rounded up); α is the heat loss coefficient of the pipe network (5%–10%); *cw* is the specific heat of water (J/(kg·K)); *tg* is production temperature (K); *th* is reinjection temperature (K); and *qp* is production rate (m3/h); and

$$m\_2 = \frac{(1+a)Q\_l}{c\_w \cdot q\_p \cdot \rho\_f \cdot \left[ (t\_\S - t\_{h1}) + \frac{\text{QCD}}{\text{CPP}-1} (t\_{h1} - t\_{h2}) \right]} \tag{11}$$

where *n*<sup>2</sup> is the number of production wells for IGDHS (rounded up); *th*<sup>1</sup> is water temperature for the first return (K); *th*<sup>2</sup> is water temperature for the second return (K); and *COP* is the coefficient of performance.

#### *3.3. Multiobjective Optimization Model for Heating Systems*

In this part, a multiobjective optimization model for geothermal heating systems is established. The project's profitability is analyzed by coupling the initial investment and operating cost of the system. The initial investment, annual operating cost, and revenue are shown in Figure 3 [35–37].

**Figure 3.** Initial investment, annual operating cost, and revenue of geothermal heating systems.

Geothermal well construction costs account for the main part of all construction costs which are based on well number, type, and construction depth. Due to limited drilling area, a vertical well is adopted for geothermal well construction, and the rest are directional wells, shown in Figure 4. The construction depth of a directional well is related to well spacing and directional angle, as in Equation (12), and the construction cost of the project can be estimated according to Equation (13):

$$h\_2 = L + (h\_1 - L) / \cos \theta \tag{12}$$

where *h*<sup>1</sup> is the depth of the vertical well (m); *h*<sup>2</sup> is the depth of the directional well (m); *L* is the depth of the first section (m); and θ is directional angle (◦); and

$$\mathcal{C}\_{z} = (h\_1 \cdot P\_1 + \mathcal{C}\_f) + (h\_2 \cdot P\_2 + \mathcal{C}\_f) \cdot (n - 1) \tag{13}$$

where *Cz* is construction cost of a geothermal well (k USD); *P*<sup>1</sup> is construction cost per meter of vertical well (k USD); *P*<sup>2</sup> is construction cost per meter of directional well (k USD); *Cf* is cost of attached ground facility (k USD); and *n* is total number of geothermal heating wells.

The payback period is the time required for the return on investment to "repay" the total investment of the geothermal heating system [38], as expressed in Equation (14), and is calculated from the beginning year of construction [39,40]:

$$\sum\_{t=0}^{Pt} (CI - CO)\_t (1 + i\_t)^{-t} = 0\tag{14}$$

where *Pt* is static payback time (year); *CI* is cash inflow; *CO* is cash outflow; and (*CI* − *CO*)*<sup>t</sup>* is net cash flow for the year, t.

**Figure 4.** Schematic diagram of vertical well and directional well.

There are two main objectives of geothermal heating system optimization: one is to have the shortest payback period, and the other is to have minimum well spacing. In order to improve the profitability of the system, its production parameters must be optimized. However, three factors restrict the optimization of a geothermal heating system. The first is the area of equipment, which must be smaller than the project area. This factor is becoming more significant due to the limited space of urban areas. Second, a geothermal system is affected by groundwater resources, which is related to the hydrogeological conditions in the specific region. Finally, the system must meet the people's requirements for thermal comfort in the built environment. The multiobjective optimization functions are shown in Equation (15).

Optimization function:

$$\min \text{Pt}(\mathbf{x}\_1, \mathbf{x}\_2, \mathbf{x}\_3 \dots \mathbf{x}\_{\mathbf{n}}) \tag{15}$$
 
$$\min \text{d}(\mathbf{x}\_1, \mathbf{x}\_2, \mathbf{x}\_3 \dots \mathbf{x}\_{\mathbf{n}})$$

Constraint conditions: Scale constraint:

*Sg* < *S*

Geological constraint:

$$q\_{\mathbf{x}} \le q\_{\text{pmax}}(\mathbf{x}\_1, \mathbf{x}\_2, \mathbf{x}\_3 \dots \mathbf{x}\_n),$$

Load constraint:

$$Q\_x \ge Q\_l(\mathbf{x}\_1, \mathbf{x}\_2, \mathbf{x}\_3 \dots \mathbf{x}\_n)$$

where *l* is well spacing (m); *Sg* is the area of equipment (m2); *S* is project area (m2); *qx* is the design value of production rate (m3/h); *q*pmax is the maximum of production rate (m3/h); *Qx* is the design value of system load (W); and x1, x2, x3 ... xn represent different working conditions.

If hydrogeological conditions in a specific region allow the use of geothermal heating systems, these multiobjective optimization functions can be used for assessment.

#### *3.4. Simulation Procedure*

According to the established thermal–hydraulic–mechanical and multiobjective optimization models, using simulation software, the calculation process is shown in Figure 5. A geothermal heating system can be optimized through the entire calculation process. First, formation parameters are input into the model. The movement of the cold water thermal front is calculated. By comparing the temperature of the production well (TP0 is the initial temperature of the production well, and TP50 is the temperature at the 50th year), the optimal spacing of wells is selected. Finally, the initial investment and operating cost are calculated. Unprofitable geothermal heating systems are simply discarded. Though the multiobjective optimization model, the optimal production parameters can be obtained.

**Figure 5.** Flowchart of simulation procedure.

#### **4. Result and Discussion**

#### *4.1. Sensitivity Analysis*

This part studies the influence mechanism of various factors on the formation of thermal breakthrough. The basic conditions of the Xinji geothermal reservoir provided by China Petroleum and Chemical Corporation are listed in Table 1. These parameters, such as well spacing, production rate, and reinjection temperature, affect the formation of thermal breakthrough. They also determine the heat extraction efficiency of a geothermal heating system. In addition, periodic heating is considered in the modeling: a geothermal heating system works for 5 months each winter. Geothermal wells are closed during the remainder of the year. Results are obtained at the final time step at the end of the 50th year.


**Table 1.** Thermal–hydraulic–mechanical parameters for the base case of the Xinji geothermal reservoir.

#### 4.1.1. Effect of Well Spacing

Well spacing is the distance between a heat production well and a reinjection well, which determines the placement of geothermal wells. In this section, different well spacings (200, 400, and 600 m) are set into the mathematical model. The remaining parameters are as follows: production rate is 100 m3/h and reinjection temperature is 311.15 K. Temporal evolutions of the water temperature and pressure at the production well are plotted in Figure 6.

**Figure 6.** Time evolution of production well in well spacing studies: (**a**) temperature at production well; (**b**) pressure at production well.

In general, large well spacing can slow down the formation of thermal breakthrough. It can be observed from Figure 6a that when the well spacing is 200 m, the temperature of the production well begins to decline at the 25th year. The reason is that the movement of cold water thermal front reaches the production wells. This temperature change would reduce the heat extraction efficiency of a geothermal heating system. Eventually, ground heating demands are not met. For the other well spacings, the temperature of production wells is constant at 333.15 K. Note that there is no thermal breakthrough.

As shown in Figure 6b, the pressure of the production well decreases with increased well spacing. This is because the reinjected water is better able to repressurize the depleted area around the production well with smaller spacing. It means that large well spacing can effectively avoid premature thermal breakthrough. The movement of the cold water thermal front takes a long time to arrive at the production well. However, it loses the significance of reinjection to maintain water pressure. As commercial, land use, economic, and policy factors should be considered in the selection of geothermal well spacing, optimizing the well spacing is a challenging task.

#### 4.1.2. Effect of Reinjection Temperature

The results of first part show that there is an apparent decline in production temperature of a geothermal heating system with well spacing of 200 m and a production rate of 100 m3/h. These parameters are helpful to visualize the influence of reinjection temperature on thermal breakthrough. In this case, four values of the parameters (288.15, 298.15, 308.15, and 318.15 K) are considered. Temporal evolution of the water temperature and pressure at a geothermal production well is plotted in Figure 7.

**Figure 7.** Time evolution of production well in reinjection temperature studies.

The temperature of reinjection water has a great influence on the formation of thermal breakthrough and the temperature of the production well. As is observed in Figure 7, thermal breakthrough is formed first when the reinjection temperature is 288.15 K. The lower the reinjection temperature is, the earlier the thermal breakthrough is formed, and the more the production water temperature drops. This is because the temperature decrease caused by the reinjected water results in the movement of the cold water thermal front originating from the reinjection well.

Although lower reinjection temperature leads to premature thermal breakthrough, it also points to a high utilization ratio of geothermal resources. Therefore, reinjection temperature must be controlled reasonably in geothermal heating projects.

#### 4.1.3. Effect of Production Rate

The production rate is the water production rate at the production well. All geothermal water must be reinjected underground, so the reinjection rate is always equal to the production rate. According to the previous analysis of well spacing and reinjection temperature, when the well spacing is 200 m and the reinjection temperature is 308.15 K, the temperature of the production well changes markedly. Therefore, those parameters help to observe the sensitivity of temperature to production rates. Different production rates (80, 100, and 120 m3/h) are input into the model. The results are shown in Figure 8.

**Figure 8.** Time evolution of production well in production rate studies: (**a**) temperature at production well; (**b**) pressure at production well.

A large production rate promotes the movement of the cold water thermal front toward the production well. As observed Figure 8a, with decreased production rate, the temperature decrease at the production well is delayed, and the drop in production temperature is not obvious. This is because the decreased well rate slows down the spread of the temperature decrease initiated by the movement of the cold water thermal front.

According to Figure 8b, the greater the production rate, the lower the pressure. This is because a lower well rate results in a smaller pressure gradient between the reinjection well and the production well.

A small production rate delays the temperature drop of production water. However, it also results in poor heating extraction of the geothermal heating system. Therefore, a reasonable range of production rates should be controlled on the basis of local test data.

#### *4.2. Numerical Optimization*

#### 4.2.1. Parameter Optimization

Taking IGDHS with a production rate of 100 m3/h and reinjection temperature of 301.15 K as an example, results are shown in Figure 9. Figure 9a shows that production temperature does not change when well spacing is over 300 m. As can be seen from Figure 9b, with increased well spacing, the pressure of the production well keeps dropping. This is similar to the first part of the simulation results, proving that the results are correct.

**Figure 9.** Time evolution of a production well for indirect geothermal district heating systems (IGDHSs) with production rate of 80 m3/h: (**a**) temperature at production well; (**b**) pressure at production well.

Temperature distribution nephograms of geothermal reservoirs under different situations are shown in Figures 10–12. The results show that the radius of the cold water thermal front is directly proportional to the production rate, and the radius of IGDHS is larger than that of DGDHS. In addition, the radius of the cold water thermal front increases with decreased reinjection temperature. As shown in Figure 11b, the radius of the cold water thermal front is about 260 m. Considering that the uneven distribution of formation and the instability of groundwater accelerate the flow of low-temperature groundwater, it is recommended that well spacing should be over 300 m in this case. The production parameters of each geothermal heating system are shown in Table 2.

**Figure 10.** Temperature distribution nephograms of direct geothermal district heating systems (DGDHSs) with different production rates in the 50th year: (**a**) radius of cold water thermal front; (**b**) production rate of 80 m3/h; (**c**) production rate of 100 m3/h; (**d**) production rate of 120 m3/h.

**Figure 11.** Temperature distribution nephograms of an IGDHS with a reinjection temperature of 301.15 K in the 50th year: (**a**) radius of cold water thermal front; (**b**) production rate of 80 m3/h; (**c**) production rate of 100 m3/h; (**d**) production rate of 120 m3/h.

**Figure 12.** Temperature distribution nephograms of a IGDHS with reinjection temperature of 288.15 K in the 50th year: (**a**) radius of cold water thermal front; (**b**) production rate of 80 m3/h; (**c**) production rate of 100 m3/h; (**d**) production rate of 120 m3/h.

**Table 2.** Optimal well spacing for different geothermal heating systems.


Projects using groundwater for heating can be optimized with the above-mentioned method. The input parameters are adjusted by the geological characteristics of each area and heating modes. The optimal well spacing of the geothermal heating system can be obtained. This methodology provides a good reference for cities developing geothermal heating technology.

#### 4.2.2. Economic Optimization of Geothermal System

In this paper, a static economic evaluation method is used to calculate the initial investment and operation cost of DGDHS and IGDHS. The basic conditions of geothermal systems provided by China Petroleum and Chemical Corporation are listed in Table 3 (1 USD = 7 RMB). Through the multiobjective optimization model, a geothermal system with high investment return and short investment payback period is selected for the case in Xinji.


**Table 3.** Economic calculation parameters for the base case.

Figure 13a shows the relationship between well price and production rate. The well price increases with increased production rate, and the same result can be achieved by lowering the reinjection temperature. The reason is that large well spacing is required for a geothermal heating system with large production rates and low reinjection temperatures to prevent thermal breakthrough.

**Figure 13.** Well price and initial construction investment for geothermal wells of different heating modes: (**a**) construction price of directional well; (**b**) construction investment of heating system.

It can be clearly seen from Figure 13b that the initial construction cost of DGDHS is the highest. Although the construction cost of DGDHS decreases with increased production rate, it is still higher than that of IGDHS. For IGDHS, increasing the production rate or reducing the reinjection temperature can reduce the initial construction cost, and the construction investment of geothermal wells will decrease by up to 30%.

There is a great difference in electricity charges and initial investments of various geothermal heating systems. As shown in Figure 14a, the electric energy consumed by IGDHS is nearly 1.5 times that of DGDHS. The power consumption of IGDHS presents an irregular change with increased well production rate. The reason is that the decreased reinjection temperature leads to increased secondary heat transfer and heat load of the heat pump, which causes increased power consumption. However, the trend in water consumption is reversed, shown in Figure 14b. Water consumption of IGDHS is much less than that of DGDHS, reduced by 10%–50%. This is because IGDHS generates heat by electricity and has better environmental adaptability than DGDHS. Even in extreme weather, IGDHS meets residential heating needs.

**Figure 14.** Annual power consumption of different geothermal heating systems: (**a**) electricity cost; (**b**) water cost.

As shown in Figure 15, DGDHS has a longer payback period than IGDHS. Actually, annual power consumption of DGDHS is lower, because there are no heat pumps, but the high cost of construction investment and equipment depreciation causes a longer payback period. On the contrary, IGDHS consumes large amounts of electricity during operation, but the initial investment for a geothermal heating system is lower, which greatly reduces the payback period. As reinjection temperature is lower, the payback period becomes longer. That is because the annual power consumption of a geothermal heating system with low reinjection temperature is significantly higher than other systems. An indirect geothermal district heating system for building heating in the case of Xinji, China, is the best heating mode. The optimal production rate, reinjection temperature, and well spacing for a 50-year heating period are 100 m3/h, 301.15 K, and 300 m, respectively.

**Figure 15.** Investment payback period of different geothermal heating systems.

In fact, the income and operating expenses of geothermal heating systems in Xinji fit nicely with economic calculations. It is shown that the calculation results by the static technical economic evaluation are correct. Due to different hydrological conditions and economic levels of cities, all items

in this paper cannot be directly used. However, the mathematical and multiobjective optimization models proposed in this paper can be used for assessment, and calculation results can be referenced for the choice of geothermal heating systems.

#### **5. Conclusions**

Renewable geothermal utilization is a significant approach for residential heating in geothermal enrichment regions. The development of geothermal reservoirs has to meet building heating demands with optimal economic benefits, and the geothermal water should be protected. In this paper, a new approach with practical engineering requirements and environmental and economic considerations is established systematically.


**Author Contributions:** Conceptualization, L.Z. and H.S.; Data curation, L.Z.; Formal analysis, L.Z., R.W., and H.S.; Funding acquisition, H.S.; Investigation, H.X., P.S., and L.D.; Methodology, R.W.; Project administration, H.S.; Resources, P.S. and L.D.; Software, L.Z.; Supervision, H.S. and H.X.; Validation, H.F.; Visualization, R.W. and H.X.; Writing—original draft, L.Z.; Writing—review & editing, H.X. and H.F.

**Funding:** This research was funded by the Beijing Nova Program, grant number Z171100001117081.

**Acknowledgments:** The work was supported by the Beijing Nova Program (grant no. Z171100001117081). The authors also wish to thank Sinopec Star Petroleum Co., Ltd. for the cooperation in this work.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).
