fipy.terms.firstOrderAdvectionTerm

Classes

FirstOrderAdvectionTerm([coeff])

The FirstOrderAdvectionTerm object constructs the b vector contribution for the advection term given by

class fipy.terms.firstOrderAdvectionTerm.FirstOrderAdvectionTerm(coeff=None)

Bases: _NonDiffusionTerm

The FirstOrderAdvectionTerm object constructs the b vector contribution for the advection term given by

\[u \abs{\nabla \phi}\]

from the advection equation given by:

\[\frac{\partial \phi}{\partial t} + u \abs{\nabla \phi} = 0\]

The construction of the gradient magnitude term requires upwinding. The formula used here is given by:

\[u_P \abs{\nabla \phi}_P = \max \left( u_P , 0 \right) \left[ \sum_A \min \left( \frac{ \phi_A - \phi_P } { d_{AP}}, 0 \right)^2 \right]^{1/2} + \min \left( u_P , 0 \right) \left[ \sum_A \max \left( \frac{ \phi_A - \phi_P } { d_{AP}}, 0 \right)^2 \right]^{1/2}\]

Here are some simple test cases for this problem:

>>> from fipy.meshes import Grid1D
>>> from fipy.solvers import *
>>> SparseMatrix = LinearLUSolver()._matrixClass
>>> mesh = Grid1D(dx = 1., nx = 3)
>>> from fipy.variables.cellVariable import CellVariable

Trivial test:

>>> var = CellVariable(value = numerix.zeros(3, 'd'), mesh = mesh)
>>> v, L, b = FirstOrderAdvectionTerm(0.)._buildMatrix(var, SparseMatrix)
>>> print(numerix.allclose(b, numerix.zeros(3, 'd'), atol = 1e-10)) 
True

Less trivial test:

>>> var = CellVariable(value = numerix.arange(3), mesh = mesh)
>>> v, L, b = FirstOrderAdvectionTerm(1.)._buildMatrix(var, SparseMatrix)
>>> print(numerix.allclose(b, numerix.array((0., -1., -1.)), atol = 1e-10)) 
True

Even less trivial

>>> var = CellVariable(value = numerix.arange(3), mesh = mesh)
>>> v, L, b = FirstOrderAdvectionTerm(-1.)._buildMatrix(var, SparseMatrix)
>>> print(numerix.allclose(b, numerix.array((1., 1., 0.)), atol = 1e-10)) 
True

Another trivial test case (more trivial than a trivial test case standing on a harpsichord singing “trivial test cases are here again”)

>>> vel = numerix.array((-1, 2, -3))
>>> var = CellVariable(value = numerix.array((4, 6, 1)), mesh = mesh)
>>> v, L, b = FirstOrderAdvectionTerm(vel)._buildMatrix(var, SparseMatrix)
>>> print(numerix.allclose(b, -vel * numerix.array((2, numerix.sqrt(5**2 + 2**2), 5)), atol = 1e-10)) 
True

Somewhat less trivial test case:

>>> from fipy.meshes import Grid2D
>>> mesh = Grid2D(dx = 1., dy = 1., nx = 2, ny = 2)
>>> vel = numerix.array((3, -5, -6, -3))
>>> var = CellVariable(value = numerix.array((3, 1, 6, 7)), mesh = mesh)
>>> v, L, b = FirstOrderAdvectionTerm(vel)._buildMatrix(var, SparseMatrix)
>>> answer = -vel * numerix.array((2, numerix.sqrt(2**2 + 6**2), 1, 0))
>>> print(numerix.allclose(b, answer, atol = 1e-10)) 
True

Create a Term.

Parameters:

coeff (float or CellVariable or FaceVariable) – Coefficient for the term. FaceVariable objects are only acceptable for diffusion or convection terms.

property RHSvector

Return the RHS vector calculated in solve() or sweep(). The cacheRHSvector() method should be called before solve() or sweep() to cache the vector.

__eq__(other)

Return self==value.

__hash__()

Return hash(self).

__mul__(other)

Multiply a term

>>> 2. * __NonDiffusionTerm(coeff=0.5)
__NonDiffusionTerm(coeff=1.0)

Test for ticket:291.

>>> from fipy import PowerLawConvectionTerm
>>> PowerLawConvectionTerm(coeff=[[1], [0]]) * 1.0
PowerLawConvectionTerm(coeff=array([[ 1.],
       [ 0.]]))
__neg__()

Negate a Term.

>>> -__NonDiffusionTerm(coeff=1.)
__NonDiffusionTerm(coeff=-1.0)
__repr__()

The representation of a Term object is given by,

>>> print(__UnaryTerm(123.456))
__UnaryTerm(coeff=123.456)
__rmul__(other)

Multiply a term

>>> 2. * __NonDiffusionTerm(coeff=0.5)
__NonDiffusionTerm(coeff=1.0)

Test for ticket:291.

>>> from fipy import PowerLawConvectionTerm
>>> PowerLawConvectionTerm(coeff=[[1], [0]]) * 1.0
PowerLawConvectionTerm(coeff=array([[ 1.],
       [ 0.]]))
cacheMatrix()

Informs solve() and sweep() to cache their matrix so that matrix can return the matrix.

cacheRHSvector()

Informs solve() and sweep() to cache their right hand side vector so that getRHSvector() can return it.

justErrorVector(var=None, solver=None, boundaryConditions=(), dt=1.0, underRelaxation=None, residualFn=None)

