Princetonlaan 63584 CB Utrecht
P.O. Box 80015
3508 TA Utrecht
The Netherlands

TNO report TNO 2019 R11549 (update 21-05-2021)

draft

Eavor Loop Audit ReportDate

21 may 2021

Author(s)J.D. van Wees

Copy no.
1
No. of copies
1
Number of pages
18
Number of appendices 1
Customer
Eavor
Project name
Eavor Audit report
Project number
060.42026

All rights reserved.
No part of this publication may be reproduced and/or published by print, photoprint,
microfilm or any other means without the previous written consent of TNO.
In case this report was drafted on instructions, the rights and obligations of contracting
parties are subject to either the General Terms and Conditions for commissions to TNO, or
the relevant agreement concluded between the contracting parties. Submitting the report for
inspection to parties who have a direct interest is permitted.
© 2021 TNO

www.tno.nl
T +31 88 866 42 56
F +31 88 866 44 75

TNO report TNO 2019 R11549 (update 21-05-2021)

Contents
1

Introduction and summary ..................................................................................... 3

2
2.1
2.2

Model benchmark .................................................................................................... 5
Robustness of Analytical and 2-dimensional axisymmetric solution ......................... 5
3d numerical simulations ........................................................................................... 83
3.1

Eavor loop Lite assessment ................................................................................... 9
Sensitivity to cement conductivity contrast, different well bore radius for main. ..... 10

4
4.1
4.2
4.3

Assesment of Eavor Loop Lite demonstration in Canada Q4 2019/Q1 2020 .. 12
Model Setup ............................................................................................................. 13
Model Results .......................................................................................................... 14
Conclusions ............................................................................................................. 14

5
5.1
5.2

Update of Assesment of Eavor Loop Lite demonstration in Canada Q4 2019/Q1
2021 ......................................................................................................................... 15
Model Results .......................................................................................................... 15
Conclusions ............................................................................................................. 16

6

Suitability of methodology for SDE+ ................................................................... 17

7

Conclusions ........................................................................................................... 18

8

References ............................................................................................................. 19

Annex 1 : python code analytical solution lateral .............................................................. 20
Annex 2 : python code analytical Eavor Lite ...................................................................... 21

2 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

1

Introduction and summary
Eavor has proposed an innovative well design for geothermal energy production, theEavor LoopTM (Figure 1). The Eavor Loop consists of a closed loop system, in whichheat is extracted from the deep subsurface through multi-lateral wells acting aseffective heat exchangers.
Eavor has proposed quantitative predictive models for the assessment of the thermal
performance and underlying uncertainty of the Eavor loop system (Price et al., 2019).
Eavor has requested in TNO to perform an audit to validate the methods used.
Furthermore, Eavor has requested TNO to evaluate the evaluation results of theEavor Loop demonstration in Canada performed in Q4 2019/Q1 2020 (Price et al.,2020), and later for this update to evaluate the Eavor Loop demonstration for a longer
period Q4 2019/Q1 2021. In addition, Eavor has requested TNO to evaluate if themethods proposed can be used in the framework of applications to the Dutch feed-in
scheme for geothermal systems (SDE+).
The finding of this audit report are the following:• The presented models for uniform operating conditions are robust and wellcapable of predicting the performance of the Eavor Loop
• Analytical approaches appear in excellent agreement with the axisymmetricmodel results for production temperatures• Results for varying operating conditions appear correct given the fact thatnumerically they are based on commercial codes and served as basis for the
excellent benchmarked results for uniform operating conditions. It is wellknown from previous studies on DBHE systems that the type of used
axisymmetric models performs well for varying operating conditions• The 3D numerical models nicely demonstrate the effects of interference andreduced performance when multilaterals are placed at insufficient spacing.
Otherwise, the axisymmetric and analytical solutions suffice.
• The Eavor Loop Lite demonstration shows production temperatures which
are consistent with analytical model prediction, marked by an excellent
correspondence of prediction and observations up to 480 days of running theloop. This close match also demonstrates that the flow rate and inlet
temperature is sufficiently stable to use an analytical model for performanceassessment.
• For SDE+ applications the analytical approach could be favoured as itprovides a transparent, fast and effective basis for assessing and judging theperformance of the Eavor Loop. Expectation curves for performance caneasily be generated by taking into account uncertainty in relevant
parameters, including thermal conductivity and geothermal gradient.

