Renewable Energy 177 (2021) 976e991

Contents lists available at ScienceDirect

Renewable Energy
journal homepage: www.elsevier.com/locate/renene

Closed-loop geothermal energy recovery from deep high enthalpy
systems
Wanju Yuan *, Zhuoheng Chen , Stephen E. Grasby , Edward Little
Geological Survey of Canada, 3303 33rd Street, NW, Calgary, T2L 2A7, Alberta, Canada

a r t i c l e i n f o

a b s t r a c t

Article history:
Received 28 October 2020
Received in revised form
4 May 2021
Accepted 4 June 2021
Available online 7 June 2021

Closed-loop geothermal energy recovery technology has advantages of being independent of reservoir
ﬂuid and permeability, experiencing less parasitic load from pumps, and being technologically ready andwidely used for heat exchange in shallow geothermal systems. Commercial application of closed-loop
geothermal technology to deep high-enthalpy systems is now feasible given advances in drilling technology. However, the technology it uses has been questioned due to differences in heat transport capacities of convective ﬂow within the wellbores and conductive ﬂux in the surrounding rock. Here wedemonstrate that closed-loop geothermal systems can provide reasonable temperature and heat duty forover 30 years using multiple laterals when installed in a suitable geological setting. Through use of two
analytical methods, our results indicate that the closed-loop geothermal system is sensitive to reservoir
thermal conductivity that controls the level of outlet temperature and interference between wells overtime. The residence time of the ﬂuid in the horizontal section, calculated as a ratio of the lateral length toﬂow rate, dictates heat transport efﬁciency. A long vertical production section could cause large drops in
ﬂuid temperature in a single lateral production system, but such heat loss can be reduced signiﬁcantly ina closed-loop system with multiple laterals.
Crown Copyright © 2021 Published by Elsevier Ltd. All rights reserved.

Keywords:
Closed-loop geothermal technology
Thermodynamics
Technology evaluation
Analytical modeling
Sensitive analysis

1. Introduction
Geothermal energy is one of a suite of renewable energies that
can help to reduce global CO2 emissions [1] and transition modernday society to a lower-carbon economy [2]. While having the
beneﬁt of being able to provide baseload power supply, geothermal
development comes with the challenge of ﬁnding thermal reservoirs that have porosity and permeability suitable for economic
production [3]. Such traditional geothermal development typicallyinvolves drilling a production well, for extraction of deep thermal
waters, along with a re-injection well to dispose of produced ﬂuids
after heat extraction. In the subsurface, heat is carried by ﬂuids
transported through porous and/or fractured media towards the
production well [4], and both heat conduction and convection
affect the overall heat production efﬁciency (Fig. 1). In contrast,
shallow geothermal-exchange systems that have been widely usedaround the world for space heating [5] make use of a closed-loopsystem that circulates a working ﬂuid to extract and store heat

* Corresponding author.
E-mail address: wanju.yuan@canada.ca (W. Yuan).
https://doi.org/10.1016/j.renene.2021.06.028
0960-1481/Crown Copyright © 2021 Published by Elsevier Ltd. All rights reserved.

