6.7.1 Coupled pore fluid diffusion and stress analysis

**Products: **ABAQUS/Standard ABAQUS/CAE

A coupled pore fluid diffusion/stress analysis:

is used to model single phase, partially or fully saturated fluid flow through porous media;

can be performed in terms of either total pore pressure or excess pore pressure by including or excluding the pore fluid weight;

requires the use of pore pressure elements with associated pore fluid flow properties defined;

can be transient or steady-state;

can be linear or nonlinear; and

can include pore pressure contact between bodies (see “Pore fluid contact properties,” Section 30.4.1).

Typical applications

Some of the more common coupled pore fluid diffusion/stress analysis problems that can be analyzed with ABAQUS/Standard are:

**Saturated flow**: Soil mechanics problems generally involve fully saturated flow, since the solid is fully saturated with ground water. Typical examples of saturated flow include consolidation of soils under foundations and excavation of tunnels in saturated soil.

**Partially saturated flow**: Partially saturated flow occurs when the wetting liquid is absorbed into or exsorbed from the medium by capillary action. Irrigation and hydrology problems generally include partially saturated flow.

**Combined flow**: Combined fully saturated and partially saturated flow occurs in problems such as seepage of water through an earth dam, where the position of the phreatic surface (the boundary between fully saturated and partially saturated soil) is of interest.

**Moisture migration**: Although not normally associated with soil mechanics, moisture migration problems can also be solved using the coupled pore fluid diffusion/stress procedure. These problems may involve partially saturated flow in polymeric materials such as paper towels and sponge-like materials; in the biomedical industry they may also involve saturated flow in hydrated soft tissues.

Flow through porous media

A porous medium is modeled in ABAQUS/Standard by a conventional approach that considers the medium as a multiphase material and adopts an effective stress principle to describe its behavior. The porous medium modeling provided considers the presence of two fluids in the medium. One is the “wetting liquid,” which is assumed to be relatively (but not entirely) incompressible. Often the other is a gas, which is relatively compressible. An example of such a system is soil containing ground water. When the medium is partially saturated, both fluids exist at a point; when it is fully saturated, the voids are completely filled with the wetting liquid. The elementary volume, , is made up of a volume of grains of solid material, ; a volume of voids, ; and a volume of wetting liquid, , that is free to move through the medium if driven. In some systems (for example, systems containing particles that absorb the wetting liquid and swell in the process) there may also be a significant volume of trapped wetting liquid, .

The porous medium is modeled by attaching the finite element mesh to the solid phase; fluid can flow through this mesh. The mechanical part of the model is based on the effective stress principle defined in “Effective stress principle for porous media,” Section 2.8.1 of the ABAQUS Theory Manual.

The model also uses a continuity equation for the mass of wetting fluid in a unit volume of the medium. This equation is described in “Continuity statement for the wetting liquid phase in a porous medium,” Section 2.8.4 of the ABAQUS Theory Manual. It is written with pore pressure (the average pressure in the wetting fluid at a point in the porous medium) as the basic variable (degree of freedom 8 at the nodes). The conjugate flux variable is the volumetric flow rate at the node, . The porous medium is partially saturated when the pore liquid pressure, , is negative.

Total and excess pore fluid pressure

The coupled pore fluid diffusion/stress analysis capability can provide solutions either in terms of total or “excess” pore fluid pressure. The excess pore fluid pressure at a point is the pore fluid pressure in excess of the hydrostatic pressure required to support the weight of pore fluid above the elevation of the material point. The difference between total and excess pore pressure is relevant only for cases in which gravitational loading is important; for example, when the loading provided by the hydrostatic pressure in the pore fluid is large or when effects like “wicking” (transient capillary suction of liquid into a dry column) are being studied. Total pore pressure solutions are provided when the gravity distributed load is used to define the gravity load on the model. Excess pore pressure solutions are provided in all other cases; for example, when gravity loading is defined with body force distributed loads.

Steady-state analysis

Steady-state coupled pore pressure/effective stress analysis assumes that there are no transient effects in the wetting liquid continuity equation; that is, the steady-state solution corresponds to constant wetting liquid velocities and constant volume of wetting liquid per unit volume in the continuum. Thus, for example, thermal expansion of the liquid phase has no effect on the steady-state solution: it is a transient effect. Therefore, the time scale chosen during steady-state analysis is relevant only to rate effects in the constitutive model used for the porous medium (excluding creep and viscoelasticity, which are disabled in steady-state analysis).

