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

revision 2252 by gross, Fri Feb 6 07:51:28 2009 UTC revision 2264 by gross, Wed Feb 11 06:48:28 2009 UTC
# Line 285  In~\cite{XXX} it has been shown that thi Line 285  In~\cite{XXX} it has been shown that thi
285  to solve the problem is coupled form, however this approach leads in some cases to a very ill-conditioned stiffness matrix in particular in the case of a very small or large permability ($\alpha\hackscore{1} \ll 1$ or $\alpha\hackscore{0} \gg 1$)    to solve the problem is coupled form, however this approach leads in some cases to a very ill-conditioned stiffness matrix in particular in the case of a very small or large permability ($\alpha\hackscore{1} \ll 1$ or $\alpha\hackscore{0} \gg 1$)
286
287  The approach we are taking is to eliminate the velocity $u$ from the problem. Assuming that $p$ is known we have  The approach we are taking is to eliminate the velocity $u$ from the problem. Assuming that $p$ is known we have
288  \begin{equation}  \begin{equation}\label{DARCY V FORM}
289  v= (I+D^*D)^{-1} (D^*f + g - Qp)  v= (I+D^*D)^{-1} (D^*f + g - Qp)
290  \end{equation}  \end{equation}
291  (notice that $(I+D^*D)$ is coercive in $W$) which is inserted into the second equation  (notice that $(I+D^*D)$ is coercive in $W$) which is inserted into the second equation
# Line 300  We use the PCG \index{linear solver!PCG} Line 300  We use the PCG \index{linear solver!PCG}
300  \begin{equation}  \begin{equation}
301  \begin{array}{rcl}  \begin{array}{rcl}
302  r & = & Q^*  \left( g -(I+D^*D)^{-1} (D^*f + g) - Qp + (I+D^*D)^{-1}Q p \right)\\  r & = & Q^*  \left( g -(I+D^*D)^{-1} (D^*f + g) - Qp + (I+D^*D)^{-1}Q p \right)\\
303  & =&  Q^* \left( - Qp - (I+D^*D)^{-1} (D^*f + g - Qp) \right) \\  & =&  Q^* \left( g- Qp - (I+D^*D)^{-1} (D^*f + g - Qp) \right) \\
304  & =&  Q^* \left( g - Qp - v \right)  & =&  Q^* \left( g - Qp - v \right)
305  \end{array}  \end{array}
306  \end{equation}  \end{equation}
307  So in a partical implementation we use the pair $(Qp,v)$ to represent the residual. This will save the  So in a partical implementation we use $\hat{r}=g-Qp-v$ to represent the residual.
308  reconstruction of the velocity $v$. In this notation the right hand side is given as  The evaluation of the iteration operator for a given $p$ is then
309  $(0,(I+D^*D)^{-1} (D^*f + g))$. The evaluation of the iteration operator for a given $p$ is then  returning $Qp+v$ where $v$ is the solution of
returning $(Qp,w)$ where $w$ is the solution of
310  \begin{equation}\label{UPDATE W}  \begin{equation}\label{UPDATE W}
311  (I+D^*D)w = Qp  (I+D^*D)v = Qp
312  \end{equation}  \end{equation}
313  We use $Q^*Q$ as a a preconditioner for the iteration operator $Q^* ( I - (I+D^*D)^{-1} ) Q$.  We use $(Q^*Q)^{-1}$ as a preconditioner for the iteration operator $Q^* ( I - (I+D^*D)^{-1} ) Q$. So the application of the preconditioner to $\hat{r}$ representing the residual is given by solving
314    implemented by solving
315  The iteration PCG \index{linear solver!PCG}\index{PCG} is terminated if  \begin{equation}\label{UPDATE P}
316  \begin{equation}\label{DARCY STOP}  Q^*Q q  = Q^*\hat{r}
\int r \cdot (Q^*Q)^{-1} r \; dx \le \mbox{ATOL}^2
317  \end{equation}  \end{equation}
318  where ATOL is a given absolute tolerance.  The residual norm used in the PCG is given as
319  The initial residual $r\hackscore{0}$ is  \begin{equation}\label{DARCY R NORM}
320  \begin{equation}\label{DARCY STOP 2}  \|r\|\hackscore{PCG}^2 = \int r \cdot (Q^*Q)^{-1} r \; dx =\int \hat{r} \cdot Q (Q^*Q)^{-1} Q^* \hat{r} \; dx \approx
321  r\hackscore{0}=Q^* \left( g - v\hackscore{ref} \right) \mbox{ with } v\hackscore{ref} = (I+D^*D)^{-1} (D^*f + g)  \|\hat{r}\|\hackscore{0}^2
322  \end{equation}  \end{equation}
323  so the  The iteration is terminated if
324  \begin{equation}\label{DARCY NORM 0}  \begin{equation}\label{DARCY STOP}
325  \int r\hackscore0 \cdot (Q^*Q)^{-1} r\hackscore0 \; dx = \int \left( g - v\hackscore{ref} \right)  \cdot  Q  p\hackscore{ref} \; dx \mbox{ with }p\hackscore{ref} = (Q^*Q)^{-1} Q^* \left( g - v\hackscore{ref} \right)  \|r\|\hackscore{PCG} \le \mbox{ATOL}
326  \end{equation}  \end{equation}
327  So we set  where we set
328  \begin{equation}\label{DARCY NORM 1}  \begin{equation}\label{DARCY ATOL DEF}
329  ATOL = atol + rtol \cdot \max(|g - v\hackscore{ref}|\hackscore{0}, |Q p\hackscore{ref} |\hackscore{0} )  \mbox{ATOL} = \mbox{atol} + \mbox{rtol} \cdot \left(\frac{1}{\|v\|\hackscore{0}} + \frac{1}{\|Qp\|\hackscore{0}} \right)^{-1}
330  \end{equation}  \end{equation}
331  where atol and rtol a given absolute and relative tolerances, respectively. The reference flux $v\hackscore{ref}$  where rtol is a given relative tolerance and $\mbox{atol}$ is a given absolute tolerance (typically $=0$).
332  and reference pressure $p\hackscore{ref}$ may be calcualated from their definition which would require to solve to  Notice that if $Qp$ and $v$ both are zero, the pair $(0,p)$ is a solution.
333  PDEs but in a practical application estimates can be used for instance solutions from previous time steps or for simplified scenarious (e.g. constant permability).  The problem is that ATOL is depending on the solution $p$ (and $v$ calculated form~\ref{DARCY V FORM}). In partcice one use the initial guess for $p$
334    to get a first value for ATOL. If the stopping crierion is met in the PCG iteration, a new $v$ is calculated from the current pressure approximation and ATOL is recalculated. If \ref{DARCY STOP} is still fullfilled the calculation is terminated and $(v,p)$ is returned. Otherwise PCG is restarted with a new ATOL.
335
336  \subsection{Functions}  \subsection{Functions}
337  \begin{classdesc}{DarcyFlow}{domain}  \begin{classdesc}{DarcyFlow}{domain}
# Line 353  The method will try to cast the given va Line 352  The method will try to cast the given va
352  \Data class objects.  \Data class objects.
353  \end{methoddesc}  \end{methoddesc}
354
355  \begin{methoddesc}[DarcyFlow]{setTolerance}{\optional{atol=0,\optional{rtol=1e-8,\optional{p_ref=None,\optional{v_ref=None}}}}}  \begin{methoddesc}[DarcyFlow]{setTolerance}{\optional{rtol=1e-4}}
356  sets the absolute tolerance ATOL according to~\ref{DARCY NORM 1}. If \var{p_ref} is not present $0$ is used.  sets the relative tolerance \mbox{rtol} in \ref{DARCY ATOL DEF}.
If \var{v_ref} is not present $0$ is used. If the final result ATOL is not positive an exception is thrown.
357  \end{methoddesc}  \end{methoddesc}
358
359    \begin{methoddesc}[DarcyFlow]{setAbsoluteTolerance}{\optional{atol=0.}}
360    sets the absolute tolerance \mbox{atol} in \ref{DARCY ATOL DEF}.
361    \end{methoddesc}
362
363    \begin{methoddesc}[DarcyFlow]{setSubProblemTolerance}{\optional{rtol=None}}
364    sets the relative tolerance used to solve the involved PDEs. If no argument is given,
365    the square of the current relative tolerance is used. The sub-problem tolerance should be choosen as large as possible to minimize the compute time. However, a too large value for the sub-problem tolerance may lead to slow convergence or even dibergence in the outer iteration.
366    \end{methoddesc}
367
368  \begin{methoddesc}[DarcyFlow]{solve}{u0,p0, \optional{max_iter=100, \optional{verbose=False \optional{sub_rtol=1.e-8}}}}  \begin{methoddesc}[DarcyFlow]{solve}{u0,p0, \optional{max_iter=100, \optional{verbose=False \optional{sub_rtol=1.e-8}}}}
369  solves the problem. and returns approximations for the flux $v$ and the pressure $p$.  solves the problem. and returns approximations for the flux $v$ and the pressure $p$.
370  \var{u0} and \var{p0} define initial guess for flux and pressure. Values marked  \var{u0} and \var{p0} define initial guess for flux and pressure. Values marked
371  by positive values \var{location_of_fixed_flux} and \var{location_of_fixed_pressure}, respectively, are kept unchanged.  by positive values \var{location_of_fixed_flux} and \var{location_of_fixed_pressure}, respectively, are kept unchanged.
\var{sub_rtol} defines the tolerance used to solve the involved PDEs. \var{sub_rtol} needs to be choosen sufficiently small to ensure convergence but users need to keep in mind that a very small value for \var{sub_rtol} will result in a long compute time. Typically  $\var{sub_rtol}=\var{rtol}^2$ is a good choice if $\var{rtol}$ is not choosen too small.
372  \end{methoddesc}  \end{methoddesc}
373
374

Legend:
 Removed from v.2252 changed lines Added in v.2264