underground (depending on the season). Advances in drilling
technology have led to consideration of such closed-loop systems
for development of deeper (higher temperature) geothermal resources (e.g., Eavor Technologies, 2021; [6]. Such closed-loop systems are similar to geothermal-exchange, but at a much larger scale
to utilize deeper high-enthalpy resources [5]. In concept, closedloop technology applies thousands of meters long horizontal
wellbores in the target thermal reservoir that are connected by two
vertical wells, one inlet well, and one production well (Fig. 2). A low
temperature working ﬂuid is injected from the inlet, heated mostly
through the horizontal interval in a geothermal reservoir withproper thermal properties, and the heat energy is harvested at theoutlet for heating and power generation purposes. In contrast to a
traditional geothermal energy system, a closed-loop system uses aworking ﬂuid that circulates with a designed ﬂow rate. In this sense
the system is independent of permeability, greatly reducing the
exploration risk, and there are no associated environmental or
scaling issues with production of geothermal ﬂuids to surface(Fig. 2). A key limitation though is that heat in a closed-loop system
can only be transported from the reservoir to the wellbore by heat
conduction [7], while within the wellbore, convection dominates.
Geothermal energy resources are widely distributed in Canada,W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991lrt

Nomenclature
A
a
C
c
D
erf
exp
h
ln
L
N
NL
Nl
q
qc
qcu
qf
Qo
Qr
ri
s
t
T
DT
Vi
v
x
y
z

cross section area of lateral wellbore, m2
factors in Matrix need to multiply heat ﬂux terms
right hand side value of the matrix calculation
speciﬁc heat capacity, J/kg/K
Inner diameter of lateral wellbore, merror function
exponential function
depth of the reservoir in case study of heat exchange
in vertical wellbore, m
the natural logarithm
length of the lateral wellbore, mnumber of segments of lateral wellbore for analyticalmodeling
numerical calculation number in Stehfest inverse
algorism
number of laterals
instantaneous conductive heat inﬂow, W
conductive heat ﬂux from reservoir to wellbore
under unit length, W/m
conductive heat ﬂux under unit length andtemperature difference, W/m/K
ﬂow rate per lateral, m3/hr
cumulative heat amount produced out by working
ﬂuids, J
cumulative heat loss from the reservoir, J
the radius of investigation in geothermal industry, mLaplace transform parameter
time, stemperature, K
temperature change caused by source/sink, K
functions in Stehfest inverse algorism
velocity of the working ﬂuid in lateral wellbore, m/s
x coordinate location, my coordinate location, mz coordinate location, mthermal conductivity of reservoir, W/m/K
density, kg/m3
instantaneous source time, s

Subscripts
aveaverage valueeboundary coordinates of the reservoir
f
ﬂuid
inlet
the variable at the inletlproperty of lateral wellbore partL
one endpoint of this segment
mid-point mid-point of each segment
outletthe variable at the outletrreservoirR
one endpoint of this segmentssurfaceseg
property of the segmentsxsource/sink function in x directionsysource/sink function in y directionszsource/sink function in z directionv
property of vertical wellbore part#1functions from sub-system #1#2
functions from sub-system #2
11,…1N,…N1,…NN location of factors in matrix a
1,2,…,Nlocation of values in matrix c
Superscriptsithe variable, or function of the ith segment
1,2,…,Nthe variable of the corresponding segment*
residence time of working ﬂuids in horizontal
wellbore‘heat source/sink location
Acronyms
NW
PTA
RTA
WCSB

northwest
Pressure Transient Analysis
Rate Transient Analysis
Western Canadian Sedimentary Basin

Greek symbols
thermal diffusivity of reservoir, m2/saRecent research work on closed-loop geothermal heat production focuses on shallow geothermal systems with four major stylesof heat exchangers (single U-tube, double U-tube, coaxial centred,
and coaxial annular) [18]. Evaluation on long-lateral U shapeclosed-loop system usually needs to apply numerical simulations
on the discretized reservoir with local ﬁne grided meshing near
wellbore [19,39]. Few works used analytical modeling of longlateral U shape closed-loop system in deep high temperature
reservoir. To evaluate the thermodynamic efﬁciency of a closedloop multilateral wellbore system, at reservoir temperatures anddepths, we developed two analytical methods to model energy
production. These two analytical methods, coupled with a transientheat transfer process, simplify nonisothermal ﬂow in a closed-loopsystem. The meshed-free modeling provides an accurate heatproduction output with less computation time. Key parameters of
the system and its sensitivity have been examined to help better
understand key controls on this geothermal energy recovery
technology. The overall thermal production and heat duty are also
discussed.
We note that a deep closed-loop system would require advanced
horizontal drilling techniques and bear additional costs. Drillingwith the highest potential areas occurring in the western andnorthwestern parts of the country [8e14]. Recently many
geothermal energy projects have been initiated in Canada (e.g.
[15e17], to demonstrate either the resource development technology or geothermal resource concepts. Notable, in 2019, Eavor
drilled the Eavor-Lite™ closed-loop project in Alberta, Canada as afull-scale prototype of the drilling technology required (but not atfull development depth), and to validate the thermosiphon effects[16].

Fig. 1. Schematic diagram of the heat transfer mechanism in conventional geothermalenergy recovery with a pair of vertical wells. Colours represent working ﬂuid temperatures ranging from cold (blue) to hot (red).
977

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

Fig. 2. Schematic diagram of a closed-loop geothermal energy recovery system and its heat transfer mechanism. Colours represent working ﬂuid temperatures ranging from cold
(blue) to hot (red). The battery symbol reﬂects end heat usage (direct heat use and/or power generation).

2. Physical and mathematical modeltechniques evolve quickly as evidenced by oil and gas drilling forshale gas production as well as drilling in the enhanced geothermal
system [20] to overcome various difﬁculties, such as long-lateral
directional drilling, multilateral completions, high-temperatureeffects, thermal expansion of casing strings, hardness [21], andlarge diameters [22]. Several recent geothermal projects in Canada
and elsewhere have successfully completed, or planned their wellsin large depth (>3000 m) [15,23,24]. And the recently completed
long lateral demonstration closed-loop system drilled in Alberta by
Eavor [16] demonstrates the state of the art of long-lateral directional drilling in geothermal production. Other new advancement
includes the sealed lateral wellbore with the absence of casing toenhance the heat transfer efﬁciency [25]. However, the technicaldifﬁculties of drilling and costs are not the intension of this study.
This study is focused on analytical methods that would adequately
characterize the closed-loop system and allow evaluation of thecritical factors that affect the heat transmission in such systems.
However, there is no publicly information available to validate
the model and evaluate the heat exchange capacity in such a system. This study examines the major factors affecting the thermaltransporting capacity for the closed-loop system and its applicability to the sedimentary basin. In Alberta, the huge amount of datafrom oil and gas exploration and production wells allows for estimation and mapping of the thermal-physical properties of thesedimentary layers. The oil and gas industry is looking for practicaltechnologies for making a transition from fossil fuel to renewableenergy. If applicable in some suitable geological settings, thisclosed-loop system would be used to help the oil and gas industry
to reduce CO2 emission by geothermal-assisted fossil fuel production, or converting oil and gas well to a geothermal well when
depleted.

Assumptions. The closed-loop geothermal recovery technology
applies a sealed horizontal wellbore along the target thermal
reservoir (Fig. 2). The working ﬂuid (e.g., water, or other ﬂuids withthermosiphon effects) circulates within the sealed wellbore fromthe vertical inlet through the horizontal leg, to the production well.
The working ﬂuid is heated by the thermal reservoir as it ﬂows to
the outlet that connects to a power generator or a heat user. The
major assumptions to model this system are as follows:
1) Heat energy resource comes from geothermal heat ﬂux from
Earth and is stored in a geothermal reservoir that has speciﬁcthermal-physical properties allowing extraction of the energy
resource through heat exchange by water circulation in tubing
efﬁciently.
2) The dimension of the geothermal reservoir is far greater thanthe extent of a horizontal wellbore in the closed-loop system sothat reservoir boundary will not affect the production
performance.
3) The initial heat reservoir is assumed to be in thermal equilibrium. In this model, the horizontal well acts as the only sink to
extract heat energy out of the system. The heat source from
Earth will not affect the heat extraction performance due to low
thermal diffusivity and large reservoir dimension.
4) Pure heat conduction in the thermal reservoir. No natural or
induced heat convection by ﬂuid ﬂow in the surrounding rockmass occurs. A heat convection mechanism dominates the heattransfer process along the ﬂuid ﬂow direction within the horizontal wellbore system. The heat transfer within cross-sectionof the working ﬂuid is neglected based on following assumption.
978

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

complexity on transient heat ﬂux calculation. In our study, twoanalytical methods were applied to tackle this dynamic heattransfer phenomenon. In general, the continuity equation needs to
be Laplace transformed as:5) The diameter of the horizontal wellbore is extremely small
compared with the dimension of the thermal reservoir. The heattransfer within the wellbore is simpliﬁed as a 1D heat convection process, and the thermal resistance by convective ﬂow isneglected in this study. The temperature proﬁle within eachcross-section of the wellbore is assumed uniform but transient.
6) The overall thermal resistance between ﬂuid and rock isneglected including the potential resistance from the steel tube
wall, casing, and cement ring. In the closed-loop system, it iscritical to enlarge the contact area between the reservoir andworking ﬂuids, and reduce the thermal resistance to improve
the thermodynamic performance. Innovations of the sealed
horizontal portion without casing have also been commercialized in long lateral closed-loop geothermal energy production[16,25]. As a result, this study will focus on the generalgeothermal production evaluation without consideration of theoverall thermal resistance within the wellbore.
7) The working ﬂuid ﬂow rate is constant throughout the thermal
production period.
8) Thermal properties of the working ﬂuid and the reservoir, suchas thermal conductivity and volumetric heat capacity are
constant.

sT ¼  v

vTðx; tÞvTðx; tÞ¼  vArf cf
þ qc ðx; tÞ
vt
vx

ðt
qc ¼
0

(1)

DT ¼ 1 ¼

ðt
0

(5)



qcu ðtÞ
1ðx  x0 Þ2 þ ðy  y0 Þ2exp dt
2pl aðt  tÞ
4aðt  tÞ(6)

and in Laplace domain, qcu can be calculated based on convolution
theory:qcu ¼sL½xðtÞ(7)ðx  x0 Þ2 þ ðy  y0 Þ2
exp 
2plat
4at

(8)

xðtÞ ¼

1“L[]” means the value of original variables after Laplace transformation. In this work, Chebyshev's Theorem [29] was applied to
calculate the numerical Laplace transformation value. And a is thethermal diffusivity of reservoir which can be deﬁned as followed:

a¼

l
rr cr

(9)

When Laplace transformation applied to Eq. (5), the conductive
heat ﬂux from reservoir to horizontal wellbore can be derived as afunction of the local temperature and the conductive heat ﬂux
under unit temperature difference:(2)

The reservoir has a closed thermal boundary in each x, y, z direction. There is no heat communication at the outer boundary.

vT¼ 0; n ¼ x; y; z
vn

vTðtÞ
$qcu ðt  tÞdtvt

In Eq. (5), qcu is the conductive heat ﬂux under unit temperature
difference calculated from the standard source and sink functions[28]:

The heat balance for each location in the horizontal wellbore isbased on the law of energy conservation. The left-hand side of Eq.
(1) is the heat amount accumulated at a certain location x at time t.
The ﬁrst term in the right-hand side is the amount of heat transferred by convection in the horizontal well. Because of the highvelocity of the convective ﬂow, heat brought by thermal conduction
is negligible. This contrasts with conventional geothermal energy
systems, where the working ﬂuid is ﬂowing at low velocity throughporous media, and both heat conduction and convection contributesigniﬁcantly to the overall heat transfer performance (Fig. 1). The
second term of the right-hand side of the equation is the heattransferred by conduction from rock forming the thermal reservoir
surrounding the wellbore. This transient heat ﬂux, qc(x,t), is
controlled by both the temperature difference between the outside
and inside of the wellbore, and the thermal conductivity of the rock.
We discuss the transient heat ﬂux derivation later.
Initial and boundary conditions. The initial temperature of thereservoir, along with the ﬂuid in the wellbore, is the reservoir
temperature, Tr, K.

Tðn; t ¼ 0Þ ¼ Tr ; n ¼ x; y; z(4)where, s is the Laplace transform parameter; symbols with “-”means parameters in Laplace domain.
Duhamel's convolution theory method. The Duhamel's
convolution theory has been widely used in Pressure Transient
Analysis (PTA) and Rate Transient Analysis (RTA) in the oil and gas
industry. Convolution theory has also been used by many researchers analyzing heat loss during steam injection for heavy oil
recovery (e.g., [26,27]. In the convolution theory, the conductive
heat ﬂux, qc, can be calculated as a convolution of temperature
derivative with standard conductive heat ﬂux under unit temperature difference:Continuity equation. A small volume of the horizontal wellbore
and its near wellbore region has been chosen (Fig. 3) to demonstrate the heat balance mathematically using the following
equation:Arf cfvTqc
þvx Arf cf

qc ¼ s $ T$qcu

(3)

(10)

Substituting Eq. (10) into Eq. (4), the outlet temperature in
Laplace domain can be expressed as:!

3. MethodologyðTr  Tinlet Þ$exp
Toutlet ¼
s

As mentioned before, determination of the conductive heal ﬂux,
qc(x,t) requires deﬁning the temperatures of the surrounding rock(outside) and working ﬂuid (inside) of the wellbore, which areunknown. The multilateral wellbore conﬁguration adds additional

Natural Coupling Method. In more complicated situations,involving multilateral wells in a heterogeneous reservoir, areservoir-scale analytical thermodynamic model is needed toevaluate the thermal production performance. In this study,
979

sL
L

$ qcu  s $Avrf cfv(11)

Fig. 3. Schematic diagram of fundamental heat transfer model of Closed-Loop geothermal energy recovery system.

Fig. 4. Schematic diagram of two sub-systems examined in the Natural Coupling Method.
980

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

Fig. 5. Schematic diagram of N segments divided from the total length of lateral usedin sub-system #1.

analytical reservoir modeling methodology has been applied to our
problem using integration of instantaneous source and sink functions [28]. This methodology has been successfully applied tosolving transient ﬂuid ﬂow problems in the oil and gas industry[30e34,40]. The original system can be seen as a combination of
two separate sub-systems (Fig. 4). Sub-system #1 is a cubic shaped
hot reservoir with a certain length line sink lying in the center ofthe reservoir. Heat is transferred into the line sink. Inﬂow heat ﬂux
varies in different part of the sink. Sub-system #2 is a cylindershaped system. The working ﬂuid is transported in one direction
at constant velocity. Along the direction of ﬂuid movement,
convective heat transfer dominates the overall performance, while
conductive heat ﬂux ﬂows laterally.
In sub-system #1, the line source can be divided into N segments(Fig. 5). In each segment, the heat ﬂux is uniform and the corresponding temperature change caused by this sink segment can be
expressed as:

DT i ¼

ðt
0Fig. 6. Schematic diagram of sub-system #2 and its N segments.

As a result of the superposition principle, the temperaturechange at any evaluation point is the summation of the contributions from all horizontal source/sink segments:

i¼1

i
i
i
$Tsx
ðt  tÞTsyðt  tÞTszðt  tÞdt

(18)

In sub-system #2, the wellbore is also divided into N segments(Fig. 6), and a uniform conducive heat ﬂux is assumed along eachhorizontal interval. Each segment has the same temperature value
at all joint points. As a result, the heat transfer governing equations
in each segment can be expressed as:

(19)

"
"
"
"
(
#
#
#
#)
∞
TRi  x þ 2nxe
TLi  x þ 2nxe
TRi þ x þ 2nxe
TLi þ x þ 2nxe
1 X
1
1
1
1
p
ﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
p
ﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
p
ﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
p
ﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
erf
 erf
þ erf
 erf
Lseg ∞ 2
22
2
2 aðt  tÞ
2 aðt  tÞ
2 aðt  tÞ
2 aðt  tÞ

(13)

(12)

i ðt  tÞ, T i ðt  tÞ, T i ðt tÞ are source/sink function in x, y, z
Tsx
sy
sz

direction, which can be expressed as following:iðt  tÞ ¼
Tsyiðt  tÞ ¼
Tsz

Lseg ¼

0

rr cr

8
>
vDT i
vDT i
qi
>
>
¼
v
þ
; 0  x  Lseg
>
>
>
vt
vx
ALseg rf cf
>
>
>
>
<
DT i ðx; t ¼ 0Þ ¼ 0; 0  x  Lseg
>
>
>


>
>
DT i ðx ¼ 0; tÞ ¼ DT i1 x ¼ Lseg ; t
>
>
>
>
>


:
DT i x ¼ Lseg ; t ¼ DT iþ1 ðx ¼ 0; tÞ

qi ðtÞ iii$T ðt  tÞTsyðt  tÞTszðt  tÞdt
rr cr sxiTsxðt  tÞ ¼

N ð i
X
q ðtÞ
t

DTðx; y; z; tÞ ¼∞
X

1
pﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃ
∞ 2 paðt  tÞ∞X

1
pﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃﬃpaðt  tÞ
2
∞

L
N

"
"
2 #
2 #)
yi  y þ 2nyeyi þ y þ 2nye
þ exp
exp
4aðt  tÞ
4aðt  tÞ

(14)

"
"
2 #
2 #)
zi  z þ 2nzezi þ z þ 2nze
þ exp4aðt  tÞ
4aðt  tÞ

(15)

(

(
exp

Sub-systems coupling. The two sub-systems share the same
conductive heat ﬂux and the same temperature value at the midpoint at each segment. In sub-system #1, the temperature change
at any location in Laplace domain can be expressed as:(16)

where, TLi , TRi , yi , zi are the location coordinates of ith segment in the
3D system. “erf” represents the error function deﬁned in Eq. (17):2erf ðxÞ ¼ pﬃﬃﬃ

pðxDTðx; y; zÞ ¼
t 2

e

dt

N
X
i¼1

(17)


qi $L

1 i
i
i
T ðtÞTsyðtÞTszðtÞ
rr cr sx
(20)

Eq. (20) indicates that the temperature change value at the midpoint of each segment can be seen as a function with N unknown

0981

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991
Table 1
Parameters of cold-water injection case after unit converted from Ramey [26]].

variables, qi ði from 1 to NÞ, with respect to the location (Eq. (21)).
The source/sink function of each segment in Laplace domain, L[],
could be numerically calculated by Chebyshev's Theorem.

iDTxi midpoint ¼ f#1
q1 ; q2 ; q3 ; …; qi ; …; qN ; i from 1 to N

(21)

In sub-system #2, the temperature value in Laplace domain ateach interval can be derived from Eq. (19):

0
1
8
s
s
qi
>
i
>
DT i ðxÞ ¼ DTx¼0$ev x þ @1  ev x A; 0  x  Lseg>>>
sLseg Arf cf
>
>
>
>
>
>
>
>
>
i
i1
<
DTx¼0
¼ DTx¼L
seg
>
>
>
>
>
>
>
>
>
>
>
>
>
>
:

(22)
As such, the temperature change at the mid-point of eachsegment is a function of location, and N unknown
variables, qi ði from 1 to NÞ.





(23)

Because the two sub-systems share the same temperaturechange at mid-point (DTxi midpoint ) and conductive heat ﬂux value as
the sink or source (qi ) of each segment, N equations with N unknown variables can be generated and solved by Gaussian elimination (Eq. (24)). The mathematical calculations of source and sink
functions, matrix computations, and following Laplace inverse
process are realized by Cþþ coding in Microsoft® Visual Studioenvironment.

2

a11
4«
aN1

/
1
/

3 2
3
C1
q1
6
6
7
7
a1N
6« 7 6« 7
5$6 qi 7 ¼ 6 Ci 7
«
6
7 6
7
aNN 4 « 5 4 « 5
CN
qN
3

Unit

21.11
14.72
0.015130.16
1828.88.732
0.5867990.654180
2.423
2.352e6



C
C

C/m
m
m
kg/s
W/m/K
kg/m3
J/K/kg
W/m/K
J/K/m3

There is no publicly available dataset for the validation of the
proposed analytical methods. We use data presented in a coldwater injection case study by Ramey [26] to validate our models
by comparing the computed results with the measured watertemperature records. The cold water was injected through a vertical
well in a reservoir with a pre-determined geothermal gradient. The
water injection rate was maintained at a constant rate, and thetemperature was measured along the wellbore after approximately
75 days. The parameters of reservoir and operation are listed inTable 1, and the results are shown in Fig. 7.
Fig. 7 presents a comparison of temperatures measured along
the vertical wellbore with temperatures computed by Ramey's
method and the two proposed analytical methods. The temperatures calculated from this study ﬁt very well with Ramey's work.
And the computed value was within 0.83  C of the measured
temperatures, which shows good application of this work on theheat transfer evaluation.

ðTr  Tinlet ÞsiDTxi midpoint ¼ f#2
q1 ; q2 ; q3 ; …; qi ; …; qN ; i from 1 to N

Value

Surface temperature
Cold water temperature at inlet
Geothermal gradient
Inner diameter
Length of the vertical well
Water injection rate
Thermal conductivity of water
Density of waterSpeciﬁc heat capacity of ﬂuid
Thermal conductivity of reservoir
Volumetric heat capacity of reservoir

4. Validation

i
iþ1
DTx¼L
¼ DTx¼0
seg

1
DTx¼0
¼

Parameter

5. Results
An application example case of a closed-loop system deﬁned bymodel parameters in Table 2 was used for developing a general

2

(24)

Stehfest inverse algorism. The solution from the matrix
calculation is done in the Laplace domain. Stehfest numericalinversion (Eq. (25)) was applied to convert the solution into a realtime domain [35].

f ðtÞ ¼



NL
ln2 X
ln2
$i
Vi F
t i¼1
t

(25)

where Vi can be calculated as following and NL is either 6, 8, or 10;
f(t) is the result in a real time domain; F() means the value inLaplace domain:  min i;N 2

Vi ¼ ð1Þ

N
þi
2

2

X

k¼iþ1
2

6
6
4

3
7
k $ð2kÞ!
7
5
N  k !k!ði  kÞ!ð2k  iÞ!
2
N
þ1
2

(26)

Fig. 7. Measured and computed temperature for a cold-water injection well.
982

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

temperature results from the heat balance between heat producedout for power generation and the heat loss from the cooled region
near the wellbore (Fig. 9).
The Natural Coupling Method has been applied to show the heattransfer equivalent. The cumulative heat amount produced by ﬂuid,
Qo, and heat loss from the reservoir, Qr, can be expressed as:

Table 2
Parameters of reservoir and the closed loop well system used in the application
example case.
Parameter

Symbol

Value

Unit

Reservoir temperature
Length lateral
Flow rate per lateral
Inner diameter lateral
Inlet temperatureDensity of ﬂuid
Speciﬁc heat capacity of ﬂuid
Thermal conductivity of rockDensity of rockSpeciﬁc heat capacity of rock
Segment length used in Natural Coupling Method

Tr
L
qf
D
Tinlet

150
4000
20
0.156
60
1000
4180
3.5
2500
1100
100


C
m
m3/hr
m

C
kg/m3
J/K/kg
W/m/K
kg/m3
J/K/kg
m

rf

cf

l
rr

cr
Lseg

ðt
Qo ¼ ðToutlet  Tinlet Þ $ rf cf $qf dt

(27)

0

Qr ¼ ∭ ðTr  Tave Þ $ rr cr dxdydz(28)

Two cumulative heat amounts have been plotted over 30 years
in Fig. 10. The two results clearly match, which shows highly accurate calculation of the proposed method in this study. The cumulative heat extracted from a thermal reservoir increases with the
projects operating time, and it will reach 9.23e5 GJ at the end of 30
years.

understanding of heat production performance and for comparison
with our two analytical methods. In this example, a 4000 m closedloop system is completed in a reservoir layer with 3.5 W/m/K
thermal conductivity and at 150  C. The temperature at the end oflateral wellbore can been seen as the outlet temperature of theworking ﬂuid, if heat exchange is neglected in the vertical wellbore.
The outlet temperatures at the end of horizontal wellbore,
calculated over 30 years by the two methods, exhibits good
agreement as shown in Fig. 8. The outlet temperature drops
dramatically from the reservoir temperature (150  C) to 108  C after
the ﬁrst year. The temperature continues to decline to 100  C over
the remaining years. As a result, this closed-loop system can provide a relatively stable energy production with ﬂuids temperature
over 100  C on a 30 year lifespan. This maintained outlet

6. Discussion
The thermal production efﬁciency and performance of a closedloop system depends on many factors. These factors come from thenature of the reservoir, such as the thermal conductivity andtemperature, and the operation variables such as the working ﬂuids
circulation rate, length of the lateral wellbore, and number andspacing of laterals. Sometimes, these factors will affect the thermalenergy production performance together.

Fig. 8. The comparison plot of results from the two methods proposed in this study.
983

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

Fig. 9. Schematic diagram of the heat balance between heat production by hot working ﬂuid cooling down and heat loss from affected reservoir region.

efﬁciency of energy conversion to electric power. As such, it is
crucial to determine the optimum production layer for a closedloop system. For sedimentary basins, thermal conductivity usually ranges from 2 to 5 W/m/K [37], and among common minerals,
quartz has the highest thermal conductivity, around 7 W/m/K. Here
we used a sensitive analysis of thermal conductivity ranging from 1
to 7 W/m/K to test the application example in this study.
Outlet temperatures, under various reservoir thermal conductivities, are plotted in Fig. 11. Results show outlet temperature overa 30 year period is clearly dependent on thermal conductivity. The
highest thermal conductivity, 7 W/m/K, yields the highest outlet
temperature of slightly over 120  C at the end of 30 years, while thelowest thermal conductivity, 1 W/m/K results in the lowest temperature of about 75  C, or half of the original reservoir temperature. The drop in outlet temperature with decreasing thermal
conductivity is nonlinear. When thermal conductivity is lower than
4 W/m/K, the difference in outlet temperature between adjacenttime series become more signiﬁcant (Fig. 11). As thermal conductivity of most sedimentary sequences is less than 4 W/m/K in
sedimentary basins, targeting of higher thermal conductivity layers
for installation of the horizontal leg is important for project success.
Initial local reservoir temperature. Target reservoir temperature is one of the most important parameters to assess if a thermal
reservoir is suitable for commercialized electricity power generation as compared to direct heat utilization. Three different local
reservoir temperature, Tr, have been studied to analyze the outlet
temperature performance and the corresponding temperature
difference between inlet and outlet temperature (Fig. 12), to help
understand heat harvesting efﬁciency.
In general, higher local temperatures will result in higher outlet
temperatures of the working ﬂuids. A reservoir rock at 200  C willproduce over 120  C ﬂuids, while a 100  C reservoir can only produce ﬂuids around 78  C after 30 years of operation. In terms of thetemperature difference between inlet and outlet temperature,higher local reservoir temperatures contribute to greater heat
harvesting, which means a higher temperature reservoir will also
give better heat recovery efﬁciency. While the local reservoir

Fig. 10. Results of cumulative produced heat amount (Qo) and heat loss from affected
reservoir region (Qr).

Thermal conductivity. Reservoir thermal conductivity is themost important parameter to evaluate the feasibility of a closedloop geothermal system to recover energy from a thermal reservoir. Because the lateral wellbore is in a closed loop and there is no
ﬂuid communication between the reservoir and wellbore, heat
conduction is the only heat transfer mechanism in a closed-loopsystem. Reservoir thermal conductivity is highly related to therock solidity (the ratio of the volume of solids to the bulk volume)and local temperature. Usually, higher rock solidity and lower
reservoir temperature will result in higher thermal conductivity[36]. Solidity and temperature also affect the drilling cost and984

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

Fig. 13. Example of the variation of thermal conductivity with depth in northeaster
British Columbia (in the NW part of the WCSB), and the related thermal gradientvariation, from Majorowicz et al. [37].

Fig. 11. Outlet temperature proﬁles over 30 years under various reservoir thermal
conductivities.

drilling and completion. The associated higher capital cost will also
increase the risk. From a geological perspective, a longer wellbore
may exceed the reservoir boundary or increase the possibility ofsigniﬁcant reservoir heterogeneity, which will all increase the uncertainty of the project. Whereas, from the engineering perspective, increasing wellbore length will also bring unexpected
problems related to wellbore integrity. On the other hand, the ﬂuid
velocity within a horizontal wellbore also affects the heat extraction performance. Higher ﬂuid velocity will reduce the contact timeof ﬂuid with the thermal reservoir, lowering the outlet temperature. However, a lower velocity may not meet the required ﬂow rate
for heat production. As such, optimized lateral length and ﬂuid
velocity in the horizontal section need to be considered carefully inthe project design phase.
It is interesting to note that the lateral length and ﬂuid velocity
can compensate for each other's effect on the outlet temperature asshown graphically in Fig. 14 and mathematically in the Duhamel's
convolution model of Eq. (11). The ratio L/v has a dimension of time
and the physical meaning is the required time of the working ﬂuid
to move through the lateral wellbore. This time, denoted as t*, isdeﬁned as residence time in the horizontal wellbore with amathematical form as:

Fig. 12. Outlet temperature and temperature difference between outlet and inlet results under three different local reservoir temperatures.

t* ¼

temperature increases with depth, thermal conductivity varies as afunction of lithology (Fig. 13). So, optimization on both local temperature and thermal conductivity is critical for closed-loop
geothermal energy production project design. A shallower andcooler layer with higher thermal conductivity may produce moreenergy than a deeper hotter layer.
Residence time. In the project design phase, the length of thehorizontal wellbore is one of the key parameters that affect theoverall heat extraction efﬁciency. Generally, a longer horizontal
wellbore will allow the working ﬂuid to have more contact areawith the thermal reservoir. However, increasing the length of thelateral part of the wellbore will also increase the difﬁculty of

L
v

(29)

A larger residence time means more time for ﬂuid contact withthe thermal reservoir, increasing heat gain. Fig. 15 displays the ﬂuid
outlet temperatures for various residence times, demonstrating
how residence time affects heat extraction performance. In Fig. 15,
residence time of the application example case in this study is
3.82 h (black). When residence time is 4.78 h (red), the outlet
temperature is higher over 30 years. When residence time issmaller (e.g., 2.87 h) the outlet temperature is lower. The difference
between the adjacent series could reach up to 10  C. However,larger residence times requires larger lateral length or lower ﬂuid
velocity. A project needs to consider carefully this value to optimise
985

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

temperature based on the principal of superposition. The time anddegree of the temperature drop depends on well spacing, reservoir
thermal conductivity, and other system parameters such as heterogeneity of the reservoir. In order to eliminate or minimize theinterference effect on a closed-loop geothermal project, wellinterference needs to considered in the design phase.
Fig. 16 compares the outlet temperature over 30 years undervarious lateral spacing of a two-lateral well conﬁguration. The
temperature proﬁle from a single lateral conﬁguration was alsoplotted in the ﬁgure for reference. The single lateral case has thehighest outlet temperature since there is no other well competingwith it. The case with 25 m spacing is the ﬁrst to deviate from thesingle lateral case, as a drop in temperature, over the 30-year
production period. The cases with 50 m and 75 m spacing deviate
from the single lateral case later in time. The smallest spacing also
displays the most signiﬁcant well interference. The difference could
be 5  C for each lateral, and the overall heat extraction efﬁciency
could be signiﬁcantly reduced by having horizontal leg spacing too
close. In this study, 75 m would provide the minimum well interference effect. However, larger spacing distance means higher
drilling costs and technical difﬁculty. Detailed simulation andsensitive analysis are needed for speciﬁc project design andoptimization.
Another important factor affecting lateral-interference is thereservoir thermal conductivity. A higher thermal conductivity
means a greater distance of thermal effect on reservoir. In the oiland gas industry, the deﬁnition of radius of investigation helps toevaluate how fast a reservoir can transfer pressure changes. Similarly, in the geothermal industry, the radius of investigation could
be deﬁned as:

Fig. 14. Outlet temperature proﬁles over 30 years under different lateral lengths andﬂow rates.

sﬃﬃﬃﬃﬃﬃﬃﬃﬃ
4lt
ri ¼
rr cr

(30)

Higher thermal conductivity will result in a greater area oftemperature change at any given time. So, a higher thermal conductivity will lead to earlier and more signiﬁcant well interference.

Fig. 15. Outlet temperature proﬁles over 30 years under various residence times.

heat extraction performance based on the fundamental geology
setting.
Multilateral well interference. Closed-loop geothermal energy
recovery technology needs multiple wells to extract enough thermal energy to meet the requirement of electricity power generation; the structure and design of multilateral wells will alsosigniﬁcantly affect heat loss reduction in the vertical production
wellbore, as discussed in the next section.
In a multilateral development, each transient heat drainage area
expands with time. When the drainage areas of adjacent wellbores
merge, the heat recovery will be affected by lowering the outlet

Fig. 16. Outlet temperature proﬁles of a two-lateral instalment with different horizontal spacing.
986

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

assumed as shown in Fig. 19. The surface temperature is assumed tobe 5  C, while the target reservoir has a temperature of 150  C at
3000 m depth. The temperature increases linearly with depth. The
working ﬂuid has a temperature of 60  C when injected at thesurface and will be heated along the 4000 m horizontal wellbore.
Heated ﬂuid will be then transported to surface by a 3000 m production well. The overall process is illustrated in Fig. 19 and theparameters of the system are listed in Table 3.
In the thermal domain, heat transfer between the inner wellbore and outside reservoir occurs where there is a temperature
difference. There are three main parts of the system having heattransfer with the outside reservoir (Fig. 20). The ﬁrst part is thevertical injection wellbore. Because the inlet temperature of theworking ﬂuid at the surface is 60  C, which is higher than thetemperature of the stratigraphic interval above a certain depth,
there will be some heat loss from the downward moving working
ﬂuid. When ﬂuids move to a depth where rock temperature ishigher than that of downward ﬂuid, it starts to gain heat fromsurrounding rocks. Through the entire lateral wellbore, the working
ﬂuid continuously gains heat. At the end of lateral wellbore, heatedﬂuids from multilaterals merge and ﬂow upwards through thevertical production well. The ﬂuid temperature reaches the
maximum at the base of the vertical production well (Fig. 20).
In this study, we applied Duhamel's Convolution Method to
calculate the heat transfer in vertical sections of the wellbore withan assumption of no well interference. While for heat transfer
modeling in the laterals, the Natural Coupling Method was applied
to deal with potential multilateral well interference and the complex geological conditions. Three modeling parts of the systemwere integrated with two connecting points (Fig. 20). While theanalytical method helped generate more accurate results.
Firstly, a simple single lateral case was studied to see how heatexchange in vertical wellbore affects the overall energy extraction
performance. Fig. 21 shows the time series of temperature variations throughout the closed-loop system over 30 years of operation.
The ﬁrst and last 3000 m of the system are the vertical injection andproduction segments. The 4000 m in the middle is the horizontal
part. In general, the temperature decreases from 60  C at the surface and increases from around 1000 m vertical depth. At the end of
the vertical injection section in the ﬁrst year, the temperaturereaches 66  C. The ﬂuid temperature at the end of lateral reaches
111  C. The temperature of the ﬂuids continues to increase in thevertical production section until the temperature reaches about thesame level of the surrounding rocks outside of the wellbore. The
temperature of working ﬂuid decreases to 99  C at the outlet. As aresult, the temperature difference between the end of lateral andoutlet of production well could be 12  C, suggesting that theworking ﬂuids may lose signiﬁcant heat in the vertical production
section in a single lateral setting.
The results from 10 multilaterals are plotted in Fig. 22. High ﬂow
rate ﬂuids (10 times in mass) are injected through the vertical injection well to the wellbore downhole. At the ﬁrst connecting point,
the ﬂuids are separated into 10 laterals and ﬂow towards the production segment. At the end of the laterals, the heated ﬂuids have atemperature of 109  C and merge together in the bottom of the
production well. There is a much larger upward convective ﬂow
mass in the vertical wellbore and conductive heat loss crossing thewall of the wellbore is limited. As a result, the temperature of theworking ﬂuid is 107  C at the outlet, with a loss of 2  C in the upward journey. For the same reason, the temperature increases only
0.5  C at the end of injection section. The multilateral well conﬁguration with much higher ﬂow rates or larger ﬂow mass in the
upward vertical section can reduce heat loss.

Fig. 17. Outlet temperature proﬁles of a two-lateral case of 25 m spacing with different
reservoir thermal conductivities.

Outlet temperature proﬁles for three different reservoir thermal
conductivities for both single-lateral and double-lateral well conﬁgurations are plotted in Fig. 17. The reservoir with 5 W/m/K conductivity has the highest outlet temperature, but with the earliest
and strongest well interference effect. When the thermal conductivity is low, 1 W/m/K, the well interference could be neglected due
to a slow heat transfer mechanism.
When the number of laterals in a well pad is more than two, thelevel of well interference for each lateral is different. To study theheating performance of each lateral, the temperature proﬁles at theend of the horizontal part of each lateral need to be determined.
Fig. 18 shows the temperature at the end of each lateral for a 10lateral multilateral well pad over 30 years under three different
reservoir thermal conductivities. In the ﬁrst several years, thetemperature at the end of each lateral shows the same value,
because there is no well interference or limited interference at thattime. However, later in time, the temperature proﬁles of the middle
laterals are lower than the wells near the edge. That is because thelaterals in the middle location experience the most signiﬁcant wellinterference effect and the heat extraction efﬁciency is reduced bysurrounding laterals. Higher reservoir thermal conductivities will
result in more signiﬁcant well interference and larger temperaturedrops. Large temperature differences between laterals may lead toless heat extraction efﬁciency and its ﬁnal efﬁciency of transition toelectricity power. Operation optimization can minimize the temperature difference of each lateral in the design phase. For example,
variable spacing distances could be applied to multilateral well pad,so that the temperature of middle laterals will be closer to theothers.
Heat exchange in vertical wellbores. The working ﬂuid transported through vertical wellbores in injection and production well
will lose or gain heat depending on the temperature variationencountered. This amount of heat could affect the overall heat
extraction performance. The application example in this study isused to examine this inﬂuence. A typical geothermal gradient is987

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

Thermal production. The most important indicator of a successful geothermal energy project is the thermal production capacity. We use the 10 lateral system to calculate the long-term
thermal energy production. The temperature difference between
the outlet and inlet multiplied by the volumetric heat capacity andﬂow rate yields the heat duty. The corresponding cumulative
thermal production with time is plotted in Fig. 23. This 10-laterals
system with application example reservoir and well parameters
generates 11 MW thermal heat energy in the early production
period, drops to 10 MW at the 6th year, and remains above 9 MW
for the rest of the period. The total heat production at 30 years willreach 2.5e9 kWh. The power consumption by pump at the injector
is limited due to the utilization of the thermosiphon effect of theworking ﬂuid. Because the temperature at the outlet decreases
gradually, the thermosiphon effect will reduce overtime. Assume
the pressure is maintained at 1300 kPa at the injector and theroughly total pressure drop of 700 kPa at the end of 30 years, thepumping power of less than 45.75 kW with 85% pump efﬁciency isneeded over the 30 years. This is just a hint of pump energy needed,
and the accurate estimation of pump power consumption needs toconsider the pipe friction, ﬂow rate, diameter changes, elbows, andpressure losses due to surface piping and equipment system [38].
As a result, this horizontal closed-loop geothermal project could
provide stable and acceptable heat production. However, the conditions to have such outcome are critical. Our base case uses 10
lateral wells of 4000 m long, and, most importantly, 3.5 W/m/K or
higher reservoir thermal conductivity and 150  C or higher reservoir temperature.

7. Conclusions
In this study, two novel analytical methods were proposed toevaluate the closed-loop system. Duhamel's convolution method
used in oil and gas industry was applied in modeling of a single
lateral closed-loop geothermal energy recovery system. For more
complicated cases such as multilateral wells, and heterogeneousand bounded reservoirs, the natural coupling analytical method
was applied to deal with the complexity of the system to help
better understand heat transfer behaviour and interaction between
the working ﬂuid and reservoir in a closed-loop system.
An application example case with typical geothermal reservoir
properties and closed-loop parameters was examined in this study.
The outlet temperature proﬁles over 30 years from two different
approaches are comparable with high accuracy. The heat balance
check between energy production and heat loss from the reservoir
helps validate the proposed methods for the closed-loop
geothermal recovery technology. The application of the methods
to the application example case in this study suggests that closedloop system provides a reasonable stable heat production.
A sensitivity study suggests reservoir thermal conductivity,
determining the efﬁciency of heat transfer from reservoir to theworking ﬂuid, is the most important parameter of the system toevaluate whether a reservoir is suitable for a closed-loop
geothermal energy recovery project. Well interference in a multilateral conﬁguration could be eliminated or reduced by increasing
well spacing. Higher reservoir thermal conductivity will also lead toearlier well interference and larger temperature drops throughtime. However, the lateral length and ﬂow rate compensate for eachother on overall heat extraction performance. For a given thermal
reservoir, optimizing the operational parameters, such as length,
Fig. 18. Temperatures after 30 years for a 10 lateral well system under 25 m wellspacing with 1 W/m/K thermal conductivity a); 3.5 W/m/K thermal conductivity b);
and 5 W/m/K thermal conductivity c).
988

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