Mechanical loads and boundary conditions can be changed gradually over the step by referring to an amplitude curve to accommodate possible geometric nonlinearities in the response.

The steady-state coupled equations are strongly unsymmetric; therefore, the unsymmetric matrix solution and storage scheme is used automatically for steady-state analysis steps (see “Procedures: overview,” Section 6.1.1).

Input File Usage: | *SOILS |

ABAQUS/CAE Usage: | Step module: |

You can specify a fixed time increment size in a coupled pore fluid diffusion/stress analysis, or ABAQUS/Standard can select the time increment size automatically. Automatic incrementation is recommended because the time increments in a typical diffusion analysis can increase by several orders of magnitude during the simulation. If you do not activate automatic incrementation, fixed time increments will be used.

Input File Usage: | Use the following option to activate automatic incrementation in steady-state analysis: |

*SOILS, UTOL= The solution does not depend on the value specified for UTOL; this value is simply a flag for automatic incrementation. |

ABAQUS/CAE Usage: | Step module: |

Transient analysis

In a transient coupled pore pressure/effective stress analysis the backward difference operator is used to integrate the continuity equation: this operator provides unconditional stability so that the only concern with respect to time integration is accuracy. You can provide the time increments, or they can be selected automatically.

The coupled partially saturated flow equations are strongly unsymmetric, so the unsymmetric solver is used automatically if you request partially saturated analysis (by including absorption/exsorption behavior in the material definition). The unsymmetric solver is also activated automatically when gravity distributed loading is used during a soils consolidation analysis.

For fully saturated flow analyses in which finite-sliding coupled pore pressure-displacement contact is modeled using contact pairs, certain contributions to the model's stiffness matrix are unsymmetric. Using the unsymmetric solver can sometimes improve convergence in such cases since ABAQUS does not automatically do so.

The integration procedure used in ABAQUS/Standard for consolidation analysis introduces a relationship between the minimum usable time increment and the element size, as shown below for fully saturated and partially saturated flows. If time increments smaller than these values are used, spurious oscillations may appear in the solution (except for partially saturated cases when linear elements or modified triangular elements are used; in these cases ABAQUS/Standard uses a special integration scheme for the wetting liquid storage term to avoid the problem). These nonphysical oscillations may cause problems if pressure-sensitive plasticity is used to model the porous medium and may lead to convergence difficulties in partially saturated analyses. If the problem requires analysis with smaller time increments than the relationships given below allow, a finer mesh is required. Generally there is no upper limit on the time step except accuracy, since the integration procedure is unconditionally stable unless nonlinearities cause convergence problems.

A simple guideline that can be used for the minimum usable time increment in the case of fully saturated flow is

where

is the time increment,

is the specific weight of the wetting liquid,

*E*

is the Young's modulus of the soil,

is the permeability of the soil (see “Permeability,” Section 20.7.2),

is the magnitude of the velocity of the pore fluid,

is the velocity coefficient in Forchheimer's flow law ( in the case of Darcy flow),

is the bulk modulus of the solid grains (see “Porous bulk moduli,” Section 20.7.3), and

is a typical element dimension.

In partially saturated flow cases the corresponding guideline for the minimum time increment is

where

*s*

is the saturation;

is the permeability-saturation relationship;

is the rate of change of saturation with respect to pore pressure (see “Sorption,” Section 20.7.4);

is the initial porosity of the material; and the other parameters are as defined for the case of fully saturated flow.

If you choose fixed time incrementation, fixed time increments equal to the size of the user-specified initial time increment, , will be used. Fixed incrementation is not generally recommended because the time increments in a typical diffusion analysis can increase over several orders of magnitude during the simulation; automatic incrementation is usually a better choice.

Input File Usage: | *SOILS, CONSOLIDATION |

ABAQUS/CAE Usage: | Step module: |

If you choose automatic time incrementation, you must specify two tolerance parameters.

The accuracy of the time integration of the flow continuity equations is governed by the maximum wetting liquid pore pressure change, , allowed in an increment. ABAQUS/Standard will restrict the time step to ensure that this value will not be exceeded at any node (except nodes with boundary conditions) during any increment in the analysis.

