/[escript]/trunk/doc/user/linearPDE.tex

Diff of /trunk/doc/user/linearPDE.tex

Parent Directory | Revision Log | View Patch Patch

-revision 659 by gross,
Fri Mar 24 07:34:34 2006 UTC
+revision 660 by gross,
Fri Mar 24 08:43:21 2006 UTC
 Line 11
  \chapter{The module \linearPDEs}
  \declaremodule{extension}{linearPDEs} \modulesynopsis{Linear partial pifferential equation handler}
  The module \linearPDEs provides an interface to define and solve linear partial
  differential equations within \escript. \linearPDEs does not provide any
  solver capabilities in itself but hands the PDE over to
  the PDE solver library defined through the \Domain of the PDE.
  The general interface is provided through the \LinearPDE class. The
  \AdvectivePDE which is derived from the \LinearPDE class
  provides an interface to PDE dominated by its advective terms. The \Poisson
  class which is also derived form the \LinearPDE class should be used
  to define the Poisson equation \index{Poisson}.
  \section{\LinearPDE Class}
  \label{SEC LinearPDE}
 Line 27 
 to define the Poisson equation \index{Po
  The \LinearPDE class is used to define a general linear, steady, second order PDE
  for an unknown function $u$ on a given $\Omega$ defined through a \Domain object.
  In the following $\Gamma$ denotes the boundary of the domain $\Omega$. $n$ denotes
  the outer normal field on $\Gamma$.
  For a single PDE with a solution with a single component the linear PDE is defined in the
  following form:
  \begin{equation}\label{LINEARPDE.SINGLE.1}
- -(A\hackscore{jl} u\hackscore{,l})\hackscore{,j}+(B\hackscore{j} u)\hackscore{,j}+C\hackscore{l} u\hackscore{,l}+D u =-X\hackscore{j,j}+Y \; .
+ -(A\hackscore{jl} u\hackscore{,l})\hackscore{,j}-(B\hackscore{j} u)\hackscore{,j}+C\hackscore{l} u\hackscore{,l}+D u =-X\hackscore{j,j}+Y \; .
  \end{equation}
  $u_{,j}$ denotes the derivative of $u$ with respect to the $j$-th spatial direction. Einstein's summation convention, ie. summation over indexes appearing twice in a term of a sum is performed, is used.
  The coefficients $A$, $B$, $C$, $D$, $X$ and $Y$ have to be specified through \Data objects in the
  \Function on the PDE or objects that can be converted into such \Data objects.
  $A$ is a \RankTwo, $B$, $C$ and $X$ are \RankOne and $D$ and $Y$ are scalar.
  The following natural
  boundary conditions are considered \index{boundary condition!natural} on $\Gamma$:
  \begin{equation}\label{LINEARPDE.SINGLE.2}
  n\hackscore{j}(A\hackscore{jl} u\hackscore{,l}+B\hackscore{j} u)+d u=n\hackscore{j}X\hackscore{j} + y  \;.
  \end{equation}
  Notice that the coefficients $A$, $B$ and $X$ are defined in the PDE. The coefficients $d$ and $y$ are
  each a \Scalar in the \FunctionOnBoundary.  Constraints \index{constraint} for the solution prescribing the value of the
  solution at certain locations in the domain. They have the form
  \begin{equation}\label{LINEARPDE.SINGLE.3}
  u=r \mbox{ where } q>0
 Line 52 
 u=r \mbox{ where } q>0
  $r$ and $q$ are each \Scalar where $q$ is the characteristic function
  \index{characteristic function} defining where the constraint is applied.
  The constraints defined by \eqn{LINEARPDE.SINGLE.3} override any other condition set by \eqn{LINEARPDE.SINGLE.1}
  or \eqn{LINEARPDE.SINGLE.2}.
  For a system of PDEs and a solution with several components the PDE has the form
  \begin{equation}\label{LINEARPDE.SYSTEM.1}
- -(A\hackscore{ijkl} u\hackscore{k,l}){,j}+(B\hackscore{ijk} u\hackscore{k})\hackscore{,j}+C\hackscore{ikl} u\hackscore{k,l}+D\hackscore{ik} u\hackscore{k} =-X\hackscore{ij,j}+Y\hackscore{i} \; .
+ -(A\hackscore{ijkl} u\hackscore{k,l})\hackscore{,j}-(B\hackscore{ijk} u\hackscore{k})\hackscore{,j}+C\hackscore{ikl} u\hackscore{k,l}+D\hackscore{ik} u\hackscore{k} =-X\hackscore{ij,j}+Y\hackscore{i} \; .
  \end{equation}
  $A$ is a \RankFour, $B$ and $C$ are each a \RankThree, $D$ and $X$ are each a \RankTwo and $Y$ is a \RankOne.
  The natural boundary conditions \index{boundary condition!natural} take the form:
  \begin{equation}\label{LINEARPDE.SYSTEM.2}
  n\hackscore{j}(A\hackscore{ijkl} u\hackscore{k,l}+B\hackscore{ijk} u\hackscore{k})+d\hackscore{ik} u\hackscore{k}=n\hackscore{j}X\hackscore{ij}+y\hackscore{i}  \;.
  \end{equation}
  The coefficient $d$ is a \RankTwo and $y$ is a
  \RankOne both in the \FunctionOnBoundary. Constraints \index{constraint} take the form
  \begin{equation}\label{LINEARPDE.SYSTEM.3}
  u\hackscore{i}=r\hackscore{i} \mbox{ where } q\hackscore{i}>0
  \end{equation}
  $r$ and $q$ are each \RankOne. Notice that not necessarily all components must
  have a constraint at all locations.
  \LinearPDE also supports solution discontinuities \index{discontinuity} over contact region $\Gamma^{contact}$
  in the domain $\Omega$. To specify the conditions across the discontinuity we are using the
  generalised flux $J$ which is in the case of a systems of PDEs and several components of the solution
  defined as
  \begin{equation}\label{LINEARPDE.SYSTEM.5}
  J\hackscore{ij}=A\hackscore{ijkl}u\hackscore{k,l}+B\hackscore{ijk}u\hackscore{k}-X\hackscore{ij}
  \end{equation}
 Line 82 
 For the case of single solution componen
  \begin{equation}\label{LINEARPDE.SINGLE.5}
  J\hackscore{j}=A\hackscore{jl}u\hackscore{,l}+B\hackscore{j}u\hackscore{k}-X\hackscore{j}
  \end{equation}
  In the context of discontinuities \index{discontinuity} $n$ denotes the normal on the
  discontinuity pointing from side 0 towards side 1. For a system of PDEs
  the contact condition takes the form
  \begin{equation}\label{LINEARPDE.SYSTEM.6}
 Line 91 
 n\hackscore{j} J^{0}\hackscore{ij}=n\hac
  where $J^{0}$ and $J^{1}$ are the fluxes on side $0$ and side $1$ of the
  discontinuity $\Gamma^{contact}$, respectively. $[u]$, which is the difference
  of the solution at side 1 and at side 0, denotes the jump of $u$ across $\Gamma^{contact}$.
  The coefficient $d^{contact}$ is a \RankTwo and $y^{contact}$ is a
  \RankOne both in the \FunctionOnContactZero or \FunctionOnContactOne.
  In case of a single PDE and a single component solution the contact condition takes the form
  \begin{equation}\label{LINEARPDE.SINGLE.6}
 Line 111 
 A\hackscore{ijkl}=A\hackscore{klij} \\
  B\hackscore{ijk}=C\hackscore{kij} \\
  D\hackscore{ik}=D\hackscore{ki} \\
  d\hackscore{ik}=d\hackscore{ki} \\
  d^{contact}\hackscore{ik}=d^{contact}\hackscore{ki}
  \end{eqnarray}
  Note that different from the scalar case~\eqn{LINEARPDE.SINGLE.4} now the coefficients $D$, $d$ abd $d^{contact}$
  have to be inspected.
  \section{\LinearPDE class}
  This is the general class to define a linear PDE in \escript. We list a selction of the most
  important methods of the class only and refer to reference guide \ReferenceGuide for a complete list.
  \begin{classdesc}{LinearPDE}{domain,numEquations=0,numSolutions=0}
  opens a linear, steady, second order PDE on the \Domain \var{domain}. \var{numEquations}
  and \var{numSolutions} gives the number of equations and the number of solutiopn components.
  If \var{numEquations} and \var{numSolutions} is non-positive, the number of equations
  and the number solutions, respctively, stay undefined until a coefficient is
  defined.
  \end{classdesc}
  \begin{methoddesc}[LinearPDE]{setValue}{
- \optional{A=Data()}\optional{, B=Data()},
+ \optional{A}\optional{, B},
- \optional{, C=Data()}\optional{, D=Data()}
+ \optional{, C}\optional{, D}
- \optional{, X=Data()}\optional{, Y=Data()}
+ \optional{, X}\optional{, Y}
- \optional{, d=Data()}\optional{, y=Data()}
+ \optional{, d}\optional{, y}
- \optional{, d_contact=Data()}\optional{, y_contact=Data()}
+ \optional{, d_contact}\optional{, y_contact}
- \optional{, q=Data()}\optional{, r=Data()}}
+ \optional{, q}\optional{, r}}
- assigns new values to coefficients.
+ assigns new values to coefficients. By dafault all values are assumed to be zero\footnote{
+ In fact it is assumed they are not present by assigning the value \code{escript.Data()}. The
+ can by used by the solver library to reduce computational costs.
+ }
  If the new coefficient value is not a \Data object, it is converted into a \Data object in the
  appropriate \FunctionSpace.
  \end{methoddesc}
  \begin{methoddesc}[LinearPDE]{getCoefficient}{name}
  return the value assigned to coefficient \var{name}. If \var{name} is not a valid name
  an exception is raised.
  \end{methoddesc}
  \begin{methoddesc}[LinearPDE]{getShapeOfCoefficient}{name}
