Skip to content

Principle of Virtual Work

Based on the equilibrium equation and boundary conditions derived in Stress and Conservation Laws, we derive the principle of virtual work, which is the weak form of the continuum-mechanics boundary-value problem. Finite element discretization starts from this weak form. This chapter presents both the current-configuration form (expressed using Cauchy stress and the linear part of Almansi strain) and the reference-configuration form (expressed using the second Piola-Kirchhoff stress and Green-Lagrange strain), demonstrates their equivalence, and then confirms the reduction to small deformation.

Equilibrium Equation and Boundary Conditions

Let \(\boldsymbol{g}\) be the body force acting on the continuum per unit mass, and consider a body occupying the domain \(\Omega\) in the current configuration. The boundary \(\Gamma\) is partitioned into the geometric boundary on which displacement is prescribed as \(\bar{\boldsymbol{u}}\), \(\Gamma_B\), and the mechanical boundary on which traction is prescribed as \(\bar{\boldsymbol{t}}\), \(\Gamma_t\), with \(\Gamma = \Gamma_B \cup \Gamma_t\) and \(\Gamma_B \cap \Gamma_t = \emptyset\). For a static problem, the equilibrium equation is obtained by omitting the inertial term from the momentum conservation law shown in Stress and Conservation Laws:

\[ \nabla_x \cdot \boldsymbol{\sigma} + \rho \boldsymbol{g} = \boldsymbol{0} \quad \text{in} \ \Omega \]

The boundary conditions are

\[ \boldsymbol{\sigma} \boldsymbol{n} = \bar{\boldsymbol{t}} \quad \text{on} \ \Gamma_t \]
\[ \boldsymbol{u} = \bar{\boldsymbol{u}} \quad \text{on} \ \Gamma_B \]

Hereafter, the principle of virtual work is derived as the weak form of the equilibrium equation and the mechanical boundary condition \(\boldsymbol{\sigma} \boldsymbol{n} = \bar{\boldsymbol{t}}\). The geometric boundary condition \(\boldsymbol{u} = \bar{\boldsymbol{u}}\) is incorporated through the choice of test function.

Weak Form in the Current Configuration

In the weak form, define the admissible space for the unknown displacement and the test-function space, respectively, as

\[ \mathcal{U} = \{ \boldsymbol{u} \in [H^1(\Omega)]^d \mid \boldsymbol{u} = \bar{\boldsymbol{u}} \ \text{on} \ \Gamma_B \} \]
\[ \mathcal{V} = \{ \delta \boldsymbol{u} \in [H^1(\Omega)]^d \mid \delta \boldsymbol{u} = \boldsymbol{0} \ \text{on} \ \Gamma_B \} \]

Here, \(d\) is the spatial dimension, \(H^1(\Omega)\) is the Sobolev space of functions square-integrable together with their first weak derivatives, and \(\delta\) denotes variation. In the current-configuration representation, \(\Omega\) is the deformed configuration; in an actual numerical solution, it is pulled back to the reference configuration or to a known intermediate configuration.

Multiplying the equilibrium equation by the weight \(\delta \boldsymbol{u} \in \mathcal{V}\) and applying the Gauss divergence theorem and the mechanical boundary condition gives the following principle of virtual work in the current configuration.

\[ \int_{\Omega} \boldsymbol{\sigma} : \delta \boldsymbol{A}_{(L)}\, dv = \int_{\Gamma_t} \delta \boldsymbol{u}^T \bar{\boldsymbol{t}}\, d\Gamma + \int_{\Omega} \delta \boldsymbol{u}^T \rho \boldsymbol{g}\, dv \]

Here, \(\boldsymbol{A}_{(L)}\) is the linear part of the Almansi strain tensor, defined by

\[ \boldsymbol{A}_{(L)} = \frac{1}{2}\left( \nabla_x \boldsymbol{u} + (\nabla_x \boldsymbol{u})^T \right), \qquad A_{(L)ij} = \frac{1}{2}\left( \frac{\partial u_i}{\partial x_j} + \frac{\partial u_j}{\partial x_i} \right) \]

Its variation is \(\delta \boldsymbol{A}_{(L)} = \tfrac{1}{2}(\nabla_x \delta \boldsymbol{u} + (\nabla_x \delta \boldsymbol{u})^T)\). Thus, we seek \(\boldsymbol{u} \in \mathcal{U}\) such that the virtual-work equation is satisfied for every \(\delta \boldsymbol{u} \in \mathcal{V}\). The left-hand side is the virtual work of internal forces, while the right-hand side is the virtual work of external forces due to prescribed traction and body force.