The accuracy of the integration of the time-dependent (creep) material behavior is governed by the maximum strain rate change allowed at any point during an increment, , as described in “Rate-dependent plasticity: creep and swelling,” Section 18.2.4.

Input File Usage: | *SOILS, CONSOLIDATION, UTOL=, CETOL= |

ABAQUS/CAE Usage: | Step module: |

Transient soils analysis can be terminated by completing a specified time period, or it can be continued until steady-state conditions are reached. By default, the analysis will end when the given time period has been completed. Alternatively, you can specify that the analysis will end when steady state is reached or the time period ends, whichever comes first. Steady state is defined by a maximum permitted rate of change of pore pressure with time: when all pore pressures are changing at less than the user-defined rate, the analysis terminates.

ABAQUS/CAE Usage: | Step module: |

You can specify that creep or viscoelastic response should be neglected during a consolidation analysis, even if creep or viscoelastic material properties have been defined.

Input File Usage: | *SOILS, CONSOLIDATION, CREEP=NONE |

ABAQUS/CAE Usage: | Step module: |

Some types of analyses may develop local instabilities, such as surface wrinkling, material instability, or local buckling. In such cases it may not be possible to obtain a quasi-static solution, even with the aid of automatic incrementation. ABAQUS/Standard offers the option to stabilize this class of problems by applying damping throughout the model in such a way that the viscous forces introduced are sufficiently large to prevent instantaneous buckling or collapse, but small enough not to affect the behavior significantly while the problem is stable. ABAQUS/Standard generates an artificial damping matrix by using the mass matrix with a unit density together with a mass-proportional damping factor, as described in “Solving nonlinear problems,” Section 7.1.1. Whenever possible, the damping factor is chosen such that, based on extrapolation of the results obtained during the first increment, the dissipated energy during the step is a small fraction of the change in strain energy during the step. This dissipated energy fraction is controlled by the user and has a default value of 2.0 × 10^{–4}. If the problem is either unstable or contains rigid body motions during the first increment, an alternative method is used to determine the damping factor; this method is based on making an averaged damping stiffness equal to the dissipated energy fraction times an averaged material stiffness.

Input File Usage: | Use the following option to activate automatic stabilization with the default dissipated energy fraction: |

*SOILS, CONSOLIDATION, STABILIZE Use the following option to specify a nondefault dissipated energy fraction: *SOILS, CONSOLIDATION, STABILIZE= Use the following option to specify the damping factor directly: *SOILS, CONSOLIDATION, STABILIZE, FACTOR= |

ABAQUS/CAE Usage: | Step module: |

Units

In coupled problems where two different fields are being solved, you must be careful when choosing the units of the problem. If the choice of units is such that the numbers generated by the equations for the two different fields differ by many orders of magnitude, the precision on some computers may be insufficient to resolve the numerical ill-conditioning of the coupled equations. Therefore, choose units that avoid badly conditioned matrices. For example, consider using units of MPascal instead of Pascal for the stress equilibrium equations to reduce the disparity between the magnitudes of the stress equilibrium equations and the pore flow continuity equations.

Initial conditions can be applied as defined in “Initial conditions,” Section 27.2.1.

Initial values of pore fluid pressures, , can be defined at the nodes.

Input File Usage: | *INITIAL CONDITIONS, TYPE=PORE PRESSURE |

ABAQUS/CAE Usage: | Initial pore fluid pressure is not supported in ABAQUS/CAE. |

Initial values of the void ratio, *e*, can be given at the nodes. The void ratio is defined as the ratio of the volume of voids to the volume of solid material (see “Effective stress principle for porous media,” Section 2.8.1 of the ABAQUS Theory Manual). The evolution of void ratio is governed by the deformation of the different phases of the material, as discussed in detail in “Constitutive behavior in a porous medium,” Section 2.8.3 of the ABAQUS Theory Manual.

Input File Usage: | *INITIAL CONDITIONS, TYPE=RATIO |

ABAQUS/CAE Usage: | Initial void ratio is not supported in ABAQUS/CAE. |

Initial saturation values, *s*, can be given at the nodes. Saturation is defined as the ratio of wetting fluid volume to void volume—see “Effective stress principle for porous media,” Section 2.8.1 of the ABAQUS Theory Manual.

Input File Usage: | *INITIAL CONDITIONS, TYPE=SATURATION |

ABAQUS/CAE Usage: | Initial saturation is not supported in ABAQUS/CAE. |