3 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

Figure 1 The Eavor LoopTM geothermal production configuration (source Eavor)

4 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

2

5 / 23Model benchmark
Eavor™ has presented three approaches of increasing complexity to predict thethermal output of Eavor-Loop™ projects; using
1. an analytical solution for a single lateral
2. 2-dimensional axisymmetric model for a single lateral
3. 3-dimensional model for multi-laterals
Each of these approaches has been described in detail, and a comparison of theresults has been shown for a benchmark case with 1-12 laterals. This benchmark
case deals with the (multi-) lateral section (thermal effects of the vertical well paths
are neglected) with a length of 5000m, inlet and reservoir temperature, andhomogeneous thermal rock properties as listed in Table 1. The flow rate is assumed
constant through time running for 25 years.

2.1

Robustness of Analytical and 2-dimensional axisymmetric solution
The results for the production temperature at the end of the lateral show an extremely
good match between the analytical solution and the 2D axisymmetric solution
reproduced in Figure 2 in this audit. The analytical solution is a so called Line Source
Model (Ramey, 1962; Carslaw and Jaeger, 1959; Kutun, 2015):
𝑇(𝑟, 𝑡) =

𝑞
2𝜋𝜆𝐿

𝑓(𝑡)

(eq. 1)𝑓(𝑡) = ln(1.4986√𝑡𝐷 ) =1
2( ln (

4𝑘𝑡
𝑟2) − 0.5772)

cf. Carlslaw& Jaeger, 1959, Ramey, 1962

𝑤𝑖𝑡ℎ
𝑘𝑡

√𝑡𝐷 = √ 𝑟2

where 𝑇(𝑟, 𝑡) is the temperature change relative to ambient temperature at the wallof the well [K], q – heat extracted at the well per unit length (negative [W]), 𝜆 – rockconductivity [J K-1 m-1], 𝑘 – rock thermal diffusivity [m2 s-1], 𝑟 – radius of well [m], L islength of the well. The derivate of 𝑇(𝑟, 𝑡) at time t [s] is a constant relating temperaturechange as a function of f(t):1𝑘𝑧=

𝑑𝑇(𝑟,𝑡)
𝑑𝑞

=

1
2𝜋𝜆𝐿

𝑓(𝑡)

(eq. 2)It has been recognized by many authors that equation (1), adopting 𝑓(𝑡) =ln(1.4986√𝑡𝐷 ) reproduces very well the characteristics of Deep Borehole Heat
Exhangers (DBHE) (e.g Kujawa and Nowak, 2000; Van Wees and Lokhorst, 2005;
Van Wees et al., 2007; Sapinska-Sliwa et al., 2015). In the analytical approximation
for wellbore temperature evolution at time ts we can write a first order differential
equation for Θ = 𝑇(𝑟, t 𝑠 ) = (𝑇 − 𝑇𝑥 ), where 𝑇 is temperature in the well and 𝑇𝑥 isambient temperature which can change linearly along the well bore as a function ofgeothermal gradient (e.g. Kujawa and Nowak, 2000). For a transformed coordinate𝑋𝑥 = , 𝑥 ∈ (0,1):
𝐿

𝑑Θ
𝑑𝑥

+ 𝑘𝑧 ∗ Θ + 𝐸 ∗ = 0(eq. 3)

TNO report TNO 2019 R11549 (update 21-05-2021)

Where𝑘𝑧 ∗ =𝑘𝑧
𝑚̇ 𝑐𝑝, 𝑇𝑥 = 𝐸 ∗ 𝑥 + 𝐹, 𝐸 ∗ = 𝐸𝐿, 𝐸 – geothermal gradient [ K m-1] , 𝐹 -surface temperature [K],1].