Fig. 20. Schematic diagram of conductive heat transfer between total wellbore andoutside reservoir in closed-loop geothermal energy recovery system with consideration of the geothermal gradient.

Fig. 19. Flow diagram of the working ﬂuid in a closed-loop geothermal energy recovery system with consideration of geothermal gradient existing (black dashed line).

Table 3
Parameters of reservoir and well system used in this case of studying heat exchange
in vertical wellbores.
Parameter

Symbol

Value

Unit

Depth
Length lateralReservoir temperature
Surface temperature
Flow rate per lateral
Number of laterals
Inner diameter lateral
Inner diameter of vertical wellbore
Inlet temperature at the surface
Density of ﬂuid
Speciﬁc heat capacity of ﬂuid
Rock thermal conductivity of lateral part
Rock thermal conductivity of vertical part
Density of rockSpeciﬁc heat capacity of rock
Segment length used in Natural Coupling Method

h
L
Tr
Ts
qf
Nl
Dl
Dv
Tinlet

3000
4000
150
5
20
10
0.156
0.21
60
1000
4180
3.5
2.0
2500
1100
100

m
m

C

C
m3/hr

rf

cf

ll
lv
rr

cr
Lseg

m
m

C
kg/m3
J/K/kg
W/m/K
W/m/K
kg/m3
J/K/kgm
Fig. 21. Working ﬂuid temperature variation in a single lateral wellbore.

