Variational forms

Functions and classes for representing and manipulating variational forms in an abstract way.

For a detailed user guide, see Assembling custom forms.

class pyiga.vform.BuiltinFuncExpr(funcname, expr)
class pyiga.vform.ConstExpr(value)

A constant scalar value.

pyiga.vform.Dx(expr, k, times=1, parametric=False)

Partial derivative of expr along the k-th coordinate axis.

class pyiga.vform.Expr

Abstract base class which all expressions derive from.

shape

the shape of the expression as a tuple, analogous to numpy.ndarray.shape. Scalar expressions have the empty tuple () as their shape.

T

For a matrix expression, return its transpose.

dot(x)

Returns the dot product of this with x; see dot() for semantics.

dx(k, times=1, parametric=False)

Compute a partial derivative. Equivalent to Dx().

find_vf()

Attempt to determine the ambient VForm from the expression tree.

is_matrix()

Returns True iff the expression is matrix-valued.

is_scalar()

Returns True iff the expression is scalar.

is_vector()

Returns True iff the expression is vector-valued.

ravel()

Ravel a matrix expression into a vector expression.

x

Return the first child expression.

y

Return the second child expression.

z

Return the third child expression.

class pyiga.vform.LiteralMatrixExpr(entries)

Matrix expression which is represented by a 2D array of individual expressions.

class pyiga.vform.LiteralVectorExpr(entries)

Vector expression which is represented by a list of individual expressions.

class pyiga.vform.MatMatExpr(A, B)
class pyiga.vform.MatVecExpr(A, x)
class pyiga.vform.NegExpr(expr)
class pyiga.vform.OuterProdExpr(x, y)
class pyiga.vform.PartialDerivExpr(basisfun, D, physical=False)

A scalar expression which refers to the value of a basis function or one of its partial derivatives.

find_vf()

Attempt to determine the ambient VForm from the expression tree.

without_derivs()

Return the underlying basis function without derivatives.

class pyiga.vform.ScalarOperExpr(oper, x, y)
class pyiga.vform.SurfaceMeasureExpr
class pyiga.vform.TensorOperExpr(oper, x, y)
class pyiga.vform.VForm(dim, geo_dim=None, arity=2, spacetime=False)

Abstract representation of a variational form.

See Assembling custom forms for a comprehensive user guide.

Parameters:
  • dim (int) – the space dimension
  • geo_dim (int) – the dimension of the image of the geometry map. By default, it is equal to dim. It should be either dim (for volume integrals) or dim + 1 (for surface integrals).
  • arity (int) – the arity of the variational form, i.e., 1 for a linear functional and 2 for a bilinear form
  • spacetime (bool) – whether the form describes a space-time discretization (deprecated)
add(expr)

Add an expression to this VForm.

all_exprs(type=None, once=True)

Deep, depth-first iteration of all expressions with dependencies.

If type is given, only exprs which are instances of that type are yielded. If once=True (default), each expr is visited only once.

basisfuns(components=(None, None), spaces=(0, 0))

Obtain expressions representing the basis functions for this vform.

Parameters:
  • components – for vector-valued problems, specify the number of components for each basis function here.
  • spaces – space indices for problems where the basis functions live in different spaces.
dependency_graph()

Compute a directed graph of the dependencies between all used variables.

finalize(do_precompute=True)

Performs standard transforms and dependency analysis.

indices_to_D(indices)

Convert a list of derivative indices into a partial derivative tuple D.

input(name, shape=(), physical=False, updatable=False)

Declare an input field with the given name and shape and return an expression representing it.

By default, input fields are assumed to be given in parametric coordinates. If the field is defined in physical geometry coordinates, pass physical=True.

If updatable is True, the generated assembler will allow updating of this field through an update(name=value) method.

let(name, expr, symmetric=False)

Define a variable with the given name which has the given expr as its value.

linearize_vars(vars)

Returns an admissible order for computing the given vars, i.e., which respects the dependency relation.

num_components()

Return number of vector components for each basis function space.

num_spaces()

Return the number of different function spaces used in this VForm.

substitute_vec_components(expr)

Given a single scalar expression in terms of vector basis functions, return a matrix of scalar expressions where each basis function has been substituted by (u,0,..,0), …, (0,…,0,u) successively (u being the underlying scalar basis function).