𝑚̇ - mass rate [kg s-1], 𝑐𝑝 - fluid heat capacity [ J kg-1 K-

The integral solution becomes:
Θ= −

𝐸∗
𝑘𝑧 ∗

+ (Θ𝑥=0 +𝐸∗∗) 𝑒 −𝑘𝑧 𝑥

𝑘𝑧 ∗(eq. 4)

With Θ𝑥=0 is the value of Θ at the entry of the fluid in the well bore. In absence athermal gradient, eq. 4 matches with the analytical formulation in Price et al., 2019
(their equation 2.18)
𝑇𝑟 −𝑇𝑓𝑜𝑢𝑡
𝑇𝑟 −𝑇𝑓𝑖𝑛

2𝜋𝜆𝐿= 𝑒𝑥𝑝 (−𝑚̇ 𝑐𝑝

𝑙𝑛(1.4986

Parameter
reservoir temperature
length lateral
Flow rate per lateral
Inner diameter lateral
Inlet temperature
Rhofluid
Cpfluid
Krock
Rhorock
Cprock√𝑘𝑡)𝑟𝑤) (eq. 5)

value

Table 1 benchmark case parameters

unit
120 C
5000 m
18 m3/h
0.156 m
60 C
1000 kg m-3
4180 J K-1 kg-1
4 W m-1 K-1
2663 kg m-31112 J K-1 kg-1

6 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

Figure 2 Comparison of 2D numerical axisymmetric model outcome from Price et al., 2019 for asingle lateral (from their figure 3) and the analytical model (equation 1, annex 1).
Parameters of the model are listed in . See text for explanation.

The results for the benchmark model in the report and reproduced here (Figure 2)
clearly show that the borehole geothermal performance can well be predicted with ananalytical model when operational conditions are constant over time (i.e. flowrate andinlet temperature), and for homogeneous reservoir conditions.
The 2D numerical models have not been numerically checked in this audit. The
approach of Price et al. (2019), using commonly adopted numerical schemes,commercial codes and high mesh resolution. They appear very robust, as is also
proven by the close match between numerical and analytical results. This is also
consistent with earlier findings of axisymmetric models applied to DBHE studies,including adopted model assumptions such as on effective heat transfer at the wellbore interface to the fluids inside the well due to turbulent flow (e.g. Kohl et al., 2002;
Lokhorst and Van Wees, 2005). Price et al., (2019) do adopt both the Joule Thomson
effect and gravitational potential energy in the numerical solution they use. Ignoring
the latter effect in adopting the Joule-Thomson effect in geothermal producton
boreholes can lead to an erroneous temperature rise of ca 2 ⁰C per km (Stauffer etal., 2014).
The comparison shows that analytical solutions suffice very well for the performanceassessment, provided operating conditions are relatively uniform in time. If flow rates
and inlet temperatures strongly vary over time, the analytical assessment can have
limitations (see further in chapter 4).

7 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

2.2

3d numerical simulations
The analytical and axisymmetric solutions are adopting infinite radius for cooling,which is fine, as the zone which is affected by cooling is limited and shows notinterference with the surface or other laterals. The cooling zone has a cylindricalshape with its axis at the borehole and a radius 0.5 x which is a function of time tand depends on the diffusivity, i.e. reproducing equation 2.4.1 from Price et al., 2019:𝑘𝑡

𝛥𝑥 = 2√𝐾(eq.5)

Where K=1.57 has been set based on best fit with numerical simulation. For thebenchmark case, marked by k=3.65x10-6 m2 s, this results cooling zones of multilaterals which do not interfere in the lifetime for the first 10-15 years when 𝛥𝑥 = 75 𝑚.
The thermal effects in 3D simulation (figure 3.) of 12 multilaterals spaced 75 m apartup to 30 years production are indeed minor and in good agreement with theaxisymmetric model. The interference for shorter spacing (25 and 50 meter ) results
in an expected reduction in production temperatures. These models nicely show thatthe multilaterals need to be spaced sufficiently far from each other, which is wellfeasible.
Given the relatively low reservoir permeability it is not expected that fluid circulation
will be triggered by the cooling of the reservoir surrounding the laterals, in such a waythat it can have sufficient impact to alter the thermal calculation results

Figure 3 Comparison of 2D (grey) and 3D numerical axisymmetric model outcome from Price et al.,2019 for 12 multilaterals with different spacing and reference solution with a single
lateral. See text for explanation.

8 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

3

9 / 23

Eavor loop Lite assessment
Price et al., (2019) makes a thermal performance assessment of the Eavor Lite loop,
which consist of a pilot in Canada with a one way loop with 2 multi-laterals, withthermal and well parameters of subvertical trajectories (denoted by main) and laterals
listed in Table 2.
The subvertical (main) path has an along hole length of ca 2565 m and consists of
two sections, which are cased and cemented. The reported inner diameter isinterpreted at the interface of the steel casing and the cement as the heat loss is
considered negligible in the steel.

Parametersurface temperature
Main Length (top+bottom section)
Main Length (top section)
Lateral length
temperature formation
Flow rate
I.D. main (top section)
O.D. main (top section)
I.D. main (bottom section)
O.D.main (bottom section)
I.D. lateral# laterals
TstartCpfluid
Krock
Rhorock
CprockKrockmain
Krockcement
Rhorockmain
Cprockmain

value

unit
5 C
2565 m
610
1920 m
75.95 C40,000 kg/h
0.16 m0.31 m
0.16
0.24
0.156 M
2
18 C
4180 J K-1 kg-1
3.5 W m-1 K-1
2663 kg m-31110 J K-1 kg-1
2 W m-1 K-1
1.2 W m-1 K-12240 kg m-3
1520 J K-1 kg-1

Table 2 Eavor loop lite input parameters

The laterals are approximately 1920 m long.
The axisymmetric model outlet temperature results for constant operating conditions
is almost in exact correspondence (< 0.1 C difference) with theTable 1 stepwise
analytical solutions for the vertical and horizontal well paths (Figure 4) reproduced inpython with annex 2, adopting parameters from Table 2.

TNO report TNO 2019 R11549 (update 21-05-2021)

Figure 4 Outlet temperature solution and intermediate temperature solutions and the start and endof the laterals for the Eavor Lite Loop (solid: the EAVOR axisymmetric numericalcalculation; dots the piece-wise analytical solution for the Eavor Lite- loop in annex 2).

Please note that the intermediate temperatures at the inlet and outlet of the laterals
are slightly different with larger difference (up to ca 1 C) between numerical andanalytical solution. These differences are primarily attributed to a) the effect of
methalpy in the numerical solution (including Joule-Thompson effect and gravitational
potential energy). In addition there may be a minor effect of the longer residence timeat relatively high depth and temperatures of the (sub) vertical main in the numericalmodel, which takes into account the full deviation survey for the numerical calculation