producing higher outlet temperatures with potential for power
generation. We have restricted our assessment to the thermodynamic viability of this technology. Demonstration projects have
proven that the drilling technology is viable. We have not assessed
the economics of this technology, which will in large part be related
to drilling costs. Overall, this technology removes the need ofpermeable zones which conventional geothermal production
methods require and avoids artiﬁcial fracturing in an enhanced
geothermal system, thereby providing an alternative solution formany high enthalpy geothermal projects. As well, a closed-loopsystem eliminates challenges related to production of geothermal
ﬂuids, including scaling and environmental risks. Further research
can greatly expand the deployment of geothermal development tosupport the transition to clean energy resources.

spacing, ﬂow rate and number of laterals in the project designphase can greatly increase efﬁciency.
The heat exchange effect in vertical injection and production
sections is minimal in a 10-multilateral closed-loop conﬁguration.
With the application example case, the closed-loop geothermalenergy recovery technology could provide 9e11 MW stable heatproduction over 30-years of operation.
Our results show that a deep closed loop system can effectively
produce energy over long time periods under suitable geologicalsettings regarding thermal gradient and thermal conductivity. Bothof these parameters are easier to determine or already known inmany petroleum basins. Our model was based on average conditions and the resultant outlet temperature may be more suitable for
district heating rather than power generation. That being said,
higher temperature target reservoirs would be capable of989

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991