An initial (effective) stress field can be specified (see “Initial conditions,” Section 27.2.1).

Most geotechnical problems begin from a geostatic state, which is a steady-state equilibrium configuration of the undisturbed soil or rock body under geostatic loading and usually includes both horizontal and vertical components. It is important to establish these initial conditions correctly so that the problem begins from an equilibrium state. The geostatic procedure can be used to verify that the user-defined initial stresses are indeed in equilibrium with the given geostatic loads and boundary conditions (see “Geostatic stress state,” Section 6.7.2).

Input File Usage: | Use one of the following options: |

*INITIAL CONDITIONS, TYPE=STRESS *INITIAL CONDITIONS, TYPE=STRESS, GEOSTATIC |

ABAQUS/CAE Usage: | Initial stress is not supported in ABAQUS/CAE. |

Boundary conditions can be applied to displacement degrees of freedom 1–6 and to pore pressure degree of freedom 8 (“Boundary conditions,” Section 27.3.1). During the analysis prescribed boundary conditions can be varied by referring to an amplitude curve (“Amplitude curves,” Section 27.1.2). If no amplitude reference is given, the default variation of a boundary condition in a coupled pore fluid diffusion/stress analysis step is as defined in “Procedures: overview,” Section 6.1.1.

If the pore pressure is prescribed with a boundary condition, fluid is assumed to enter and leave through the node as needed to maintain the prescribed pressure.

The following loading types can be prescribed in a coupled pore fluid diffusion/stress analysis:

Concentrated nodal forces can be applied to the displacement degrees of freedom (1–6); see “Concentrated loads,” Section 27.4.2.

Distributed pressure forces or body forces can be applied; see “Distributed loads,” Section 27.4.3. The distributed load types available with particular elements are described in Part VI, “Elements.” The magnitude and direction of gravitational loading are usually defined by using the gravity distributed load type.

Pore fluid flow is controlled as described in “Pore fluid flow,” Section 27.4.6.

The following predefined fields can be prescribed, as described in “Predefined fields,” Section 27.6.1:

Although temperature is not a degree of freedom in coupled pore fluid diffusion/stress analysis, nodal temperatures can be specified. Any difference between the applied and initial temperatures will cause thermal strain if a thermal expansion coefficient is given for the material (“Thermal expansion,” Section 20.1.2). The specified temperature also affects temperature-dependent material properties, if any.

The values of user-defined field variables can be specified; these values affect only field-variable-dependent material properties, if any.

Any of the mechanical constitutive models available in ABAQUS/Standard can be used to model the porous material.

In problems formulated in terms of total pore pressure, you must include the density of the dry material in the material definition (see “Density,” Section 16.2.1).

You can use a permeability material property to define the specific weight of the wetting liquid, ; the permeability, , and its dependence on the void ratio, *e*, and saturation, ; and the flow velocity, (see “Permeability,” Section 20.7.2).

You can define the compressibility of the solid grains and of the permeating fluid in both fully and partially saturated flow problems (see “Elastic behavior of porous materials,” Section 17.3.1). If you do not specify the porous bulk moduli, the materials are assumed to be fully incompressible.

For partially saturated flow you must define the porous medium's absorption/exsorption behavior (see “Sorption,” Section 20.7.4).

Gel swelling (“Swelling gel,” Section 20.7.5) and volumetric moisture swelling of the solid skeleton (“Moisture swelling,” Section 20.7.6) can be included in partially saturated cases. These effects are usually associated with modeling of moisture migration in polymeric systems rather than with geotechnical systems.

Thermal expansion can be defined separately for the solid material and for the permeating fluid. In such a case you should repeat the expansion material property, with the necessary parameters, to define the different thermal expansion effects (see “Thermal expansion,” Section 20.1.2). Thermal expansion will be active only in a consolidation (transient) analysis.

Input File Usage: | To define the thermal expansion of the permeating fluid: |

*EXPANSION, TYPE=ISO, PORE FLUID To define the thermal expansion of the solid material: *EXPANSION, TYPE=ISO or ORTHO or ANISO |

ABAQUS/CAE Usage: | To define the thermal expansion of the permeating fluid: |

Property module: material editor: To define the thermal expansion of the solid material: Property module: material editor: |