3.1

Sensitivity to cement conductivity contrast, different well bore radius formain.

The axisymmetic numerical solution does not allow for variations in conductivity inradial direction. We did check with the analytical solution what could be the effect of
adopting a realistic cement conductivity in between the inner and outer diameter inthe borehole. For the subvertical sections the calculation of 𝑘𝑧 is (cf. Lokhorst and
Van Wees, 2004):𝑟ln( 𝑟𝑜 )

1
𝑘𝑧

=

𝑓(𝑡)

( 𝜆 𝑖 + 𝜆 )𝑑𝑇(𝑟𝑖 ,𝑡)
𝑐
=
𝑑𝑞
2𝜋𝐿

(eq. 6)

10 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

Where 𝑟𝑜 and 𝑟𝑖 are outer and inner radius of the vertical sections and 𝜆𝑐 is the thermal
conductivity of the cement. 𝑓(𝑡) is calculated taking 𝑟 = 𝑟𝑜 . When the inner and outer
diameter are equal eq. 6 reduces to eq. 2.
The results for a cement conductivity of 𝜆𝑐 = 1.2 do show a very small difference(<<0.1 C) with the reference case presented in Figure 4. Similarly, increase of theinner diameter by 10% of the vertical main has hardly influence on the results.

11 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

4

Assesment of Eavor Loop Lite demonstration in
Canada Q4 2019/Q1 2020
Eavor has performed Q4 2019/Q1 2020 a demonstration project of the Eavor loop in
Canada (Eavor Lite). The setup of the Eavor Loop demonstration is given in Figure 5
and Eavor’s detailed evaluation of the results of the demonstration is given in Priceet al., 2020. In this section we perform an independent assessment of thedemonstration results, in particular in view of thermodynamics analysis.
The geometric setup of the demonstration site and geological conditions are inaccordance with section 3, Table 2, except for flow rate, inlet and outlet temperaturesand pressures, and adopted thermal conductivities.