of Energy Research and Development (OERD) provided funding forthis study. The authors appreciate discussions with Eavor on thebasic principals, however the modeling work conducted here wascompletely independent and any errors or omissions are the authors alone. The authors also thank our internal reviewer, Heather
King, and four anonymous external reviewers for their constructive
comments and suggestions on this paper. This paper has been
assigned NRCan contribution number: 20200775.

References[1] International Energy Agency (IEA), Renewables 2018: Market Analysis andForecast from 2018 to 2023, 2018. https://www.iea.org/renewables2018/.
[2] D. Gielen, F. Boshell, D. Saygin, M.D. Bazilian, N. Wagner, R. Gorini, The role ofrenewable energy in the global energy transformation, Energy Strategy Reviews 24 (2019) 38e50, https://doi.org/10.1016/j.esr.2019.01.006.
[3] S. Schumacher, R. Pierau, W. Wirth, Probability of success studies forgeothermal projects in clastic reservoirs: from subsurface data to geologicalrisk analysis, Geothermics 83 (2020) 375e6505, https://doi.org/10.1016/
j.geothermics.2019.101725.
[4] M. Abbasi, N. Khazali, M. Shari, Analytical model for convection-conduction
heat transfer during water injection in fractured geothermal reservoirs withvariable rock matrix block size, Geothermics 69 (April) (2017) 1e14, https://doi.org/10.1016/j.geothermics.2017.04.002.
[5] Y. Zhang, C. Yu, G. Li, X. Guo, G. Wang, Y. Shi, C. Peng, Y. Tan, Performance
analysis of a downhole coaxial heat exchanger geothermal system withvarious working ﬂuids, Appl. Therm. Eng. 163 (September) (2019), 114317,
https://doi.org/10.1016/j.applthermaleng.2019.114317.
[6] Greenﬁre Energy, GreenLoopTM technology overview, 2021. Available on theinternet: https://www.greenﬁreenergy.com/greenloop-technology/.
[7] P.B. Umwelt, Investigations on the Improvement of the Energy Output of a
Closed Loop Geothermal System (CLGS), Technical University of Berlin, 2008.
[8] GeothermamEx Inc, Report on the South Meager Geothermal Resource British
Columbia, Canada, for Western Geopower Corp, 2004, p. 156p.
[9] S.E. Grasby, D.M. Allen, S. Bell, Z. Chen, G. Ferguson, A. Jessop, M. Kelman,
M. Ko, J. Majorowicz, M. Moore, J. Raymond, R. Therrien, Geothermal energy
resource potential of Canada, Geological Survey of Canada, Open File 6914
(2011) 322.
[10] J. Majorowicz, S.E. Grasby, High potential regions for enhanced geothermal
systems in Canada, Nat. Resour. Res. 19 (2010) 177e188.
[11] J. Majorowicz, S.E. Grasby, Geothermal energy for northern Canada e is it
economic? Nat. Resour. Res. 23 (2014) 159e173.
[12] J. Majorowicz, S.E. Grasby, Deep geothermal energy in Canadian sedimentarybasins vs. fossil based energy we try to replace e exergy [KJ/Kg] compared,
Renew. Energy 141 (2019) 259e277.
[13] S. Weides, J. Majorowicz, Implications of spatial variability in heat ﬂow for
geothermal resource evaluation in large foreland basins: the case of thewestern Canada sedimentary basin, Energy 7 (2014) 2573e2594, 2014.
[14] S.B. Mahbaz, A.R. Dehghani-sanij, M.B. Dusseault, J.S. Nathwani, Enhanced andintegrated geothermal systems for sustainable development of Canada’ s
northern communities, Sustainable Energy Technologies and Assessments 37(2020), 100565, https://doi.org/10.1016/j.seta.2019.100565.
[15] Deep Earth Energy Production, About DEEP, 2021. Available on the Internet:
https://deepcorp.ca/.
[16] Eavor Technologies Inc, Available on the Internet: https://eavor.com/about/technology, 2021.
[17] C. Hickson, M. Kumataka, P. Akto, D. Cotterill, D. Benoit, R. Eccles, K. Huang,
M. Colombina, S. Collins, Alberta # 1: the Province's First Electrical
Geothermal Project, vols. 1e13, 2020.
[18] H. Shao, P. Hein, A. Sachse, O. Kolditz, Geoenergy Modeling II, Springer International Publishing, 2016, https://doi.org/10.1007/978-3-319-45057-5.
[19] Z. Hu, T. Xu, B. Feng, Y. Yuan, F. Li, G. Feng, Z. Jiang, Thermal and ﬂuid processes in a closed-loop geothermal system using CO2 as a working ﬂuid,
Renew.
Energy154(2020)
351e367,https://doi.org/10.1016/
j.renene.2020.02.096, 2020.
[20] J.W. Tester, B.J. Anderson, A.S. Batchelor, D.D. Blackwell, R. DiPippo, The Futureof Geothermal Energy - Impact of Enhanced Geothermal Systems (EGS) on the
United States in the 21st Century, 0-615-13438-6, Massachusetts Institute ofTechnology, 2006. Available from: the Internet: www1.eere.energy.gov/
geothermal/pdfs/future_geo_energy.pdf.
[21] M. Akhter, J. Mallams, X. Tang, A. Kazi, Y. Kao, S. Kumar, B. Tai, D. Antao,
A. Palazollo, D. Staack, Shockwave and plasma accelerated rock cracking(SPARC) for hard rock drilling, in: Bulletin of the American Physical Society,73rd Annual Gaseous Electronics Virtual Conference, 2020. October 5e9,2020.
[22] J. Finger, D. Blankenship, Handbook of Best Practices for Geothermal Drilling,
2010. Sandia Report (SAND2010-6048), http://artikel-software.com/ﬁle/geothermal drilling handbook.pdf.
[23] Geothermal Engineering ltd., 2021. Available on the Internet: https://geothermalengineering.co.uk/united-downs/#.
[24] K. Huang, C. Hickson, D. Cotterill, Y. Champollion, Geothermal assessment ofFig. 22. Temperature proﬁles along total wellbore length of 10-laterals case.

