 # Diff of /trunk/doc/user/stokessolver.tex

revision 2843 by gross, Thu Jan 14 05:51:28 2010 UTC revision 2851 by gross, Fri Jan 15 07:37:47 2010 UTC
# Line 6  We want to calculate the velocity \index Line 6  We want to calculate the velocity \index
6  \begin{equation}\label{Stokes 1}  \begin{equation}\label{Stokes 1}
7  -\left(\eta(v\hackscore{i,j}+ v\hackscore{j,i})\right)\hackscore{,j}+p\hackscore{,i}=f\hackscore{i}-\sigma\hackscore{ij,j}  -\left(\eta(v\hackscore{i,j}+ v\hackscore{j,i})\right)\hackscore{,j}+p\hackscore{,i}=f\hackscore{i}-\sigma\hackscore{ij,j}
8  \end{equation}  \end{equation}
9  where  $f\hackscore{i}$ defines an internal force \index{force, internal} and $\sigma\hackscore{ij}$ is an intial stress \index{stress, initial}. The viscosity $\eta$ may weakly depend on pressure and velocity. If relevant we will use the notation $\eta(v,p)$ to express this dependency.  where  $f\hackscore{i}$ defines an internal force \index{force, internal} and $\sigma\hackscore{ij}$ is an initial stress \index{stress, initial}. The viscosity $\eta$ may weakly depend on pressure and velocity. If relevant we will use the notation $\eta(v,p)$ to express this dependency.
10
11  We assume an incompressible media:  We assume an incompressible media:
12  \begin{equation}\label{Stokes 2}  \begin{equation}\label{Stokes 2}
# Line 62  for the pressure and an additional corre Line 62  for the pressure and an additional corre
62   \label{STOKES ITER STEP 2}   \label{STOKES ITER STEP 2}
63   \end{equation}   \end{equation}
64  The new velocity and pressure are then given by $v\hackscore{2}=v\hackscore{1}-dv\hackscore{2}$ and  The new velocity and pressure are then given by $v\hackscore{2}=v\hackscore{1}-dv\hackscore{2}$ and
65  $p\hackscore{2}=p\hackscore{0}+dp\hackscore{2}$ which will fullfill the block system~\ref{STOKES}.  $p\hackscore{2}=p\hackscore{0}+dp\hackscore{2}$ which will fulfill the block system~\ref{STOKES}.
66  This solution strategy is called the Uzawa scheme \index{Uzawa scheme}.  This solution strategy is called the Uzawa scheme \index{Uzawa scheme}.
67
68  There is a problem with this scheme: In practice we will use an iterative scheme  There is a problem with this scheme: In practice we will use an iterative scheme
# Line 80  In the following we will use the two nor Line 80  In the following we will use the two nor
80  \|p\|\hackscore{0}^2= \int\hackscore{\Omega} p^2 \; dx.  \|p\|\hackscore{0}^2= \int\hackscore{\Omega} p^2 \; dx.
81  \label{STOKES STOP}  \label{STOKES STOP}
82  \end{equation}  \end{equation}
83  for velicity $v$ and pressure $p$. The iteration is terminated if the stopping criterium  for velocity $v$ and pressure $p$. The iteration is terminated if the stopping criteria
84   \begin{equation} \label{STOKES STOPPING CRITERIA}   \begin{equation} \label{STOKES STOPPING CRITERIA}
85  \max(\|Bv\hackscore{1}\|\hackscore{0},\|v\hackscore{2}-v\hackscore{0}\|\hackscore{1}) \le \tau \cdot \|v\hackscore{2}\|\hackscore{1}  \max(\|Bv\hackscore{1}\|\hackscore{0},\|v\hackscore{2}-v\hackscore{0}\|\hackscore{1}) \le \tau \cdot \|v\hackscore{2}\|\hackscore{1}
86   \end{equation}   \end{equation}
# Line 90  norm of $B A^{-1} B^{*} dp\hackscore{2}$ Line 90  norm of $B A^{-1} B^{*} dp\hackscore{2}$
90
91  We want to optimize the tolerance choice for solving~\ref{STOKES ITER STEP 1}  We want to optimize the tolerance choice for solving~\ref{STOKES ITER STEP 1}
92  and~\ref{STOKES ITER STEP 2}. To do this we write the iteration scheme as a fixed point problem. Here  and~\ref{STOKES ITER STEP 2}. To do this we write the iteration scheme as a fixed point problem. Here
93  we consider the errors produced by the iterative solvers beeing used.  we consider the errors produced by the iterative solvers being used.
94  From Equation~\ref{STOKES ITER STEP 1} we have  From Equation~\ref{STOKES ITER STEP 1} we have
95  \begin{equation} \label{STOKES total V1}  \begin{equation} \label{STOKES total V1}
96  v\hackscore{1} = e\hackscore{1} + v\hackscore{0} + A^{-1} ( G - Av\hackscore{0} - B^{*} p\hackscore{0} )  v\hackscore{1} = e\hackscore{1} + v\hackscore{0} + A^{-1} ( G - Av\hackscore{0} - B^{*} p\hackscore{0} )
# Line 123  the problem Line 123  the problem
123  \begin{equation} \label{STOKES P OPERATOR}  \begin{equation} \label{STOKES P OPERATOR}
124  A w = B^{*} s  A w = B^{*} s
125  \end{equation}  \end{equation}
126  with sufficient accuracy to return $q=Bw$. We assume that the choicen tolerance is  with sufficient accuracy to return $q=Bw$. We assume that the desired tolerance is
127  sufficiently small, for instane one can take $\tau\hackscore{2}^2$  sufficiently small, for instance one can take $\tau\hackscore{2}^2$
128  where $\tau\hackscore{2}$ is the tolerance for~\ref{STOKES ITER STEP 2}.  where $\tau\hackscore{2}$ is the tolerance for~\ref{STOKES ITER STEP 2}.
129
130  In an implementation we use the fact that the residual $r$ is given as  In an implementation we use the fact that the residual $r$ is given as
# Line 155  setting the errors to zero we can write Line 155  setting the errors to zero we can write
155  \begin{equation}  \begin{equation}
156  v = \Phi(v,p) \mbox{ and } p = \Psi(u,p)  v = \Phi(v,p) \mbox{ and } p = \Psi(u,p)
157  \end{equation}  \end{equation}
158  with suitable functions $\Phi(v,p)$ and $\Psi(v,p)$ representing the itartion operator without  with suitable functions $\Phi(v,p)$ and $\Psi(v,p)$ representing the iteration operator without
159  errors. In fact for a linear problem,  $\Phi$ and $\Psi$ are constant. With this notation we can write  errors. In fact for a linear problem,  $\Phi$ and $\Psi$ are constant. With this notation we can write
160  the update step in the form $p\hackscore{2}= \delta p + \Psi(v\hackscore{0},p\hackscore{0})$ and  the update step in the form $p\hackscore{2}= \delta p + \Psi(v\hackscore{0},p\hackscore{0})$ and
161  $v\hackscore{2}= \delta v + \Phi(v\hackscore{0},p\hackscore{0})$ where  $v\hackscore{2}= \delta v + \Phi(v\hackscore{0},p\hackscore{0})$ where
# Line 201  So we choose the value for $\tau\hacksco Line 201 So we choose the value for$\tau\hacksco
201   M \tau\hackscore{2} \|B v\hackscore{1}\|\hackscore{0}  \le (\chi^{-})^2 \epsilon^{-}   M \tau\hackscore{2} \|B v\hackscore{1}\|\hackscore{0}  \le (\chi^{-})^2 \epsilon^{-}
202  \end{equation}  \end{equation}
203  in order to make the perturbation for the termination of the pressure iteration a second order effect. We use a  in order to make the perturbation for the termination of the pressure iteration a second order effect. We use a
204  similar argument for the velcity:  similar argument for the velocity:
205  \begin{equation}\label{STOKES EST 2}  \begin{equation}\label{STOKES EST 2}
206  \begin{array}{rcl}  \begin{array}{rcl}
207  \|v\hackscore{2}-v\|\hackscore{1} & \le & \|v\hackscore{2}-\delta v-v\|\hackscore{1} + \| \delta v\|\hackscore{1} \\  \|v\hackscore{2}-v\|\hackscore{1} & \le & \|v\hackscore{2}-\delta v-v\|\hackscore{1} + \| \delta v\|\hackscore{1} \\
# Line 216  So we choose the value for $\tau\hacksco Line 216 So we choose the value for$\tau\hacksco
216  \begin{equation} \label{STOKES TOL1}  \begin{equation} \label{STOKES TOL1}
217  K \tau\hackscore{1} \le \chi^{-}  K \tau\hackscore{1} \le \chi^{-}
218  \end{equation}  \end{equation}
219  Assuming we have etsimates for $M$ and $K$ (if none is available, we use the value $1$.)  Assuming we have estimates for $M$ and $K$ (if none is available, we use the value $1$.)
220  we can use~\ref{STOKES TOL1} and~\ref{STOKES TOL2} to get appropriate values for the tolerances. After  we can use~\ref{STOKES TOL1} and~\ref{STOKES TOL2} to get appropriate values for the tolerances. After
221  the step has been completed we can calculate a new convergence rate $\chi =\frac{\epsilon}{\epsilon^{-}}$.  the step has been completed we can calculate a new convergence rate $\chi =\frac{\epsilon}{\epsilon^{-}}$.
222  For partical reasons we restrict $\chi$ to be less or equal a given maximum value $\chi\hackscore{max}\le 1$.  For partial reasons we restrict $\chi$ to be less or equal a given maximum value $\chi\hackscore{max}\le 1$.
223  If we see $\chi \le \chi^{-} (1+\chi^{-})$ our choices for the tolerances was suitable. Otherwise, we need to adjust the values for $K$ and $M$. From the estimates~\ref{STOKES EST 1} and~\ref{STOKES EST 2} we establish  If we see $\chi \le \chi^{-} (1+\chi^{-})$ our choices for the tolerances was suitable. Otherwise, we need to adjust the values for $K$ and $M$. From the estimates~\ref{STOKES EST 1} and~\ref{STOKES EST 2} we establish
224  \begin{equation}\label{STOKES EST 3}  \begin{equation}\label{STOKES EST 3}
225  \chi \le ( 1 + \max(M \frac{\tau\hackscore{2} \|B v\hackscore{1}\|\hackscore{0}}{\chi^{-} \epsilon^{-}},K \tau\hackscore{1}  ) ) \cdot \chi^{-}  \chi \le ( 1 + \max(M \frac{\tau\hackscore{2} \|B v\hackscore{1}\|\hackscore{0}}{\chi^{-} \epsilon^{-}},K \tau\hackscore{1}  ) ) \cdot \chi^{-}
# Line 230  $M^{+}$ and $K^{+}$ then we get Line 230  $M^{+}$ and $K^{+}$ then we get
230  \chi =  ( 1 + \max(M^{+} \frac{\chi^{-}}{M},K^{+} \frac{\chi^{-}}{K}) ) \cdot \chi^{-}  \chi =  ( 1 + \max(M^{+} \frac{\chi^{-}}{M},K^{+} \frac{\chi^{-}}{K}) ) \cdot \chi^{-}
231  \end{equation}  \end{equation}
232  From this equation we set for  From this equation we set for
233  than our choice for $K$ was not good anough. In this case we can calaculate a new value  than our choice for $K$ was not good enough. In this case we can calculate a new value
234   \begin{equation}   \begin{equation}
235  K^{+} =  \frac{\chi-\chi^{-}}{(\chi^{-})^2} K  K^{+} =  \frac{\chi-\chi^{-}}{(\chi^{-})^2} K
236  \end{equation}  \end{equation}
# Line 274  the solution procedure is implemented as Line 274  the solution procedure is implemented as
274    \begin{enumerate}    \begin{enumerate}
275    \item update $p\hackscore{2}\leftarrow p$ and $v\hackscore{2}\leftarrow v\hackscore{1}$    \item update $p\hackscore{2}\leftarrow p$ and $v\hackscore{2}\leftarrow v\hackscore{1}$
276     \end{enumerate}     \end{enumerate}
277     \item calculate convergence messure $\epsilon$ and convergence rate $\chi$     \item calculate convergence measure $\epsilon$ and convergence rate $\chi$
278  \item if stopping criterium~\ref{STOKES STOPPING CRITERIA} holds:  \item if stopping criteria~\ref{STOKES STOPPING CRITERIA} holds:
279   \begin{enumerate}   \begin{enumerate}
280   \item return $v\hackscore{2}$ and $p\hackscore{2}$   \item return $v\hackscore{2}$ and $p\hackscore{2}$
281   \end{enumerate}   \end{enumerate}
# Line 291  the solution procedure is implemented as Line 291  the solution procedure is implemented as
291
292  \begin{classdesc}{StokesProblemCartesian}{domain}  \begin{classdesc}{StokesProblemCartesian}{domain}
293  opens the Stokes problem\index{Stokes problem} on the \Domain domain. The domain  opens the Stokes problem\index{Stokes problem} on the \Domain domain. The domain
294  needs to support LBB compliant elements for the Stokes problem, see~\cite{LBB} for detials~\index{LBB condition}.  needs to support LBB compliant elements for the Stokes problem, see~\cite{LBB} for details~\index{LBB condition}.
295  For instance one can use second order polynomials for velocity and  For instance one can use second order polynomials for velocity and
296  first order polynomials for the pressure on the same element. Alternativly, one can use  first order polynomials for the pressure on the same element. Alternatively, one can use
297  macro elements\index{macro elements} using linear polynomial for both pressure and velocity bu with a subdivided  macro elements\index{macro elements} using linear polynomial for both pressure and velocity    with a subdivided
298  element for the velocity. Typically, the macro element is more cost effective. The fact that pressure and velocity are represented in different way is expressed by  element for the velocity. Typically, the macro element is more cost effective. The fact that pressure and velocity are represented in different way is expressed by
299  \begin{python}  \begin{python}
300  velocity=Vector(0.0, Solution(mesh))  velocity=Vector(0.0, Solution(mesh))
# Line 307  restoration_factor=0}}}}}} Line 307  restoration_factor=0}}}}}}
307  assigns values to the model parameters. In any call all values must be set.  assigns values to the model parameters. In any call all values must be set.
308  \var{f} defines the external force $f$, \var{eta} the viscosity $\eta$,  \var{f} defines the external force $f$, \var{eta} the viscosity $\eta$,
309  \var{surface_stress} the surface stress $s$ and \var{stress} the initial stress $\sigma$.  \var{surface_stress} the surface stress $s$ and \var{stress} the initial stress $\sigma$.
310  The locations and compontents where the velocity is fixed are set by  The locations and components where the velocity is fixed are set by
311  the values of \var{fixed_u_mask}. \var{restoration_factor} defines the restoring force factor $\alpha$.  the values of \var{fixed_u_mask}. \var{restoration_factor} defines the restoring force factor $\alpha$.
312  The method will try to cast the given values to appropriate  The method will try to cast the given values to appropriate
313  \Data class objects.  \Data class objects.
# Line 338  returns the current relative tolerance. Line 338  returns the current relative tolerance.
338  \end{methoddesc}  \end{methoddesc}
339  \begin{methoddesc}[StokesProblemCartesian]{setAbsoluteTolerance}{\optional{tolerance=0.}}  \begin{methoddesc}[StokesProblemCartesian]{setAbsoluteTolerance}{\optional{tolerance=0.}}
340  sets the absolute tolerance for the error in the relevant norm. The tolerance must be non-negative. Typically the  sets the absolute tolerance for the error in the relevant norm. The tolerance must be non-negative. Typically the
341  absolute talerance is set to 0.  absolute tolerance is set to 0.
342  \end{methoddesc}  \end{methoddesc}
343
344  \begin{methoddesc}[StokesProblemCartesian]{getAbsoluteTolerance}{}  \begin{methoddesc}[StokesProblemCartesian]{getAbsoluteTolerance}{}
345  sreturns the current absolute tolerance.  returns the current absolute tolerance.
346  \end{methoddesc}  \end{methoddesc}
347
348  \begin{methoddesc}[StokesProblemCartesian]{getSolverOptionsVelocity}{}  \begin{methoddesc}[StokesProblemCartesian]{getSolverOptionsVelocity}{}

Legend:
 Removed from v.2843 changed lines Added in v.2851