distinguish two cases.
The continuous and convex function $\mathcal{F}$ denotes the von Mises flow function
-$$\mathcal{F}(\tau) = \vert\tau^D\vert - \sigma_0¸\quad\text{with}\quad \tau^D
-= \tau - \dfrac{1}{3}tr(\tau)I,$$
-$\sigma_0$ as yield stress and $\vert .\vert$ as the frobenius norm. If there
+@f{gather*}\mathcal{F}(\tau) = \vert\tau^D\vert - \sigma_0¸\quad\text{with}\quad \tau^D
+= \tau - \dfrac{1}{3}tr(\tau)I,@f}
+$\sigma_0$ as yield stress and $\vert .\vert$ as the Frobenius norm. If there
are no plastic deformations in a particular point - that is $\lambda=0$ - this yields $\vert\sigma^D\vert <
\sigma_0$ and otherwise if $\lambda > 0$ it follows that $\vert\sigma^D\vert = \sigma_0$.
That means if the stress is smaller than the yield stress there are only elastic
The boundary of $\Omega$ separates as follows $\Gamma=\Gamma_D\bigcup\Gamma_C$ and $\Gamma_D\bigcap\Gamma_C=\emptyset$.
At the boundary $\Gamma_D$ we have zero Dirichlet conditions. $\Gamma_C$ denotes the potential contact boundary.\\
The last two lines decribe the so-called Signorini contact conditions. If there is no contact the normal stress
-$$ \sigma_n = \sigma n\cdot n$$
+@f{gather*} \sigma_n = \sigma n\cdot n@f}
is zero with the outward normal $n$. If there is contact ($u_n = g$) the tangential stress $\sigma_t = \sigma\cdot n - \sigma_n n$
vanishes, because we consider a frictionless situation and the normal stress is
negative. The gap $g$ comes with the start configuration of the obstacle and the
As a starting point to derive the equations above, let us imagine that we want
to minimise an energy functional:
-$$E(\tau) := \dfrac{1}{2}\int\limits_{\Omega}\tau A \tau d\tau,\quad \tau\in \Pi W^{\textrm{div}}$$
+@f{gather*}E(\tau) := \dfrac{1}{2}\int\limits_{\Omega}\tau A \tau d\tau,\quad \tau\in \Pi W^{\textrm{div}}@f}
with
-$$W^{\textrm{div}}:=\lbrace \tau\in
-L^2(\Omega,\mathbb{R}^{\textrm{dim}\times\textrm{dim}}_{\textrm{sym}}):\textrm{div}(\tau)\in L^2(\Omega,\mathbb{R}^{\textrm{dim}})\rbrace$$ and
-$$\Pi \Sigma:=\lbrace \tau\in \Sigma, \mathcal{F}(\tau)\leq 0\rbrace$$
+@f{gather*}W^{\textrm{div}}:=\lbrace \tau\in
+L^2(\Omega,\mathbb{R}^{\textrm{dim}\times\textrm{dim}}_{\textrm{sym}}):\textrm{div}(\tau)\in L^2(\Omega,\mathbb{R}^{\textrm{dim}})\rbrace@f} and
+@f{gather*}\Pi \Sigma:=\lbrace \tau\in \Sigma, \mathcal{F}(\tau)\leq 0\rbrace@f}
as the set of admissible stresses which is defined
by a continious, convex flow function $\mathcal{F}$.
With the goal of deriving the dual formulation of the minimisation
problem, we define a lagrange function:
-$$L(\tau,\varphi) := E(\tau) + (\varphi, \textrm{div}(\tau)),\quad \lbrace\tau,\varphi\rbrace\in\Pi W^{\textrm{div}}\times V^+$$
+@f{gather*}L(\tau,\varphi) := E(\tau) + (\varphi, \textrm{div}(\tau)),\quad \lbrace\tau,\varphi\rbrace\in\Pi W^{\textrm{div}}\times V^+@f}
with
-$$V^+ := \lbrace u\in V: u_n\leq g \text{ on } \Gamma_C \rbrace$$
-$$V:=\left[ H_0^1 \right]^{\textrm{dim}}:=\lbrace u\in \left[H^1(\Omega)\right]^{\textrm{dim}}: u
-= 0 \text{ on } \Gamma_D\rbrace$$
+@f{gather*}V^+ := \lbrace u\in V: u_n\leq g \text{ on } \Gamma_C \rbrace@f}
+@f{gather*}V:=\left[ H_0^1 \right]^{\textrm{dim}}:=\lbrace u\in \left[H^1(\Omega)\right]^{\textrm{dim}}: u
+= 0 \text{ on } \Gamma_D\rbrace@f}
By building the Fréchet derivatives of $L$ for both components we obtain the
-dual formulation for the stationary case which is known as \textbf{Hencky-Type-Model}:\\
+dual formulation for the stationary case which is known as <i>Hencky-Type-Model</i>:\\
Find a pair $\lbrace\sigma,u\rbrace\in \Pi W\times V^+$ with
-$$\left(A\sigma,\tau - \sigma\right) + \left(u, \textrm{div}(\tau) - \textrm{div}(\sigma)\right) \geq 0,\quad \forall \tau\in \Pi W^{\textrm{div}}$$
-$$-\left(\textrm{div}(\sigma),\varphi - u\right) \geq 0,\quad \forall \varphi\in V^+.$$
+@f{gather*}\left(A\sigma,\tau - \sigma\right) + \left(u, \textrm{div}(\tau) - \textrm{div}(\sigma)\right) \geq 0,\quad \forall \tau\in \Pi W^{\textrm{div}}@f}
+@f{gather*}-\left(\textrm{div}(\sigma),\varphi - u\right) \geq 0,\quad \forall \varphi\in V^+.@f}
By integrating by parts and multiplying the first inequality by the elastic
tensor $C=A^{-1}$ we achieve the primal-mixed version of our problem:
Find a pair $\lbrace\sigma,u\rbrace\in \Pi W\times V^+$ with
-$$\left(\sigma,\tau - \sigma\right) - \left(C\varepsilon(u), \tau - \sigma\right) \geq 0,\quad \forall \tau\in \Pi W$$
-$$\left(\sigma,\varepsilon(\varphi) - \varepsilon(u)\right) \geq 0,\quad \forall \varphi\in V^+.$$
+@f{gather*}\left(\sigma,\tau - \sigma\right) - \left(C\varepsilon(u), \tau - \sigma\right) \geq 0,\quad \forall \tau\in \Pi W@f}
+@f{gather*}\left(\sigma,\varepsilon(\varphi) - \varepsilon(u)\right) \geq 0,\quad \forall \varphi\in V^+.@f}
Therein $\varepsilon$ denotes the linearised deformation tensor with $\varepsilon(u) := \dfrac{1}{2}\left(\nabla u + \nabla u^T\right)$ for small deformations.\\
Most materials - especially metals - have the property that they show some hardening effects during the forming process.
There are different constitutive laws to describe those material behaviors. The
A\sigma\vert$.
It can be considered by establishing an additional term in our primal-mixed formulation:\\
Find a pair $\lbrace(\sigma,\xi),u\rbrace\in \Pi (W\times L^2(\Omega,\mathbb{R}))\times V^+$ with
-$$\left(\sigma,\tau - \sigma\right) - \left(C\varepsilon(u), \tau - \sigma\right) + \gamma\left( \xi, \eta - \xi\right) \geq 0,\quad \forall (\tau,\eta)\in \Pi (W,L^2(\Omega,\mathbb{R}))$$
-$$\left(\sigma,\varepsilon(\varphi) - \varepsilon(u)\right) \geq 0,\quad \forall \varphi\in V^+,$$
+@f{gather*}\left(\sigma,\tau - \sigma\right) - \left(C\varepsilon(u), \tau - \sigma\right) + \gamma\left( \xi, \eta - \xi\right) \geq 0,\quad \forall (\tau,\eta)\in \Pi (W,L^2(\Omega,\mathbb{R}))@f}
+@f{gather*}\left(\sigma,\varepsilon(\varphi) - \varepsilon(u)\right) \geq 0,\quad \forall \varphi\in V^+,@f}
with the hardening parameter $\gamma > 0$.
Now we want to derive a primal problem which only depends on the displacement $u$. For that purpose we
Equations, Springer-Verlag Berlin Heidelberg, 2007 and Frohne: FEM-Simulation
der Umformtechnik metallischer Oberflächen im Mikrokosmos, Ph.D. thesis,
University of Siegen, Germany, 2011) on
-$$\left(\sigma - C\varepsilon(u), \tau - \sigma\right) \geq 0,\quad \forall \tau\in \Pi W,$$
+@f{gather*}\left(\sigma - C\varepsilon(u), \tau - \sigma\right) \geq 0,\quad \forall \tau\in \Pi W,@f}
which yields with the second inequality:\\
Find the displacement $u\in V^+$ with
-$$\left(P_{\Pi}(C\varepsilon(u)),\varepsilon(\varphi) - \varepsilon(u)\right) \geq 0,\quad \forall \varphi\in V^+,$$
+@f{gather*}\left(P_{\Pi}(C\varepsilon(u)),\varepsilon(\varphi) - \varepsilon(u)\right) \geq 0,\quad \forall \varphi\in V^+,@f}
with the projection:
-$$P_{\Pi}(\tau):=@f{cases}
+@f{gather*}P_{\Pi}(\tau):=@f{cases}
\tau, & \text{if }\vert\tau^D\vert \leq \sigma_0 + \gamma\xi,\\
\hat\alpha\dfrac{\tau^D}{\vert\tau^D\vert} + \dfrac{1}{3}tr(\tau), & \text{if }\vert\tau^D\vert > \sigma_0 + \gamma\xi,
- @f}$$
+ @f}@f}
with the radius
-$$\hat\alpha := \sigma_0 + \gamma\xi .$$
+@f{gather*}\hat\alpha := \sigma_0 + \gamma\xi .@f}
With the relation $\xi = \vert\varepsilon(u) - A\sigma\vert$ it is possible to eliminate $\xi$ inside the projection $P_{\Pi}$:\\
-$$P_{\Pi}(\tau):=@f{cases}
+@f{gather*}P_{\Pi}(\tau):=@f{cases}
\tau, & \text{if }\vert\tau^D\vert \leq \sigma_0,\\
\alpha\dfrac{\tau^D}{\vert\tau^D\vert} + \dfrac{1}{3}tr(\tau), & \text{if }\vert\tau^D\vert > \sigma_0,
- @f}$$
-$$\alpha := \sigma_0 + \dfrac{\gamma}{2\mu+\gamma}\left(\vert\tau^D\vert - \sigma_0\right) ,$$
+ @f}@f}
+@f{gather*}\alpha := \sigma_0 + \dfrac{\gamma}{2\mu+\gamma}\left(\vert\tau^D\vert - \sigma_0\right) ,@f}
with a further material parameter $\mu>0$ called shear modulus. We refer that
this only possible for isotropic plasticity.
So what we do is to calculate the stresses by using Hooke's law for linear elastic, isotropic materials
-$$\sigma = C \varepsilon(u) = 2\mu \varepsilon^D(u) + \kappa tr(\varepsilon(u))I = \left[2\mu\left(\mathbb{I} -\dfrac{1}{3} I\otimes I\right) + \kappa I\otimes I\right]\varepsilon(u)$$
+@f{gather*}\sigma = C \varepsilon(u) = 2\mu \varepsilon^D(u) + \kappa tr(\varepsilon(u))I = \left[2\mu\left(\mathbb{I} -\dfrac{1}{3} I\otimes I\right) + \kappa I\otimes I\right]\varepsilon(u)@f}
with the material parameter $\kappa>0$ (bulk modulus). The variables $I$ and
$\mathbb{I}$ denote the identity tensors of second and forth order. In that
notation $2\mu \varepsilon^D(u)$ is the deviatoric part and $\kappa
<h3>Linearisation of the constitutive law for the Newton method</h3>
For the Newton method we have to linearise the following semi-linearform
-$$a(\psi;\varphi) := \left(P_{\Pi}(C\varepsilon(\varphi)),\varepsilon(\varphi)\right).$$
+@f{gather*}a(\psi;\varphi) := \left(P_{\Pi}(C\varepsilon(\varphi)),\varepsilon(\varphi)\right).@f}
Because we have to find the solution $u$ in the convex set $V^+$, we have to
apply an SQP-method (SQP: sequential quadratic programming). That means we have
to solve a minimisation problem for a known $u^i$ in every SQP-step of the form
&\rightarrow& \textrm{min},\quad u^{i+1}\in V^+.
@f}
Neglecting the constant terms $ a(u^i;u^i)$ and $ a'(u^i;u^i,u^i)$ we obtain the
-following minimisation problem $$\dfrac{1}{2} a'(u^i;u^{i+1},u^{i+1}) - F(u^i)\rightarrow \textrm{min},\quad u^{i+1}\in V^+$$ with
-$$F(\varphi) := \left(a'(\varphi;\varphi,u^{i+1}) - a(\varphi;u^{i+1}) \right).$$
+following minimisation problem @f{gather*}\dfrac{1}{2} a'(u^i;u^{i+1},u^{i+1}) - F(u^i)\rightarrow \textrm{min},\quad u^{i+1}\in V^+@f} with
+@f{gather*}F(\varphi) := \left(a'(\varphi;\varphi,u^{i+1}) - a(\varphi;u^{i+1}) \right).@f}
In the case of our constitutive law the Fréchet derivative of the
semi-linearform $a(.;.)$ at the point $u^i$ is
-$$a'(u^i;\psi,\varphi) =
-(I(x)\varepsilon(\psi),\varepsilon(\varphi)),\quad x\in\Omega,$$ $$
+@f{gather*}a'(u^i;\psi,\varphi) =
+(I(x)\varepsilon(\psi),\varepsilon(\varphi)),\quad x\in\Omega,@f} @f{gather*}
I(x) := @f{cases}
2\mu\left(\mathbb{I} - \dfrac{1}{3} I\otimes I\right) + \kappa I\otimes I, &
\quad \vert \tau^D \vert \leq \sigma_0\\
- \dfrac{\tau^D\otimes\tau^D}{\vert\tau^D\vert}\right) + \kappa I\otimes I,
&\quad \vert \tau^D \vert > \sigma_0
@f}
-$$
+@f}
with
-$$\tau^D := C\varepsilon^D(u^i).$$
+@f{gather*}\tau^D := C\varepsilon^D(u^i).@f}
Remark that $a(.;.)$ is not differentiable in the common sense but it is
slantly differentiable like the function for the contact problem and again we refer to
Hintermueller, Ito, Kunisch: The primal-dual active set strategy as a semismooth newton method, SIAM J. OPTIM., 2003, Vol. 13, No. 3, pp. 865-888.
We introduce a Lagrange multiplier $\lambda$ and the convex cone $K\subset W'$,
$W'$ dual space of the trace space $W$ of $V$ restricted to $\Gamma_C$,
-$$K:=\{\mu\in W':\mu_T = 0,\quad\langle\mu n,v\rangle_{\Gamma_C}\geq 0,\quad
-\forall v\in W, v \ge 0\text{ on }\Gamma_C \}$$
+@f{gather*}K:=\{\mu\in W':\mu_T = 0,\quad\langle\mu n,v\rangle_{\Gamma_C}\geq 0,\quad
+\forall v\in W, v \ge 0\text{ on }\Gamma_C \}@f}
of Lagrange multipliers, where $\langle\cdot,\cdot\rangle$
denotes the duality pairing, i.e. a boundary integral, between $W'$ and $W$.
Intuitively, $K$ is the cone of all "non-positive functions", except that $ K\subset
Now we describe an algorithm that combines the Newton-method, which we use for
the nonlinear constitutive law, with the semismooth Newton method for the contact. It
sums up the results of the sections before and works as follows:
-@f{itemize}
- \item[(0)] Initialize $\mathcal{A}_k$ and $\mathcal{F}_k$, such that $\mathcal{S} = \mathcal{A}_k \cup \mathcal{F}_k$ and $\mathcal{A}_k \cap \mathcal{F}_k = \emptyset$ and set $k = 1$.
- \item[(1)] Assemble the Newton matrix $a'(U^k;\varphi_i,\varphi_j)$ and the right-hand-side $F(U^k)$.
- \item[(2)] Find the primal-dual pair $(U^k,\Lambda^k)$ that satisfies
+<ol>
+ <li> Initialize $\mathcal{A}_k$ and $\mathcal{F}_k$, such that $\mathcal{S} = \mathcal{A}_k \cup \mathcal{F}_k$ and $\mathcal{A}_k \cap \mathcal{F}_k = \emptyset$ and set $k = 1$.
+ <li> Assemble the Newton matrix $a'(U^k;\varphi_i,\varphi_j)$ and the right-hand-side $F(U^k)$.
+ <li> Find the primal-dual pair $(U^k,\Lambda^k)$ that satisfies
@f{align*}
AU^k + B\Lambda^k & = F, &\\
\left[B^TU^k\right]_i & = G_i & \forall i\in\mathcal{A}_k,\\
\Lambda^k_i & = 0 & \forall i\in\mathcal{F}_k.
@f}
- \item[(3)] Define the new active and inactive sets by
- $$\mathcal{A}_{k+1}:=\lbrace i\in\mathcal{S}:\Lambda^k_i +
- c\left(\left[B^TU^k\right]_i - G_i\right) > 0\rbrace,$$
- $$\mathcal{F}_{k+1}:=\lbrace i\in\mathcal{S}:\Lambda^k_i +
- c\left(\left[B^TU^k\right]_i - G_i\right) \leq 0\rbrace.$$
- \item[(4)] If $\mathcal{A}_{k+1} = \mathcal{A}_k$ and $\vert
+ <li> Define the new active and inactive sets by
+ @f{gather*}\mathcal{A}_{k+1}:=\lbrace i\in\mathcal{S}:\Lambda^k_i +
+ c\left(\left[B^TU^k\right]_i - G_i\right) > 0\rbrace,@f}
+ @f{gather*}\mathcal{F}_{k+1}:=\lbrace i\in\mathcal{S}:\Lambda^k_i +
+ c\left(\left[B^TU^k\right]_i - G_i\right) \leq 0\rbrace.@f}
+ <li> If $\mathcal{A}_{k+1} = \mathcal{A}_k$ and $\vert
F\left(U^{k+1}\right) \vert < \delta$ then stop, else set $k=k+1$ and go to
- step (1).
-@f}
+ step (2).
+</ol>
The mass matrix $B\in\mathbb{R}^{n\times m}$, $n>m$, is quadratic in our
situation since $\Lambda^k$ is only defined on $\Gamma_C$:
-$$B_{ij} = @f{cases}
+@f{gather*}B_{ij} = \begin{cases}
\int\limits_{\Gamma_C}\varphi_i^2(x)dx, & \text{if}\quad i=j\\
0, & \text{if}\quad i\neq j.
-@f}$$
+\end{cases}@f}
So $m$ denotes the size of $\Lambda^k$ and $i$ a degree of freedom. In our
programm we use the structure of a quadratic sparse for $B\in\mathbb{R}^{n\times
n}$ and the length of $\Lambda^k$ is $n$ with $\Lambda^k_i = 0$ for $i>m$.
The vector $G$ is defined by a suitable approximation $g_h$ of the gap $g$
-$$G_i = @f{cases}
+@f{gather*}G_i = \begin{cases}
\int\limits_{\Gamma_C}g_h(x)\varphi_i(x)dx, & \text{if}\quad i\leq m\\
0, & \text{if}\quad i>m.
-@f}$$
+\end{cases}@f}
Compared to step-41, step (1) is added but it should be clear
from the sections above that we only linearize the problem. In step (2) we have to solve a linear
<h3>Implementation</h3>
-This tutorial is essentailly a mixture of step-40 and step-41 but instead of
+This tutorial is essentially a mixture of step-40 and step-41 but instead of
PETSc we let the Trilinos library deal with parallelizing the linear algebra
(like in step-32). Since we are trying to solve a similar problem like in
step-41 we will use the same methods but now in parallel.
motivation in Chinese. If your audience is Japanese, please see the other entry
for motivation. This is a word in Japanese and Korean, but it means "motive
power" or "kinetic energy" (without the motivation meaning that you are
-probably looking for)").
\ No newline at end of file
+probably looking for)").
+
projects / dealii-svn.git / commitdiff