Fig. 23. Heat duty and cumulative thermal production proﬁles over 30 years of 10laterals case.

CRediT authorship contribution statement
Wanju Yuan: Methodology, Formal analysis, Software, Writinge original draft. Zhuoheng Chen: Conceptualization, Supervision,
Writing e review & editing, Funding acquisition. Stephen E.
Grasby: Conceptualization, Writing e review & editing. Edward
Little: Project administration, Writing e review & editing.
Declaration of competing interest
The authors declare that they have no known competing
ﬁnancial interests or personal relationships that could haveappeared to inﬂuence the work reported in this paper.
Acknowledgements
This work originated as an independent evaluation of deep
closed-loop system and is an output of the Geoscience for New
Energy Supply (GNES) Program of Natural Resources Canada. Ofﬁce
990

W. Yuan, Z. Chen, S.E. Grasby et al.

Renewable Energy 177 (2021) 976e991
bounded reservoir with multistage fractured horizontal wells, J. Nat. Gas Sci.
Eng. 45 (2017b) 272e290, https://doi.org/10.1016/j.jngse.2017.05.003, 2017.
[34] C. Su, Semi-analytical Modeling of Fluid Flow in and Formation Evaluation of
Unconventional Reservoir Using Boundary Integration Strategies, PhD
Dissertation, Faculty of Graduate Study and Research, University of Regina,
2018. July, 2018.
[35] H. Stehfest, Numerical inversion of Laplace transforms, Commun. ACM 13 (1)(1970) 47, 1970.
[36] E.C. Robertson, Thermal Properties of Rocks, vol. 88, US Department of theInterior: Geological Survey, 1988.
[37] J.A. Majorowicz, F.W. Jones, A.M. Jessop, Preliminary geothermics of thesedimentary basins in the Yukon and Northwest territories (60-70 degrees
north) - estimates from petroleum bottom-hole temperature data, Bull. Can.
Petrol. Geol. 36 (1988) 39e51.
[38] S. Schulz, Investigations on the Improvement of the Energy Output of a Closed
Loop Geothermal System (CLGS), PhD Dissertation, Faculty of Planning,
Building, Environment, Technical University of Berlin, 2008, 2008.
[39] C.M. Oldenburg, L. Pan, M.P. Muir, A.D. Eastman, B.S. Huggins, Numerical
Simulation of Critical Factors Controlling Heat Extraction from Geothermal
Systems Using a Closed-Loop Heat Exchange Method, Lawrence Berkeley
National Laboratory Recent Work, 2019. Permalink: https://escholarship.org/
uc/item/4nj3s9ws.
[40] G. Zhao, Reservoir modeling method.” U. S. Patent provisional application No.
61/226, 197, ﬁled on july 16, 2009. Pub. No.: US 2011/0015909 A1. January 20,
2011, US Patent No 8 (2012) 275, 593, granted on Sept. 25.

