Introduction
Large deformation refers to a state in which a material or structure undergoes changes in shape or size that are too significant to be ignored, due to external loads, temperature variations, or internal stresses.
In such cases, displacement gradients and rotations are no longer small, and linear approximations cease to be valid. Therefore, the traditional small strain assumption cannot accurately describe the behavior, and geometric nonlinearity must be taken into account.
Theory
Deformation Gradient
F (Deformation Gradient) represents how the body is deformed. It contains all information: stretching, compression, rotation, and shear.
C (Right Cauchy–Green Tensor) is obtained from F after removing the rotational part. It represents only the change of shape, such as changes in lengths and angles.
\mathbf{F} = \mathbf{I} + \nabla \mathbf{u},
\quad
\mathbf{C} = \mathbf{F}^\mathsf{T} \mathbf{F},
\quad
J = \sqrt{\det \mathbf{C}}Green–Lagrange Strain
E (Green–Lagrange Strain) is a strain measure derived from C. It correctly describes strain even under large deformation or rotation. For small deformations, it becomes the same as the usual linear strain.
\mathbf{E}
= \frac{1}{2}\bigl(\mathbf{C} - \mathbf{I}\bigr)
Variation of the Green–Lagrange Strain
δE (Variation of the Green–Lagrange Strain) represents how the strain E changes when the displacement field is perturbed slightly. This small change is introduced by the virtual displacement δu.
\delta \mathbf{E}(\delta\mathbf{u},\mathbf{u})
= \frac{1}{2}\Bigl(
(\nabla \delta\mathbf{u})\,\mathbf{F}
+ \bigl((\nabla \delta\mathbf{u})\,\mathbf{F}\bigr)^\mathsf{T}
\Bigr)
Second Piola–Kirchhoff Stress
S(u) (Second Piola–Kirchhoff Stress) is a stress measure defined in the reference configuration. It describes how much the material is stretched or compressed, based on the strain E. The parameters μ and λ are material constants (Lamé parameters). The formula uses the volume change J and the inverse tensor C^{-1} to represent nonlinear stress behavior.
\mathbf{S}(\mathbf{u})
= \mu\bigl(\mathbf{I} - \mathbf{C}^{-1}\bigr)
+ \lambda \,\ln J\, \mathbf{C}^{-1}Internal Virtual Work
In nonlinear continuum mechanics, the internal virtual work is given by the variation of the strain energy density:
\delta W = \int_{\Omega_0} \mathbf{S} : \delta \mathbf{E} dVBut, to represent this virtual work as weak from, test function should be introduced to this equation. Therefore, additional stress measure below is used.
Stress Measure Used in the Weak Form
Although the weak form in the previous section is written in terms of the Second Piola–Kirchhoff stress S\mathbf{S} and the variation of the Green–Lagrange strain δE\delta \mathbf{E}, it is important to note that this representation is primarily energy-based and not the most natural form for finite element implementation. In the weak form, which is actual implementation of the code, First Piola–Kirchhoff Stress ( \mathbf{P(u)} )is used.
\mathbf{P}(u) = \mathbf{F}\mathbf{S}(u)Using this First Piola–Kirchhoff stress, the weak form of the equilibrium equations can be written in a more standard and implementation-friendly form:
R(\mathbf{u}; \mathbf{v})
=
\int_{\Omega_0}
\mathbf{P}(\mathbf{u}) : \nabla \mathbf{v}
dV - l_f\int_{\Gamma_N}
\mathbf{t} \cdot \mathbf{v}
d\GammaThe l_f here is scaling factor that is used as technique to make equation converged.
First and Second PIola–Kirchhoff
Second Piola–Kirchhoff stress \mathbf{S} is primarily used as an intermediate stress measure for constitutive modeling, whereas the First Piola-Kirchhoff stress \mathbf{P} directly enters the weak form and numerical implementation.
Although both formulations are mathematically equivalent, the formulation based on \mathbf{P} aligns more naturally with standard finite element discretizations.
Practical Implementation
In general, textbooks state that nonlinear structural problems are solved using the Newton–Raphson method, and the explanation usually ends there. However, in practice this is often not sufficient. In fact, even a seemingly simple cantilever example failed to converge.
Why does this happen?
It is likely because nonlinear equations frequently fall into local solutions or otherwise fail to converge depending on the initial guess and loading conditions.
To address this, I adopted a method in which the Neumann boundary load is increased gradually. This is so-called load-scaling. Although this is a fairly heuristic approach, it worked quite effectively. There are certainly more sophisticated methods available, arc-length for example, but in terms of ease of implementation, this strategy offers significant value.
Result
The result of cantilever is shown below, numerically converged enough as well.
