Solid infinite elements

Abaqus uses infinite elements to solve boundary value problems defined in unbounded domains or problems in which the region of interest is small in size compared to the surrounding medium.

The following topics are discussed:

Static analysis
Dynamic response

Static analysis

The solution in the far field is assumed to be linear, so only linear behavior is provided in the infinite elements. The static behavior of the infinite elements is based on modeling the basic solution variable, u (in stress analysis u is a displacement component) with respect to spatial distance r measured from a “pole” of the solution, so that $u \to 0$ as $r \to \infty$ , and $u \to \infty$ as $r \to 0$ . The interpolation provides terms of order $1 / r$ , $1 / r^{2}$ , and, when the solution variable is a stress-like variable (such as the pore liquid pressure in the analysis of flow through a porous medium), $1 / r^{3}$ as $r \to 0$ . The far-field behavior of many common cases, such as a point load on a half-space, is thereby included. This modeling is achieved by using standard quadratic or cubic interpolation for $u (s)$ in $- 1 \leq s \leq 1$ , where s is a mapped coordinate that is chosen such that the mapping $r (s)$ causes $r \to \infty$ as $s \to 1$ . We obtain two- and three-dimensional models of domains that reach to infinity by combining this interpolation in the s-direction in a product form with standard linear or quadratic interpolation in orthogonal directions in the mapped space.

In using infinite elements for static analysis, the pole must be located so as to provide a reasonable far-field solution for the particular problem being modeled. The infinite elements in Abaqus are written with nodes on the interface between the finite and infinite elements and, on each edge that stretches to infinity, a node that must be placed in the infinite direction such that the straight line from that node through the corresponding interface node passes through the pole for that ray at a distance on the other side of the interface from the infinite element equal to the distance between these nodes (Figure 1).

Figure 1. Pole node location for an infinite element.

The one-dimensional concept is, thus, based on a node (node 1) on the interface between the finite and infinite elements, distance $r_{1} = a$ from the pole and at $s = - 1$ in the mapped space, and node 2, at $r_{2} = 2 a$ from the pole (the pole is at $r = 0$ ) and at $s = 0$ is the mapped space. The $r (s)$ mapping is chosen as

r = - \frac{2 s}{1 - s} r_{1} + \frac{1 + s}{1 - s} r_{2},

so that

r = \frac{2 a}{1 - s},

which inverts to give

s = 1 - \frac{2 a}{r} .

When an element with $1 / r$ and $1 / r^{2}$ behavior is required, we combine this geometric mapping with standard quadratic interpolation of u with respect to s, written in terms of its values at node 1 and at node 2:

u = \frac{1}{2} s (s - 1) u_{1} + (1 - s^{2}) u_{2}

(this gives $u = 0$ at $s = 1$ , where $r \to \infty$ ). Using the inverted geometric mapping to define $u (r)$ then gives

u = (- u_{1} + 4 u_{2}) \frac{a}{r} + (2 u_{1} - 4 u_{2}) {(\frac{a}{r})}^{2},

which provides the desired behavior. Likewise, when $1 / r^{3}$ behavior is also required, we use cubic interpolation of u with respect to s, written in terms of its values at nodes 1 and 2 and at a third node, which we choose to place at $s = \frac{1}{2}$ :

u = - \frac{1}{3} s (s - 1) (s - \frac{1}{2}) u_{1} - 2 (1 - s^{2}) (s - \frac{1}{2}) u_{2} + \frac{8}{3} s (1 - s^{2}) u_{3} .

The inverted geometric mapping then provides

u = (\frac{1}{3} u_{1} - 4 u_{2} + \frac{32}{3} u_{3}) \frac{a}{r} + (- 2 u_{1} + 20 u_{2} - 32 u_{3}) {(\frac{a}{r})}^{2} + (\frac{8}{3} u_{1} - 16 u_{2} + \frac{64}{3} u_{3}) {(\frac{a}{r})}^{3} .

The infinite elements in Abaqus consist of two- and three-dimensional elements for uncoupled stress analysis that use quadratic interpolation for displacement components and two- and three-dimensional elements for coupled stress-pore liquid pressure elements, in which the displacements use quadratic interpolation and the pore liquid pressure uses cubic interpolation in the infinite direction. This higher-order interpolation is used for the pore liquid pressure for compatibility: since the displacement varies as $1 / r^{2}$ , the strain (and, therefore, the stress) may vary as $1 / r^{3}$ .

Dynamic response

The dynamic response of the infinite elements is based on consideration of plane body waves traveling orthogonally to the boundary. Again, we assume the response adjacent to the boundary is of small enough amplitude so that the medium responds in a linear elastic fashion.

The equilibrium equation is

- ρ \ddot{u} + \frac{\partial}{\partial x} \cdot σ = 0,