Figure 5 The schematic layout of the Eavor loop Lite demonstration (Price et al., 2020).

Eavor constructed the loop successfully mid 2019 and operated the loop for threemonths starting December 4 2019, up to March 4 2020. During that time, Eavor
measured inlet and outlet temperatures and pressures, flow rate (Figure 6), as well
as fluid losses. The latter were minor (average 0.72 m3/day), and are unlikely to affectthermodynamics.

12 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

13 / 23

Figure 6 Measured inlet and outlet temperatures, flow rate of the Eavor Lite Loop (data provided by
Eavor), starting date is December 4 2019, end date is march 4 2020

4.1

Model Setup
TNO applied its independent analytical model (section 2) to evaluate the Eavor Loopdemonstration outcomes. It has been assumed that the analytical model can be usedas the flow rate conditions and inlet temperature have been relatively constant (Figure6). The aim was to determine if the measured outlet temperature given the flow rateand inlet temperature is consistent with the Eavor Lite loop well layout, ambient rockproperties and reservoir temperature. The model adopts parameters as listed inTable 2, with modified flow rate, inlet temperature and thermal conductivities inagreement with the yellow marked rows in Table 3. The incorporated rockconductivities have been calibrated by Eavor based on history matching with theirnumerical performance model, and deviate little (<3%) from pre-spud estimates forthe horizontal laterals (Price et al., 2020).

Parameter
Value Units
Flow Rate
23.0 m3/h
Inlet Temperature
23.4 ºC
Inlet Pressure
31.0 kPag
Outlet Temperature
52.4 ºC
Outlet Pressure
304.9 kPag
Horizontal Thermal Conductivity
4.64 W/mK
Vertical Thermal Conductivity
2.25 W/mK

Comments
Average over time span
Average over time span
Average over time span
Last value
Last value
History Matched Dec 4 - Feb 2
History Matched Dec 4 - Feb 2

Table 3 Parameters in the Eavor Loop Lite demonstration (modified from Price et al., 2020). Thosemarked in yellow have been adopted in the model.

TNO report TNO 2019 R11549 (update 21-05-2021)

4.2

Model Results

Figure 7 shows the predicted temperatures of the Eavor Loop Lite demonstration.
The predicted outlet temperatures show an excellent match with the observed outlet
temperatures after 30 days of operation.

Figure 7 Eavor Loop Lite predicted temperatures (dotted lines) at the subsurface inlet of the lateralsegment, end of the lateral, and the outlet temperature.

4.3

Conclusions
The Eavor Loop Lite demonstration shows production temperatures which areconsistent with analytical model prediction, marked by an excellent correspondenceof prediction and observations after 30 days. This close match also demonstrates
that the flow rate and inlet temperature is sufficiently stable to use an analytical model
for performance assessment.

14 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

5

Update of Assesment of Eavor Loop Litedemonstration in Canada Q4 2019/Q1 2021
The data on the Eavor Loop Lite demonstration have been updated based on data
running up to Q1 2021, covering 480 days of operation

5.1

Model Results
Figure 8 shows the predicted temperatures of the Eavor Loop Lite demonstration,with the same parameter settings as in section 4. The long term predicted outlet
temperatures show an excellent match with the observed outlet temperatures up to
480 days of operation. The short lived deviations from the long term trend are aconsequence of temporary modifications (or arrest) of the flow rate. These short termchanges in flow rate cannot be incorporated in the semi-analytical model we use,
however we expect that the effects of the short term deviations in flow rate are limited
for the long term production temperature. Furthermore, we note that the inlet
temperature in the second half of the run time is about 1⁰C higher than in the first
half. In the semi-analytical model we cannot vary the inlet temperature over time.
However, if we run the model with 1⁰C higher inlet temperature for the whole duration
the effects on the resulting production temperature are minor.

