Table of contents
The partial differential equation for the static linear elasticity problem
− ∇ ⋅ σ = f in Ω σ = 2 μ ε + λ I ∇ ⋅ u in Ω u = g on ∂ Ω D σ ⋅ n = t on ∂ Ω N
\def\arraystretch{1.5}
\begin{array}{lc}
-\nabla \cdot \bm \sigma =f & \text{in }\Omega \\
\bm \sigma = 2\mu \bm \varepsilon + \lambda \bm I \nabla \cdot \bm u & \text{in }\Omega \\
\bm u = \bm g & \text{on }\partial \Omega_D \\
\bm \sigma \cdot \bm n = \bm t & \text{on }\partial \Omega_N
\end{array}
− ∇ ⋅ σ = f σ = 2 μ ε + λ I ∇ ⋅ u u = g σ ⋅ n = t in Ω in Ω on ∂ Ω D on ∂ Ω N We cannot directly tackle the strong form of a differential equation to create the FEM. The first step in the FEM is to derive the weak form (or variational form) of the strong form.
Following the recipe for deriving the weak form we multiply the PDE by v v v v = 0 v=0 v = 0 ∂ Ω D \partial \Omega_D ∂ Ω D Ω \Omega Ω
∫ Ω ( − ∇ ⋅ σ ) ⋅ v d Ω = ∫ Ω f ⋅ v d Ω \int_{\Omega } \left(-\nabla \cdot \bm \sigma \right)\cdot \bm v d\Omega =\int_{\Omega } \bm f \cdot \bm v d\Omega ∫ Ω ( − ∇ ⋅ σ ) ⋅ v d Ω = ∫ Ω f ⋅ v d Ω
1: The 2D domain.
In Figure 1 , ∂ Ω D \partial \Omega_D ∂ Ω D Dirichlet boundary conditions (also known as the essential boundary conditions), where the solution u \bm u u ∂ Ω N \partial \Omega_N ∂ Ω N Neumann boundary conditions (also known as the natural boundary conditions), where the force is applied.
Recall that (in 2D)
σ = [ σ x τ x y τ x y σ y ] \bm \sigma =\left\lbrack \begin{array}{cc}
\sigma_x & \tau_{xy} \\
\tau_{xy} & \sigma_y
\end{array}\right\rbrack σ = [ σ x τ x y τ x y σ y ] so the divergence of the stress is
∇ ⋅ σ = [ ∂ σ x ∂ x + ∂ τ x y ∂ y ∂ τ x y ∂ x + ∂ σ y ∂ y ]
\def\arraystretch{2.0}
\nabla \cdot \bm \sigma =\left\lbrack \begin{array}{c}
\dfrac{\partial \sigma_x }{\partial x}+\dfrac{\partial \tau_{xy} }{\partial y}\\
\dfrac{\partial \tau_{xy} }{\partial x}+\dfrac{\partial \sigma_y }{\partial y}
\end{array}\right\rbrack
∇ ⋅ σ = ⎣ ⎡ ∂ x ∂ σ x + ∂ y ∂ τ x y ∂ x ∂ τ x y + ∂ y ∂ σ y ⎦ ⎤ Let
q 1 = [ σ x τ x y ] and q 2 = [ τ x y σ y ] \bm q_1 =\left\lbrack \begin{array}{c}
\sigma_x \\
\tau_{xy}
\end{array}\right\rbrack \text{ and } \bm q_2 =\left\lbrack \begin{array}{c}
\tau_{xy} \\
\sigma_y
\end{array}\right\rbrack q 1 = [ σ x τ x y ] and q 2 = [ τ x y σ y ] where q 1 \bm q_1 q 1 x x x q 2 \bm q_2 q 2 y y y
we can write
( ∇ ⋅ σ ) ⋅ v = ∑ i , j ∂ σ i j ∂ x j v j = ∇ ⋅ q 1 v x + ∇ ⋅ q 2 v y \left(\nabla \cdot \bm \sigma \right)\cdot \bm v=\sum_{i,j} \dfrac{\partial \sigma_{ij} }{\partial x_j }v_j =\nabla \cdot \bm q_1 v_x +\nabla \cdot \bm q_2 v_y ( ∇ ⋅ σ ) ⋅ v = i , j ∑ ∂ x j ∂ σ ij v j = ∇ ⋅ q 1 v x + ∇ ⋅ q 2 v y then we can formulate the weak form as such:
∫ Ω ( − ∇ ⋅ σ ) ⋅ v d Ω = − ∫ Ω ( ∇ ⋅ q 1 v x + ∇ ⋅ q 2 v y ) d Ω \int_{\Omega } \left(-\nabla \cdot \bm \sigma \right)\cdot \bm v d\Omega = -\int_{\Omega } \left(\nabla \cdot q_1 v_x +\nabla \cdot q_2 v_y \right)d\Omega ∫ Ω ( − ∇ ⋅ σ ) ⋅ v d Ω = − ∫ Ω ( ∇ ⋅ q 1 v x + ∇ ⋅ q 2 v y ) d Ω We can expand ∇ ⋅ ( v q ) \nabla \cdot \left(vq\right) ∇ ⋅ ( v q )
∇ ⋅ ( v q ) = ∂ ∂ x ( v q x ) + ∂ ∂ y ( v q y ) = product rule ∂ q x ∂ x v + q x ∂ v ∂ x + ∂ y ∂ y v + q y ∂ v ∂ y \nabla \cdot \left(vq\right)=\dfrac{\partial }{\partial x}\left(vq_x \right)+\dfrac{\partial }{\partial y}\left(vq_y \right)\overset{\text{product rule} }{=} \dfrac{\partial q_x }{\partial x}v+q_x \dfrac{\partial v}{\partial x}+\dfrac{\partial y}{\partial y}v+q_y \dfrac{\partial v}{\partial y} ∇ ⋅ ( v q ) = ∂ x ∂ ( v q x ) + ∂ y ∂ ( v q y ) = product rule ∂ x ∂ q x v + q x ∂ x ∂ v + ∂ y ∂ y v + q y ∂ y ∂ v ∇ ⋅ ( v q ) = q ⋅ ∇ v + v ∇ ⋅ q \nabla \cdot \left(vq\right)=q\cdot \nabla v+v\nabla \cdot q ∇ ⋅ ( v q ) = q ⋅ ∇ v + v ∇ ⋅ q So
∫ Ω ∇ ⋅ ( v q ) d Ω = ∫ Ω q ⋅ ∇ v d Ω + ∫ Ω v ∇ ⋅ q d Ω \int_{\Omega } \nabla \cdot \left(vq\right)d\Omega =\int_{\Omega } q\cdot \nabla vd\Omega +\int_{\Omega } v\nabla \cdot qd\Omega ∫ Ω ∇ ⋅ ( v q ) d Ω = ∫ Ω q ⋅ ∇ v d Ω + ∫ Ω v ∇ ⋅ q d Ω Recall the Gauss divergence theorem
∫ Ω ∇ ⋅ ( v q ) d Ω = ∫ ∂ Ω N v q ⋅ n d s \int_{\Omega } \nabla \cdot \left(vq\right)d\Omega =\int_{\partial \Omega_N } vq\cdot nds ∫ Ω ∇ ⋅ ( v q ) d Ω = ∫ ∂ Ω N v q ⋅ n d s which can be seen as the multivariate integration by parts formula. The theorem also describes conservation of mass in classical mechanics. What ever mass is generated (or lost) inside the body needs to be transferred across its boundary.
We then have
∫ Ω ∇ ⋅ ( v q ) d Ω = ∫ Ω q ⋅ ∇ v d Ω + ∫ Ω v ∇ ⋅ q ⏟ t w o d i f f e r e n t i a l s d Ω = ∫ ∂ Ω N v q ⋅ n d s \int_{\Omega } \nabla \cdot \left(vq\right)d\Omega =\int_{\Omega } q\cdot \nabla vd\Omega +\int_{\Omega } v\underset{\mathrm{two}\mathrm{differentials} }{\underbrace{\nabla \cdot q} } d\Omega =\int_{\partial \Omega_N } vq\cdot nds ∫ Ω ∇ ⋅ ( v q ) d Ω = ∫ Ω q ⋅ ∇ v d Ω + ∫ Ω v two differentials ∇ ⋅ q d Ω = ∫ ∂ Ω N v q ⋅ n d s We move the term with the higher order differential to the left hand side
∫ Ω v ∇ ⋅ q d Ω = ∫ ∂ Ω N v q ⋅ n d s − ∫ Ω q ⋅ ∇ v d Ω \int_{\Omega } v\nabla \cdot qd\Omega =\int_{\partial \Omega_N } vq\cdot nds-\int_{\Omega } q\cdot \nabla vd\Omega ∫ Ω v ∇ ⋅ q d Ω = ∫ ∂ Ω N v q ⋅ n d s − ∫ Ω q ⋅ ∇ v d Ω So we bring it all together
− ∫ Ω v x ∇ ⋅ q 1 + v y ∇ ⋅ q 2 d Ω = ∫ Ω q 1 ⋅ ∇ v x + q 2 ⋅ ∇ v y d Ω − ∫ ∂ Ω N v x q 1 ⋅ n + v y q 2 ⋅ n d s -\int_{\Omega } v_x \nabla \cdot \bm q_1 +v_y \nabla \cdot \bm q_2 d\Omega =\int_{\Omega } \bm q_1 \cdot \nabla v_x + \bm q_2 \cdot \nabla v_y d\Omega -\int_{\partial \Omega_N } v_x \bm q_1 \cdot n+ v_y \bm q_2 \cdot nds − ∫ Ω v x ∇ ⋅ q 1 + v y ∇ ⋅ q 2 d Ω = ∫ Ω q 1 ⋅ ∇ v x + q 2 ⋅ ∇ v y d Ω − ∫ ∂ Ω N v x q 1 ⋅ n + v y q 2 ⋅ n d s We substitute back q 1 \bm q_1 q 1 q 2 \bm q_2 q 2
− ∫ ∂ Ω N v x q 1 ⋅ n + v y q 2 ⋅ n d s = ∫ ∂ Ω N ( σ x n x + τ x y n y ) v x + ( τ x y n x + σ x n y ) v y d s = ∫ ∂ Ω N ( σ ⋅ n ) ⋅ v d s = ∫ ∂ Ω N t ⋅ v d s
\def\arraystretch{2.5}
\begin{align*}
& -\int_{\partial\Omega_N}v_{x}\bm{q}_{1}\cdot\bm{n}+v_{y}\bm{q}_{2}\cdot\bm{n}ds=\int_{\partial\Omega_N}\left(\sigma_{x}n_{x}+\tau_{xy}n_{y}\right)v_{x}+\left(\tau_{xy}n_{x}+\sigma_{x}n_{y}\right)v_{y}ds\\
& =\int_{\partial\Omega_N}\left(\bm{\sigma}\cdot\bm{n}\right)\cdot\bm{v}ds=\int_{\partial\Omega_N}\bm{t}\cdot\bm{v}ds
\end{align*}
− ∫ ∂ Ω N v x q 1 ⋅ n + v y q 2 ⋅ n d s = ∫ ∂ Ω N ( σ x n x + τ x y n y ) v x + ( τ x y n x + σ x n y ) v y d s = ∫ ∂ Ω N ( σ ⋅ n ) ⋅ v d s = ∫ ∂ Ω N t ⋅ v d s Furthermore
q 1 ⋅ ∇ v x + q 2 ⋅ ∇ v y = σ x ∂ v x ∂ x + τ x y ∂ v x ∂ y + τ x y ∂ v y ∂ x + σ x ∂ v y ∂ y
q_1 \cdot \nabla v_x +q_2 \cdot \nabla v_y =\sigma_x \dfrac{\partial v_x }{\partial x}+\tau_{xy} \dfrac{\partial v_x }{\partial y}+\tau_{xy} \dfrac{\partial v_y }{\partial x}+\sigma_x \dfrac{\partial v_y }{\partial y}
q 1 ⋅ ∇ v x + q 2 ⋅ ∇ v y = σ x ∂ x ∂ v x + τ x y ∂ y ∂ v x + τ x y ∂ x ∂ v y + σ x ∂ y ∂ v y σ x ∂ v x ∂ x + 2 τ x y 1 2 ( ∂ v x ∂ y + ∂ v y ∂ x ) + σ x ∂ v y ∂ y = σ : ε ( v )
\sigma_x \dfrac{\partial v_x }{\partial x}+2\tau_{xy} \dfrac{1}{2}\left(\dfrac{\partial v_x }{\partial y}+\dfrac{\partial v_y }{\partial x}\right)+\sigma_x \dfrac{\partial v_y }{\partial y}= \bm \sigma : \bm \varepsilon \left(v\right)
σ x ∂ x ∂ v x + 2 τ x y 2 1 ( ∂ y ∂ v x + ∂ x ∂ v y ) + σ x ∂ y ∂ v y = σ : ε ( v ) σ : ε ( v ) = ∑ i , j σ i j ε i j \bm \sigma : \bm \varepsilon \left(v\right)=\sum_{i,j} \sigma_{ij} \varepsilon_{ij} σ : ε ( v ) = i , j ∑ σ ij ε ij here we identified the double dot product between the stress and virtual strain. A similar derivation is used as a proof of equivalence between the principle of virtual work and the equilibrium equation. See here .
Thus, from the starting expression
∫ Ω ( − ∇ ⋅ σ ( u ) ) ⋅ v d Ω = ∫ Ω f ⋅ v d Ω \int_{\Omega } \left(-\nabla \cdot \bm \sigma (\bm u) \right)\cdot \bm v d\Omega =\int_{\Omega } \bm f \cdot \bm v d\Omega ∫ Ω ( − ∇ ⋅ σ ( u ) ) ⋅ v d Ω = ∫ Ω f ⋅ v d Ω we formulate our final weak form: Find the displacement field u \bm u u u = g \bm u = \bm g u = g ∂ Ω D \partial \Omega_D ∂ Ω D
∫ Ω σ ( u ) : ε ( v ) d Ω = ∫ Ω f ⋅ v d Ω + ∫ ∂ Ω N t ⋅ v d s ∀ v that are 0 on ∂ Ω D \boxed{ \int_{\Omega } \bm \sigma (\bm u) : \bm \varepsilon \left( \bm v\right) d\Omega =\int_{\Omega } \bm f\cdot \bm vd\Omega +\int_{\partial \Omega_N } \bm{t} \cdot \bm v ds \quad \forall v \text{ that are } 0 \text{ on } \partial \Omega_D} ∫ Ω σ ( u ) : ε ( v ) d Ω = ∫ Ω f ⋅ v d Ω + ∫ ∂ Ω N t ⋅ v d s ∀ v that are 0 on ∂ Ω D More info on the weak form can be found in the section on Variational Formulation - The Weak Form
