Cased ThermoElastic Wellbore Problem
Problem description
This example uses the thermal option of the SinglePhasePoromechanics
solver to handle a cased wellbore problem subject to a uniform temperature change on the inner surface of the casing. The wellbore is composed of a steel casing, a cement sheath and rock formation. Isotropic linear thermoelastic behavior is assumed for all three materials. No separation or thermal barrier is allowed for the casing-cement and cement-rock contact interfaces. Plane strain condition is assumed.
Solution to this axisymmetric problem can be obtained in the cylindrical coordinate system by using an implicit 1D finite difference method (Jane and Lee 1999). Results of such analysis will be considered as reference solutions to validate GEOS results.
Input file
This benchmark example uses no external input files and everything required is contained within two GEOS XML files located at:
inputFiles/wellbore/CasedThermoElasticWellbore_base.xml
and
inputFiles/wellbore/CasedThermoElasticWellbore_benchmark.xml
The corresponding integrated test is
inputFiles/wellbore/CasedThermoElasticWellbore_smoke.xml
Geometry and mesh
The internal wellbore mesh generator InternalWellbore
is employed to create the mesh of this wellbore problem. The radii of the casing cylinder, the cement sheath cylinder and the far-field boundary of the surrounding rock formation are defined by a vector radius
. In the tangent direction, theta
angle is specified from 0 to 90 degrees to simulate the problem on a quarter of the wellbore geometry. The problem is under plane strain condition and therefore we only consider radial thermal diffusion on a single horizontal layer. The trajectory of the well is defined by trajectory
, which is vertical in this case. The autoSpaceRadialElems
parameters allow for optimally increasing the element size from the wellbore to the far-field zone. In this example, the auto spacing option is only applied to the rock formation. The useCartesianOuterBoundary
with a value 3 specified for the rock layer transforms the far-field boundary to a circular shape. The cellBlockNames
and elementTypes
define the regions and related element types associated to casing, cement sheath, and rock.
<Mesh>
<InternalWellbore
name="mesh1"
elementTypes="{ C3D8, C3D8, C3D8 }"
radius="{ 0.15707, 0.17780, 0.21272, 1.5707 }"
theta="{ 0, 90 }"
zCoords="{ 0, 0.1 }"
nr="{ 5, 5, 5 }"
nt="{ 10 }"
nz="{ 1 }"
trajectory="{ { 0.0, 0.0, 0.0 },
{ 0.0, 0.0, 0.1 } }"
autoSpaceRadialElems="{ 0, 0, 1 }"
cellBlockNames="{ casing, cement, rock }"
/>
</Mesh>
Material properties
The bulk and shear drained elastic moduli of the materials as well as its drained linear thermal expansion coefficient relating stress change to temperature change are defined within the Constitutive
tag as follows:
<ElasticIsotropic
name="casingSolid"
defaultDensity="7500"
defaultBulkModulus="159.4202899e9"
defaultShearModulus="86.61417323e9"
defaultDrainedLinearTEC="1.2e-5"/>
<ElasticIsotropic
name="cementSolid"
defaultDensity="2700"
defaultBulkModulus="2.298850575e9"
defaultShearModulus="1.652892562e9"
defaultDrainedLinearTEC="2.0e-5"/>
<ElasticIsotropic
name="rockSolid"
defaultDensity="2700"
defaultBulkModulus="5.535714286e9"
defaultShearModulus="3.81147541e9"
defaultDrainedLinearTEC="2.0e-5"/>
Here the solid density is also defined but it is not used because the gravitational effect is ignored in this example. To mimic a thermoelastic coupling without fluid flow, a negligible porosity and a zero Biot coefficient are defined as:
<BiotPorosity
name="casingPorosity"
defaultReferencePorosity="1e-6"
defaultGrainBulkModulus="159.4202899e9"/>
<BiotPorosity
name="cementPorosity"
defaultReferencePorosity="1e-6"
defaultGrainBulkModulus="2.298850575e9"/>
<BiotPorosity
name="rockPorosity"
defaultReferencePorosity="1e-6"
defaultGrainBulkModulus="5.535714286e9"/>
In this XML block, the Biot coefficient is defined using the elastic bulk modulus of the solid skeleton as . In this example, we define a skeleton bulk modulus that is identical to the drained bulk modulus defined above to enforce the Biot coefficient to zero.
The thermal conductivities and the volumetric heat capacities of casing, cement, and rock are defined by following XML blocks:
<SinglePhaseThermalConductivity
name="casingThermalCond"
defaultThermalConductivityComponents="{ 15, 15, 15 }"/>
<SinglePhaseThermalConductivity
name="cementThermalCond"
defaultThermalConductivityComponents="{ 1.0, 1.0, 1.0 }"/>
<SinglePhaseThermalConductivity
name="rockThermalCond"
defaultThermalConductivityComponents="{ 1.66, 1.66, 1.66 }"/>
and
<SolidInternalEnergy
name="casingInternalEnergy"
referenceVolumetricHeatCapacity="1.375e6"
referenceTemperature="0"
referenceInternalEnergy="0"/>
<SolidInternalEnergy
name="cementInternalEnergy"
referenceVolumetricHeatCapacity="4.2e6"
referenceTemperature="0"
referenceInternalEnergy="0"/>
<SolidInternalEnergy
name="rockInternalEnergy"
referenceVolumetricHeatCapacity="4.56e6"
referenceTemperature="0"
referenceInternalEnergy="0"/>
An ultra-low permeability is defined for the three layers to simulate a thermoelastic problem without the impact of fluid flow.
<ConstantPermeability
name="casingPerm"
permeabilityComponents="{ 1.0e-100, 1.0e-100, 1.0e-100 }"/>
<ConstantPermeability
name="cementPerm"
permeabilityComponents="{ 1.0e-100, 1.0e-100, 1.0e-100 }"/>
<ConstantPermeability
name="rockPerm"
permeabilityComponents="{ 1.0e-100, 1.0e-100, 1.0e-100 }"/>
Also, a negligible volumetric heat capacity is defined for the fluid to completely ignore the thermal convection effect such that only thermal transfers via the diffusion phenomenon are considered.
<ThermalCompressibleSinglePhaseFluid
name="fluid"
defaultDensity="1000"
defaultViscosity="1e-3"
referencePressure="0.0"
referenceTemperature="20.0"
compressibility="5e-10"
thermalExpansionCoeff="1e-10"
viscosibility="0.0"
specificHeatCapacity="1"
referenceInternalEnergy="1"/>
Other fluid properties such as viscosity, thermal expansion coefficient, etc. are not relevant to this example because fluid flow is ignored and pore pressure is zero everywhere.
Boundary conditions
The mechanical boundary conditions are applied to ensure the axisymmetric plane strain conditions such as:
<FieldSpecification
name="tNegConstraint"
objectPath="nodeManager"
fieldName="totalDisplacement"
component="1"
scale="0.0"
setNames="{ tneg }"/>
<FieldSpecification
name="tPosConstraint"
objectPath="nodeManager"
fieldName="totalDisplacement"
component="0"
scale="0.0"
setNames="{ tpos }"/>
<FieldSpecification
name="zconstraint"
objectPath="nodeManager"
fieldName="totalDisplacement"
component="2"
scale="0.0"
setNames="{ zneg, zpos }"/>
Besides, the far-field boundary is assumed to be fixed because the local changes on the wellbore must have negligible effect on the far-field boundary.
<FieldSpecification
name="rPosConstraint_x"
objectPath="nodeManager"
fieldName="totalDisplacement"
component="0"
scale="0.0"
setNames="{ rpos }"/>
<FieldSpecification
name="rPosConstraint_y"
objectPath="nodeManager"
fieldName="totalDisplacement"
component="1"
scale="0.0"
setNames="{ rpos }"/>
The traction-free condition on the inner surface of the casing is defined by:
<Traction
name="innerPressure"
objectPath="faceManager"
tractionType="normal"
scale="0.0e6"
setNames="{ rneg }"/>
The initial reservoir temperature (that is also the far-field boundary temperature) and the temperature of a cold fluid applied on the inner surface of the casing are defined as
<FieldSpecification
name="initialTemperature"
initialCondition="1"
setNames="{ all }"
objectPath="ElementRegions"
fieldName="temperature"
scale="100"/>
<FieldSpecification
name="farfieldTemperature"
setNames="{ rpos }"
objectPath="faceManager"
fieldName="temperature"
scale="100"/>
<FieldSpecification
name="innerTemperature"
setNames="{ rneg }"
objectPath="faceManager"
fieldName="temperature"
scale="-20.0"/>
It is important to remark that the initial effective stress of each layers must be set with accordance to the initial temperature: where is the initial effective principal stress, is the initial temperature change, is the drained bulk modulus and is the drained linear thermal expansion coefficient of the materials.
<FieldSpecification
name="initialSigma_x_casing"
initialCondition="1"
setNames="{ all }"
objectPath="ElementRegions/casing/casing"
fieldName="casingSolid_stress"
component="0"
scale="573913043.5"/>
<FieldSpecification
name="initialSigma_y_casing"
initialCondition="1"
setNames="{ all }"
objectPath="ElementRegions/casing/casing"
fieldName="casingSolid_stress"
component="1"
scale="573913043.5"/>
<FieldSpecification
name="initialSigma_z_casing"
initialCondition="1"
setNames="{ all }"
objectPath="ElementRegions/casing/casing"
fieldName="casingSolid_stress"
component="2"
scale="573913043.5"/>
Zero pore pressure is set everywhere to simulate a thermoelastic problem in which fluid flow is ignored:
<FieldSpecification
name="zeroPressure"
setNames="{ all }"
objectPath="ElementRegions"
fieldName="pressure"
scale="0e6"/>
<FieldSpecification
name="sourcePressure"
setNames="{ rneg }"
objectPath="faceManager"
fieldName="pressure"
scale="0"/>
<FieldSpecification
name="sinkPressure"
setNames="{ rpos }"
objectPath="faceManager"
fieldName="pressure"
scale="0"/>
Collecting output data
It is convenient to collect data in hdf5 format that can be easily post-processed using Python. To collect the temperature field in the three layers for all the time steps, the following XML blocks need to be defined:
<PackCollection
name="temperatureCollection_casing"
objectPath="ElementRegions/casing/casing"
fieldName="temperature"/>
<PackCollection
name="temperatureCollection_cement"
objectPath="ElementRegions/cement/cement"
fieldName="temperature"/>
<PackCollection
name="temperatureCollection_rock"
objectPath="ElementRegions/rock/rock"
fieldName="temperature"/>
<TimeHistory
name="temperatureHistoryOutput_casing"
sources="{ /Tasks/temperatureCollection_casing }"
filename="temperatureHistory_casing"/>
<TimeHistory
name="temperatureHistoryOutput_cement"
sources="{ /Tasks/temperatureCollection_cement }"
filename="temperatureHistory_cement"/>
<TimeHistory
name="temperatureHistoryOutput_rock"
sources="{ /Tasks/temperatureCollection_rock }"
filename="temperatureHistory_rock"/>
Similarly, the following blocks are needed to collect the solid stress:
<PackCollection
name="stressCollection_casing"
objectPath="ElementRegions/casing/casing"
fieldName="casingSolid_stress"/>
<PackCollection
name="stressCollection_cement"
objectPath="ElementRegions/cement/cement"
fieldName="cementSolid_stress"/>
<PackCollection
name="stressCollection_rock"
objectPath="ElementRegions/rock/rock"
fieldName="rockSolid_stress"/>
<TimeHistory
name="stressHistoryOutput_casing"
sources="{ /Tasks/stressCollection_casing }"
filename="stressHistory_casing"/>
<TimeHistory
name="stressHistoryOutput_cement"
sources="{ /Tasks/stressCollection_cement }"
filename="stressHistory_cement"/>
<TimeHistory
name="stressHistoryOutput_rock"
sources="{ /Tasks/stressCollection_rock }"
filename="stressHistory_rock"/>
The displacement field can be collected for the whole domain using nodeManager
as follows
<PackCollection
name="displacementCollection"
objectPath="nodeManager"
fieldName="totalDisplacement"/>
<TimeHistory
name="displacementHistoryOutput"
sources="{ /Tasks/displacementCollection }"
filename="displacementHistory"/>
Also, periodic events are required to trigger the collection of this data on the mesh. For example, the periodic events for collecting the displacement field are defined as:
<PeriodicEvent
name="displacementHistoryCollection"
endTime="1e5"
forceDt="1e4"
target="/Tasks/displacementCollection"/>
<PeriodicEvent
name="displacementTimeHistoryOutput"
endTime="1e5"
forceDt="1e4"
target="/Outputs/displacementHistoryOutput"/>
Results and benchmark
A good agreement between the GEOS results and analytical results for temperature distribution around the cased wellbore is shown in the figures below:
and the validation for the radial displacement around the cased wellbore is shown below:
The validations of the total radial and hoop stress (tangent stress) components computed by GEOS against reference results are shown in the figure below:
To go further
Feedback on this example
This concludes the cased wellbore example. For any feedback on this example, please submit a GitHub issue on the project’s GitHub page.