Figure 8 (top) Eavor Loop Lite predicted temperatures (dotted lines) at the subsurface inlet of thelateral segment, end of the lateral, and the outlet temperature (adopting parametersettings from section 4); (bottom) with 1⁰C higher inlet temperature.

15 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

5.2

Conclusions
The Eavor Loop Lite demonstration shows production temperatures which areconsistent with analytical model prediction, marked by a very good agreement ofprediction and observations up to 480 days of running the loop.

16 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

6

Suitability of methodology for SDE+
The Eavor Loop is a novel geothermal production concept. In this document I
assume that RVO can accept the Eavor concept in the 2020 geothermal <4000m
category for SDE+, but this needs to be discussed and confirmed from RVO.
If RVO accepts it in the current category, the SDE application requires an expectation
curve for the produced geothermal power. For the P50 estimate, I think theapproach of Price et al. (2019) for predicting an expectation curve of production
temperature at given flow rate is appropriate. The resulting P50 power value isprimarily a function of the flowrate [kg/h], thermal properties (which are varied),and I think this may also need to include uncertainty in other subsurface parameters
beyond engineering control if this affects the power prediction (e.g. the thermal
gradient).
To correctly forecast thermal output of the Eavor loop, thermal conductivity andspecific heat of the rock must be correctly estimated. Please note that the adopted
values for the benchmark and the Eavor loop lite are specific to that geologicalformation.
For an estimate of performance for SDE+ the analytical models can be used,provided laterals are sufficiently spaced (see 2.2), and under the assumption thatduring the SDE load hours the system runs at a flow rate higher or equal that theremainder of the year. The analytical (or axisymmetric) models need to incorporate
the effects of heating and cooling of the vertical trajectories of the Eavor loop as
these have an effect on the thermal outcome (Figure 4). Interestingly, the effects of
the Joule-Thomson and gravitational potential energy effects appear to be of minor
influence for the predicted outlet temperatures.

17 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

7

Conclusions

The finding of this audit report are the following with respect to the assessment ofappropriate predictive models for the hydrothermal performance of the Eavor Loop:• The presented models for uniform operating conditions are robust and wellcapable of predicting the performance of the Eavor Loop
• Analytical approaches appear in excellent agreement with the axisymmetricmodel results for production temperatures• Results for varying operating conditions appear correct given the fact thatnumerically they are based on commercial codes and served as basis for the
excellent benchmarked results for uniform operating conditions. It is wellknown from previous studies on DBHE systems that the type of used
axisymmetric models performs well for varying operating conditions• The 3D numerical models nicely demonstrate the effects of interference andreduced performance when multilaterals are placed at insufficient spacing.
Otherwise, the axisymmetric and analytical solutions suffice• The Eavor Loop Lite demonstration shows production temperatures which
are consistent with analytical model prediction, marked by an excellent
correspondence of prediction and observations up to 480 days of running theloop. This close match also demonstrates that the flow rate and inlet
temperature is sufficiently stable to use an analytical model for performanceassessment.
• For SDE+ applications the analytical approach could be favoured as itprovides a transparent, fast and effective basis for assessing and judging theperformance of the Eavor Loop. Expectation curves for performance caneasily be generated by taking into account uncertainty in relevant
parameters, including thermal conductivity and geothermal gradient.

18 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

8

