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.

../../../../../../_images/sketch8.png

Fig. 56 Sketch of a cased thermoelastic wellbore

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>
../../../../../../_images/mesh11.png

Fig. 57 An optimized mesh for the cased wellbore.

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"
      grainBulkModulus="159.4202899e9"/>

    <BiotPorosity
      name="cementPorosity"
      defaultReferencePorosity="1e-6"
      grainBulkModulus="2.298850575e9"/>

    <BiotPorosity
      name="rockPorosity"
      defaultReferencePorosity="1e-6"
      grainBulkModulus="5.535714286e9"/>

In this XML block, the Biot coefficient is defined using the elastic bulk modulus K_{s} of the solid skeleton as b_{Biot} = 1 - K/K_{s}. In this example, we define a skeleton bulk modulus that is identical to the drained bulk modulus K 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:

    <SinglePhaseConstantThermalConductivity
      name="casingThermalCond"
      thermalConductivityComponents="{ 15, 15, 15 }"/>

    <SinglePhaseConstantThermalConductivity
      name="cementThermalCond"
      thermalConductivityComponents="{ 1.0, 1.0, 1.0 }"/>

    <SinglePhaseConstantThermalConductivity
      name="rockThermalCond"
      thermalConductivityComponents="{ 1.66, 1.66, 1.66 }"/>

and

    <SolidInternalEnergy
      name="casingInternalEnergy"
      volumetricHeatCapacity="1.375e6"
      referenceTemperature="0"
      referenceInternalEnergy="0"/>

    <SolidInternalEnergy
      name="cementInternalEnergy"
      volumetricHeatCapacity="4.2e6"
      referenceTemperature="0"
      referenceInternalEnergy="0"/>

    <SolidInternalEnergy
      name="rockInternalEnergy"
      volumetricHeatCapacity="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: \sigma_{0} = 3K\alpha \delta T_{0} where \sigma_{0} is the initial effective principal stress, \delta T_{0} is the initial temperature change, K is the drained bulk modulus and \alpha is the drained linear thermal expansion coefficient of the materials.

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:

(Source code)

../../../../../../_images/temperature1.png

Fig. 58 Validation of the temperature.

and the validation for the radial displacement around the cased wellbore is shown below:

(Source code)

../../../../../../_images/displacement1.png

Fig. 59 Validation of the displacement.

The validations of the total radial and hoop stress (tangent stress) components computed by GEOS against reference results are shown in the figure below:

(Source code)

../../../../../../_images/stress1.png

Fig. 60 Validation of the stresses.

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.