-Line 168 
 Otherwise \False is returned.
+Line 171 
 Otherwise \False is returned.
  \begin{methoddesc}[LinearPDE]{setSolverMethod}{\optional{solver=LinearPDE.DEFAULT}\optional{, preconditioner=LinearPDE.DEFAULT}}
  sets the solver method and preconditioner to be used. It is pointed out that a PDE solver library
  may not know the specified solver method but may choose a similar method and preconditioner.
  \end{methoddesc}
  \begin{methoddesc}[LinearPDE]{getSolverMethodName}{}
-Line 180 
 returns the solver method and preconditi
+Line 183 
 returns the solver method and preconditi
  \end{methoddesc}
  \begin{methoddesc}[LinearPDE]{setSolverPackage}{\optional{package=LinearPDE.DEFAULT}}
  Set the solver package to be used by PDE library to solve the linear systems of equations. The
  specified package may not be supported by the PDE solver library. In this case, dependng on
  the PDE solver, the default solver is used or an exeption is thrown.
  If \var{package} is not specified, the default package of the PDE solver library is used.
  \end{methoddesc}
  \begin{methoddesc}[LinearPDE]{getSolverPackage}{}
-Line 193 
 returns the linear solver package curren
+Line 196 
 returns the linear solver package curren
  \begin{methoddesc}[LinearPDE]{setTolerance}{\optional{tol=1.e-8}}:
  resets the tolerance for solution. The actually meaning of tolerance is
  depending on the underlying PDE library. In most cases, the tolerance
  will only consider the error from solving the discerete problem but will
  not consider any discretization error.
  \end{methoddesc}