target formations using recorded temperature measurements for the Alberta
No. 1 geothermal project, Appl. Sci. 11 (2) (2021) 1e10, https://doi.org/
10.3390/app11020608.
[25] M. Toews, “Method and Apparatus for Power Production” U. S. Patent
Application No. 16/181,492, Filed on Nov. 6. 2018, 2019. Pub. No.: US 2019/0154101 A1. May 23.
[26] H.J. Ramey, Wellbore heat transmission, J. Petrol. Technol. (1962) 427e435.
[27] J. Hagoort, Ramey's Wellbore heat transmission revisited, SPE J. 9 (4) (2004)
465e474, https://doi.org/10.2118/87305-PA e441.
[28] H.S. Carslaw, J.C. Jaeger, Conduction of Heat in Solids, second ed., University
Press, Oxford, 1959, p. 339. Oxford at the Clarendon Press.
[29] W.H. Press, S.A. Teukolsky, W.T. Vetterling, B.P. Flannery, Numerical recipes in
C. The Art of Scientiﬁc Computing, second ed., Cambridge University Press,
USA, 1992, ISBN 978-0-521-43108-8.
[30] G. Zhao, Modeling Complex Natural Fracture Network in Heterogeneous Tight
Formations Using Semi-analytical Strategy, Society of Petroleum Engineers,
2013, https://doi.org/10.2118/167127-MS.
[31] G. Zhao, L. Xiao, C. Su, Z. Chen, K. Hu, Model-based type curves and theirapplications for horizontal wells with multi-staged hydraulic fractures, Journal of Canadian Energy Technology Innovation (CETI) 2 (3) (2016) 29e43.
January 2016.
[32] L. Xiao, G. Zhao, H. Qing, A compatible boundary element approach withgeologic modeling techniques to model transient ﬂuid ﬂow in heterogeneous
systems, J. Petrol. Sci. Eng. 151 (2017a) 318e329, https://doi.org/10.1016/
j.petrol.2017.01.011, 2017.
[33] L. Xiao, G. Zhao, H. Qing, Depletion-induced stress change in a hydraulically991