where $ρ$ is the material's density, $\ddot{u}$ is the material particle acceleration, $σ$ is the stress, and $x$ is position.

We assume the material's response is isotropic, linear elastic, and—thus—can be written as

σ = λ I I : ε + 2 G ε,

where $ε$ is the strain and

λ = \frac{E ν}{(1 + ν) (1 - 2 ν)}

and

G = \frac{E}{2 (1 + ν)}

are Lamé's constants (E is Young's modulus and $ν$ is Poisson's ratio). Introducing this material response in the equilibrium equation, and assuming small strain:

ε = \frac{1}{2} {\frac{\partial u}{\partial x} + {[\frac{\partial u}{\partial x}]}^{T}},

provides the governing equation for the motion

ρ {\ddot{u}}_{i} = G \frac{\partial^{2} u_{i}}{\partial x_{j} \partial x_{j}} + (λ + G) \frac{\partial^{2} u_{j}}{\partial x_{i} \partial x_{j}},

where index notation has been used for simplicity.

We consider plane waves traveling along the x-axis. Two body wave solutions of this form exist for this equation. One describes plane, longitudinal (“push”) waves, which have the form

u_{x} = f (x \pm c_{p} t), u_{y} = u_{z} = 0,

where, by substitution in the governing equation above, we find that the wave speed, $c_{p}$ , is

c_{p} = \sqrt{\frac{λ + 2 G}{ρ}} .

The other solution of this form is the “shear” wave solution

u_{y} = f (x \pm c_{s} t), u_{x} = u_{z} = 0

u_{z} = f (x \pm c_{s} t), u_{x} = u_{y} = 0,

where—again by substitution in the governing equation—we obtain

c_{s} = \sqrt{\frac{G}{ρ}} .

In each case the solution $f (x - c t)$ represents waves moving in the direction of increasing x, while $f (x + c t)$ represents waves moving in the direction of decreasing x.

Now consider a boundary at $x = L$ of a medium modeled by finite elements in $x < L$ . We introduce distributed damping on this boundary, such that

σ_{x x} = - d_{p} {\dot{u}}_{x}

and

\begin{matrix} σ_{x y} & = - d_{s} {\dot{u}}_{y} \\ σ_{x z} & = - d_{s} {\dot{u}}_{z}, \end{matrix}

where we will now choose the damping constants $d_{p}$ and $d_{s}$ to avoid reflection of longitudinal and shear wave energy back into the medium in $x < L$ . Plane, longitudinal waves approaching the boundary have the form $u_{x} = f_{1} (x - c_{p} t)$ , $u_{y} = u_{z} = 0$ . If they are reflected at all as plane, longitudinal waves, their reflection will travel away from the boundary in some form $u_{x} = f_{2} (x + c_{p} t)$ , $u_{y} = u_{z} = 0$ . Since the problem is linear, superposition provides the total displacement $f_{1} + f_{2}$ , with corresponding stresses $σ_{x x} = (λ + 2 G) (f_{1}^{'} + f_{2}^{'})$ , all other $σ_{i j} = 0$ , and velocity ${\dot{u}}_{x} = - c_{p} (f_{1}^{'} - f_{2}^{'})$ . For this solution to satisfy the damping behavior introduced on the boundary at $x = L$ requires

(λ + 2 G - d_{p} c_{p}) f_{1}^{'} + (λ + 2 G + d_{p} c_{p}) f_{2}^{'} = 0.

We can, therefore, ensure that $f_{2} = 0$ (so that $f_{2}^{'} = 0$ ) for any $f_{1}$ by choosing

d_{p} = \frac{λ + 2 G}{c_{p}} = ρ c_{p} .

A similar argument for shear waves provides

d_{s} = ρ c_{s} .

These values of boundary damping are built into the infinite elements in Abaqus. From the above discussion we see that they transmit all normally impinging plane body waves exactly (provided that the material behavior close to the boundary is linear elastic). General problems involve nonplane body waves that do not impinge on the boundary from an orthogonal direction and may also involve Rayleigh surface waves and Love waves. Nevertheless, these “quiet” boundaries work quite well even for such general cases, provided that they are arranged so that the dominant direction of wave propagation is orthogonal to the boundary or, at free surfaces and interfaces where Rayleigh or Love waves are of concern, they are orthogonal to the surface (see, for example, Cohen and Jennings, 1983). As the boundaries are “quiet” rather than silent (perfect transmitters of all waveforms), and because the boundaries rely on the solution adjacent to them being linear elastic, they should be placed some reasonable distance from the region of main interest.

During dynamic response analysis following static preload (as is common in geotechnical applications), the traction provided by the infinite elements to the boundary of the finite element mesh consists of the constant stress obtained from the static response with the quiet boundary damping stress added. Since the elements have no stiffness during dynamic analysis, they allow a net rigid body motion to occur, which is usually not a significant effect.