Because this equation is written over the deformed (current-configuration) domain, an actual solution procedure reselects the initial configuration \(\Omega_0\) (reference configuration) or a known intermediate configuration as the reference configuration, rewrites the equation in incremental form, and then solves it. See Framework for Incremental Analysis for the specific choice of reference configuration (Total Lagrange / Updated Lagrange) and the incremental decomposition.

Weak Form in the Initial Configuration

Consider a body occupying the domain \(\Omega_0\) in the reference configuration, with boundary \(\Gamma_0\) partitioned into \(\Gamma_{0B} \cup \Gamma_{0t}\). Pulling the current-configuration representation back to the reference configuration gives the conjugate stress-strain pair of the second Piola-Kirchhoff stress \(\boldsymbol{S}\) and Green-Lagrange strain \(\boldsymbol{E}\). The principle of virtual work in the initial configuration is then

\[ \int_{\Omega_0} \boldsymbol{S} : \delta \boldsymbol{E}\, dV = \int_{\Gamma_{0t}} \delta \boldsymbol{u}^T \bar{\boldsymbol{t}}\, d\Gamma_0 + \int_{\Omega_0} \delta \boldsymbol{u}^T \rho_0 \boldsymbol{g}\, dV \]

Here, \(\rho_0\) is the mass density in the reference configuration, and by the mass-conservation relation \(\rho_0 = J\rho\), this is equivalent to the body-force representation in the current configuration.

Equivalence of the Current- and Initial-Configuration Representations

The virtual work of internal forces in the two representations coincides through the transformation by the deformation gradient \(\boldsymbol{F}\) and the volume ratio \(J = \det \boldsymbol{F}\); namely,

\[ \int_{\Omega_0} \boldsymbol{S} : \delta \boldsymbol{E}\, dV = \int_{\Omega} \boldsymbol{\sigma} : \delta \boldsymbol{A}_{(L)}\, dv \]

The external-force terms are also equivalent by mass conservation and traction transformation. Therefore, the virtual-work equation in the current configuration and that in the initial configuration express the same principle in different configurations. A solution method referring to the reference configuration corresponds to the Total Lagrange method, whereas a solution method referring to the current configuration (the immediately preceding converged configuration) corresponds to the Updated Lagrange method.

Reduction to Small Deformation

Under the small-deformation assumptions \(\boldsymbol{F} \approx \boldsymbol{I}\) and \(J \approx 1\), the distinction between the current and reference configurations disappears, the second PK stress coincides with the Cauchy stress (\(\boldsymbol{S} \to \boldsymbol{\sigma}\)), and both the Green-Lagrange strain and the linear part of the Almansi strain reduce to the infinitesimal strain \(\boldsymbol{\varepsilon}\).

\[ \boldsymbol{\varepsilon} = \nabla_S \boldsymbol{u} = \frac{1}{2}\left( \nabla \boldsymbol{u} + (\nabla \boldsymbol{u})^T \right), \qquad \varepsilon_{ij} = \frac{1}{2}\left( \frac{\partial u_i}{\partial x_j} + \frac{\partial u_j}{\partial x_i} \right) \]

The principle of virtual work then reduces to the weak form expressed in terms of the Cauchy stress \(\boldsymbol{\sigma}\) and infinitesimal strain \(\boldsymbol{\varepsilon}\):

\[ \int_{\Omega} \boldsymbol{\sigma} : \delta \boldsymbol{\varepsilon}\, dV = \int_{\Gamma_t} \delta \boldsymbol{u}^T \bar{\boldsymbol{t}}\, d\Gamma + \int_{\Omega} \delta \boldsymbol{u}^T \rho \boldsymbol{g}\, dV \]
\[ \delta \boldsymbol{u} = \boldsymbol{0} \quad \text{on} \ \Gamma_B \]

This is the weak form used directly for discretization in small-deformation linear-elastic static analysis (Linear Elastic Static Analysis (Introduction and Appendix) starts from this form and shows the construction of the element stiffness \(\boldsymbol{K}^e\) through assembly of the global equation \(\boldsymbol{K}\boldsymbol{U} = \boldsymbol{F}\)).

Substituting the linear-elastic constitutive law \(\boldsymbol{\sigma} = \boldsymbol{\mathsf{C}} : \boldsymbol{\varepsilon}\) and writing \(\hat{\sigma} = D\, \hat{\varepsilon}\) in Voigt notation, the weak form becomes

\[ \int_{\Omega} \delta \hat{\varepsilon}^T D\, \hat{\varepsilon}\, dV = \int_{\Gamma_t} \delta \boldsymbol{u}^T \bar{\boldsymbol{t}}\, d\Gamma + \int_{\Omega} \delta \boldsymbol{u}^T \rho \boldsymbol{g}\, dV \]

which has this form.