1 |
|
2 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
3 |
% Copyright (c) 2003-2018 by The University of Queensland |
4 |
% http://www.uq.edu.au |
5 |
% |
6 |
% Primary Business: Queensland, Australia |
7 |
% Licensed under the Apache License, version 2.0 |
8 |
% http://www.apache.org/licenses/LICENSE-2.0 |
9 |
% |
10 |
% Development until 2012 by Earth Systems Science Computational Center (ESSCC) |
11 |
% Development 2012-2013 by School of Earth Sciences |
12 |
% Development from 2014 by Centre for Geoscience Computing (GeoComp) |
13 |
% |
14 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
15 |
|
16 |
\chapter{The \linearPDEs Module} |
17 |
|
18 |
\section{Linear Partial Differential Equations} |
19 |
\label{SEC LinearPDE} |
20 |
|
21 |
The \LinearPDE class is used to define a general linear, steady, second order |
22 |
PDE for an unknown function $u$ on a given $\Omega$ defined through a \Domain object. |
23 |
In the following $\Gamma$ denotes the boundary of the domain $\Omega$ and $n$ |
24 |
denotes the outer normal field on $\Gamma$. |
25 |
|
26 |
For a single PDE with a solution that has a single component the linear PDE is |
27 |
defined in the following form \footnote{This PDE system can be written in the equivalent form, \begin{align*} \centering \, -\nabla \cdot \left( \textbf{A} \nabla u \right) - \nabla \cdot \left( \textbf{B} \, u \right) + \textbf{C} \cdot \nabla u + \textbf{D} \, u + \sum_p \delta(p) u(p) &= - \nabla \cdot \textbf{X} + \textbf{Y} + \sum_p \delta(p), \\ |
28 |
\hat{\textbf{n}} \cdot \left( \textbf{A} \nabla u + \textbf{B} \, u \right) + \textbf{d} \, u &= \hat{\textbf{n}} \cdot \textbf{X} + \textbf{y}, |
29 |
\end{align*} where bold font denotes a data object and $\hat{\textbf{n}}$ denotes a unit vector normal. The weak form for this is |
30 |
\begin{equation*}\int_{\Omega} (\nabla v)\cdot \left( \textbf{A} \nabla u \right) +( \nabla v) \cdot \left( \textbf{B} \, u \right) + v\textbf{C} \cdot \nabla u + \textbf{D} \,v u + \sum_p \delta(p) u(p)v \;d\Omega = \int_{\Omega}(\nabla v) \cdot \textbf{X} + \textbf{Y}v + \sum_p \delta(p)v\; d\Omega. |
31 |
\end{equation*} |
32 |
} |
33 |
\begin{equation}\label{LINEARPDE.SINGLE.1} |
34 |
-(A_{jl} u_{,l})_{,j}-(B_{j} u)_{,j}+C_{l} u_{,l}+D u + \sum_p d^{dirac}(p) \; u(p) =-X_{j,j}+Y + \sum_p y^{dirac}(p) \; . |
35 |
\end{equation} |
36 |
$u_{,j}$ denotes the derivative of $u$ with respect to the $j$-th spatial direction. |
37 |
Einstein's summation convention, i.e. summation over indexes appearing twice |
38 |
in a term of a sum, is used in this chapter. $y^{dirac}(p)$ represent a nodal source term |
39 |
at point $p$, cf. $y^{dirac}(p) $ and similar $d^{dirac}(p)$ define Dirac delta-function |
40 |
terms. |
41 |
The coefficients $A$, $B$, $C$, $D$, $X$ and $Y$ have to be specified through |
42 |
\Data objects in the \Function on the PDE or objects that can be converted |
43 |
into such \Data objects. $d^{dirac}$ |
44 |
$A$ is a \RankTwo, $B$, $C$ and $X$ are each a \RankOne and $D$ and $Y$ are |
45 |
scalars. $y^{dirac}$ and $d^{dirac}$ are each scalars in the \DiracDeltaFunctions. |
46 |
The following natural boundary conditions are considered\index{boundary condition!natural} on $\Gamma$: |
47 |
\begin{equation}\label{LINEARPDE.SINGLE.2} |
48 |
n_{j}(A_{jl} u_{,l}+B_{j} u)+d u=n_{j}X_{j} + y \;. |
49 |
\end{equation} |
50 |
Notice that the coefficients $A$, $B$ and $X$ are defined in the PDE. |
51 |
The coefficients $d$ and $y$ are each a \Scalar in the \FunctionOnBoundary. |
52 |
Constraints\index{constraint} for the solution prescribe the value of the |
53 |
solution at certain locations in the domain. They have the form |
54 |
\begin{equation}\label{LINEARPDE.SINGLE.3} |
55 |
u=r \mbox{ where } q>0 |
56 |
\end{equation} |
57 |
$r$ and $q$ are each a \Scalar where $q$ is the characteristic function\index{characteristic function} |
58 |
defining where the constraint is applied. |
59 |
The constraints defined by \eqn{LINEARPDE.SINGLE.3} override any other |
60 |
condition set by \eqn{LINEARPDE.SINGLE.1} or \eqn{LINEARPDE.SINGLE.2}. |
61 |
|
62 |
For a system of PDEs and a solution with several components the PDE has the form |
63 |
\begin{equation}\label{LINEARPDE.SYSTEM.1} |
64 |
-(A_{ijkl} u_{k,l})_{,j}-(B_{ijk} u_{k})_{,j}+C_{ikl} u_{k,l}+D_{ik} u_{k}+ |
65 |
\sum_p d^{dirac}_{ik}(p) \; u_i(p) |
66 |
=-X_{ij,j}+Y_{i} +\sum_p \; y^{dirac}_i(p) \; . |
67 |
\end{equation} |
68 |
$A$ is a \RankFour, $B$ and $C$ are each a \RankThree, $D$, $d^{dirac}$ and $X$ are each a \RankTwo and $Y$ and $y^{dirac}$ is a \RankOne. |
69 |
The natural boundary conditions\index{boundary condition!natural} take the form: |
70 |
\begin{equation}\label{LINEARPDE.SYSTEM.2} |
71 |
n_{j}(A_{ijkl} u_{k,l}+B_{ijk} u_{k})+d_{ik} u_{k}=n_{j}X_{ij}+y_{i} \;. |
72 |
\end{equation} |
73 |
The coefficient $d$ is a \RankTwo and $y$ is a \RankOne both in the |
74 |
\FunctionOnBoundary. Constraints\index{constraint} take the form |
75 |
\begin{equation}\label{LINEARPDE.SYSTEM.3} |
76 |
u_{i}=r_{i} \mbox{ where } q_{i}>0 |
77 |
\end{equation} |
78 |
$r$ and $q$ are each a \RankOne. Notice that not necessarily all components |
79 |
must have a constraint at all locations. An example for a system of PDEs is shown in the elastic deformation example in section \Refe{ELASTIC CHAP} |
80 |
|
81 |
\LinearPDE also supports solution discontinuities\index{discontinuity} over a |
82 |
contact region $\Gamma^{contact}$ in the domain $\Omega$. |
83 |
To specify the conditions across the discontinuity we are using the |
84 |
generalised flux $J$\footnote{In some applications the definition of flux used |
85 |
here can be different from the commonly used definition. |
86 |
For instance, if $T$ is a temperature field the heat flux $q$ is defined as |
87 |
$q_{,i}=-\kappa T_{,i}$ ($\kappa$ is the diffusivity) which differs from the |
88 |
definition used here by the sign. This needs to be kept in mind when defining |
89 |
natural boundary conditions\index{boundary condition!natural}.} which in |
90 |
the case of a system of PDEs and several components of the solution, is |
91 |
defined as |
92 |
\begin{equation}\label{LINEARPDE.SYSTEM.5} |
93 |
J_{ij}=A_{ijkl}u_{k,l}+B_{ijk}u_{k}-X_{ij} |
94 |
\end{equation} |
95 |
For the case of single solution component and single PDE, $J$ is defined as |
96 |
\begin{equation}\label{LINEARPDE.SINGLE.5} |
97 |
J_{j}=A_{jl}u_{,l}+B_{j}u_{k}-X_{j} |
98 |
\end{equation} |
99 |
In the context of discontinuities\index{discontinuity} $n$ denotes the normal |
100 |
on the discontinuity pointing from side 0 towards side 1. |
101 |
For a system of PDEs the contact condition takes the form |
102 |
\begin{equation}\label{LINEARPDE.SYSTEM.6} |
103 |
n_{j} J^{0}_{ij}=n_{j} J^{1}_{ij}=y^{contact}_{i} - d^{contact}_{ik} [u]_{k} \; . |
104 |
\end{equation} |
105 |
where $J^{0}$ and $J^{1}$ are the fluxes on side $0$ and side $1$ of the |
106 |
discontinuity $\Gamma^{contact}$, respectively. $[u]$, which is the difference |
107 |
of the solution at side 1 and at side 0, denotes the jump of $u$ across $\Gamma^{contact}$. |
108 |
The coefficient $d^{contact}$ is a \RankTwo and $y^{contact}$ is a |
109 |
\RankOne both in the \FunctionOnContactZero or \FunctionOnContactOne. |
110 |
In the case of a single PDE and a single component solution the contact |
111 |
condition takes the form |
112 |
\begin{equation}\label{LINEARPDE.SINGLE.6} |
113 |
n_{j} J^{0}_{j}=n_{j} J^{1}_{j}=y^{contact} - d^{contact}[u] |
114 |
\end{equation} |
115 |
In this case the coefficient $d^{contact}$ and $y^{contact}$ are each |
116 |
\Scalars both in either\\ |
117 |
the \FunctionOnContactZero or the \FunctionOnContactOne. |
118 |
|
119 |
The PDE is symmetrical\index{symmetrical} if |
120 |
\begin{equation}\label{LINEARPDE.SINGLE.4} |
121 |
A_{jl}=A_{lj} \mbox{ and } B_{j}=C_{j} |
122 |
\end{equation} |
123 |
The system of PDEs is symmetrical\index{symmetrical} if |
124 |
\begin{eqnarray} |
125 |
\label{LINEARPDE.SYSTEM.4} |
126 |
A_{ijkl}&=&A_{klij} \\ |
127 |
B_{ijk}&=&C_{kij} \\ |
128 |
D_{ik}&=&D_{ki} \\ |
129 |
d_{ik}&=&d_{ki} \\ |
130 |
d^{contact}_{ik}&=&d^{contact}_{ki} |
131 |
\end{eqnarray} |
132 |
Note that in contrast to the scalar case \eqn{LINEARPDE.SINGLE.4} now the |
133 |
coefficients $D$, $d$ and $d^{contact}$ have to be inspected. |
134 |
|
135 |
The following example illustrates a typical usage of the \LinearPDE class: |
136 |
\begin{python} |
137 |
from esys.escript import * |
138 |
from esys.escript.linearPDEs import LinearPDE |
139 |
from esys.finley import Rectangle |
140 |
mydomain = Rectangle(l0=1., l1=1., n0=40, n1=20) |
141 |
mypde=LinearPDE(mydomain) |
142 |
mypde.setSymmetryOn() |
143 |
mypde.setValue(A=kappa*kronecker(mydomain), D=1, Y=1) |
144 |
u=mypde.getSolution() |
145 |
\end{python} |
146 |
We refer to \Chap{CHAP: Tutorial} for more details. |
147 |
|
148 |
An instance of the \SolverOptions class is attached to the \LinearPDE class |
149 |
object. It holds options for the solver that may be set before solving the PDE. |
150 |
In the following example the \method{getSolverOptions} method is used to |
151 |
access the \SolverOptions object attached to \var{mypde}: |
152 |
\begin{python} |
153 |
from esys.escript import * |
154 |
from esys.escript.linearPDEs import LinearPDE, SolverOptions |
155 |
from esys.finley import Rectangle |
156 |
mydomain = Rectangle(l0=1., l1=1., n0=40, n1=20) |
157 |
mypde=LinearPDE(mydomain) |
158 |
mypde.setValue(A=kappa*kronecker(mydomain), D=1, Y=1) |
159 |
mypde.getSolverOptions().setVerbosityOn() |
160 |
mypde.getSolverOptions().setSolverMethod(SolverOptions.PCG) |
161 |
mypde.getSolverOptions().setPreconditioner(SolverOptions.AMG) |
162 |
mypde.getSolverOptions().setTolerance(1e-8) |
163 |
mypde.getSolverOptions().setIterMax(1000) |
164 |
u=mypde.getSolution() |
165 |
\end{python} |
166 |
In this example, the preconditioned conjugate gradient method \PCG is used |
167 |
with preconditioner \AMG. The relative tolerance is set to $10^{-8}$ and |
168 |
the maximum number of iteration steps to $1000$. |
169 |
After a completed call to \method{getSolution()}, the attached \SolverOptions |
170 |
object gives access to diagnostic information: |
171 |
\begin{python} |
172 |
u=mypde.getSolution() |
173 |
print("Number of iteration steps =", mypde.getDiagnostics("num_iter")) |
174 |
print("Total solution time =", mypde.getDiagnostics("time")) |
175 |
print("Set-up time =", mypde.getDiagnostics("set_up_time")) |
176 |
print("Net time =", mypde.getDiagnostics("net_time")) |
177 |
print("Residual norm of returned solution =", |
178 |
mypde.getDiagnostics('residual_norm')) |
179 |
\end{python} |
180 |
Typically, a negative value for a diagnostic variable indicates that it is |
181 |
undefined. For more details on \SolverOptions and Trilinos see chapter \Refe{TRILINOS}. |
182 |
|
183 |
\subsection{Classes} |
184 |
%\declaremodule{extension}{esys.escript.linearPDEs} |
185 |
%\modulesynopsis{Linear partial differential equation handler} |
186 |
The module \linearPDEs provides an interface to define and solve linear partial |
187 |
differential equations within \escript. The module \linearPDEs does not |
188 |
provide any solver capabilities in itself but hands the PDE over to the PDE |
189 |
solver library defined through the \Domain of the PDE, e.g. \finley. |
190 |
The general interface is provided through the \LinearPDE class. The \Poisson |
191 |
class which is also derived form the \LinearPDE class can be used to define |
192 |
the Poisson equation\index{Poisson}. |
193 |
|
194 |
\subsection{\LinearPDE class} |
195 |
This is the general class to define a linear PDE in \escript. |
196 |
We list a selection of the most important methods of the class. |
197 |
For a complete list, see the reference at \ReferenceGuide. |
198 |
|
199 |
\begin{classdesc}{LinearPDE}{domain,numEquations=0,numSolutions=0} |
200 |
opens a linear, steady, second order PDE on the \Domain \var{domain}. |
201 |
The parameters \var{numEquations} and \var{numSolutions} give the number of |
202 |
equations and the number of solution components. |
203 |
If \var{numEquations} and \var{numSolutions} are non-positive, then the number |
204 |
of equations and the number of solutions, respectively, stay undefined until a |
205 |
coefficient is defined. |
206 |
\end{classdesc} |
207 |
|
208 |
\subsubsection{\LinearPDE methods} |
209 |
\begin{methoddesc}[LinearPDE]{setValue}{ |
210 |
\optional{A}\optional{, B}, |
211 |
\optional{, C}\optional{, D} |
212 |
\optional{, X}\optional{, Y} |
213 |
\optional{, d}\optional{, y} |
214 |
\optional{, d_contact}\optional{, y_contact} |
215 |
\optional{, d_dirac}\optional{, y_dirac} |
216 |
\optional{, q}\optional{, r}, } |
217 |
assigns new values to coefficients. By default all values are assumed to be |
218 |
zero\footnote{In fact, it is assumed they are not present by assigning the |
219 |
value \code{escript.Data()}. This can be used by the solver library to reduce |
220 |
computational costs.}. |
221 |
If the new coefficient value is not a \Data object, it is converted into a |
222 |
\Data object in the appropriate \FunctionSpace. |
223 |
\end{methoddesc} |
224 |
|
225 |
\begin{methoddesc}[LinearPDE]{getCoefficient}{name} |
226 |
returns the value assigned to coefficient \var{name}. If \var{name} is not a |
227 |
valid name an exception is raised. |
228 |
\end{methoddesc} |
229 |
|
230 |
\begin{methoddesc}[LinearPDE]{getShapeOfCoefficient}{name} |
231 |
returns the shape of the coefficient \var{name} even if no value has been |
232 |
assigned to it. |
233 |
\end{methoddesc} |
234 |
|
235 |
\begin{methoddesc}[LinearPDE]{getFunctionSpaceForCoefficient}{name} |
236 |
returns the \FunctionSpace of the coefficient \var{name} even if no value has |
237 |
been assigned to it. |
238 |
\end{methoddesc} |
239 |
|
240 |
\begin{methoddesc}[LinearPDE]{setDebugOn}{} |
241 |
switches on debug mode so more diagnostic messages will be printed. |
242 |
\end{methoddesc} |
243 |
|
244 |
\begin{methoddesc}[LinearPDE]{setDebugOff}{} |
245 |
switches off debug mode. |
246 |
\end{methoddesc} |
247 |
|
248 |
\begin{methoddesc}[LinearPDE]{getSolverOptions}{} |
249 |
returns the solver options for solving the PDE. In fact, the method returns |
250 |
a \SolverOptions class object which can be used to modify the tolerance, |
251 |
the solver or the preconditioner, see \Sec{SEC Solver Options} for details. |
252 |
\end{methoddesc} |
253 |
|
254 |
\begin{methoddesc}[LinearPDE]{setSolverOptions}{\optional{options=None}} |
255 |
sets the solver options for solving the PDE. If argument \var{options} is |
256 |
present it must be a \SolverOptions class object, see \Sec{SEC Solver Options} |
257 |
for details. Otherwise the solver options are reset to the default. |
258 |
\end{methoddesc} |
259 |
|
260 |
\begin{methoddesc}[LinearPDE]{isUsingLumping}{} |
261 |
returns \True if matrix lumping is set as the solver for the system of linear |
262 |
equations, \False otherwise. |
263 |
\end{methoddesc} |
264 |
|
265 |
\begin{methoddesc}[LinearPDE]{getDomain}{} |
266 |
returns the \Domain of the PDE. |
267 |
\end{methoddesc} |
268 |
|
269 |
\begin{methoddesc}[LinearPDE]{getDim}{} |
270 |
returns the number of spatial dimensions of the PDE. |
271 |
\end{methoddesc} |
272 |
|
273 |
\begin{methoddesc}[LinearPDE]{getNumEquations}{} |
274 |
returns the number of equations. |
275 |
\end{methoddesc} |
276 |
|
277 |
\begin{methoddesc}[LinearPDE]{getNumSolutions}{} |
278 |
returns the number of components of the solution. |
279 |
\end{methoddesc} |
280 |
|
281 |
\begin{methoddesc}[LinearPDE]{checkSymmetry}{verbose=\False} |
282 |
returns \True if the PDE is symmetric, \False otherwise. |
283 |
The method is very computationally expensive and should only be called for |
284 |
testing purposes. The symmetry flag is not altered. |
285 |
If \var{verbose=True} information about where symmetry is violated is printed. |
286 |
\end{methoddesc} |
287 |
|
288 |
\begin{methoddesc}[LinearPDE]{getFlux}{u} |
289 |
returns the flux $J_{ij}$\index{flux} for given solution \var{u} defined by |
290 |
\eqn{LINEARPDE.SYSTEM.5} and \eqn{LINEARPDE.SINGLE.5}. |
291 |
\end{methoddesc} |
292 |
|
293 |
\begin{methoddesc}[LinearPDE]{isSymmetric}{} |
294 |
returns \True if the PDE has been indicated to be symmetric, \False otherwise. |
295 |
\end{methoddesc} |
296 |
|
297 |
\begin{methoddesc}[LinearPDE]{setSymmetryOn}{} |
298 |
indicates that the PDE is symmetric which enables the use of certain solvers |
299 |
and can potentially speed up the solver. |
300 |
\end{methoddesc} |
301 |
|
302 |
\begin{methoddesc}[LinearPDE]{setSymmetryOff}{} |
303 |
indicates that the PDE is not symmetric. |
304 |
\end{methoddesc} |
305 |
|
306 |
\begin{methoddesc}[LinearPDE]{setReducedOrderOn}{} |
307 |
enables the reduction of polynomial order for the solution and equation |
308 |
evaluation even if a quadratic or higher interpolation order is defined in the |
309 |
\Domain. This feature may not be supported by all PDE libraries. |
310 |
\end{methoddesc} |
311 |
|
312 |
\begin{methoddesc}[LinearPDE]{setReducedOrderOff}{} |
313 |
disables the reduction of polynomial order for the solution and equation evaluation. |
314 |
\end{methoddesc} |
315 |
|
316 |
\begin{methoddesc}[LinearPDE]{getOperator}{} |
317 |
returns the \Operator of the PDE. |
318 |
\end{methoddesc} |
319 |
|
320 |
\begin{methoddesc}[LinearPDE]{getRightHandSide}{} |
321 |
returns the right hand side of the PDE as a \Data object. |
322 |
\end{methoddesc} |
323 |
|
324 |
\begin{methoddesc}[LinearPDE]{getSystem}{} |
325 |
returns the \Operator and right hand side of the PDE as a tuple. |
326 |
\end{methoddesc} |
327 |
|
328 |
\begin{methoddesc}[LinearPDE]{getSolution}{} |
329 |
returns (an approximation of) the solution of the PDE. This call will invoke |
330 |
the discretization of the PDE and the solution of the resulting system of |
331 |
linear equations. Keep in mind that this call is typically computationally |
332 |
expensive and -- depending on the PDE and the discretization -- can take a |
333 |
long time to complete. |
334 |
\end{methoddesc} |
335 |
|
336 |
\subsection{The \Poisson Class} |
337 |
The \Poisson class provides an easy way to define and solve the Poisson |
338 |
equation |
339 |
\begin{equation}\label{POISSON.1} |
340 |
-u_{,ii}=f |
341 |
\end{equation} |
342 |
with homogeneous boundary conditions |
343 |
\begin{equation}\label{POISSON.2} |
344 |
n_{i}u_{,i}=0 |
345 |
\end{equation} |
346 |
and homogeneous constraints |
347 |
\begin{equation}\label{POISSON.3} |
348 |
u=0 \mbox{ where } q>0 . |
349 |
\end{equation} |
350 |
$f$ has to be a \Scalar in the \Function and $q$ must be a \Scalar in the \SolutionFS. |
351 |
|
352 |
\begin{classdesc}{Poisson}{domain} |
353 |
opens a Poisson equation on the \Domain domain. \Poisson is derived from \LinearPDE. |
354 |
\end{classdesc} |
355 |
\begin{methoddesc}[Poisson]{setValue}{f=escript.Data(),q=escript.Data()} |
356 |
assigns new values to \var{f} and \var{q}. |
357 |
\end{methoddesc} |
358 |
|
359 |
\subsection{The \Helmholtz Class} |
360 |
The \Helmholtz class defines the Helmholtz problem |
361 |
\begin{equation}\label{HZ.1} |
362 |
\omega \; u - (k\; u_{,j})_{,j} = f |
363 |
\end{equation} |
364 |
with natural boundary conditions |
365 |
\begin{equation}\label{HZ.2} |
366 |
k\; u_{,j} n_{,j} = g- \alpha \; u |
367 |
\end{equation} |
368 |
and constraints |
369 |
\begin{equation}\label{HZ.3} |
370 |
u=r \mbox{ where } q>0 . |
371 |
\end{equation} |
372 |
$\omega$, $k$, and $f$ each have to be a \Scalar in the \Function, $g$ and |
373 |
$\alpha$ must be a \Scalar in the \FunctionOnBoundary, and $q$ and $r$ must be |
374 |
a \Scalar in the \SolutionFS or must be mapped or interpolated into the |
375 |
particular \FunctionSpace. |
376 |
|
377 |
\begin{classdesc}{Helmholtz}{domain} |
378 |
opens a Helmholtz equation on the \Domain domain. \Helmholtz is derived from \LinearPDE. |
379 |
\end{classdesc} |
380 |
\begin{methoddesc}[Helmholtz]{setValue}{ \optional{omega} \optional{, k} \optional{, f} \optional{, alpha} \optional{, g} \optional{, r} \optional{, q}} |
381 |
assigns new values to \var{omega}, \var{k}, \var{f}, \var{alpha}, \var{g}, |
382 |
\var{r}, and \var{q}. By default all values are set to zero. |
383 |
\end{methoddesc} |
384 |
|
385 |
\subsection{The \Lame Class} |
386 |
The \Lame class defines a Lam\'e equation problem |
387 |
\begin{equation}\label{LE.1} |
388 |
-(\mu (u_{i,j}+u_{j,i})+\lambda u_{k,k}\delta_{ij})_{j} = F_{i}-\sigma_{ij,j} |
389 |
\end{equation} |
390 |
with natural boundary conditions |
391 |
\begin{equation}\label{LE.2} |
392 |
n_{j}(\mu \; (u_{i,j}+u_{j,i})+\lambda u_{k,k}\delta_{ij}) = f_{i}+n_{j}\sigma_{ij} |
393 |
\end{equation} |
394 |
and constraint |
395 |
\begin{equation}\label{LE.3} |
396 |
u_{i}=r_{i} \mbox{ where } q_{i}>0 . |
397 |
\end{equation} |
398 |
$\mu$, $\lambda$ have to be a \Scalar in the \Function, $F$ has to be a |
399 |
\Vector in the \Function, $\sigma$ has to be a \Tensor in the \Function, |
400 |
$f$ must be a \Vector in the \FunctionOnBoundary, and $q$ and $r$ must be a |
401 |
\Vector in the \SolutionFS or must be mapped or interpolated into the |
402 |
particular \FunctionSpace. |
403 |
|
404 |
\begin{classdesc}{Lame}{domain} |
405 |
opens a Lam\'e equation on the \Domain domain. \Lame is derived from \LinearPDE. |
406 |
\end{classdesc} |
407 |
\begin{methoddesc}[Lame]{setValue}{ \optional{lame_lambda} \optional{, lame_mu} \optional{, F} \optional{, sigma} \optional{, f} \optional{, r} \optional{, q}} |
408 |
assigns new values to \var{lame_lambda}, \var{lame_mu}, \var{F}, \var{sigma}, |
409 |
\var{f}, \var{r}, and \var{q}. By default all values are set to zero. |
410 |
\end{methoddesc} |
411 |
|
412 |
\section{Projection} |
413 |
%\declaremodule{extension}{esys.escript.pdetools} |
414 |
\label{SEC Projection} |
415 |
|
416 |
Using the \LinearPDE class provides an option to change the \FunctionSpace |
417 |
attribute in addition to the standard interpolation mechanism\index{interpolation} |
418 |
as discussed in \Chap{ESCRIPT CHAP}. If you consider the stripped-down version |
419 |
\begin{equation}\label{PROJ.1} |
420 |
u = Y |
421 |
\end{equation} |
422 |
of the general scalar PDE~\Refe{LINEARPDE.SINGLE.1} (boundary conditions are |
423 |
irrelevant), you can see the solution $u$ of this PDE as a projection of the |
424 |
input function $Y$ which has the \Function attribute to a function with the |
425 |
\SolutionFS or \ReducedSolutionFS attribute. |
426 |
In fact, the solution maps values defined at element centers representing a |
427 |
possibly discontinuous function onto a continuous function represented by its |
428 |
values at the nodes of the FEM mesh. |
429 |
This mapping is called a projection\index{projection}. Projection can be a |
430 |
useful tool but needs to be applied with some care due to the possibility of |
431 |
projecting a potentially discontinuous function onto a continuous function, |
432 |
although this may also be a desirable effect, for instance to smooth a function. |
433 |
The projection of the gradient of a function typically calculated on the |
434 |
element center to the nodes of a FEM mesh can be evaluated on the domain |
435 |
boundary and so projection provides a tool to extrapolate the gradient from |
436 |
the internal to the boundary. This is only a reasonable procedure in the |
437 |
absence of singularities at the boundary. |
438 |
|
439 |
As projection is often used in simulations \escript provides an easy to use |
440 |
class \class{Projector} which is part of the \pdetools module. |
441 |
The following script demonstrates the usage of the class to project the |
442 |
piecewise constant function ($=1$ for $x_{0}\ge 0.5$ and $=-1$ for $x_{0}<0.5$) |
443 |
to a function with the \ReducedSolutionFS attribute (default target): |
444 |
\begin{python} |
445 |
from esys.escript.pdetools import Projector |
446 |
proj=Projector(domain) |
447 |
x0=domain.getX()[0] |
448 |
jmp=1.-2.*wherePositive(x0-0.5) |
449 |
u=proj.getValue(jmp) |
450 |
# alternative call: |
451 |
u=proj(jmp) |
452 |
\end{python} |
453 |
By default the class uses lumping to solve the PDE~\Refe{PROJ.1}. |
454 |
This technique is faster than using the standard solver techniques of PDEs. |
455 |
In essence it leads to using the average of neighbour element values to |
456 |
calculate the value at each FEM node. |
457 |
|
458 |
The following script illustrates how to evaluate the normal stress on the |
459 |
boundary from a given displacement field \var{u}: |
460 |
\begin{python} |
461 |
from esys.escript.pdetools import Projector |
462 |
u=... |
463 |
proj=Projector(u.getDomain()) |
464 |
e=symmetric(grad(u)) |
465 |
stress = G*e+ (K-2./3.*G)*trace(e)*kronecker(u.getDomain()) |
466 |
normal_stress = inner(u.getDomain().getNormal(), proj(stress)) |
467 |
\end{python} |
468 |
|
469 |
\begin{classdesc}{Projector}{domain\optional{, reduce=\True \optional{, fast=\True}}} |
470 |
This class defines a projector on the domain \var{domain}. |
471 |
If \var{reduce} is set to \True the projection will be returned as a |
472 |
\ReducedSolutionFS \Data object. |
473 |
Otherwise the \SolutionFS representation is returned. |
474 |
If \var{reduce} is set to \True lumping is used when the \eqn{PROJ.1} is |
475 |
solved, otherwise the standard PDE solver is used. |
476 |
Notice, that lumping requires significantly less computation time and memory. |
477 |
The class is callable. |
478 |
\end{classdesc} |
479 |
|
480 |
\begin{methoddesc}[Projector]{getSolverOptions}{} |
481 |
returns the solver options for solving the PDE. In fact, the method returns |
482 |
a \SolverOptions class object which can be used to modify the tolerance, |
483 |
the solver or the preconditioner, see \Sec{SEC Solver Options} for details. |
484 |
\end{methoddesc} |
485 |
|
486 |
\begin{methoddesc}[Projector]{getValue}{input_data} |
487 |
projects the \var{input_data}. This method is equivalent to call an instance |
488 |
of the class with argument \var{input_data} |
489 |
\end{methoddesc} |
490 |
|
491 |
% \section{Transport Problems} |
492 |
% \label{SEC Transport} |
493 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
494 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
495 |
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
496 |
\section{Solver Options} |
497 |
\label{SEC Solver Options} |
498 |
|
499 |
\begin{classdesc}{SolverOptions}{} |
500 |
This class defines the solver options for a linear or non-linear solver. |
501 |
The option also supports the handling of diagnostic information. |
502 |
\end{classdesc} |
503 |
|
504 |
\begin{methoddesc}[SolverOptions]{getSummary}{} |
505 |
returns a string reporting the current settings. |
506 |
\end{methoddesc} |
507 |
|
508 |
\begin{methoddesc}[SolverOptions]{getName}{key} |
509 |
returns the name as a string of a given key. |
510 |
\end{methoddesc} |
511 |
|
512 |
\begin{methoddesc}[SolverOptions]{setSolverMethod}{method} |
513 |
sets the solver method to be used. |
514 |
Use \var{method}=\member{SolverOptions.DIRECT} to indicate that a direct |
515 |
rather than an iterative solver should be used and use |
516 |
\var{method}=\member{SolverOptions.ITERATIVE} to indicate that an iterative |
517 |
rather than a direct solver should be used. Note that SolverOptions needs to be |
518 |
imported from linearPDEs and is not the same as the object returned by pde.getSolverOptions(). |
519 |
The value of \var{method} must be one of the constants:\\ |
520 |
\member{SolverOptions.DEFAULT} -- use default solver depending on other options\\ |
521 |
\member{SolverOptions.BICGSTAB} -- Biconjugate Gradient Stabilized iterative method\\ |
522 |
\member{SolverOptions.CGLS} -- Conjugate Gradient with Least Squares method\\ |
523 |
\member{SolverOptions.CGS} -- Conjugate Gradient Square method\\ |
524 |
\member{SolverOptions.CHOLEVSKY} -- Direct solver based on LDLT factorization\\ |
525 |
\member{SolverOptions.CR} -- Conjugate Residual method\\ |
526 |
\member{SolverOptions.DIRECT} -- use a direct solver if available\\ |
527 |
%\member{Solver.Options.DIRECT_MUMPS}--runtime error\\ |
528 |
%\member{SolverOptions.DIRECT_PARDISO}--runtime error\\ |
529 |
\member{SolverOptions.DIRECT_SUPERLU} -- use a direct SUPERLU solver\\ |
530 |
\member{SolverOptions.DIRECT_TRILINOS} -- use default TRILINOS default solver (KLU2) (chapter \Refe{TRILINOS})\\ |
531 |
\member{SolverOptions.GMRES} -- Gram-Schmidt minimum residual method\\ |
532 |
\member{SolverOptions.HRZ_LUMPING} -- Matrix lumping using the HRZ approach (section \Refe{WAVE CHAP})\\ |
533 |
\member{SolverOptions.ITERATIVE} -- use a suitable iterative solver\\ |
534 |
\member{SolverOptions.LSQR} -- Least squares based LSQR solver\\ |
535 |
\member{SolverOptions.LUMPING} -- Matrix lumping\\ |
536 |
\member{SolverOptions.MINRES} -- Minimum Residual method\\ |
537 |
\member{SolverOptions.NONLINEAR_GMRES} -- restarted GMRES for nonlinear systems\\ |
538 |
\member{SolverOptions.PCG} -- Preconditioned Conjugate Gradient method\\ |
539 |
\member{SolverOptions.PRES20} -- GMRES with restart after 20 steps and truncations after 5 residuals\\ |
540 |
\member{SolverOptions.ROWSUM_LUMPING} -- Matrix lumping using row sum\\ |
541 |
\member{SolverOptions.TFQMR} -- Transpose Free Quasi Minimum Residual method.\\ |
542 |
Not all packages support all solvers. It can be assumed that a package makes a |
543 |
reasonable choice if it encounters an unknown solver. |
544 |
See Table~\Refe{TAB FINLEY SOLVER OPTIONS 1} for the solvers supported by |
545 |
\finley. |
546 |
\end{methoddesc} |
547 |
|
548 |
\begin{methoddesc}[SolverOptions]{getSolverMethod}{} |
549 |
returns the key of the solver method to be used. |
550 |
\end{methoddesc} |
551 |
|
552 |
\begin{methoddesc}[SolverOptions]{setPreconditioner}{preconditioner} |
553 |
sets the preconditioner to be used. |
554 |
The value of \var{preconditioner} must be one of the constants:\\ |
555 |
\member{SolverOptions.AMG} -- Algebraic Multi Grid\\ |
556 |
%\member{SolverOptions.AMLI} -- Algebraic Multi Level Iteration\\ |
557 |
\member{SolverOptions.GAUSS_SEIDEL} -- Gauss-Seidel\\ |
558 |
\member{SolverOptions.ILU0} -- Incomplete LU-factorization with no fill-in\\ |
559 |
\member{SolverOptions.ILUT} -- Incomplete LU-factorization with fill-in\\ |
560 |
\member{SolverOptions.JACOBI} -- Jacobi preconditioner\\ |
561 |
\member{SolverOptions.NO_PRECONDITIONER} -- do not apply a preconditioner\\ |
562 |
%\member{SolverOptions.REC_ILU} -- recursive ILU0\\ |
563 |
\member{SolverOptions.RILU} -- relaxed ILU0.\\ |
564 |
Not all packages support all preconditioners. It can be assumed that a package |
565 |
makes a reasonable choice if it encounters an unknown preconditioner. |
566 |
See Table~\Refe{TAB FINLEY SOLVER OPTIONS 2} for the preconditioners supported |
567 |
by \finley. |
568 |
\end{methoddesc} |
569 |
|
570 |
\begin{methoddesc}[SolverOptions]{getPreconditioner}{} |
571 |
returns the key of the preconditioner to be used. |
572 |
\end{methoddesc} |
573 |
|
574 |
\begin{methoddesc}[SolverOptions]{setPackage}{package} |
575 |
sets the solver package to be used as a solver. |
576 |
The value of \var{package} must be one of the constants:\\ |
577 |
\member{SolverOptions.DEFAULT} -- choose a default depending on other options\\ |
578 |
\member{SolverOptions.CUSP} -- CUDA sparse linear algebra package\\ |
579 |
\member{SolverOptions.MKL} -- Intel MKL direct solver\\ |
580 |
\member{SolverOptions.PASO} -- built-in PASO solver library\\ |
581 |
\member{SolverOptions.TRILINOS} -- Trilinos solver package (chapter \Refe{TRILINOS})\\ |
582 |
\member{SolverOptions.UMFPACK} -- direct solver from the UMFPACK library.\\ |
583 |
Not all packages are supported on all implementations. An exception may be |
584 |
thrown on some platforms if a particular package is requested. |
585 |
Currently \finley supports \member{SolverOptions.PASO} (as default) and, if |
586 |
available, \member{SolverOptions.MKL}\footnote{If the stiffness matrix is |
587 |
non-regular \MKL may return without returning a proper error code. |
588 |
If you observe suspicious solutions when using MKL, this may be caused by a |
589 |
non-invertible operator.} and \member{SolverOptions.UMFPACK}. |
590 |
\end{methoddesc} |
591 |
|
592 |
\begin{methoddesc}[SolverOptions]{getPackage}{} |
593 |
returns the solver package key. |
594 |
\end{methoddesc} |
595 |
|
596 |
\begin{methoddesc}[SolverOptions]{setTrilinosParameter}{\optional{ name , value}} |
597 |
sets Trilinos parameters for use by MueLu and Trilinos iterative and direct solvers. It is important that the exact space between words is included in the names. |
598 |
The following keywords are supported:\\ |
599 |
\var{"verbosity"}: controls level of AMG output, |
600 |
\var{"none"}, \var{"low"}, \var{"medium"}, \var{"high"}, and \var{"extreme"}\\ |
601 |
%\var{"print initial parameters"}: \var{True} |
602 |
\var{"problem:type"}: \var{"unknown"}, \var{"Poisson-2D"}, \var{"Poisson-3D"}, \var{"Elasticity-2D"},\\ \hspace{3.2cm}\var{"Elasticity-3D"}, \var{"Poisson-2D-complex"}, \var{"Poisson-3D-complex"},\\ \hspace{3.2cm}\var{"Elasticity-2D-complex"}, \var{"Elasticity-3D-complex"},\\ |
603 |
\hspace{3.2cm}\var{"ConvectionDiffusion"} and \var{"MHD"}. |
604 |
\\ |
605 |
\var{"number~of~equations"}: number of equations\\ |
606 |
\var{"multigrid~algorithm"}: \var{ "sa"}, \var{"unsmoothed"}, \var{ "pg"}, \var{"interp"}, \var{"emin"},\\ \hspace{4.2cm}\var{"semicoarsen"}\\ |
607 |
\var{"sa:~damping factor"}: for smoothed aggregation AMG alorithm, default =1.3\\ |
608 |
\var{"sa:~use~filtered~matrix"}: Booleon\\ |
609 |
\var{"filtered~matrix:~use~lumping"}: Boolean\\ |
610 |
\var{"filtered~matrix:~reuse~eigenvalue"}: Boolean\\ |
611 |
\var{"interp:~interpolation~order"}: 0,~1\\ |
612 |
\var{"interp:~build~coarse~coordinates}: Boolean\\ |
613 |
\var{"emin:~iterative~method"}: \var{"cg"}, \var{"gmres"} or \var{"sd"}\\ |
614 |
\var{"emin:~num~iterations"}: integer\\ |
615 |
\var{"emin:~num~reuse~iterations"}: integer\\ |
616 |
%options.setTrilinosParameter("emin: pattern", "AkPtent") |
617 |
%options.setTrilinosParameter("emin: pattern order", 1) |
618 |
%# semicoarsen |
619 |
%options.setTrilinosParameter("semicoarsen: coarsen rate", 3) |
620 |
%# |
621 |
%options.setTrilinosParameter("transpose: use implicit", False) |
622 |
%options.setTrilinosParameter("max levels", 10) |
623 |
%options.setTrilinosParameter("coarse: max size", 2000) |
624 |
%options.setTrilinosParameter("coarse: type", "SuperLU") |
625 |
%options.setTrilinosParameter("aggregation: type", "structured") |
626 |
%options.setTrilinosParameter("aggregation: ordering", "natural") |
627 |
% # "natural", "graph", "random" |
628 |
%options.setTrilinosParameter("aggregation: drop scheme", "classical") |
629 |
% # "classical", "distance laplacian" |
630 |
%options.setTrilinosParameter("aggregation: drop tol", 0.0) |
631 |
%options.setTrilinosParameter("aggregation: min agg size", 2) |
632 |
%options.setTrilinosParameter("aggregation: max agg size", -1) |
633 |
% # -1 means unlimited |
634 |
%options.setTrilinosParameter("aggregation: Dirichlet threshold", 1e-5) |
635 |
%options.setTrilinosParameter("problem: symmetric", True) |
636 |
%options.setTrilinosParameter("smoother: pre or post", "both") |
637 |
%options.setTrilinosParameter("smoother: type", "RELAXATION") |
638 |
%options.setTrilinosParameter("smoother: pre type", "CHEBYSHEV") |
639 |
%options.setTrilinosParameter("smoother: post type", "RELAXATION") |
640 |
\var{"cycle~type"}: \var{"V"} or \var{"W"}\\ |
641 |
%options.setTrilinosParameter("reuse: type", "full") |
642 |
%options.setTrilinosParameter("repartition: enable", False) |
643 |
%options.setTrilinosParameter("repartition: start level", 2) |
644 |
%options.setTrilinosParameter("repartition: min rows per proc", 800) |
645 |
%options.setTrilinosParameter("repartition: max imbalance", 1.2) |
646 |
%options.setTrilinosParameter("repartition: remap parts", True) |
647 |
%options.setTrilinosParameter("repartition: rebalance P and R", False) |
648 |
\var{"xml~parameter~file"}: name of XML file containing chosen XML parameters\\ |
649 |
\end{methoddesc} |
650 |
|
651 |
\begin{methoddesc}[SolverOptions]{resetDiagnostics}{\optional{all=False}} |
652 |
resets the diagnostics. If \var{all} is \True all diagnostics, including |
653 |
accumulative counters, are reset. |
654 |
\end{methoddesc} |
655 |
|
656 |
\begin{methoddesc}[SolverOptions]{getDiagnostics}{\optional{ name}} |
657 |
returns the diagnostic information \var{name}. The following keywords are |
658 |
supported:\\ |
659 |
\var{"num_iter"}: number of iteration steps\\ |
660 |
\var{"cum_num_iter"}: cumulative number of iteration steps\\ |
661 |
\var{"num_level"}: number of levels in the multi level solver\\ |
662 |
\var{"num_inner_iter"}: number of inner iteration steps\\ |
663 |
\var{"cum_num_inner_iter"}: cumulative number of inner iteration steps\\ |
664 |
\var{"time"}: execution time\\ |
665 |
\var{"cum_time"}: cumulative execution time\\ |
666 |
\var{"set_up_time"}: time to set up the solver, typically this includes |
667 |
factorization and reordering\\ |
668 |
\var{"cum_set_up_time"}: cumulative time to set up the solver\\ |
669 |
\var{"net_time"}: net execution time, excluding setup time for the solver |
670 |
and execution time for preconditioner\\ |
671 |
\var{"cum_net_time"}: cumulative net execution time\\ |
672 |
\var{"residual_norm"}: norm of the final residual\\ |
673 |
\var{"converged"}: status of convergence\\ |
674 |
\var{"preconditioner_size"}: size of preconditioner in MBytes \\ |
675 |
\var{"time_step_backtracking_used"}: whether the time step size was reduced after convergence failure \\ |
676 |
\var{"coarse_level_sparsity"}: the sparsity at coarse level (AMG only) \\ |
677 |
\var{"num_coarse_unknowns"}: number of unknowns at coarse level (AMG only) . |
678 |
|
679 |
\end{methoddesc} |
680 |
|
681 |
\begin{methoddesc}[SolverOptions]{hasConverged}{} |
682 |
returns \True if the last solver call has been finalized successfully. |
683 |
If an exception has been thrown by the solver the status of this flag is undefined. |
684 |
\end{methoddesc} |
685 |
|
686 |
|
687 |
\begin{methoddesc}[SolverOptions]{setReordering}{ordering} |
688 |
sets the key of the reordering method to be applied if supported by the solver. |
689 |
Some direct solvers support reordering to optimize compute time and storage |
690 |
use during elimination. The value of \var{ordering} must be one of the |
691 |
constants:\\ |
692 |
\member{SolverOptions.DEFAULT_REORDERING} -- as recommended by the solver\\ |
693 |
\member{SolverOptions.MINIMUM_FILL_IN} -- reorder matrix to reduce fill-in during factorization\\ |
694 |
\member{SolverOptions.NESTED_DISSECTION} -- reorder matrix to improve load balancing during factorization\\ |
695 |
\member{SolverOptions.NO_REORDERING} -- no matrix reordering applied. |
696 |
\end{methoddesc} |
697 |
|
698 |
\begin{methoddesc}[SolverOptions]{getReordering}{} |
699 |
returns the key of the reordering method to be applied if supported by the solver. |
700 |
\end{methoddesc} |
701 |
|
702 |
\begin{methoddesc}[SolverOptions]{setRestart}{\optional{restart=None}} |
703 |
sets the number of iterations steps after which \GMRES is to perform a restart. |
704 |
If \var{restart} is equal to \var{None} no restart is performed. |
705 |
\end{methoddesc} |
706 |
|
707 |
\begin{methoddesc}[SolverOptions]{getRestart}{} |
708 |
returns the number of iterations steps after which \GMRES performs a restart. |
709 |
\end{methoddesc} |
710 |
|
711 |
\begin{methoddesc}[SolverOptions]{setTruncation}{\optional{truncation=20}} |
712 |
sets the number of residuals in \GMRES to be stored for orthogonalization. |
713 |
The more residuals are stored the faster \GMRES converges but the higher the storage needs are and the more expensive |
714 |
a single iteration step becomes. |
715 |
\end{methoddesc} |
716 |
|
717 |
\begin{methoddesc}[SolverOptions]{getTruncation}{} |
718 |
returns the number of residuals in \GMRES to be stored for orthogonalization. |
719 |
\end{methoddesc} |
720 |
|
721 |
|
722 |
\begin{methoddesc}[SolverOptions]{setIterMax}{\optional{iter_max=10000}} |
723 |
sets the maximum number of iteration steps. |
724 |
\end{methoddesc} |
725 |
|
726 |
\begin{methoddesc}[SolverOptions]{getIterMax}{} |
727 |
returns maximum number of iteration steps. |
728 |
\end{methoddesc} |
729 |
|
730 |
%\begin{methoddesc}[SolverOptions]{setLevelMax}{\optional{level_max=10}} |
731 |
%sets the maximum number of coarsening levels to be used in the \AMG solver or |
732 |
%preconditioner. |
733 |
%\end{methoddesc} |
734 |
|
735 |
%\begin{methoddesc}[SolverOptions]{getLevelMax}{} |
736 |
%returns the maximum number of coarsening levels to be used in an algebraic |
737 |
%multi level solver or preconditioner. |
738 |
%\end{methoddesc} |
739 |
|
740 |
%\begin{methoddesc}[SolverOptions]{setCoarseningThreshold}{\optional{theta=0.25}} |
741 |
%sets the threshold for coarsening in the \AMG solver or preconditioner. |
742 |
%\end{methoddesc} |
743 |
|
744 |
%\begin{methoddesc}[SolverOptions]{getCoarseningThreshold}{} |
745 |
%returns the threshold for coarsening in the \AMG solver or preconditioner. |
746 |
%\end{methoddesc} |
747 |
|
748 |
%\begin{methoddesc}[SolverOptions]{setDiagonalDominanceThreshold}{\optional{value=0.5}} |
749 |
%sets the threshold for diagonal dominant rows which are eliminated during \AMG coarsening. |
750 |
%\end{methoddesc} |
751 |
|
752 |
%\begin{methoddesc}[SolverOptions]{getDiagonalDominanceThreshold}{} |
753 |
%returns the threshold for diagonal dominant rows which are eliminated during \AMG coarsening. |
754 |
%\end{methoddesc} |
755 |
|
756 |
%\begin{methoddesc}[SolverOptions]{setMinCoarseMatrixSize}{\optional{size=500}} |
757 |
%sets the minimum size of the coarsest level matrix in \AMG. |
758 |
%\end{methoddesc} |
759 |
|
760 |
%\begin{methoddesc}[SolverOptions]{getMinCoarseMatrixSize}{} |
761 |
%returns the minimum size of the coarsest level matrix in \AMG. |
762 |
%\end{methoddesc} |
763 |
|
764 |
%\begin{methoddesc}[SolverOptions]{setSmoother}{\optional{smoother=\GAUSSSEIDEL}} |
765 |
%sets the \JACOBI or \GAUSSSEIDEL smoother to be used with \AMG. |
766 |
%\end{methoddesc} |
767 |
|
768 |
%\begin{methoddesc}[SolverOptions]{getSmoother}{} |
769 |
%returns the key of the smoother used in \AMG. |
770 |
%\end{methoddesc} |
771 |
|
772 |
%\begin{methoddesc}[SolverOptions]{setAMGInterpolation}{\optional{method=\var{None}}} |
773 |
%sets interpolation method for \AMG to |
774 |
%\member{CLASSIC_INTERPOLATION_WITH_FF_COUPLING}, |
775 |
%\member{CLASSIC_INTERPOLATION}, or |
776 |
%\member{DIRECT_INTERPOLATION}. |
777 |
%\end{methoddesc} |
778 |
|
779 |
%\begin{methoddesc}[SolverOptions]{getAMGInterpolation}{} |
780 |
%returns the key of the interpolation method for \AMG. |
781 |
%\end{methoddesc} |
782 |
|
783 |
%\begin{methoddesc}[SolverOptions]{setNumSweeps}{\optional{sweeps=2}} |
784 |
%sets the number of sweeps in a \JACOBI or \GAUSSSEIDEL preconditioner. |
785 |
%\end{methoddesc} |
786 |
|
787 |
%\begin{methoddesc}[SolverOptions]{getNumSweeps}{} |
788 |
%returns the number of sweeps in a \JACOBI or \GAUSSSEIDEL preconditioner. |
789 |
%\end{methoddesc} |
790 |
|
791 |
%\begin{methoddesc}[SolverOptions]{setNumPreSweeps}{\optional{sweeps=2}} |
792 |
%sets the number of sweeps in the pre-smoothing step of \AMG. |
793 |
%\end{methoddesc} |
794 |
|
795 |
%\begin{methoddesc}[SolverOptions]{getNumPreSweeps}{} |
796 |
%returns the number of sweeps in the pre-smoothing step of \AMG. |
797 |
%\end{methoddesc} |
798 |
|
799 |
%\begin{methoddesc}[SolverOptions]{setNumPostSweeps}{\optional{sweeps=2}} |
800 |
%sets the number of sweeps in the post-smoothing step of \AMG. |
801 |
%\end{methoddesc} |
802 |
|
803 |
%\begin{methoddesc}[SolverOptions]{getNumPostSweeps}{} |
804 |
%returns he number of sweeps in the post-smoothing step of \AMG. |
805 |
%\end{methoddesc} |
806 |
|
807 |
\begin{methoddesc}[SolverOptions]{setTolerance}{\optional{rtol=1.e-8}} |
808 |
sets the relative tolerance for the solver. The actual meaning of tolerance |
809 |
depends on the underlying PDE library. In most cases, the tolerance |
810 |
will only consider the error from solving the discrete problem but will |
811 |
not consider any discretization error. |
812 |
\end{methoddesc} |
813 |
|
814 |
\begin{methoddesc}[SolverOptions]{getTolerance}{} |
815 |
returns the relative tolerance for the solver. |
816 |
\end{methoddesc} |
817 |
|
818 |
\begin{methoddesc}[SolverOptions]{setAbsoluteTolerance}{\optional{atol=0.}} |
819 |
sets the absolute tolerance for the solver. The actual meaning of tolerance |
820 |
depends on the underlying PDE library. In most cases, the tolerance |
821 |
will only consider the error from solving the discrete problem but will |
822 |
not consider any discretization error. |
823 |
\end{methoddesc} |
824 |
|
825 |
\begin{methoddesc}[SolverOptions]{getAbsoluteTolerance}{} |
826 |
returns the absolute tolerance for the solver. |
827 |
\end{methoddesc} |
828 |
|
829 |
\begin{methoddesc}[SolverOptions]{setInnerTolerance}{\optional{rtol=0.9}} |
830 |
sets the relative tolerance for an inner iteration scheme, for instance |
831 |
on the coarsest level in a multi-level scheme. |
832 |
\end{methoddesc} |
833 |
|
834 |
\begin{methoddesc}[SolverOptions]{getInnerTolerance}{} |
835 |
returns the relative tolerance for an inner iteration scheme. |
836 |
\end{methoddesc} |
837 |
|
838 |
\begin{methoddesc}[SolverOptions]{setRelaxationFactor}{\optional{factor=0.3}} |
839 |
sets the relaxation factor used to add dropped elements in \RILU to the main diagonal. |
840 |
\end{methoddesc} |
841 |
|
842 |
\begin{methoddesc}[SolverOptions]{getRelaxationFactor}{} |
843 |
returns the relaxation factor used to add dropped elements in \RILU to the main diagonal. |
844 |
\end{methoddesc} |
845 |
|
846 |
\begin{methoddesc}[SolverOptions]{isSymmetric}{} |
847 |
returns \True if the discrete system is indicated as symmetric. |
848 |
\end{methoddesc} |
849 |
|
850 |
\begin{methoddesc}[SolverOptions]{setSymmetryOn}{} |
851 |
sets the symmetry flag to indicate that the coefficient matrix is symmetric. |
852 |
\end{methoddesc} |
853 |
|
854 |
\begin{methoddesc}[SolverOptions]{setSymmetryOff}{} |
855 |
clears the symmetry flag for the coefficient matrix. |
856 |
\end{methoddesc} |
857 |
|
858 |
\begin{methoddesc}[SolverOptions]{isVerbose}{} |
859 |
returns \True if the solver is expected to be verbose. |
860 |
\end{methoddesc} |
861 |
|
862 |
\begin{methoddesc}[SolverOptions]{setVerbosityOn}{} |
863 |
switches the verbosity of the solver on. |
864 |
\end{methoddesc} |
865 |
|
866 |
\begin{methoddesc}[SolverOptions]{setVerbosityOff}{} |
867 |
switches the verbosity of the solver off. |
868 |
\end{methoddesc} |
869 |
|
870 |
\begin{methoddesc}[SolverOptions]{adaptInnerTolerance}{} |
871 |
returns \True if the tolerance of the inner solver is selected automatically. |
872 |
Otherwise the inner tolerance set by \member{setInnerTolerance} is used. |
873 |
\end{methoddesc} |
874 |
|
875 |
\begin{methoddesc}[SolverOptions]{setInnerToleranceAdaptionOn}{} |
876 |
switches the automatic selection of inner tolerance on. |
877 |
\end{methoddesc} |
878 |
|
879 |
\begin{methoddesc}[SolverOptions]{setInnerToleranceAdaptionOff}{} |
880 |
switches the automatic selection of inner tolerance off. |
881 |
\end{methoddesc} |
882 |
|
883 |
\begin{methoddesc}[SolverOptions]{setInnerIterMax}{\optional{iter_max=10}} |
884 |
sets the maximum number of iteration steps for the inner iteration. |
885 |
\end{methoddesc} |
886 |
|
887 |
\begin{methoddesc}[SolverOptions]{getInnerIterMax}{} |
888 |
returns the maximum number of inner iteration steps. |
889 |
\end{methoddesc} |
890 |
|
891 |
\begin{methoddesc}[SolverOptions]{acceptConvergenceFailure}{} |
892 |
returns \True if a failure to meet the stopping criteria within the given |
893 |
number of iteration steps is not raising in exception. This is useful |
894 |
if a solver is used in a non-linear context where the non-linear solver can |
895 |
continue even if the returned solution does not necessarily meet the stopping |
896 |
criteria. One can use the \member{hasConverged} method to check if the |
897 |
last call to the solver was successful. |
898 |
\end{methoddesc} |
899 |
|
900 |
\begin{methoddesc}[SolverOptions]{setAcceptanceConvergenceFailureOn}{} |
901 |
switches the acceptance of a failure of convergence on. |
902 |
\end{methoddesc} |
903 |
|
904 |
\begin{methoddesc}[SolverOptions]{setAcceptanceConvergenceFailureOff}{} |
905 |
switches the acceptance of a failure of convergence off. |
906 |
\end{methoddesc} |
907 |
|
908 |
\begin{memberdesc}[SolverOptions]{DEFAULT} |
909 |
default method, preconditioner or package to be used to solve the PDE. |
910 |
An appropriate method should be chosen by the used PDE solver library. |
911 |
\end{memberdesc} |
912 |
|
913 |
\begin{memberdesc}[SolverOptions]{MKL} |
914 |
the \MKL library by Intel, \Refe{MKL}\footnote{The \MKL library will only be |
915 |
available when the Intel compilation environment was used to build \escript.}. |
916 |
\end{memberdesc} |
917 |
|
918 |
\begin{memberdesc}[SolverOptions]{UMFPACK} |
919 |
the \UMFPACK library, \Refe{UMFPACK}. Note that \UMFPACK is not parallelized. |
920 |
\end{memberdesc} |
921 |
|
922 |
\begin{memberdesc}[SolverOptions]{PASO} |
923 |
\PASO is the default solver library of \finley, see \Sec{chap:finley}. |
924 |
\end{memberdesc} |
925 |
|
926 |
\begin{memberdesc}[SolverOptions]{ITERATIVE} |
927 |
the default iterative method and preconditioner. The actual method used |
928 |
depends on the PDE solver library and the chosen solver package. |
929 |
Typically, \PCG is used for symmetric PDEs and \BiCGStab otherwise, both with |
930 |
\JACOBI preconditioner. |
931 |
\end{memberdesc} |
932 |
|
933 |
\begin{memberdesc}[SolverOptions]{DIRECT} |
934 |
the default direct linear solver. |
935 |
\end{memberdesc} |
936 |
|
937 |
\begin{memberdesc}[SolverOptions]{CHOLEVSKY} |
938 |
direct solver based on Cholevsky factorization (or similar), see \Refe{Saad}. |
939 |
The solver requires a symmetric PDE. |
940 |
\end{memberdesc} |
941 |
|
942 |
\begin{memberdesc}[SolverOptions]{PCG} |
943 |
preconditioned conjugate gradient method, see \Refe{WEISS}\index{linear solver!PCG}\index{PCG}. |
944 |
The solver requires a symmetric PDE. |
945 |
\end{memberdesc} |
946 |
|
947 |
\begin{memberdesc}[SolverOptions]{TFQMR} |
948 |
transpose-free quasi-minimal residual method, see \Refe{WEISS}\index{linear solver!TFQMR}\index{TFQMR}. |
949 |
\end{memberdesc} |
950 |
|
951 |
\begin{memberdesc}[SolverOptions]{GMRES} |
952 |
the GMRES method, see \Refe{WEISS}\index{linear solver!GMRES}\index{GMRES}. |
953 |
Truncation and restart are controlled by the \var{truncation} and \var{restart} |
954 |
parameters of \method{getSolution}. |
955 |
\end{memberdesc} |
956 |
|
957 |
\begin{memberdesc}[SolverOptions]{MINRES} |
958 |
minimal residual method\index{linear solver!MINRES}\index{MINRES} |
959 |
\end{memberdesc} |
960 |
|
961 |
\begin{memberdesc}[SolverOptions]{ROWSUM_LUMPING} |
962 |
row sum lumping of the stiffness matrix, see Section~\Refe{LUMPING} for details\index{linear solver!row sum lumping}\index{row sum lumping}. |
963 |
Lumping does not use the linear system solver library. |
964 |
\end{memberdesc} |
965 |
|
966 |
\begin{memberdesc}[SolverOptions]{HRZ_LUMPING} |
967 |
HRZ lumping of the stiffness matrix, see Section~\Refe{LUMPING} for details\index{linear solver!HRZ lumping}\index{HRZ lumping}. |
968 |
Lumping does not use the linear system solver library. |
969 |
\end{memberdesc} |
970 |
|
971 |
\begin{memberdesc}[SolverOptions]{PRES20} |
972 |
the GMRES method with truncation after five residuals and restart after |
973 |
20 steps, see \Refe{WEISS}. |
974 |
\end{memberdesc} |
975 |
|
976 |
\begin{memberdesc}[SolverOptions]{CGS} |
977 |
conjugate gradient squared method, see \Refe{WEISS}. |
978 |
\end{memberdesc} |
979 |
|
980 |
\begin{memberdesc}[SolverOptions]{BICGSTAB} |
981 |
stabilized bi-conjugate gradients methods, see \Refe{WEISS}. |
982 |
\end{memberdesc} |
983 |
|
984 |
\begin{memberdesc}[SolverOptions]{SSOR} |
985 |
symmetric successive over-relaxation method, see \Refe{WEISS}. |
986 |
Typically used as preconditioner but some linear solver libraries support |
987 |
this as a solver. |
988 |
\end{memberdesc} |
989 |
|
990 |
\begin{memberdesc}[SolverOptions]{ILU0} |
991 |
the incomplete LU factorization preconditioner with no fill-in, see \Refe{Saad}. |
992 |
\end{memberdesc} |
993 |
|
994 |
\begin{memberdesc}[SolverOptions]{JACOBI} |
995 |
the Jacobi preconditioner, see \Refe{Saad}. |
996 |
\end{memberdesc} |
997 |
|
998 |
\begin{memberdesc}[SolverOptions]{AMG} |
999 |
the algebraic multi grid method, see \Refe{AMG}. This method can be used as |
1000 |
linear solver method but is more robust when used as a preconditioner. |
1001 |
\end{memberdesc} |
1002 |
|
1003 |
\begin{memberdesc}[SolverOptions]{GAUSS_SEIDEL} |
1004 |
the symmetric Gauss-Seidel preconditioner, see \Refe{Saad}. |
1005 |
\member{getNumSweeps()} is the number of sweeps used. |
1006 |
\end{memberdesc} |
1007 |
|
1008 |
%\begin{memberdesc}[SolverOptions]{RILU} |
1009 |
%relaxed incomplete LU factorization preconditioner, see \Refe{RELAXILU}. |
1010 |
%This method is similar to the one used for \ILU but dropped elements are added |
1011 |
%to the main diagonal with the relaxation factor \member{getRelaxationFactor}. |
1012 |
%\end{memberdesc} |
1013 |
|
1014 |
\begin{memberdesc}[SolverOptions]{REC_ILU} |
1015 |
recursive incomplete LU factorization preconditioner, see \Refe{RILU}. |
1016 |
This method is similar to the one used for \ILU but applies reordering during |
1017 |
the factorization. |
1018 |
\end{memberdesc} |
1019 |
|
1020 |
\begin{memberdesc}[SolverOptions]{NO_REORDERING} |
1021 |
no reordering is used during factorization. |
1022 |
\end{memberdesc} |
1023 |
|
1024 |
\begin{memberdesc}[SolverOptions]{DEFAULT_REORDERING} |
1025 |
the default reordering method during factorization. |
1026 |
\end{memberdesc} |
1027 |
|
1028 |
\begin{memberdesc}[SolverOptions]{MINIMUM_FILL_IN} |
1029 |
applies reordering before factorization using a fill-in minimization strategy. |
1030 |
You have to check with the particular solver library or linear solver package |
1031 |
if this is supported. In any case, it is advisable to apply reordering on the |
1032 |
mesh to minimize fill-in. |
1033 |
\end{memberdesc} |
1034 |
|
1035 |
\begin{memberdesc}[SolverOptions]{NESTED_DISSECTION} |
1036 |
applies reordering before factorization using a nested dissection strategy. |
1037 |
You have to check with the particular solver library or linear solver package |
1038 |
if this is supported. In any case, it is advisable to apply reordering on the |
1039 |
mesh to minimize fill-in. |
1040 |
\end{memberdesc} |
1041 |
|
1042 |
\begin{memberdesc}[SolverOptions]{TRILINOS} |
1043 |
the Trilinos library~\cite{Trilinos} is used as a solver. |
1044 |
\end{memberdesc} |
1045 |
|
1046 |
\begin{memberdesc}[SolverOptions]{NO_PRECONDITIONER} |
1047 |
no preconditioner is applied. |
1048 |
\end{memberdesc} |
1049 |
|
1050 |
\begin{memberdesc}[SolverOptions]{DIRECT_INTERPOLATION} |
1051 |
direct interpolation in \AMG, see \cite{AMG} |
1052 |
\end{memberdesc} |
1053 |
\begin{memberdesc}[SolverOptions]{CLASSIC_INTERPOLATION} |
1054 |
classic interpolation in \AMG, see \cite{AMG} |
1055 |
\end{memberdesc} |
1056 |
\begin{memberdesc}[SolverOptions]{CLASSIC_INTERPOLATION_WITH_FF_COUPLING} |
1057 |
classic interpolation with enforced FF coupling in \AMG, see \cite{AMG} |
1058 |
\end{memberdesc} |
1059 |
|
1060 |
\input{lumping} |