-Line 219 
 returns the number of components of the
+Line 222 
 returns the number of components of the
  \end{methoddesc}
  \begin{methoddesc}[LinearPDE]{checkSymmetry}{verbose=\False}
  returns \True if the PDE is symmetric and \False otherwise.
  The method is very computational expensive and should only be
  called for testing purposes. The symmetry flag is not altered.
  If \var{verbose}=\True information about where symmetry is violated
  are printed.
-Line 246 
 indicates that the PDE is not symmetric.
+Line 249 
 indicates that the PDE is not symmetric.
  \end{methoddesc}
  \begin{methoddesc}[LinearPDE]{setReducedOrderOn}{}
  switches on the reduction of polynomial order for the solution and equation evaluation even if
  a quadratic or higher interpolation order is defined in the \Domain. This feature may not
  be supported by all PDE libraries.
  \end{methoddesc}
  \begin{methoddesc}[LinearPDE]{setReducedOrderOff}{}
  switches off the reduction of polynomial order for the solution and
  equation evaluation.
  \end{methoddesc}
-Line 280 
 returns the \Operator and right hand sid
+Line 283 
 returns the \Operator and right hand sid
  \optional{, restart=-1}
  }
  returns (an approximation of) the solution of the PDE. If \code{verbose=\True} some information during the solution process printed.
  \var{reordering} selects a reordering methods that is applied before or during the solution process
  (=\NOREORDERING ,\MINIMUMFILLIN ,\NESTEDDESCTION).
  \var{iter_max} specifies the maximum number of iteration steps that are allowed to reach the specified tolerance.
  \var{drop_tolerance} specifies a relative tolerance for small elements to be dropped when building a preconditioner
  (eg. in \ILUT). \var{drop_storage} limits the extra storage allowed when building a preconditioner
  (eg. in \ILUT). The extra storage is given relative to the size of the stiffness matrix, eg.
  \var{drop_storage=1.2} will allow the preconditioner to use the $1.2$ fold storage space than used
  for the stiffness matrix. \var{truncation} defines the truncation.
  \end{methoddesc}
  \begin{memberdesc}[LinearPDE]{DEFAULT}
  default method, preconditioner or package to be used to solve the PDE. An appropriate method should be
  chosen by the used PDE solver library.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{SCSL}
  the SCSL library by SGI,~\Ref{SCSL}\footnote{The SCSL library will only be available on SGI systems}
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{MKL}
-Line 313 
 the solver library of \finley, see \Sec{
+Line 316 
 the solver library of \finley, see \Sec{
  \begin{memberdesc}[LinearPDE]{ITERATIVE}
  the default iterative method and preconditioner. The actually used method depends on the
  PDE solver library and the solver package been choosen. Typically, \PCG is used for symmetric PDEs
  and \BiCGStab otherwise, both with \JACOBI preconditioner.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{DIRECT}
  the default direct linear solver.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{CHOLEVSKY}
  direct solver based on Cholevsky factorization (or similar), see~\Ref{Saad}. The solver will require a symmetric PDE.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{PCG}
-Line 335 
 the GMRES method, see~\Ref{WEISS}\index{
+Line 338 
 the GMRES method, see~\Ref{WEISS}\index{
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{LUMPING}
  uses lumping to solve the system of linear equations~\index{linear solver!lumping}\index{lumping}. This solver technique
  condenses the stiffness matrix to a diagonal matrix so the solution of the linear systems becomes very cheap. It can be used when
  only \var{D} is present but in any case has to applied with care. The difference in the solutions with and without lumping can be significant
  but is expect to converge to zero when the mesh gets finer.
  Lumping does not use the linear system solver library.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{PRES20}
-Line 352 
 conjugate gradient squared method, see~\
+Line 355 
 conjugate gradient squared method, see~\
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{BICGSTAB}
  stabilized bi-conjugate gradients methods, see~\Ref{WEISS}.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{SSOR}
  symmetric successive over-relaxation method, see~\Ref{WEISS}. Typically used as preconditioner but some linear solver libraries support
  this as a solver.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{ILU0}
  the incomplete LU factorization preconditioner with no fill-in, see~\Ref{Saad}.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{ILUT}
  the incomplete LU factorization preconditioner with fill-in, see~\Ref{Saad}. During the  LU-factorization element with
  relative size less then \var{drop_tolerance} are dropped. Moreover, the size of the LU-factorization is restricted to the
  \var{drop_storage}-fold of the stiffness matrix. \var{drop_tolerance} and \var{drop_storage} are both set in the
  \method{getSolution} call.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{JACOBI}
-Line 375 
 the Jacobi preconditioner, see~\Ref{Saad
+Line 378 
 the Jacobi preconditioner, see~\Ref{Saad
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{AMG}
  the algebraic--multi grid method, see~\Ref{AMG}. This method can be used as linear solver method but is more robust when used
  in a preconditioner.
  \end{memberdesc}
-Line 383 
 in a preconditioner.
+Line 386 
 in a preconditioner.
  recursive incomplete LU factorization preconditioner, see~\Ref{RILU}. This method is similar to \ILUT but uses smoothing
  between levels. During the  LU-factorization element with
  relative size less then \var{drop_tolerance} are dropped. Moreover, the size of the LU-factorization is restricted to the
  \var{drop_storage}-fold of the stiffness matrix. \var{drop_tolerance} and \var{drop_storage} are both set in the
  \method{getSolution} call.
  \end{memberdesc}
  \begin{memberdesc}[LinearPDE]{NO_REORDERING}
-Line 416 
 and homogeneous constraints
+Line 419 
 and homogeneous constraints
  u=0 \mbox{ where } q>0
  \end{equation}
  $f$ has to be a \Scalar in the \Function and $q$ must be
  a \Scalar in  the \SolutionFS.
  \begin{classdesc}{Poisson}{domain}
  opens a Poisson equation on the \Domain domain. \Poisson is derived from \LinearPDE.
-Line 426 
 assigns new values to \var{f} and \var{q
+Line 429 
 assigns new values to \var{f} and \var{q
  \end{methoddesc}
  \section{The \Helmholtz Class}
+ The \Helmholtz class defines the Helmholtz problem
+ \begin{equation}\label{HZ.1}
+ \omega \; u - (k\; u\hackscore{,j})\hackscore{,j} = f
+ \end{equation}
+  with natural boundary conditons
+ \begin{equation}\label{HZ.2}
+ k\; u\hackscore{,j} n\hackscore{,j} = g- \alpha \; u
+ \end{equation}
+ and constraints:
+ \begin{equation}\label{HZ.3}
+ u=r \mbox{ where } q>0
+ \end{equation}
+ $\omega$, $k$, $f$ have to be a \Scalar in the \Function,
+ $g$ and $\alpha$ must be a \Scalar in  the \FunctionOnBoundary,
+ and $q$ and $r$ must be a \Scalar in  the \SolutionFS or must be mapped or interpolated into the particular \FunctionSpace.
+ \begin{classdesc}{Helmholtz}{domain}
+ opens a Helmholtz equation on the \Domain domain. \Helmholtz is derived from \LinearPDE.
+ \end{classdesc}
+ \begin{methoddesc}[Helmholtz]{setValue}{ \optional{omega} \optional{, k} \optional{, f} \optional{, alpha} \optional{, g} \optional{, r} \optional{, q}}
+ assigns new values to \var{omega}, \var{k}, \var{f}, \var{alpha}, \var{g}, \var{r}, \var{q}. By default all values are set to be zero.
+ \end{methoddesc}
  \section{The \Lame Class}
+ The \Lame class defines a Lame equation problem:
+ \begin{equation}\label{LE.1}
+ -\mu (u\hackscore{i,j}+u\hackscore{j,i})+\lambda u\hackscore{k,k})\hackscore{j} = F\hackscore{i}-\sigma\hackscore{ij,j}
+ \end{equation}
+ with natural boundary conditons:
+ \begin{equation}\label{LE.2}
+ n\hackscore{j}(\mu \; (u\hackscore{i,j}+u\hackscore{j,i})+\lambda*u\hackscore{k,k}) = f\hackscore{i}+n\hackscore{j}\sigma\hackscore{ij}
+ \end{equation}
+ and constraint
+ \begin{equation}\label{LE.3}
+ u\hackscore{i}=r\hackscore{i} \mbox{ where } q\hackscore{i}>0
+ \end{equation}
+ $\mu$, $\lambda$ have to be a \Scalar in the \Function,
+ $F$ has to be a \Vector in the \Function,
+ $\sigma$ has to be a \Tensor in the \Function,
+ $f$ must be a \Vector in  the \FunctionOnBoundary,
+ and $q$ and $r$ must be a \Vector in  the \SolutionFS or must be mapped or interpolated into the particular \FunctionSpace.
+ \begin{classdesc}{Lame}{domain}
+ opens a Lame equation on the \Domain domain. \Lame is derived from \LinearPDE.
+ \end{classdesc}
+ \begin{methoddesc}[Lame]{setValue}{ \optional{lame_lambda} \optional{, lame_mu} \optional{, F} \optional{, sigma} \optional{, f} \optional{, r} \optional{, q}}
+ assigns new values to
+ \var{lame_lambda},
+ \var{lame_mu},
+ \var{F},
+ \var{sigma},
+ \var{f},
+ \var{r} and
+ \var{q}
+ By default all values are set to be zero.
+ \end{methoddesc}

 Legend:



Removed from v.659
 


changed lines


 
Added in v.660
 Legend:



Removed from v.659
 


changed lines


 
Added in v.660
-Removed from v.659
+Added in v.660

	ViewVC Help
Powered by ViewVC 1.1.26