| 16 |
|
|
| 17 |
The following sections give a breif overview of the model classes and their corresponding methods. |
The following sections give a breif overview of the model classes and their corresponding methods. |
| 18 |
|
|
| 19 |
\section{Stokes Cartesian (Saddle Point Problem)} |
\section{Stokes Problem} |
| 20 |
|
The velocity field $v$ and pressure $p$ of an incompressible fluid is given as the solution of their |
| 21 |
\subsection{Description} |
Stokes problem |
| 22 |
|
\begin{equation}\label{Stokes 1} |
| 23 |
|
-\left(\eta(v\hackscore{i,j}+ v\hackscore{i,j})\right)\hackscore{,j}+p\hackscore{,i}=F\hackscore{i} |
| 24 |
|
\end{equation} |
| 25 |
|
where $\eta$ is the viscosity and $F_i$ defines an internal force. We assume an incompressible media: |
| 26 |
|
\begin{equation}\label{Stokes 2} |
| 27 |
|
-v\hackscore{i,i}=0 |
| 28 |
|
\end{equation} |
| 29 |
|
Natural boundary conditions are taken in the form |
| 30 |
|
\begin{equation}\label{Stokes Boundary} |
| 31 |
|
\left(\eta(v\hackscore{i,j}+ v\hackscore{i,j})\right)n\hackscore{j}-n\hackscore{i}p=f |
| 32 |
|
\end{equation} |
| 33 |
|
which can be overwritten by a constraint |
| 34 |
|
\begin{equation}\label{Stokes Boundary0} |
| 35 |
|
v\hackscore{i}(x)=0 |
| 36 |
|
\end{equation} |
| 37 |
|
where the index $i$ may depend on the location $x$ on the bondary. |
| 38 |
|
|
|
Saddle point type problems emerge in a number of applications throughout physics and engineering. Finite element discretisation of the Navier-Stokes (momentum) equations for incompressible flow leads to equations of a saddle point type, which can be formulated as a solution of the following operator problem for $u \in V$ and $p \in Q$ with suitable Hilbert spaces $V$ and $Q$: |
|
| 39 |
|
|
| 40 |
|
\subsection{Solution Method \label{STOKES SOLVE}} |
| 41 |
|
In block form equation equations~\ref{Stokes 1} and~\ref{Stokes 2} takes the form of a saddle point problem |
| 42 |
\begin{equation} |
\begin{equation} |
| 43 |
\left[ \begin{array}{cc} |
\left[ \begin{array}{cc} |
| 44 |
A & B \\ |
A & B^{*} \\ |
| 45 |
b^{*} & 0 \\ |
B & 0 \\ |
| 46 |
\end{array} \right] |
\end{array} \right] |
| 47 |
\left[ \begin{array}{c} |
\left[ \begin{array}{c} |
| 48 |
u \\ |
u \\ |
| 49 |
p \\ |
p \\ |
| 50 |
\end{array} \right] |
\end{array} \right] |
| 51 |
=\left[ \begin{array}{c} |
=\left[ \begin{array}{c} |
| 52 |
f \\ |
F \\ |
| 53 |
g \\ |
0 \\ |
| 54 |
\end{array} \right] |
\end{array} \right] |
| 55 |
\label{SADDLEPOINT} |
\label{SADDLEPOINT} |
| 56 |
\end{equation} |
\end{equation} |
| 57 |
|
where $A$ is coercive, self-adjoint linear operator in $V$, $B$ is a divergence operator and $B^{*}$ is it adjoint operator (=gradient operator)). For more details on the mathematics see references \cite{AAMIRBERKYAN2008,MBENZI2005}. |
| 58 |
|
We apply the preconditioner matrix |
| 59 |
|
\begin{equation} |
| 60 |
|
\left[ \begin{array}{cc} |
| 61 |
|
A^{-1} & 0 \\ |
| 62 |
|
S^{-1} B A^{-1} & -S^{-1} \\ |
| 63 |
|
\end{array} \right] |
| 64 |
|
\label{SADDLEPOINT PRECODITIONER} |
| 65 |
|
\end{equation} |
| 66 |
|
with the Schur complement $S=BA^{-1}B^{*}$ to solve the problem iteratively. The updates $dv$ and $dp$ for |
| 67 |
|
given velocity $v$ and pressure $p$ is given as |
| 68 |
|
\begin{equation} |
| 69 |
|
\begin{array}{rcl} |
| 70 |
|
A dv & = & F-Av-B^{*}p \\ |
| 71 |
|
S dp & = & B(v+dv) \\ |
| 72 |
|
\end{array} |
| 73 |
|
\label{SADDLEPOINT iteration} |
| 74 |
|
\end{equation} |
| 75 |
|
This scheme is called the Uzawa scheme. |
| 76 |
|
|
| 77 |
|
In the case of the Stokes problem it can be shown that $S$ can be approximated by $\frac{1}{\eta}$. With this the iteration scheme can be implemented as |
| 78 |
|
\begin{equation} |
| 79 |
|
\begin{array}{rcl} |
| 80 |
|
-\left(\eta(dv\hackscore{i,j}+ dv\hackscore{i,j})\right)\hackscore{,j} & = & F\hackscore{i}+\left(\eta(v\hackscore{i,j}+ v\hackscore{i,j})\right)\hackscore{,j}-p\hackscore{,i} \\ |
| 81 |
|
\frac{1}{\eta} dp & = & - \left(v\hackscore{,i}+dv\hackscore{,i} \right)\hackscore{,i} \\ |
| 82 |
|
\end{array} |
| 83 |
|
\label{SADDLEPOINT iteration 2} |
| 84 |
|
\end{equation} |
| 85 |
|
To accelerate the convergence we are using the restarted $GMRES$ method using the norm |
| 86 |
|
\begin{equation} |
| 87 |
|
\|(v,p)\|^2= \int\hackscore{\Omega} \eta^2 v\hackscore{i,j}^2 + p^2 \; dx |
| 88 |
|
\label{SADDLEPOINT iteration 3} |
| 89 |
|
\end{equation} |
| 90 |
|
or alternatively the $PCG$ method on the pressure only using the norm |
| 91 |
|
\begin{equation} |
| 92 |
|
\|p\|^2= \int\hackscore{\Omega} \frac{1}{\eta^2} p^2 \; dx |
| 93 |
|
\label{SADDLEPOINT iteration 4} |
| 94 |
|
\end{equation} |
| 95 |
|
|
| 96 |
|
|
| 97 |
|
|
|
where $A$ is coercive, self-adjoint linear operator in $V$, $B$ is a linear operator from $Q$ into $V$ and $B^{*}$ is the adjoint operator of $B$. $f$ and $g$ are given elements from $V$ and $Q$ respectivitly. For more details on the mathematics see references \cite{AAMIRBERKYAN2008,MBENZI2005}. |
|
| 98 |
|
|
|
The Uzawa scheme scheme is used to solve the momentum equation with the secondary condition of incompressibility \cite{GROSS2006,AAMIRBERKYAN2008}. |
|
| 99 |
|
|
| 100 |
\begin{classdesc}{StokesProblemCartesian}{domain,debug} |
\begin{classdesc}{StokesProblemCartesian}{domain,debug} |
| 101 |
opens the stokes equations on the \Domain domain. Setting debug=True switches the debug mode to on. |
opens the stokes equations on the \Domain domain. Setting debug=True switches the debug mode to on. |
| 170 |
\end{equation} |
\end{equation} |
| 171 |
with the elastic strain given as |
with the elastic strain given as |
| 172 |
\begin{equation}\label{IKM-EQU-3} |
\begin{equation}\label{IKM-EQU-3} |
| 173 |
D\hackscore{ij}'^{el}=\frac{1}{2 \mu} \sigma\hackscore{ij}' |
D\hackscore{ij}'^{el}=\frac{1}{2 \mu} \dot{\sigma}\hackscore{ij}' |
| 174 |
\end{equation} |
\end{equation} |
| 175 |
where $\sigma'\hackscore{ij}$ is the deviatoric stress (Notice that $\sigma'\hackscore{ii}=0$). |
where $\sigma'\hackscore{ij}$ is the deviatoric stress (Notice that $\sigma'\hackscore{ii}=0$). |
| 176 |
If the material is composed by materials $q$ the visco-plastic strain can be decomposed as |
If the material is composed by materials $q$ the visco-plastic strain can be decomposed as |
| 213 |
\begin{equation}\label{IKM-EQU-6b} |
\begin{equation}\label{IKM-EQU-6b} |
| 214 |
\frac{1}{\eta^{vp}}=\max(\sum\hackscore{q} \frac{1}{\eta^{q}}, \frac{\dot{\gamma}^{vp}} {\tau\hackscore{Y} + \beta \; p}) |
\frac{1}{\eta^{vp}}=\max(\sum\hackscore{q} \frac{1}{\eta^{q}}, \frac{\dot{\gamma}^{vp}} {\tau\hackscore{Y} + \beta \; p}) |
| 215 |
\end{equation} |
\end{equation} |
|
Now we can combine equations~\ref{IKM-EQU-2}, \ref{IKM-EQU-3} and~\ref{IKM-EQU-6b} to get |
|
|
\begin{equation}\label{IKM-EQU-10} |
|
|
D\hackscore{ij}'=\frac{1}{2 \eta\hackscore{eff}} \sigma\hackscore{ij}' \mbox{ with } |
|
|
\frac{1}{\eta\hackscore{eff}}=\frac{1}{\mu}+\frac{1}{\eta^{vp}} |
|
|
\end{equation} |
|
| 216 |
The deviatoric stress needs to fullfill the equilibrion equation |
The deviatoric stress needs to fullfill the equilibrion equation |
| 217 |
\begin{equation}\label{IKM-EQU-1} |
\begin{equation}\label{IKM-EQU-1} |
| 218 |
-\sigma'\hackscore{ij,j}+p\hackscore{,i}=F\hackscore{i} |
-\sigma'\hackscore{ij,j}+p\hackscore{,i}=F\hackscore{i} |
| 221 |
\begin{equation}\label{IKM-EQU-2} |
\begin{equation}\label{IKM-EQU-2} |
| 222 |
-v\hackscore{i,i}=0 |
-v\hackscore{i,i}=0 |
| 223 |
\end{equation} |
\end{equation} |
| 224 |
|
Natural boundary conditions are taken in the form |
| 225 |
|
\begin{equation}\label{IKM-EQU-Boundary} |
| 226 |
|
\sigma'\hackscore{ij}n\hackscore{j}-n\hackscore{i}p=f |
| 227 |
|
\end{equation} |
| 228 |
|
which can be overwritten by a constraint |
| 229 |
|
\begin{equation}\label{IKM-EQU-Boundary0} |
| 230 |
|
v\hackscore{i}(x)=0 |
| 231 |
|
\end{equation} |
| 232 |
|
where the index $i$ may depend on the location $x$ on the bondary. |
| 233 |
|
|
| 234 |
|
\subsection{Solution Method \label{IKM-SOLVE}} |
| 235 |
|
By using a first order finite difference approximation wit step size $dt>0$~\ref{IKM-EQU-3} get the form |
| 236 |
|
\begin{equation}\label{IKM-EQU-3b} |
| 237 |
|
D\hackscore{ij}'^{el}=\frac{1}{2 \mu dt } \left( \sigma\hackscore{ij}' - \sigma\hackscore{ij}^{'-} \right) |
| 238 |
|
\end{equation} |
| 239 |
|
where $\sigma\hackscore{ij}^{'-}$ is the deviatoric stress at the precious time step. |
| 240 |
|
Now we can combine equations~\ref{IKM-EQU-2}, \ref{IKM-EQU-3b} and~\ref{IKM-EQU-6b} to get |
| 241 |
|
\begin{equation}\label{IKM-EQU-10} |
| 242 |
|
\sigma\hackscore{ij}' = 2 \eta\hackscore{eff} D\hackscore{ij}' + |
| 243 |
|
\frac{\eta\hackscore{eff}}{\mu \; dt} |
| 244 |
|
\sigma\hackscore{ij}^{'-} \mbox{ with } |
| 245 |
|
\frac{1}{\eta\hackscore{eff}}=\frac{1}{\mu \; dt}+\frac{1}{\eta^{vp}} |
| 246 |
|
\end{equation} |
| 247 |
|
Notice that $\eta\hackscore{eff}$ is a function of stress. |
| 248 |
After inserting~\ref{IKM-EQU-10} into~\ref{IKM-EQU-1} we get |
After inserting~\ref{IKM-EQU-10} into~\ref{IKM-EQU-1} we get |
| 249 |
\begin{equation}\label{IKM-EQU-1ib} |
\begin{equation}\label{IKM-EQU-1ib} |
| 250 |
-\left(\eta\hackscore{eff} (v\hackscore{i,j}+ v\hackscore{i,j})\right)\hackscore{,j}+p\hackscore{,i}=F\hackscore{i} |
-\left(\eta\hackscore{eff} (v\hackscore{i,j}+ v\hackscore{i,j})\right)\hackscore{,j}+p\hackscore{,i}=F\hackscore{i}+ |
| 251 |
|
\frac{\eta\hackscore{eff}}{\mu dt } \sigma\hackscore{ij,j}^{'-} |
| 252 |
|
\end{equation} |
| 253 |
|
Together with the incomressibilty condition~\ref{IKM-EQU-2} we need to solve a problem with a form almost identical |
| 254 |
|
to the Stokes problem discussed in section~\label{STOKES SOLVE} but with the difference that $\eta\hackscore{eff}$ is depending on the solution. Analog to the iteration scheme~\ref{SADDLEPOINT iteration 2} we can run |
| 255 |
|
\begin{equation} |
| 256 |
|
\begin{array}{rcl} |
| 257 |
|
-\left(\eta\hackscore{eff}(dv\hackscore{i,j}+ dv\hackscore{i,j})\right)\hackscore{,j} & = & F\hackscore{i}+ \sigma\hackscore{ij,j}' \\ |
| 258 |
|
\frac{1}{\eta\hackscore{eff}} dp & = & - v\hackscore{i,i}^+ \\ |
| 259 |
|
d\sigma\hackscore{ij}' & = & 2 \eta\hackscore{eff} D\hackscore{ij}^{+'} + \frac{\eta\hackscore{eff}}{\mu \; dt} \sigma\hackscore{ij}' - \sigma\hackscore{ij}\\ |
| 260 |
|
\end{array} |
| 261 |
|
\label{IKM iteration 2} |
| 262 |
|
\end{equation} |
| 263 |
|
where $v^+=v+dv$. As this problem is non-linear the Jacobi-free Newton-GMRES method is used with the norm |
| 264 |
|
\begin{equation} |
| 265 |
|
\|(\sigma, p)\|^2= \int\hackscore{\Omega} \sigma\hackscore{i,j}^2 + p^2 \; dx |
| 266 |
|
\label{IKM iteration 3} |
| 267 |
\end{equation} |
\end{equation} |
| 268 |
|
|
| 269 |
|
|
Parent Directory
| 
Revision Log
| 
Patch