Package esys :: Package escript :: Module pdetools
[hide private]
[frames] | no frames]

Module pdetools


Author: Lutz Gross, l.gross@uq.edu.au

Copyright: Copyright (c) 2003-2009 by University of Queensland Earth Systems Science Computational Center (ESSCC) http://www.uq.edu.au/esscc Primary Business: Queensland, Australia

License: Licensed under the Open Software License version 3.0 http://www.opensource.org/licenses/osl-3.0.php

Classes [hide private]
  ArithmeticTuple
Tuple supporting inplace update x+=y and scaling x=a*y where x,y is an ArithmeticTuple and a is a float.
  CorrectionFailed
Exception thrown if no convergence has been achieved in the solution correction scheme.
  Defect
Defines a non-linear defect F(x) of a variable x.
  HomogeneousSaddlePointProblem
This class provides a framework for solving linear homogeneous saddle point problems of the form:
  IndefinitePreconditioner
Exception thrown if the preconditioner is not positive definite.
  IterationBreakDown
Exception thrown if the iteration scheme encountered an incurable breakdown.
  Locator
Locator provides access to the values of data objects at a given spatial coordinate x.
  MaxIterReached
Exception thrown if the maximum number of iteration steps is reached.
  NegativeNorm
Exception thrown if a norm calculation returns a negative norm.
  NoPDE
Solves the following problem for u:
  Projector
The Projector is a factory which projects a discontinuous function onto a continuous function on a given domain.
  SolverSchemeException
This is a generic exception thrown by solvers.
  TimeIntegrationManager
A simple mechanism to manage time dependend values.
Functions [hide private]
tuple
GMRES(r, Aprod, x, bilinearform, atol=0, rtol=1e-08, iter_max=100, iter_restart=20, verbose=False, P_R=None)
Solver for
tuple
MINRES(r, Aprod, x, Msolve, bilinearform, atol=0, rtol=1e-08, iter_max=100)
Solver for
escript.Data of rank 0
MaskFromBoundaryTag(domain, *tags)
Creates a mask on the Solution(domain) function space where the value is one for samples that touch the boundary tagged by tags.
escript.Data of rank 0
MaskFromTag(domain, *tags)
Creates a mask on the Solution(domain) function space where the value is one for samples that touch regions tagged by tags.
same type as the initial guess
NewtonGMRES(defect, x, iter_max=100, sub_iter_max=20, atol=0, rtol=0.0001, sub_tol_max=0.5, gamma=0.9, verbose=False)
Solves a non-linear problem F(x)=0 for unknown x using the stopping criterion:
tuple
PCG(r, Aprod, x, Msolve, bilinearform, atol=0, rtol=1e-08, iter_max=100, initial_guess=True, verbose=False)
Solver for
tuple
TFQMR(r, Aprod, x, bilinearform, atol=0, rtol=1e-08, iter_max=100)
Solver for
 
_GMRESm(r, Aprod, x, bilinearform, atol, iter_max=100, iter_restart=20, verbose=False, P_R=None)
 
__FDGMRES(F0, defect, x0, atol, iter_max=100, iter_restart=20)
 
__givapp(c, s, vin)
Applies a sequence of Givens rotations (c,s) recursively to the vector vin
Variables [hide private]
  __url__ = 'https://launchpad.net/escript-finley'
Function Details [hide private]

GMRES(r, Aprod, x, bilinearform, atol=0, rtol=1e-08, iter_max=100, iter_restart=20, verbose=False, P_R=None)

 

Solver for

Ax=b

with a general operator A (more details required!). It uses the generalized minimum residual method (GMRES).

The iteration is terminated if

|r| <= atol+rtol*|r0|

where r0 is the initial residual and |.| is the energy norm. In fact

|r| = sqrt( bilinearform(r,r))

Parameters:
  • r (any object supporting inplace add (x+=y) and scaling (x=scalar*y)) - initial residual r=b-Ax. r is altered.
  • x (same like r) - an initial guess for the solution
  • Aprod (function Aprod(x) where x is of the same object like argument x. The returned object needs to be of the same type like argument r.) - returns the value Ax
  • bilinearform (function bilinearform(x,r) where x is of the same type like argument x and r. The returned value is a float.) - inner product <x,r>
  • atol (non-negative float) - absolute tolerance
  • rtol (non-negative float) - relative tolerance
  • iter_max (int) - maximum number of iteration steps
  • iter_restart (int) - in order to save memory the orthogonalization process is terminated after iter_restart steps and the iteration is restarted.
Returns: tuple
the solution approximation and the corresponding residual

Warning: r and x are altered.

MINRES(r, Aprod, x, Msolve, bilinearform, atol=0, rtol=1e-08, iter_max=100)

 

Solver for

Ax=b

with a symmetric and positive definite operator A (more details required!). It uses the minimum residual method (MINRES) with preconditioner M providing an approximation of A.

The iteration is terminated if

|r| <= atol+rtol*|r0|

where r0 is the initial residual and |.| is the energy norm. In fact

|r| = sqrt( bilinearform(Msolve(r),r))

For details on the preconditioned conjugate gradient method see the book:

Templates for the Solution of Linear Systems by R. Barrett, M. Berry, T.F. Chan, J. Demmel, J. Donato, J. Dongarra, V. Eijkhout, R. Pozo, C. Romine, and H. van der Vorst.

