ThermoPoroElastic Wellbore Problem

Problem description

This example uses the thermal option of the SinglePhasePoromechanics solver to handle an open wellbore problem subjected to a uniform temperature change on its inner surface. Isotropic linear thermoporoelastic behavior is considered for the rock formation around the wellbore. Plane strain and axisymmetric conditions are assumed.

../../../../../../_images/sketch9.png — Fig. 61 Sketch of a thermoporoelastic wellbore

Analytical solutions to this problem were first derived by (Wang and Papamichos 1994) using a one-way coupling simplification. They are also reformulated for the full coupling assumption in the book of (Cheng 2016). These solutions will be considered to validate GEOS results.

Input file

This benchmark example uses no external input files and everything required is contained within two GEOS xml files that are located at:

inputFiles/wellbore/ThermoPoroElasticWellbore_base.xml

and

inputFiles/wellbore/ThermoPoroElasticWellbore_benchmark.xml

The corresponding integrated test is

inputFiles/wellbore/ThermoPoroElasticWellbore_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 open wellbore 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 thermal diffusion along the radial direction within 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 near wellbore zone to the far-field one.

  <Mesh>
    <InternalWellbore
      name="mesh1"
      elementTypes="{ C3D8 }"
      radius="{ 0.1, 5.0 }"
      theta="{ 0, 90 }"
      zCoords="{ 0, 0.1 }"
      nr="{ 100 }"
      nt="{ 40 }"
      nz="{ 1 }"
      trajectory="{ { 0.0, 0.0, 0.0 },
                    { 0.0, 0.0, 0.1 } }"
      autoSpaceRadialElems="{ 1 }"
      cellBlockNames="{ rock }"
      />
  </Mesh>

../../../../../../_images/mesh15.png — Fig. 62 An optimized mesh for the open wellbore.

Material properties

The bulk and shear drained elastic moduli of rock as well as its drained linear thermal expansion coefficient relating stress change to temperature variation are defined within the Constitutive tag as follows:

    <ElasticIsotropic
      name="rockSolid"
      defaultDensity="2700"
      defaultBulkModulus="20.7e9"
      defaultShearModulus="12.4e9"
      defaultDrainedLinearTEC="4e-5"/>

Here the solid density is also defined, but it is not used as the gravitational effect is ignored in this example. The porosity and the elastic bulk modulus $K_{s}$ of the solid skeleton are defined as:

    <BiotPorosity
      name="rockPorosity"
      defaultReferencePorosity="0.001"
      defaultGrainBulkModulus="23.5e9"
      defaultPorosityTEC="4e-5"/>

The thermal conductivities and the volumetric heat capacities of rock are defined by following XML blocks:

    <SinglePhaseThermalConductivity
      name="rockThermalCond"
      defaultThermalConductivityComponents="{ 6.6, 6.6, 6.6 }"/>

and

    <SolidInternalEnergy
      name="rockInternalEnergy"
      referenceVolumetricHeatCapacity="1.89e6"
      referenceTemperature="0"
      referenceInternalEnergy="0"/>

The permeability of rock is defined by:

    <ConstantPermeability
      name="rockPerm"
      permeabilityComponents="{ 1.0e-21, 1.0e-21, 1.0e-21 }"/>

Fluid properties such as viscosity, thermal expansion coefficient, etc. are defined by the XML block below. A negligible volumetric heat capacity is defined for fluid to ignore the thermal convection effect. This way, only thermal transfer via the diffusion phenomenon is considered.

    <ThermalCompressibleSinglePhaseFluid
      name="fluid"
      defaultDensity="1000"
      defaultViscosity="1e-3"
      referencePressure="0.0"
      referenceTemperature="20.0"
      compressibility="5e-10"
      thermalExpansionCoeff="3e-4"
      viscosibility="0.0"
      specificHeatCapacity="1"
      referenceInternalEnergy="1"/>

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 wellbore is defined by:

    <Traction
      name="innerTraction"
      objectPath="faceManager"
      tractionType="normal"
      scale="0.0e6"
      setNames="{ rneg }"/>

