Finite Differences

class Diff(argument, target=-1, superscript=-1)

Sympy Node representing a derivative.

The difference to sympy’s built in differential is:
  • shortened latex representation

  • all simplifications have to be done manually

  • optional marker displayed as superscript


Returns the argument the derivative acts on, for nested derivatives the inner argument is returned


Returns a Diff node with the given ‘new_arg’ instead of the current argument. For nested derivatives a new nested derivative is returned where the inner Diff has the ‘new_arg’


Applies linearity property of Diff: i.e. ‘Diff(c*a+b)’ is transformed to ‘c * Diff(a) + Diff(b)’ The parameter functions is a list of all symbols that are considered functions, not constants. For the example above: functions=[a, b]

property arg

Expression the derivative acts on

property target

Subscript, usually the variable the Diff is w.r.t.

property superscript

Superscript, for example used as the Chapman-Enskog order index

interpolated_access(offset, **kwargs)

Represents an interpolated access on a spatially differentiated field


offset (Tuple[sympy.Expr]) – Absolute position to determine the value of the spatial derivative

diff(expr, *args)

Shortcut function to create nested derivatives

>>> f = sp.Symbol("f")
>>> diff(f, 0, 0, 1) == Diff(Diff( Diff(f, 1), 0), 0)
class DiffOperator(target=-1, superscript=-1)

Un-applied differential, i.e. differential operator

  • target – the differential is w.r.t to this variable. This target is mainly for display purposes (its the subscript) and to distinguish DiffOperators If the target is ‘-1’ no subscript is displayed

  • superscript – optional marker displayed as superscript is not displayed if set to ‘-1’

The DiffOperator behaves much like a variable with special name. Its main use is to be applied later, using the DiffOperator.apply(expr, arg) which transforms ‘DiffOperator’s to applied ‘Diff’s

static apply(expr, argument, apply_to_constants=True)

Returns a new expression where each ‘DiffOperator’ is replaced by a ‘Diff’ node. Multiplications of ‘DiffOperator’s are interpreted as nested application of differentiation: i.e. DiffOperator(‘x’)*DiffOperator(‘x’) is a second derivative replaced by Diff(Diff(arg, x), t)


Returns set of all derivatives in an expression.

This function yields different results than ‘expr.atoms(Diff)’ when nested derivatives are in the expression, since this function only returns the outer derivatives


>>> x, y = sp.symbols("x, y")
>>> diff_terms( diff(x, 0, 0) )
{Diff(Diff(x, 0, -1), 0, -1)}
>>> diff_terms( diff(x, 0, 0) + y )
{Diff(Diff(x, 0, -1), 0, -1)}

Rewrites expression into a sum of distinct derivatives with pre-factors

zero_diffs(expr, label)

Replaces all differentials with the given target by 0


>>> x, y, f = sp.symbols("x y f")
>>> expression = Diff(f, x) + Diff(f, y) + Diff(Diff(f, y), x) + 7
>>> zero_diffs(expression, x)
Diff(f, y, -1) + 7
evaluate_diffs(expr, var=None)

Replaces pystencils diff objects by sympy diff objects and evaluates them.

Replaces Diff nodes by sp.diff , the free variable is either the target (if var=None) otherwise the specified var

normalize_diff_order(expression, functions=None, constants=None, sort_key=<function _default_diff_sort_key>)

Assumes order of differentiation can be exchanged. Changes the order of nested Diffs to a standard order defined by the sorting key ‘sort_key’ such that the derivative terms can be further simplified

expand_diff_linear(expr, functions=None, constants=None)

Expands all derivative nodes by applying Diff.split_linear

  • expr – expression containing derivatives

  • functions – sequence of symbols that are considered functions and can not be pulled before the derivative. if None, all symbols are viewed as functions

  • constants – sequence of symbols which are considered constants and can be pulled before the derivative


Fully expands all derivatives by applying product rule


Inverse product rule

functional_derivative(functional, v)

Computes functional derivative of functional with respect to v using Euler-Lagrange equation

\[\frac{\delta F}{\delta v} = \frac{\partial F}{\partial v} - \nabla \cdot \frac{\partial F}{\partial \nabla v}\]
  • assumes that gradients are represented by Diff() node

  • Diff(Diff(r)) represents the divergence of r

  • the constants parameter is a list with symbols not affected by the derivative. This is used for simplification of the derivative terms.

advection(advected_scalar, velocity_field, idx=None)

Advection term ∇·(velocity_field · advected_scalar)

Term that describes the advection of a scalar quantity in a velocity field.

diffusion(scalar, diffusion_coeff, idx=None)

Diffusion term ∇·( diffusion_coeff · ∇(scalar))


>>> f = Field.create_generic('f', spatial_dimensions=2)
>>> d = sp.Symbol("d")
>>> dx = sp.Symbol("dx")
>>> diffusion_term = diffusion(scalar=f, diffusion_coeff=d)
>>> discretization = Discretization2ndOrder()
>>> expected_output = ((f[-1, 0] + f[0, -1] - 4 * f[0, 0] + f[0, 1] + f[1, 0]) * d) / dx**2
>>> sp.simplify(discretization(diffusion_term) - expected_output)
transient(scalar, idx=None)

Transient term ∂_t(scalar)

class FVM1stOrder(field, flux=0, source=0)

Finite-volume discretization

  • field (Field) – the field with the quantity to calculate, e.g. a concentration

  • flux – a list of sympy expressions that specify the flux, one for each cartesian direction

  • source – a list of sympy expressions that specify the source


Return a list of assignments for the discrete fluxes


flux_field (Field) – a staggered field to which the fluxes should be assigned


Return a list of assignments for the discrete source term


Return a list of assignments for the continuity equation, which includes the source term


flux_field (Field) – a staggered field from which the fluxes are taken

VOF(j, v, ρ)

Volume-of-fluid discretization of advection

  • j (Field) – the staggered field to write the fluxes to. Should have a D2Q9/D3Q27 stencil. Other stencils work too, but incur a small error (D2Q5/D3Q7: v^2, D3Q19: v^3).

  • v (Field) – the flow velocity field

  • ρ (Field) – the quantity to advect