Parameters:
  • r (any object supporting inplace add (x+=y) and scaling (x=scalar*y)) - initial residual r=b-Ax. r is altered.
  • x (any object supporting inplace add (x+=y) and scaling (x=scalar*y)) - an initial guess for the solution
  • Aprod (function Aprod(x) where x is of the same object like argument x. The returned object needs to be of the same type like argument r.) - returns the value Ax
  • Msolve (function Msolve(r) where r is of the same type like argument r. The returned object needs to be of the same type like argument x.) - solves Mx=r
  • bilinearform (function bilinearform(x,r) where x is of the same type like argument x and r is. The returned value is a float.) - inner product <x,r>
  • atol (non-negative float) - absolute tolerance
  • rtol (non-negative float) - relative tolerance
  • iter_max (int) - maximum number of iteration steps
Returns: tuple
the solution approximation and the corresponding residual

Warning: r and x are altered.

MaskFromBoundaryTag(domain, *tags)

 

Creates a mask on the Solution(domain) function space where the value is one for samples that touch the boundary tagged by tags.

Usage: m=MaskFromBoundaryTag(domain, "left", "right")

Parameters:
  • domain (escript.Domain) - domain to be used
  • tags (str) - boundary tags
Returns: escript.Data of rank 0
a mask which marks samples that are touching the boundary tagged by any of the given tags

MaskFromTag(domain, *tags)

 

Creates a mask on the Solution(domain) function space where the value is one for samples that touch regions tagged by tags.

Usage: m=MaskFromTag(domain, "ham")

Parameters:
  • domain (escript.Domain) - domain to be used
  • tags (str) - boundary tags
Returns: escript.Data of rank 0
a mask which marks samples that are touching the boundary tagged by any of the given tags

NewtonGMRES(defect, x, iter_max=100, sub_iter_max=20, atol=0, rtol=0.0001, sub_tol_max=0.5, gamma=0.9, verbose=False)

 

Solves a non-linear problem F(x)=0 for unknown x using the stopping criterion:

norm(F(x) <= atol + rtol * norm(F(x0)

where x0 is the initial guess.

Parameters:
  • defect (Defect) - object defining the function F. defect.norm defines the norm used in the stopping criterion.
  • x (any object type allowing basic operations such as numpy.ndarray, Data) - initial guess for the solution, x is altered.
  • iter_max (positive int) - maximum number of iteration steps
  • sub_iter_max (positive int) - maximum number of inner iteration steps
  • atol (positive float) - absolute tolerance for the solution
  • rtol (positive float) - relative tolerance for the solution
  • gamma (positive float, less than 1) - tolerance safety factor for inner iteration
  • sub_tol_max (positive float, less than 1) - upper bound for inner tolerance
Returns: same type as the initial guess
an approximation of the solution with the desired accuracy

PCG(r, Aprod, x, Msolve, bilinearform, atol=0, rtol=1e-08, iter_max=100, initial_guess=True, verbose=False)

 

Solver for

Ax=b

with a symmetric and positive definite operator A (more details required!). It uses the conjugate gradient method with preconditioner M providing an approximation of A.

The iteration is terminated if

|r| <= atol+rtol*|r0|

where r0 is the initial residual and |.| is the energy norm. In fact

|r| = sqrt( bilinearform(Msolve(r),r))

For details on the preconditioned conjugate gradient method see the book:

Templates for the Solution of Linear Systems by R. Barrett, M. Berry, T.F. Chan, J. Demmel, J. Donato, J. Dongarra, V. Eijkhout, R. Pozo, C. Romine, and H. van der Vorst.

Parameters:
  • r (any object supporting inplace add (x+=y) and scaling (x=scalar*y)) - initial residual r=b-Ax. r is altered.
  • x (any object supporting inplace add (x+=y) and scaling (x=scalar*y)) - an initial guess for the solution
  • Aprod (function Aprod(x) where x is of the same object like argument x. The returned object needs to be of the same type like argument r.) - returns the value Ax
  • Msolve (function Msolve(r) where r is of the same type like argument r. The returned object needs to be of the same type like argument x.) - solves Mx=r
  • bilinearform (function bilinearform(x,r) where x is of the same type like argument x and r is. The returned value is a float.) - inner product <x,r>
  • atol (non-negative float) - absolute tolerance
  • rtol (non-negative float) - relative tolerance
  • iter_max (int) - maximum number of iteration steps
Returns: tuple
the solution approximation and the corresponding residual

Warning: r and x are altered.

TFQMR(r, Aprod, x, bilinearform, atol=0, rtol=1e-08, iter_max=100)

 

Solver for

Ax=b

with a general operator A (more details required!). It uses the Transpose-Free Quasi-Minimal Residual method (TFQMR).

The iteration is terminated if

|r| <= atol+rtol*|r0|

where r0 is the initial residual and |.| is the energy norm. In fact

|r| = sqrt( bilinearform(r,r))

Parameters:
  • r (any object supporting inplace add (x+=y) and scaling (x=scalar*y)) - initial residual r=b-Ax. r is altered.
  • x (same like r) - an initial guess for the solution
  • Aprod (function Aprod(x) where x is of the same object like argument x. The returned object needs to be of the same type like argument r.) - returns the value Ax
  • bilinearform (function bilinearform(x,r) where x is of the same type like argument x and r. The returned value is a float.) - inner product <x,r>
  • atol (non-negative float) - absolute tolerance
  • rtol (non-negative float) - relative tolerance
  • iter_max (int) - maximum number of iteration steps
Returns: tuple

Warning: r and x are altered.

__givapp(c, s, vin)

 

Applies a sequence of Givens rotations (c,s) recursively to the vector vin

Warning: vin is altered.