The initial temperature (that is also the far-field boundary temperature) and the temperature applied on the inner surface of the wellbore are defined as

    <FieldSpecification
      name="initialTemperature"
      initialCondition="1"
      setNames="{ all }"
      objectPath="ElementRegions"
      fieldName="temperature"
      scale="0"/>

    <FieldSpecification
      name="farfieldTemperature"
      setNames="{ rpos }"       
      objectPath="faceManager"
      fieldName="temperature"
      scale="0"/>

    <FieldSpecification
      name="innerTemperature"
      setNames="{ rneg }"     
      objectPath="faceManager"
      fieldName="temperature"
      scale="100.0"/>

It is important to remark that the initial effective stress of rock must be set with accordance to the initial temperature change: $\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. In this example, the initial effective stresses are set to zero because the initial temperature change is set to zero.

    <FieldSpecification
      name="initialSigma_x_rock"
      initialCondition="1"
      setNames="{ all }"
      objectPath="ElementRegions/rock/rock"
      fieldName="rockSolid_stress"
      component="0"
      scale="0"/>
    <FieldSpecification
      name="initialSigma_y_rock"
      initialCondition="1"
      setNames="{ all }"
      objectPath="ElementRegions/rock/rock"
      fieldName="rockSolid_stress"
      component="1"
      scale="0"/>
    <FieldSpecification
      name="initialSigma_z_rock"
      initialCondition="1"
      setNames="{ all }"
      objectPath="ElementRegions/rock/rock"
      fieldName="rockSolid_stress"
      component="2"
      scale="0"/>

The initial and boundary conditions for pore pressure are defined in the block below:

    <FieldSpecification
      name="initialPressure"
      initialCondition="1"
      setNames="{ all }"
      objectPath="ElementRegions"
      fieldName="pressure"
      scale="0e6"/>

    <FieldSpecification
      name="innerPressure"
      setNames="{ rneg }"     
      objectPath="faceManager"
      fieldName="pressure"
      scale="0e6"/>

    <FieldSpecification
      name="farfieldPressure"
      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 for all the time steps, the following XML blocks need to be defined:

    <PackCollection
      name="temperatureCollection_rock"
      objectPath="ElementRegions/rock/rock"
      fieldName="temperature"/>

    <TimeHistory
      name="temperatureHistoryOutput_rock"
      sources="{ /Tasks/temperatureCollection_rock }"
      filename="temperatureHistory_rock"/>

Similarly, the following blocks are needed to collect the effective stress field across the domain:

    <PackCollection
      name="stressCollection_rock"
      objectPath="ElementRegions/rock/rock"
      fieldName="rockSolid_stress"/>

    <TimeHistory
      name="stressHistoryOutput_rock"
      sources="{ /Tasks/stressCollection_rock }"
      filename="stressHistory_rock"/>

The displacement field can be collected 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 during the entire simulation. For example, the periodic events for collecting the displacement field are defined as:

    <PeriodicEvent
      name="displacementHistoryCollection"
      beginTime="0"
      endTime="360"
      forceDt="60"
      target="/Tasks/displacementCollection"/>
    <PeriodicEvent
      name="displacementTimeHistoryOutput_1"
      beginTime="0"
      endTime="360"
      forceDt="60"
      target="/Outputs/displacementHistoryOutput"/>
    <PeriodicEvent
      name="displacementHistoryCollection_2"
      beginTime="360"
      endTime="3700"
      forceDt="360"
      target="/Tasks/displacementCollection"/>
    <PeriodicEvent
      name="displacementTimeHistoryOutput_2"
      beginTime="360"
      endTime="3700"
      forceDt="360"
      target="/Outputs/displacementHistoryOutput"/>

Results and benchmark

A good agreement between the GEOS results and analytical results for temperature and pore pressure distribution around the wellbore is shown in the figures below:

../../../../../../_images/temperature_pressure.png — Fig. 63 Validation of temperature and pore pressure.

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

../../../../../../_images/displacement2.png — Fig. 64 Validation of the radial 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:

../../../../../../_images/stress2.png — Fig. 65 Validation of the radial and tangent 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.