The analysis of flow through porous media in ABAQUS/Standard is available for plane strain, axisymmetric, and three-dimensional problems. Continuum pore pressure elements are provided for modeling fluid flow through a deforming porous medium in a coupled pore fluid diffusion/stress analysis. These elements have pore pressure degree of freedom 8 in addition to displacement degrees of freedom 1–3. Stress/displacement elements can be used in parts of the model without pore fluid flow. See “Choosing the appropriate element for an analysis type,” Section 21.1.3, for more information.

The element output available for a coupled pore fluid diffusion/stress analysis includes the usual mechanical quantities such as (effective) stress; strain; energies; and the values of state, field, and user-defined variables. In addition, the following quantities associated with pore fluid flow are available:

VOIDR | Void ratio, |

POR | Pore pressure, . |

SAT | Saturation, |

GELVR | Gel volume ratio, . |

FLUVR | Total fluid volume ratio, . |

FLVEL | Magnitude and components of the pore fluid effective velocity vector, . |

FLVELM | Magnitude, , of the pore fluid effective velocity vector. |

FLVEL n | Component |

The nodal output available includes the usual mechanical quantities such as displacements, reaction forces, and coordinates. In addition, the following quantities associated with pore fluid flow are available:

POR | Pore pressure at a node. |

RVF | Reaction fluid volume flux due to prescribed pressure. This flux is the rate at which fluid volume is entering or leaving the model through the node to maintain the prescribed pressure boundary condition. A positive value of RVF indicates that fluid is entering the model. |

RVT | Reaction total fluid volume (computed only in a transient analysis). This value is the time integrated value of RVF. |

All of the output variable identifiers are outlined in “ABAQUS/Standard output variable identifiers,” Section 4.2.1.

*HEADING … *********************************** ** ** Material definition ** *********************************** *MATERIAL, NAME=soilData lines to define mechanical properties of the solid material… *EXPANSIONData lines to define the thermal expansion coefficient of the solid grains*EXPANSION, TYPE=ISO, PORE FLUIDData lines to define the thermal expansion coefficient of the permeating fluid*PERMEABILITY, SPECIFIC=Data lines to define permeability, , as a function of the void ratio,*PERMEABILITY, TYPE=SATURATIONeData lines to define the dependence of permeability on saturation,*PERMEABILITY, TYPE=VELOCITYData lines to define the velocity coefficient,*POROUS BULK MODULIData line to define the bulk moduli of the solid grains and the permeating fluid*SORPTION, TYPE=ABSORPTIONData lines to define absorption behavior*SORPTION, TYPE=EXSORPTIONData lines to define exsorption behavior*SORPTION, TYPE=SCANNINGData lines to define scanning behavior (between absorption and exsorption)*GELData line to define gel behavior in partially saturated flow*MOISTURE SWELLINGData lines to define moisture swelling strain as a function of saturation in partially saturated flow… *********************************** ** ** Boundary conditions and initial conditions ** *********************************** *BOUNDARYData lines to specify zero-valued boundary conditions*INITIAL CONDITIONS, TYPE=STRESS, GEOSTATICData lines to specify initial stresses*INITIAL CONDITIONS, TYPE=PORE PRESSUREData lines to define initial values of pore fluid pressures*INITIAL CONDITIONS, TYPE=RATIOData lines to define initial values of the void ratio*INITIAL CONDITIONS, TYPE=SATURATIONData lines to define initial saturation*AMPLITUDE, NAME=nameData lines to define amplitude variations*********************************** ** ** Step 1: Optional step to ensure an equilibrium ** geostatic stress field ** *********************************** *STEP *GEOSTATIC *CLOAD and/or *DLOAD and/or *TEMPERATURE and/or *FIELDData lines to specify mechanical loading*FLOW and/or *SFLOW and/or *DFLOW and/or *DSFLOWData lines to specify pore fluid flow*BOUNDARYData lines to specify displacements or pore pressures*END STEP *********************************** ** ** Step 2: Coupled pore diffusion/stress analysis step ** *********************************** *STEP (,NLGEOM) ** Use NLGEOM to include geometric nonlinearities *SOILSData line to define incrementation*CLOAD and/or *DLOAD and/or *DSLOADData lines to specify mechanical loading*FLOW and/or *SFLOW and/or *DFLOW and/or *DSFLOWData lines to specify pore fluid flow*BOUNDARYData lines to specify displacements or pore pressures*END STEP