Builds the Term’s linear system once.

This method also recalculates and returns the error as well as applying under-relaxation.

justErrorVector returns the overlapping local value in parallel (not the non-overlapping value).

>>> from fipy.solvers import DummySolver
>>> from fipy import *
>>> m = Grid1D(nx=10)
>>> v = CellVariable(mesh=m)
>>> len(DiffusionTerm().justErrorVector(v, solver=DummySolver())) == m.numberOfCells
True
Parameters:
  • var (CellVariable) – Variable to be solved for. Provides the initial condition, the old value and holds the solution on completion.

  • solver (Solver) – Iterative solver to be used to solve the linear system of equations. The default sovler depends on the solver package selected.

  • boundaryConditions (tuple of BoundaryCondition) –

  • dt (float) – Timestep size.

  • underRelaxation (float) – Usually a value between 0 and 1 or None in the case of no under-relaxation

  • residualFn (function) – Takes var, matrix, and RHSvector arguments, used to customize the residual calculation.

Returns:

error – The residual vector \(\vec{e}=\mathsf{L}\vec{x}_\text{old} - \vec{b}\)

Return type:

CellVariable

justResidualVector(var=None, solver=None, boundaryConditions=(), dt=None, underRelaxation=None, residualFn=None)

Builds the Term’s linear system once.

This method also recalculates and returns the residual as well as applying under-relaxation.

justResidualVector returns the overlapping local value in parallel (not the non-overlapping value).

>>> from fipy import *
>>> m = Grid1D(nx=10)
>>> v = CellVariable(mesh=m)
>>> len(numerix.asarray(DiffusionTerm().justResidualVector(v))) == m.numberOfCells
True
Parameters:
  • var (CellVariable) – Variable to be solved for. Provides the initial condition, the old value and holds the solution on completion.

  • solver (Solver) – Iterative solver to be used to solve the linear system of equations. The default sovler depends on the solver package selected.

  • boundaryConditions (tuple of BoundaryCondition) –

  • dt (float) – Timestep size.

  • underRelaxation (float) – Usually a value between 0 and 1 or None in the case of no under-relaxation

  • residualFn (function) – Takes var, matrix, and RHSvector arguments, used to customize the residual calculation.

Returns:

residual – The residual vector \(\vec{r}=\mathsf{L}\vec{x} - \vec{b}\)

Return type:

CellVariable

property matrix

Return the matrix calculated in solve() or sweep(). The cacheMatrix() method should be called before solve() or sweep() to cache the matrix.

residualVectorAndNorm(var=None, solver=None, boundaryConditions=(), dt=None, underRelaxation=None, residualFn=None)

Builds the Term’s linear system once.

This method also recalculates and returns the residual as well as applying under-relaxation.

Parameters:
  • var (CellVariable) – Variable to be solved for. Provides the initial condition, the old value and holds the solution on completion.

  • solver (Solver) – Iterative solver to be used to solve the linear system of equations. The default sovler depends on the solver package selected.

  • boundaryConditions (tuple of BoundaryCondition) –

  • dt (float) – Timestep size.

  • underRelaxation (float) – Usually a value between 0 and 1 or None in the case of no under-relaxation

  • residualFn (function) – Takes var, matrix, and RHSvector arguments, used to customize the residual calculation.

Returns:

  • residual (~fipy.variables.cellVariable.CellVariable) – The residual vector \(\vec{r}=\mathsf{L}\vec{x} - \vec{b}\)

  • norm (float) – The L2 norm of residual, \(\|\vec{r}\|_2\)

solve(var=None, solver=None, boundaryConditions=(), dt=None)

Builds and solves the Term’s linear system once. This method does not return the residual. It should be used when the residual is not required.

Parameters:
  • var (CellVariable) – Variable to be solved for. Provides the initial condition, the old value and holds the solution on completion.

  • solver (Solver) – Iterative solver to be used to solve the linear system of equations. The default sovler depends on the solver package selected.

  • boundaryConditions (tuple of BoundaryCondition) –

  • dt (float) – Timestep size.

sweep(var=None, solver=None, boundaryConditions=(), dt=None, underRelaxation=None, residualFn=None, cacheResidual=False, cacheError=False)

Builds and solves the Term’s linear system once. This method also recalculates and returns the residual as well as applying under-relaxation.

Parameters:
  • var (CellVariable) – Variable to be solved for. Provides the initial condition, the old value and holds the solution on completion.

  • solver (Solver) – Iterative solver to be used to solve the linear system of equations. The default sovler depends on the solver package selected.

  • boundaryConditions (tuple of BoundaryCondition) –

  • dt (float) – Timestep size.

  • underRelaxation (float) – Usually a value between 0 and 1 or None in the case of no under-relaxation

  • residualFn (function) – Takes var, matrix, and RHSvector arguments, used to customize the residual calculation.

  • cacheResidual (bool) – If True, calculate and store the residual vector \(\vec{r}=\mathsf{L}\vec{x} - \vec{b}\) in the residualVector member of Term

  • cacheError (bool) – If True, use the residual vector \(\vec{r}\) to solve \(\mathsf{L}\vec{e}=\vec{r}\) for the error vector \(\vec{e}\) and store it in the errorVector member of Term

Returns:

residual – The residual vector \(\vec{r}=\mathsf{L}\vec{x} - \vec{b}\)

Return type:

CellVariable

Last updated on Jun 26, 2024. Created using Sphinx 7.1.2.