References
Carslaw, H.S., and Jaeger, J.C.: Conduction of Heat in Solids, Oxford at theClarendon Press, Second Edition (1959), University Press, Oxford, p. 339
Kohl T., Brenni R., and Eugster W.: System performance of a deep borehole heatexchanger, Geothermics, 31, (2002), 687-708
Price, G., Toews, M., Holmes, M., Schwarz, B., Vany, J., 2019. Eavor-Loop™
Thermodynamic Modelling Report. Draft Version received august 28 andcommented draft sep 20.
Price, G., Toews, M., Holmes, M., Schwarz, B., Vany, J., 2020. Eavor-Loop™
Demonstration Project Preliminary Technical Report March 9, 2020
Kutun, K., Tureyen, O.I., Satman, A., 2015. Analysis of Wellhead Production
Temperature Derivatives. PROCEEDINGS, Fortieth Workshop on Geothermal
Reservoir Engineering. Stanford University, Stanford, California, January 26-28,
2015 SGP-TR-204
Kujawa, T., Nowak, W. 2000. Shallow and Deep Vertical Geothermal heatexchangers as low temperature sources for heat pumps. World Geothermal
Congress extended abstract.
Lokhorst, A and J.D van Wees, 2005. Technical and economical feasibility of a deep
borehole heat exchanger (DBHE) in LTG-01. TNO, NITG 05-068-B0510. Nlog.nl
Ramey, H.J.: Wellbore Heat Transmission, Journal of Petroleum Technology, April
1962, 427-435.
Sapinska-Sliwa, A, and others, 2015. Deep Borehole Heat exchangers – aconceptual review. Proceedings World Geothermal Congress 2015, Melbourne,
Australia
Stauffer, P.H., Lewis, K.C., Stein, J.S., Bryan, S. ,Travis, J., Lichtner, P., Zyvoloski,G., 2014. Joule –Thomson Effects on the Flow of Liquid Water. Transp Porous
Med. DOI 10.1007/s11242-014-0379-3
Van Wees, J.D., Lokhorst, A., Zoethout, J., 2007. Re-use of E&P-boreholes forgeothermal energy production. EAGE conference extended abstract.
https://www.tno.nl/media/1662/370beno2.pdf

19 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

Annex 1 : python code analytical solution lateralimport numpy as np
import matplotlib.pyplot as plt
import os.path
x = np.linspace(0, 25, 26)
x[0] =0.1
tstart = 60
tres = 120
theta0 = tstart-tres
Tinlet = x*0 +tstart

secinyear = 365.25*24*3600 # s
krock = 4.0 # W m-1 K-1
rhorock = 2663 # kg m-3
cprock = 1112 #J/kg K-1
radius = 0.078 # m
L = 5000 # m
mrate = 18e3 # kg/h
cpfluid = 4180 # J/kg K-1
wdot=
mrate *cpfluid /3600
radius2 = radius*radius# K formula for C&J
Kfac = 1.4986
overkzK = np.log( Kfac *np.sqrt((krock/(rhorock*cprock))*x*secinyear/
radius2))/(2*np.pi*krock*L)
kzstar = 1/ (overkzK*wdot)
theta = theta0*np.exp(-kzstar)
Toutlet4 = tres + theta

fig, ax = plt.subplots()
dir = "C:/Users/weesjdamv/review_eavor/"
data = np.loadtxt(os.path.join(dir,'figure3_ana2d.txt'))
line1, = plt.plot(data[:,0], data[:,1], color='grey', linewidth=4,
label='Toutlet figure3')
line2, = ax.plot(x, Toutlet4, "ro",markersize=2, label='Toutlet')
plt.xlim(0.1,25)
plt.xlabel('time (years)')
plt.ylabel('Temperature (C)')#plt.semilogx()
plt.grid(True)ax.legend()
plt.show()

20 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

Annex 2 : python code analytical Eavor Lite
import numpy as np
import matplotlib.pyplot as plt
import os.path
x = np.linspace(0, 30, 31)
x[0] =0.1
tsurface = 5.0 # C
tres = 75.95 # 75.95 # C
depth =
2565.0 #2560.0 # 2365.0 # 2570 # m
depthsurface = 610.0 # surface to first section depth
Estar = tres-tsurface
depthratio = depthsurface/depth
Estar1 = Estar * depthratio
Estar2 = Estar - Estar1
tsection = tsurface + depthratio*Estar
secinyear = 365.25*24*3600# mass rate of Eavor lite system
mrate = 40e3 # kg/h
# inlet temperature of Eavor systemtstart = 18

# main systemkrockmain = 2.0 # W m-1 K-1
krockcement = 1.2 # W m-1 K-1
rhorockmain = 2240 #2240.0 # kg m-3
cprockmain = 1520.0 #J/kg K-1
radiussurface = 0.08 # adopt 0.155 to take into account cement
radiussurfaceinner = 0.08 #0.08 # 0.12
radiussurface2 = radiussurface*radiussurface
radiusmain = 0.08 # adopt 0.111
to take into account cement
radiusmain2 = radiusmain*radiusmain
radiusmaininner = 0.08 #0.08 #0.088 # inner radius# lateral properties
krock = 3.5 # W m-1 K-1
rhorock = 2663.0 # kg m-3
cprock = 1110.0 #J/kg K-1
radius = 0.078
L = 1920.0 # m
nlateral = 2 # number of laterals
cpfluid = 4180.0 # J/kg K-1
wdot=
mrate *cpfluid /3600.0
Tinlet = x*0 + tstart
Kfacvertical = 1.4986

