Chris Nkinthorn
Finite Element Methods
Problem 1
Let $\Omega$ be a region in $\mathbb{R}^2$ and let its boundary $\Gamma = \overline{\Gamma_1\cup\Gamma_2\cup\Gamma_3\cup\Gamma_4}$ be composed of non overlapping subregions of $\Gamma$. Let $n$ be the unit outward normal vector to $\Gamma$ such that $s$ and $n$ form a right-hand rule basis. Consider the following BVP in classical linear elastostatics
Given:
f : Ī© ā R ; q i : Ī 1 ā R ; h i : Ī 2 ā R ; q n & h s : Ī 4 ā R ; q s & h n : Ī 4 ā R ; f:\Omega\rightarrow \mathbb{R}; \quad
q_i:\Gamma_1\rightarrow \mathbb{R}; \quad
h_i:\Gamma_2\rightarrow \mathbb{R}; \quad
q_n \& h_s : \Gamma_4\rightarrow \mathbb{R};\quad
q_s \& h_n : \Gamma_4\rightarrow \mathbb{R};\quad f : Ī© ā R ; q i ā : Ī 1 ā ā R ; h i ā : Ī 2 ā ā R ; q n ā & h s ā : Ī 4 ā ā R ; q s ā & h n ā : Ī 4 ā ā R ; find $u_i : \overline{\Omega}\rightarrow\mathbb{R}$ such that
Ļ i j , j + f i = 0 u i = q i Ļ i j n j = h i \begin{align}
\sigma_{ij,j} + f_i = 0
\\
u_i = q_i
\\
\sigma_{ij}n_j = h_i
\end{align} Ļ ij , j ā + f i ā = 0 u i ā = q i ā Ļ ij ā n j ā = h i ā ā ā where $\sigma{ij} = c {ijkl}u_{(k,l)}$. Establish a weak formulation for this problem in which all ā$q$-typeā boundary conditions are essential and all ā$h$-type ā boundary conditions are natural. State all requirements on the spaces $\delta$ and $\mathcal{V}$. Hint $w_i = w_n n_i+ w_s s_i$.
Method of Weighted Residuals
\begin{align}
\sigma_{ij,j}+f_i = 0
\rightarrow
0 &=
\int_\Omega w_i
\left(
\sigma_{ij,j}+f_i
\right)
d\Omega
\\
&=
\overbrace{\int_\Omega w_i
\sigma_{ij,j}
d\Omega
}^\text{Integrate by Parts}+
\int_\Omega w_if_i
d\Omega
\\
&=
\int_\Omega w_if_i
d\Omega
-
\int_\Omega
w_{(i,j)}\sigma_{ij}
d\Omega
+
\int_{\Gamma}
\left(
w_i \sigma_{ij}n_j
\right)
d\Gamma
\\
&=
\int_\Omega w_if_i
d\Omega
-
\int_\Omega
w_{(i,j)}\sigma_{ij}
d\Omega
+
\overbrace{
\cancelto{0}{
\int_{\Gamma_{g_i}}
\left(
w_i \sigma_{ij}n_j
\right)
d\Gamma}
+
\int_{\Gamma_{h_i}}
\left(
w_i \sigma_{ij}n_j
\right)
d\Gamma
}^\text{from IBP}
\\
&=
\int_\Omega w_if_i
d\Omega
-
\int_\Omega
w_{(i,j)}\sigma_{ij}
d\Omega
+
\int_{\Gamma_{h_i}}
\left(
w_i \sigma_{ij}n_j
\right)
d\Gamma
\end{align}
The boundary $\Gamma_{h_i}$ is the union of the 4 non overlapping boundaries per $\Gamma = \overline{\Gamma_1\cup\Gamma_2\cup\Gamma_3\cup\Gamma_4}$. Each correspond to one of the given kinds of boundary conditions, which we are asked to generalize. Investigation of the boundary contour integral requires one integral.
Function Space Assumptions
Ī“ = { u ā£ u ā H 1 , u i = g i Ā onĀ Ī g i } V = { w ā£ w ā H 1 , w i = 0 Ā onĀ Ī g i } \delta = \{{\bf u}|{\bf u}\in \mathcal{H}^1,u_i=g_i \text{ on } \Gamma_{g_i}\}
\qquad
\mathcal{V} = \{{\bf w}|{\bf w}\in \mathcal{H}^1,w_i=0 \text{ on } \Gamma_{g_i}\} Ī“ = { u ā£ u ā H 1 , u i ā = g i ā Ā onĀ Ī g i ā ā } V = { w ā£ w ā H 1 , w i ā = 0 Ā onĀ Ī g i ā ā } The usual decomposition of these function spaces such that our finite spaces are $\delta^h \in \delta \text{ and } \mathcal{V}^h \in \mathcal{V}$.
ā« Ī h i ( w i Ļ i j n j ) d Ī = ā« Ī h 1 ( w i Ļ i j n j ) d Ī ā B C 1 : g i + ā« Ī h 2 ( w i Ļ i j n j ) ā w i h i d Ī ā B C 2 + ā« Ī h 3 ( w i Ļ i j n j ) d Ī ā B C 3 + ā« Ī h 4 ( w i Ļ i j n j ) d Ī ā B C 4 ā applyĀ hint:Ā w i = w n n i + w s s i = g i + ā« Ī h 2 w i h i d Ī + ā« Ī h 3 ( w n n i + w s s i ) Ļ i j n j d Ī + ā« Ī h 4 ( w n n i + w s s i ) Ļ i j n j d Ī = g i + ā« Ī h 2 w i h i d Ī + ā« Ī h 3 w n n i Ļ i j n j d Ī ā g n + ā« Ī h 3 w s s i Ļ i j n j ā h s d Ī + ā« Ī h 4 w s s i Ļ i j n j d Ī ā g s + ā« Ī h 4 w n n i Ļ i j n j ā h n d Ī \int_{\Gamma_{h_i}}
\left(
w_i \sigma_{ij}n_j
\right)
d\Gamma
=
\overbrace{\int_{\Gamma_{h_1}}
\left(
w_i \sigma_{ij}n_j
\right)
d\Gamma
}^{{BC}_1:g_i}
+
\overbrace{\int_{\Gamma_{h_2}}
\underbrace{\left(
w_i \sigma_{ij}n_j
\right)}_{w_ih_i}
d\Gamma
}^{BC_2}
+
\underbrace{\overbrace{\int_{\Gamma_{h_3}}
\left(
w_i \sigma_{ij}n_j
\right)
d\Gamma
}^{BC_3}
+
\overbrace{\int_{\Gamma_{h_4}}
\left(
w_i \sigma_{ij}n_j
\right)
d\Gamma
}^{BC_4}}_{\text{apply hint: }w_i = w^nn_i+w^ss_i}
\\
=
g_i
+
\int\limits_{\Gamma_{h_2}}
w_i h_i
d\Gamma
+
\int_{\Gamma_{h_3}}
\left(w^nn_i+w^ss_i\right) \sigma_{ij}n_j
d\Gamma
+
\int_{\Gamma_{h_4}}
\left(w^nn_i+w^ss_i\right) \sigma_{ij}n_j
d\Gamma
\\
=
g_i
+
\int\limits_{\Gamma_{h_2}}
w_i h_i
d\Gamma
+
\overbrace{\int_{\Gamma_{h_3}}
w^nn_i \sigma_{ij}n_j
d\Gamma}^{g_n}
+
\int_{\Gamma_{h_3}}
w^s
\overbrace{s_i \sigma_{ij}n_j}^{h_s}
d\Gamma
\\
+
\overbrace{
\int_{\Gamma_{h_4}}
w^ss_i \sigma_{ij}n_j
d\Gamma}^{g_s}
+
\int_{\Gamma_{h_4}}
\overbrace{w^nn_i \sigma_{ij}n_j}^{h_n}
d\Gamma ā« Ī h i ā ā ā ( w i ā Ļ ij ā n j ā ) d Ī = ā« Ī h 1 ā ā ā ( w i ā Ļ ij ā n j ā ) d Ī ā BC 1 ā : g i ā ā + ā« Ī h 2 ā ā ā w i ā h i ā ( w i ā Ļ ij ā n j ā ) ā ā d Ī ā B C 2 ā ā + applyĀ hint:Ā w i ā = w n n i ā + w s s i ā ā« Ī h 3 ā ā ā ( w i ā Ļ ij ā n j ā ) d Ī ā B C 3 ā ā + ā« Ī h 4 ā ā ā ( w i ā Ļ ij ā n j ā ) d Ī ā B C 4 ā ā ā ā = g i ā + Ī h 2 ā ā ā« ā w i ā h i ā d Ī + ā« Ī h 3 ā ā ā ( w n n i ā + w s s i ā ) Ļ ij ā n j ā d Ī + ā« Ī h 4 ā ā ā ( w n n i ā + w s s i ā ) Ļ ij ā n j ā d Ī = g i ā + Ī h 2 ā ā ā« ā w i ā h i ā d Ī + ā« Ī h 3 ā ā ā w n n i ā Ļ ij ā n j ā d Ī ā g n ā ā + ā« Ī h 3 ā ā ā w s s i ā Ļ ij ā n j ā ā h s ā ā d Ī + ā« Ī h 4 ā ā ā w s s i ā Ļ ij ā n j ā d Ī ā g s ā ā + ā« Ī h 4 ā ā ā w n n i ā Ļ ij ā n j ā ā h n ā ā d Ī Finally, we can express the simplified form:
0 = g i ā£ Ī h 1 + ā« Ī h 2 w i h i d Ī + g n ā£ Ī h 3 + ā« Ī h 3 w s h s d Ī + g s ā£ Ī h 4 + ā« Ī h 4 h n d Ī \boxed{0=g_i\Big|_{\Gamma_{h_1}}
+
\int\limits_{\Gamma_{h_2}}
w_i h_i
d\Gamma
+
g_n\Big|_{\Gamma_{h_3}}
+
\int_{\Gamma_{h_3}}
w^s
h_s
d\Gamma
+
g_s\Big|_{\Gamma_{h_4}}
+
\int_{\Gamma_{h_4}}
h_n
d\Gamma
} 0 = g i ā ā Ī h 1 ā ā ā + Ī h 2 ā ā ā« ā w i ā h i ā d Ī + g n ā ā Ī h 3 ā ā ā + ā« Ī h 3 ā ā ā w s h s ā d Ī + g s ā ā Ī h 4 ā ā ā + ā« Ī h 4 ā ā ā h n ā d Ī ā Problem 2
In practice it is often useful to generalize the constitutive equation of classical elasticity to the form
Ļ i j = c i j k l ( Ļµ k l ā Ļµ k l 0 ) + Ļ i j 0 \sigma_{i j}=c_{i j k l}\left(\epsilon_{k l}-\epsilon_{k l}^{0}\right)+\sigma_{i j}^{0} Ļ ij ā = c ijk l ā ( Ļµ k l ā ā Ļµ k l 0 ā ) + Ļ ij 0 ā Where $\epsilon{ij}^0$ and $\sigma {ij}^0$
are the initial strain and initial stress both given functions of $x$. The initial strain term may be used to represent thermal expansion effects by way of
Ļµ k l 0 = ā Īø c k l \epsilon_{k l}^{0}=-\theta c_{k l} Ļµ k l 0 ā = ā Īø c k l ā where $\theta$ is the temperature and the $c_{kl}$ās are the thermal expansion coefficients (both given functions). Clearly [const eqn] will in no way change the stiffness matrix . However there will be additional contributions to $f_p^e$ Generalize the definition of $f_p^e$ to account for these additional terms.
ā« Ī© w ( i , j ) Ļ i j d Ī© = ā« Ī© w i f i d Ī© + ā i = 1 n s d ( ā« Ī h i w i h i d Ī ) \int_{\Omega} w_{(i, j)} \sigma_{i j} d \Omega=\int_{\Omega} w_{i} f_{i} d \Omega+\sum_{i=1}^{n_{s d}}\left(\int_{\Gamma_{h_{i}}} w_{i} h_{i} d \Gamma\right) ā« Ī© ā w ( i , j ) ā Ļ ij ā d Ī© = ā« Ī© ā w i ā f i ā d Ī© + i = 1 ā n s d ā ā ( ā« Ī h i ā ā ā w i ā h i ā d Ī ) Substitution of first given to the second adds to the new contributions of the
ā« Ī© w ( i , j ) c i j k l Ļµ k l d Ī© = ā« Ī© w i ā i d Ī© + ā i = 1 n s d ( ā« Ī k i w i h i d Ī ) + ā« Ī© w ( i , j ) c i j k l Ļµ k l 0 d Ī© ā ā« Ī© w ( i , j ) Ļ i j 0 d Ī© \int_{\Omega} w_{(i, j)} c_{i j k l} \epsilon_{k l} d \Omega=\int_{\Omega} w_{i} \ell_{i} d \Omega+\sum_{i=1}^{n_{s d}}\left(\int_{\Gamma_{k_{i}}} w_{i} h_{i} d \Gamma\right)
% the rest of it
+\int_{\Omega} w_{(i, j)} c_{i j k l} \epsilon_{k l}^{0} d \Omega
-\int_{\Omega} w_{(i, j)} \sigma_{i j}^{0} d \Omega ā« Ī© ā w ( i , j ) ā c ijk l ā Ļµ k l ā d Ī© = ā« Ī© ā w i ā ā i ā d Ī© + i = 1 ā n s d ā ā ( ā« Ī k i ā ā ā w i ā h i ā d Ī ) + ā« Ī© ā w ( i , j ) ā c ijk l ā Ļµ k l 0 ā d Ī© ā ā« Ī© ā w ( i , j ) ā Ļ ij 0 ā d Ī© So, the new component is
f p e = āÆ + e i T ā« Ī© e B a T D Īø ā newĀ scalar c d Ī© ā e i T ā« Ī© e B a T Ļ 0 d Ī© \boxed{f_{p}^{e}=\cdots+e_{i}^{T} \int_{\Omega^{e}} \boldsymbol{B}_{a}^{T} \boldsymbol{D} \overbrace{\theta}^\text{new scalar}\mathbf{c} d \Omega-e_{i}^{T} \int_{\Omega^{e}} \boldsymbol{B}_{a}^{T} \boldsymbol{\sigma}^{0} d \Omega} f p e ā = āÆ + e i T ā ā« Ī© e ā B a T ā D Īø newĀ scalar c d Ī© ā e i T ā ā« Ī© e ā B a T ā Ļ 0 d Ī© ā The elasticity and stress tensors were reduced by 3 orders and 1 order, respectively, leveraging their symmetry.
c = { c 11 c 22 2 c 12 } ; Ļ 0 = { Ļ 11 0 Ļ 22 0 Ļ 12 0 } {\bf c}=\left\{\begin{array}{c}{c_{11}} \\ {c_{22}} \\ {2 c_{12}}\end{array}\right\} ; \quad \boldsymbol{\sigma}^{0}=\left\{\begin{array}{c}{\boldsymbol{\sigma}_{11}^{0}} \\ {\boldsymbol{\sigma}_{22}^{0}} \\ {\boldsymbol{\sigma}_{12}^{0}}\end{array}\right\} c = ā© āØ ā§ ā c 11 ā c 22 ā 2 c 12 ā ā ā ā¬ ā« ā ; Ļ 0 = ā© āØ ā§ ā Ļ 11 0 ā Ļ 22 0 ā Ļ 12 0 ā ā ā ā¬ ā« ā Problem 3
The mesh and boundary conditions given below are for a 2-D elastostatics problem. Following the method given in the class handout, determine the ID matrix, IEN vectors and LM matrix. Also indicate the locations in the stiffness matrix, or force vector the ācontributionsā for the terms in the element stiffness matrix for element 7 will go. Essential BC information (be sure to account for it): The vertical edge (nodes 1,4,7,8 are on it) can not move in the in either direction The horizontal edge (nodes 15, 3,13 on it) has a prescribed vertical displacement of 0.05 units. Note ā Nodes are marked with a black dot (there are corner nodes and three edge nodes). Element numbers have the circles around them.
Let $n$ be $null$, $f$ be a degree of freedom, and $g$ be a nonzero essential boundary condition. As $n_\text{sd} =2$, the first entry corresponds to the $x$ direction and the second with $y$.
Element 7
$\left{\begin{matrix}3:\text{dof}&13:\text{dof}&6:\text{dof}&12:\text{dof}\1:\text{dog}&2:\text{dog}&7:\text{dof}&14:\text{dof}\end{matrix}\right}$
Stiffness Matrix
k 7 = [ k 3 , 3 f 3 k 3 , 13 f 3 k 3 , 6 k 3 , 7 k 3 , 12 k 3 , 14 ā ā ā ā ā ā ā ā k 13 , 3 f 4 k 3 , 13 f 3 k 3 , 6 k 3 , 7 k 3 , 12 k 3 , 14 ā ā ā ā ā ā ā ā k 6 , 3 f 6 k 6 , 13 f 6 k 6 , 6 k 6 , 7 k 6 , 12 k 6 , 14 k 7 , 3 f 7 k 7 , 13 f 7 k 7 , 6 k 7 , 7 k 7 , 12 k 7 , 14 k 12 , 3 f 12 k 12 , 13 f 12 k 3 , 6 k 3 , 7 k 3 , 12 k 3 , 14 k 14 , 3 f 14 k 14 , 13 f 3 k 14 , 6 k 14 , 7 k 14 , 12 k 14 , 14 ] {\bf k}^7
=
\begin{bmatrix}
k_{3,3}&f_3&k_{3,13}&f_3&k_{3,6}&k_{3,7}&k_{3,12}&k_{3,14}\\
- & - & - & - & - & - & - & - \\
k_{13,3}&f_4&k_{3,13}&f_3&k_{3,6}&k_{3,7}&k_{3,12}&k_{3,14}\\
- & - & - & - & - & - & - & - \\
k_{6,3}&f_{6}&k_{6,13}&f_{6}&k_{6,6}&k_{6,7}&k_{6,12}&k_{6,14}\\
k_{7,3}&f_7&k_{7,13}&f_7&k_{7,6}&k_{7,7}&k_{7,12}&k_{7,14}\\
k_{12,3}&f_{12}&k_{12,13}&f_{12}&k_{3,6}&k_{3,7}&k_{3,12}&k_{3,14}\\
k_{14,3}&f_{14}&k_{14,13}&f_3&k_{14,6}&k_{14,7}&k_{14,12}&k_{14,14}\\
\end{bmatrix} k 7 = ā k 3 , 3 ā ā k 13 , 3 ā ā k 6 , 3 ā k 7 , 3 ā k 12 , 3 ā k 14 , 3 ā ā f 3 ā ā f 4 ā ā f 6 ā f 7 ā f 12 ā f 14 ā ā k 3 , 13 ā ā k 3 , 13 ā ā k 6 , 13 ā k 7 , 13 ā k 12 , 13 ā k 14 , 13 ā ā f 3 ā ā f 3 ā ā f 6 ā f 7 ā f 12 ā f 3 ā ā k 3 , 6 ā ā k 3 , 6 ā ā k 6 , 6 ā k 7 , 6 ā k 3 , 6 ā k 14 , 6 ā ā k 3 , 7 ā ā k 3 , 7 ā ā k 6 , 7 ā k 7 , 7 ā k 3 , 7 ā k 14 , 7 ā ā k 3 , 12 ā ā k 3 , 12 ā ā k 6 , 12 ā k 7 , 12 ā k 3 , 12 ā k 14 , 12 ā ā k 3 , 14 ā ā k 3 , 14 ā ā k 6 , 14 ā k 7 , 14 ā k 3 , 14 ā k 14 , 14 ā ā ā