esys.escript.pdetools Package¶
Classes¶
- class esys.escript.pdetools.ArithmeticTuple(*args)¶
Tuple supporting inplace update x+=y and scaling x=a*y where
x,y
is an ArithmeticTuple anda
is a float.Example of usage:
from esys.escript import Data from numpy import array a=eData(...) b=array([1.,4.]) x=ArithmeticTuple(a,b) y=5.*x
- __init__(*args)¶
Initializes object with elements
args
.- Parameters
args – tuple of objects that support inplace add (x+=y) and scaling (x=a*y)
- class esys.escript.pdetools.CorrectionFailed¶
Exception thrown if no convergence has been achieved in the solution correction scheme.
- __init__(*args, **kwargs)¶
- class esys.escript.pdetools.Defect¶
Defines a non-linear defect F(x) of a variable x. This class includes two functions (bilinearform and eval) that must be overridden by subclassing before use.
- __init__()¶
Initializes defect.
- bilinearform(x0, x1)¶
Returns the inner product of x0 and x1
NOTE: MUST BE OVERRIDDEN BY A SUBCLASS
- Parameters
x0 – value for x0
x1 – value for x1
- Returns
the inner product of x0 and x1
- Return type
float
- derivative(F0, x0, v, v_is_normalised=True)¶
Returns the directional derivative at
x0
in the direction ofv
.- Parameters
F0 – value of this defect at x0
x0 – value at which derivative is calculated
v – direction
v_is_normalised – True to indicate that
v
is nomalized (self.norm(v)=0)
- Returns
derivative of this defect at x0 in the direction of
v
- Note
by default numerical evaluation (self.eval(x0+eps*v)-F0)/eps is used but this method maybe overwritten to use exact evaluation.
- eval(x)¶
Returns the value F of a given
x
.NOTE: MUST BE OVERRIDDEN BY A SUBCLASS
- Parameters
x – value for which the defect
F
is evaluated- Returns
value of the defect at
x
- getDerivativeIncrementLength()¶
Returns the relative increment length used to approximate the derivative of the defect. :return: value of the defect at
x
:rtype: positivefloat
- norm(x)¶
Returns the norm of argument
x
.- Parameters
x – a value
- Returns
norm of argument x
- Return type
float
- Note
by default
sqrt(self.bilinearform(x,x)
is returned.
- setDerivativeIncrementLength(inc=1.4901161193847656e-05)¶
Sets the relative length of the increment used to approximate the derivative of the defect. The increment is inc*norm(x)/norm(v)*v in the direction of v with x as a starting point.
- Parameters
inc (positive
float
) – relative increment length
- class esys.escript.pdetools.HomogeneousSaddlePointProblem(**kwargs)¶
This class provides a framework for solving linear homogeneous saddle point problems of the form:
*Av+B^*p=f* *Bv =0*
for the unknowns v and p and given operators A and B and given right hand side f. B^* is the adjoint operator of B. A may depend weakly on v and p.
- __init__(**kwargs)¶
initializes the saddle point problem
- Bv(v, tol)¶
Returns Bv with accuracy
tol
(overwrite)- Return type
equal to the type of p
- Note
boundary conditions on p should be zero!
- getAbsoluteTolerance()¶
Returns the absolute tolerance.
- Returns
absolute tolerance
- Return type
float
- getDV(p, v, tol)¶
return a correction to the value for a given v and a given p with accuracy
tol
(overwrite)- Parameters
p – pressure
v – pressure
- Returns
dv given as dv= A^{-1} (f-A v-B^*p)
- Note
Only A may depend on v and p
- getTolerance()¶
Returns the relative tolerance.
- Returns
relative tolerance
- Return type
float
- inner_p(p0, p1)¶
Returns inner product of p0 and p1 (overwrite).
- Parameters
p0 – a pressure
p1 – a pressure
- Returns
inner product of p0 and p1
- Return type
float
- inner_pBv(p, Bv)¶
Returns inner product of element p and Bv (overwrite).
- Parameters
p – a pressure increment
Bv – a residual
- Returns
inner product of element p and Bv
- Return type
float
- Note
used if PCG is applied.
- norm_Bv(Bv)¶
Returns the norm of Bv (overwrite).
- Return type
equal to the type of p
- Note
boundary conditions on p should be zero!
- norm_p(p)¶
calculates the norm of
p
- Parameters
p – a pressure
- Returns
the norm of
p
using the inner product for pressure- Return type
float
- norm_v(v)¶
Returns the norm of v (overwrite).
- Parameters
v – a velovity
- Returns
norm of v
- Return type
non-negative
float
- resetControlParameters(K_p=1.0, K_v=1.0, rtol_max=0.01, rtol_min=1e-07, chi_max=0.5, reduction_factor=0.3, theta=0.1)¶
sets a control parameter
- Parameters
K_p (
float
) – initial value for constant to adjust pressure toleranceK_v (
float
) – initial value for constant to adjust velocity tolerancertol_max (
float
) – maximuim relative tolerance used to calculate presssure and velocity increment.chi_max (
float
) – maximum tolerable converegence rate.reduction_factor (
float
) – reduction factor for adjustment factors.
- setAbsoluteTolerance(tolerance=0.0)¶
Sets the absolute tolerance.
- Parameters
tolerance (non-negative
float
) – tolerance to be used
- setControlParameter(K_p=None, K_v=None, rtol_max=None, rtol_min=None, chi_max=None, reduction_factor=None, theta=None)¶
sets a control parameter
- Parameters
K_p (
float
) – initial value for constant to adjust pressure toleranceK_v (
float
) – initial value for constant to adjust velocity tolerancertol_max (
float
) – maximuim relative tolerance used to calculate presssure and velocity increment.chi_max (
float
) – maximum tolerable converegence rate.
- setTolerance(tolerance=0.0001)¶
Sets the relative tolerance for (v,p).
- Parameters
tolerance (non-negative
float
) – tolerance to be used
- solve(v, p, max_iter=20, verbose=False, usePCG=True, iter_restart=20, max_correction_steps=10)¶
Solves the saddle point problem using initial guesses v and p.
- Parameters
v (
Data
) – initial guess for velocityp (
Data
) – initial guess for pressureusePCG (
bool
) – indicates the usage of the PCG rather than GMRES scheme.max_iter (
int
) – maximum number of iteration steps per correction attemptverbose (
bool
) – if True, shows information on the progress of the saddlepoint problem solver.iter_restart (
int
) – restart the iteration afteriter_restart
steps (only used if useUzaw=False)
- Return type
tuple
ofData
objects- Note
typically this method is overwritten by a subclass. It provides a wrapper for the
_solve
method.
- solve_AinvBt(dp, tol)¶
Solves A dv=B^*dp with accuracy
tol
- Parameters
dp – a pressure increment
- Returns
the solution of A dv=B^*dp
- Note
boundary conditions on dv should be zero! A is the operator used in
getDV
and must not be altered.
- solve_prec(Bv, tol)¶
Provides a preconditioner for (BA^{-1}B^ * ) applied to Bv with accuracy
tol
- Return type
equal to the type of p
- Note
boundary conditions on p should be zero!
- class esys.escript.pdetools.IndefinitePreconditioner¶
Exception thrown if the preconditioner is not positive definite.
- __init__(*args, **kwargs)¶
- class esys.escript.pdetools.IterationBreakDown¶
Exception thrown if the iteration scheme encountered an incurable breakdown.
- __init__(*args, **kwargs)¶
- class esys.escript.pdetools.Locator(where, x=array([0., 0., 0.]))¶
Locator provides access to the values of data objects at a given spatial coordinate x.
In fact, a Locator object finds the sample in the set of samples of a given function space or domain which is closest to the given point x.
- __init__(where, x=array([0., 0., 0.]))¶
Initializes a Locator to access values in Data objects on the Doamin or FunctionSpace for the sample point which is closest to the given point x.
- Parameters
where (
escript.FunctionSpace
) – function spacex (
numpy.ndarray
orlist
ofnumpy.ndarray
) – location(s) of the Locator
- getFunctionSpace()¶
Returns the function space of the Locator.
- getId(item=None)¶
Returns the identifier of the location.
- getValue(data)¶
Returns the value of
data
at the Locator ifdata
is aData
object otherwise the object is returned.
- getX()¶
Returns the exact coordinates of the Locator.
- setValue(data, v)¶
Sets the value of the
data
at the Locator.
- class esys.escript.pdetools.MaxIterReached¶
Exception thrown if the maximum number of iteration steps is reached.
- __init__(*args, **kwargs)¶
- class esys.escript.pdetools.NegativeNorm¶
Exception thrown if a norm calculation returns a negative norm.
- __init__(*args, **kwargs)¶
- class esys.escript.pdetools.NoPDE(domain, D=None, Y=None, q=None, r=None)¶
Solves the following problem for u:
kronecker[i,j]*D[j]*u[j]=Y[i]
with constraint
u[j]=r[j] where q[j]>0
where D, Y, r and q are given functions of rank 1.
In the case of scalars this takes the form
D*u=Y
with constraint
u=r where q>0
where D, Y, r and q are given scalar functions.
The constraint overwrites any other condition.
- Note
This class is similar to the
linearPDEs.LinearPDE
class with A=B=C=X=0 but has the intention that all input parameters are given inSolution
orReducedSolution
.
- __init__(domain, D=None, Y=None, q=None, r=None)¶
Initializes the problem.
- Parameters
domain (
Domain
) – domain of the PDED (
float
,int
,numpy.ndarray
,Data
) – coefficient of the solutionY (
float
,int
,numpy.ndarray
,Data
) – right hand sideq (
float
,int
,numpy.ndarray
,Data
) – location of constraintsr (
float
,int
,numpy.ndarray
,Data
) – value of solution at locations of constraints
- getSolution()¶
Returns the solution.
- Returns
the solution of the problem
- Return type
Data
object in theFunctionSpace
Solution
orReducedSolution
- setReducedOff()¶
Sets the
FunctionSpace
of the solution toSolution
.
- setReducedOn()¶
Sets the
FunctionSpace
of the solution toReducedSolution
.
- setValue(D=None, Y=None, q=None, r=None)¶
Assigns values to the parameters.
- Parameters
D (
float
,int
,numpy.ndarray
,Data
) – coefficient of the solutionY (
float
,int
,numpy.ndarray
,Data
) – right hand sideq (
float
,int
,numpy.ndarray
,Data
) – location of constraintsr (
float
,int
,numpy.ndarray
,Data
) – value of solution at locations of constraints
- class esys.escript.pdetools.Projector(domain, reduce=True, fast=True)¶
The Projector is a factory which projects a discontinuous function onto a continuous function on a given domain.
- __init__(domain, reduce=True, fast=True)¶
Creates a continuous function space projector for a domain.
- Parameters
domain – Domain of the projection.
reduce – Flag to reduce projection order
fast – Flag to use a fast method based on matrix lumping
- getSolverOptions()¶
Returns the solver options of the PDE solver.
- Return type
linearPDEs.SolverOptions
- getValue(input_data)¶
Projects
input_data
onto a continuous function.- Parameters
input_data – the data to be projected
- class esys.escript.pdetools.SolverSchemeException¶
This is a generic exception thrown by solvers.
- __init__(*args, **kwargs)¶
- class esys.escript.pdetools.TimeIntegrationManager(*inital_values, **kwargs)¶
A simple mechanism to manage time dependend values.
Typical usage is:
dt=0.1 # time increment tm=TimeIntegrationManager(inital_value,p=1) while t<1. v_guess=tm.extrapolate(dt) # extrapolate to t+dt v=... tm.checkin(dt,v) t+=dt
- Note
currently only p=1 is supported.
- __init__(*inital_values, **kwargs)¶
Sets up the value manager where
inital_values
are the initial values and p is the order used for extrapolation.
- checkin(dt, *values)¶
Adds new values to the manager. The p+1 last values are lost.
- extrapolate(dt)¶
Extrapolates to
dt
forward in time.
- getTime()¶
- getValue()¶
Functions¶
- esys.escript.pdetools.BoundaryValuesFromVolumeTag(domain, **values)¶
Creates a mask on the Solution(domain) function space where the value is one for samples that touch regions tagged by tags.
Usage: m=BoundaryValuesFromVolumeTag(domain, ham=1, f=6)
- Parameters
domain (
escript.Domain
) – domain to be used- Returns
a mask which marks samples that are touching the boundary tagged by any of the given tags
- Return type
escript.Data
of rank 0
- esys.escript.pdetools.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 solutionAprod (function
Aprod(x)
wherex
is of the same object like argumentx
. The returned object needs to be of the same type like argumentr
.) – returns the value Axbilinearform (function
bilinearform(x,r)
wherex
is of the same type like argumentx
andr
. The returned value is afloat
.) – inner product<x,r>
atol (non-negative
float
) – absolute tolerancertol (non-negative
float
) – relative toleranceiter_max (
int
) – maximum number of iteration stepsiter_restart (
int
) – in order to save memory the orthogonalization process is terminated afteriter_restart
steps and the iteration is restarted.
- Returns
the solution approximation and the corresponding residual
- Return type
tuple
- Warning
r
andx
are altered.
- esys.escript.pdetools.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)
wherex
is of the same object like argumentx
. The returned object needs to be of the same type like argumentr
.) – returns the value AxMsolve (function
Msolve(r)
wherer
is of the same type like argumentr
. The returned object needs to be of the same type like argumentx
.) – solves Mx=rbilinearform (function
bilinearform(x,r)
wherex
is of the same type like argumentx
andr
is. The returned value is afloat
.) – inner product<x,r>
atol (non-negative
float
) – absolute tolerancertol (non-negative
float
) – relative toleranceiter_max (
int
) – maximum number of iteration steps
- Returns
the solution approximation and the corresponding residual
- Return type
tuple
- Warning
r
andx
are altered.
- esys.escript.pdetools.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 usedtags (
str
) – boundary tags
- Returns
a mask which marks samples that are touching the boundary tagged by any of the given tags
- Return type
escript.Data
of rank 0
- esys.escript.pdetools.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 usedtags (
str
) – boundary tags
- Returns
a mask which marks samples that are touching the boundary tagged by any of the given tags
- Return type
escript.Data
of rank 0
- esys.escript.pdetools.NewtonGMRES(defect, x, iter_max=100, sub_iter_max=20, atol=0, rtol=0.0001, subtol_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 stepssub_iter_max (positive
int
) – maximum number of inner iteration stepsatol (positive
float
) – absolute tolerance for the solutionrtol (positive
float
) – relative tolerance for the solutiongamma (positive
float
, less than 1) – tolerance safety factor for inner iterationsubtol_max (positive
float
, less than 1) – upper bound for inner tolerance
- Returns
an approximation of the solution with the desired accuracy
- Return type
same type as the initial guess
- esys.escript.pdetools.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)
wherex
is of the same object like argumentx
. The returned object needs to be of the same type like argumentr
.) – returns the value AxMsolve (function
Msolve(r)
wherer
is of the same type like argumentr
. The returned object needs to be of the same type like argumentx
.) – solves Mx=rbilinearform (function
bilinearform(x,r)
wherex
is of the same type like argumentx
andr
is. The returned value is afloat
.) – inner product<x,r>
atol (non-negative
float
) – absolute tolerancertol (non-negative
float
) – relative toleranceiter_max (
int
) – maximum number of iteration steps
- Returns
the solution approximation and the corresponding residual
- Return type
tuple
- Warning
r
andx
are altered.
- esys.escript.pdetools.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 solutionAprod (function
Aprod(x)
wherex
is of the same object like argumentx
. The returned object needs to be of the same type like argumentr
.) – returns the value Axbilinearform (function
bilinearform(x,r)
wherex
is of the same type like argumentx
andr
. The returned value is afloat
.) – inner product<x,r>
atol (non-negative
float
) – absolute tolerancertol (non-negative
float
) – relative toleranceiter_max (
int
) – maximum number of iteration steps
- Return type
tuple
- Warning
r
andx
are altered.
- esys.escript.pdetools.getInfLocator(arg)¶
Return a Locator for a point with the inf value over all arg.
- esys.escript.pdetools.getSupLocator(arg)¶
Return a Locator for a point with the sup value over all arg.