transform(fun, type=None, deep=True)

Apply fun to all exprs (or all exprs of the given type). If fun returns an expr, replace the old expr by this new one.

transitive_closure(dep_graph, vars, exclude={})

Linearized transitive closure (w.r.t. dependency) of the given vars.

transitive_deps(dep_graph, vars)

Return all vars on which the given vars depend directly or indirectly, in linearized order.

vars_without_dep_on(dep_graph, exclude)

Return a linearized list of all expr vars which do not depend on the given vars.

class pyiga.vform.VarRefExpr(var, I, D=None, parametric=False)

A scalar expression which refers to an entry of a variable or a derivative thereof.

find_vf()

Attempt to determine the ambient VForm from the expression tree.

without_derivs()

Return a reference to the underlying variable without derivatives.

class pyiga.vform.VectorCrossExpr(x, y)
class pyiga.vform.VolumeMeasureExpr
pyiga.vform.as_expr(x)

Interpret input as an expression; useful for constants.

pyiga.vform.as_matrix(x)

Convert a sequence of sequence of expressions to a matrix expression.

pyiga.vform.as_vector(x)

Convert a sequence of expressions to a vector expression.

pyiga.vform.cos(x)

Cosine of an expression.

pyiga.vform.cross(x, y)

Cross product of two 3D vectors.

pyiga.vform.curl(expr)

The curl (or rot) of a 3D vector expression.

pyiga.vform.det(A)

Determinant of a matrix.

pyiga.vform.div(expr, parametric=False)

The divergence of a vector-valued expressions, resulting in a scalar.

pyiga.vform.dot(a, b)

The dot product of two expressions.

Depending on the shapes of the arguments a and b, its semantics differ:

  • vector, vector: inner product (see inner())
  • matrix, vector: matrix-vector product
  • matrix, matrix: matrix-matrix product
pyiga.vform.ds = <pyiga.vform.SurfaceMeasureExpr object>

A symbolic expression representing the integration weight stemming from the geometry map for surface integrals.

pyiga.vform.dx = <pyiga.vform.VolumeMeasureExpr object>

A symbolic expression representing the integration weight stemming from the geometry map.

pyiga.vform.exp(x)

Exponential of an expression.

pyiga.vform.grad(expr, dims=None, parametric=False)

Gradient of an expression.

If expr is scalar, results in a vector of all partial derivatives.

If expr is a vector, results in the Jacobian matrix whose rows are the gradients of the individual components.

If dims is specified, it is a tuple of dimensions along which to take the derivative. By default, all space dimensions are used.

If parametric is true, the gradient with respect to the coordinates in the parameter domain is computed. By default, the gradient is computed in physical coordinates (transformed by the geometry map).

pyiga.vform.hess(expr, parametric=False)

Hessian matrix of a scalar expression.

If parametric is true, the Hessian with respect to the coordinates in the parameter domain is computed. By default, the Hessian is computed in physical coordinates (transformed by the geometry map).

pyiga.vform.inner(x, y)

The inner product of two vector or matrix expressions.

pyiga.vform.inv(A)

Inverse of a matrix.

pyiga.vform.iterexprs(exprs, deep=False, type=None, once=True)

Iterate through all subexpressions of the list of expressions exprs in depth-first order.

If deep=True, follow variable references. If type is given, only exprs which are instances of that type are yielded. If once=True (default), each expr is visited only once.

pyiga.vform.log(x)

Natural logarithm of an expression.

pyiga.vform.make_var_expr(vf, var)

Create an expression of the proper shape which refers to the variable var.

pyiga.vform.mapexprs(exprs, fun, deep=False)

Replace each expr e in a list of expr trees by fun(e), depth first.

If deep=True, follow variable references.

pyiga.vform.norm(x)

Euclidean norm of a vector.

pyiga.vform.outer(x, y)

Outer product of two vectors, resulting in a matrix.

pyiga.vform.sin(x)

Sine of an expression.

pyiga.vform.sqrt(x)

Square root of an expression.

pyiga.vform.tan(x)

Tangent of an expression.

pyiga.vform.tr(A)

Trace of a matrix.