21 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)# K factor for mains and lateral
Kfac = Kfacvertical # 1.57 # 1.4986
lsm = np.log(Kfac *np.sqrt((krockmain/(rhorockmain*cprockmain))*x*secinyear/
radiussurface2))/ krockmain
cement = np.log(radiussurface/radiussurfaceinner)/krockcement
overkzKsurface = (cement + lsm) /(2*np.pi*(depthsurface))kzstarsurface = 1/ (overkzKsurface *wdot)

lsm = np.log(Kfac *np.sqrt((krockmain/(rhorockmain*cprockmain))*x*secinyear/
radiusmain2))/ krockmain
cement = np.log(radiusmain/radiusmaininner)/krockcement
overkzKmain = (cement + lsm) /(2*np.pi*(depth-depthsurface))kzstarmain = 1/ (overkzKmain *wdot)
theta0 = Tinlet - tsurface
theta = -Estar1/kzstarsurface + (theta0+Estar1/kzstarsurface) *np.exp(-kzstarsurface)
Tsectionbreak = tsurface + theta+Estar1

theta0 = Tsectionbreak - tsection
theta = -Estar2/kzstarmain + (theta0+Estar2/kzstarmain) * np.exp(kzstarmain)
Tinletlateral = tsection + theta+Estar2
Kfac = 1.4986 # 1.57 # 1.4986
theta0 = Tinletlateral - tresradius2 = radius*radius
overkzK = np.log(Kfac *np.sqrt((krock/(rhorock*cprock))*x*secinyear/
radius2))/(2*np.pi*krock*L)
kzstar = 1/ (overkzK *wdot/nlateral)theta = theta0*np.exp(-kzstar)
Tendlateral = tres + theta
Kfac = Kfacvertical # 1.57 # 1.4986
theta0 = Tendlateral - tres
Estar2 = -Estar2
theta = -Estar2/kzstarmain + (theta0+Estar2/kzstarmain) * np.exp(kzstarmain)
Tsectionbreak2 = tres + theta+Estar2
theta0 = Tsectionbreak2 - tsection
Estar1 = -Estar1
theta = -Estar1/kzstarsurface + (theta0+Estar1/kzstarsurface) *np.exp(-kzstarsurface)
Toutlet = tsection + theta+Estar1
fig, ax = plt.subplots()
dir = "C:/Users/weesjdamv/OneDrive TNO/wees/geothermie/39.TKIwaylandeavor/review_eavor/"
data =
np.loadtxt(os.path.join(dir,'eavorlite_axicalc_20_nocement.txt'))

22 / 23

TNO report TNO 2019 R11549 (update 21-05-2021)

line1, = plt.plot(data[:,0], data[:,1], color='blue', linewidth=2,
label='Tinletlateral EAVOR')
line2, = plt.plot(data[:,0], data[:,2], color='green', linewidth=2,
label='Tendlateral EAVOR')
line3, = plt.plot(data[:,0], data[:,3], color='brown', linewidth=2,
label='Toutlet EAVOR')# Using set_dashes() to modify dashing of an existing line
line1a, = ax.plot(x, Tinletlateral, "bo",markersize=2,
label='Tinletlateral')line2a, = ax.plot(x, Tendlateral, "go",markersize=2,
label='Tendlateral')line3a, = ax.plot(x, Toutlet, "ro",markersize=2, label='Toutlet')#line4, = ax.plot(x, Tsectionbreak2, "bo",markersize=2,
label='Tsectionbreak2')

plt.xlim(0.1,30)
plt.ylim(10, 50)
plt.xlabel('time (years)')
plt.ylabel('Temperature (C)')#plt.semilogx()
plt.grid(True)ax.legend()
plt.show()

23 / 23