Updated Lagrangian formulation: plane stress state, isoparametric 4-node quadrilateral, Hencky strains

01

Formulation and hypotheses

The isoparametric 4-node quadrilateral element is described in Section 9.4.

In an Updated Lagrangian formulation, the displacements and stresses at load step $i+1$ are evaluated with respect to the known deformed configuration at the preceding load step $i$. That known configuration becomes the new reference configuration. Consequently, the derivatives and integrals used for load step $i+1$ are evaluated with respect to the configuration at load step $i$; see also Chapter 15.

Cauchy stresses are used together with either the Hencky or the Biot strain tensor. Both strain measures are implemented in the MATLAB program and compared in the numerical example. For small strains, the Biot strain tensor represents the first-order approximation of the Hencky strain tensor.

02

Strain energy and notation

Strain energy at load step $i+1$, written with respect to the configuration at load step $i$ (see also Sections 15.1 and 15.2):

Working in the natural coordinates $s,t$, where $\boldsymbol{\varepsilon}$, $\boldsymbol{\sigma}$, and $\boldsymbol{J}$ are functions of $s,t$:

Here, $h$ is the initial element thickness. The current thickness is approximated by accounting for transverse contraction as in Section 15.1:

This thickness relation is a small-strain approximation; it is not a general exact relation for finite elastic strains.

The strains and stresses at load step $i+1$ are:

This additive strain update is an approximation consistent with the present small-strain incremental formulation; Hencky and Biot strain tensors are not generally additive under finite successive deformations.

For a linear elastic material in plane stress:

The finite-element area and Jacobian are evaluated in the known configuration at load step $i$. The Jacobian is:

The coordinates $x_1,\ldots,x_4$ and $y_1,\ldots,y_4$ are the nodal coordinates of the current finite element at load step $i$. After load step $i+1$ is completed, the nodal displacements, strains, and corresponding Cauchy stresses at that load step are known and define the reference state for the next increment.

03

Virtual work and strain increment

The virtual-work equation is evaluated using a $2\times2$ Gauss quadrature:

Here, $n_G=2$ is the number of Gauss points in each natural-coordinate direction. Thus the Q4 element uses $2\times2=4$ Gauss points in total.

For any geometrically admissible virtual displacement, the strain variation is written as:

The element nodal-displacement vector is:

The Hencky strain tensor for the increment is:

For small strains, $\boldsymbol{C}$ remains close to the identity tensor. Hence, to first order:

This approximation is used in the following derivation. Since $\boldsymbol{C}=\boldsymbol{F}^{T}\boldsymbol{F}$ and $\boldsymbol{C}^{-1}=\boldsymbol{U}^{-2}$:

The deformation gradient of the increment, evaluated with respect to the configuration at load step $i$, is:

04

Construction of the strain-displacement matrix

The displacement derivatives follow from the same Q4 interpolation used in Section 9.4:

Here:

The increment of the deformation gradient is produced by the eight element nodal-displacement increments:

For each node $a=1,2,3,4$:

Consequently:

where the $2\times2$ matrices are:

Define the $2\times16$ matrix:

The engineering strain-displacement matrix has three rows and eight columns:

The matrix $\boldsymbol{B}$ is evaluated for every finite element and at every Gauss point in B_matrix.m, which is called by stiff.m.

05

Internal force and approximate tangent

The element internal nodal force vector is:

To keep the program simple, the tangent stiffness used in this section is approximated by:

This expression neglects the dependence of $\boldsymbol{B}$ on the nodal displacements. It is therefore an approximate tangent stiffness matrix. The equilibrium state is still governed by the element internal force vector, but the approximate tangent gives poorer iterative convergence.

For this reason, the underrelaxation factor used in main.m is retained:

%*** main ***
dS=K\F*sr;
S=S-dS;

A value of $s_r=0.6$ provides satisfactory convergence for this example.

An improved approximation of the tangent stiffness matrix is presented in Section 17.3.

06

Hencky and Biot implementation

Both Hencky and Biot strain tensors are implemented in the same MATLAB program.

%*** stiff ***
istrain=1;   % 1 - Hencky strain
             % 2 - Biot strain

At the beginning of the stiff subroutine, the variable istrain selects the strain measure used in the analysis: istrain=1 for Hencky strain and istrain=2 for Biot strain.

if istrain==1
    ee=logm(U);
    Uic=Ui^2;       % Hencky strain tensor
elseif istrain==2
    ee=U-eye(2);
    Uic=Ui;         % Biot strain tensor
else
    error('Unknown strain tensor option')
end

For small strains, the Biot strain tensor may be regarded as the first-order approximation of the Hencky strain tensor:

07

Numerical example

Example 1. Consider a cantilever beam of length $L=400\ \mathrm{mm}$ and rectangular cross-section $5\ \mathrm{mm}\times20\ \mathrm{mm}$. The material properties are $E=1000\ \mathrm{MPa}$ and $\nu=0.3$. A force $F=100\ \mathrm{N}$ is applied at the free end. The mesh contains $60\times8=480$ Q4 elements and 549 nodes. The load is applied in 20 steps.

The displacement and normal-stress results are compared below. The Hencky and Biot values are reported in that order.

The Hencky- and Biot-based results are extremely close for this example, as expected from their first-order relation at small strain. The comparison also shows the corresponding Euler–Almansi and ANSYS results without changing the stress measure: all reported stresses are Cauchy stresses.

08

Programs and downloads

The MATLAB package used for the example and the accompanying ANSYS macro can be downloaded below.