Symbolic operator representation

Based on pymbolic.

Basic objects

Object types

Based on the mathematical quantity being represented, the following types of objects occur as part of a symbolic operator representation:

pyopencl.array.Array instances do not occur on the symbolic of pymbolic at all. Those hold per-node degrees of freedom (and only those), which is not visible as an array axis in symbolic code. (They’re visible only once evaluated.)

Placeholders

pytential.symbolic.primitives.var

alias of pymbolic.primitives.Variable

pytential.symbolic.primitives.make_sym_vector(name, components, var_factory=None)[source]

Return an object array of components subscripted Variable (or subclass) instances.

Parameters:
  • components – The number of components in the vector.
  • var_factory – The Variable subclass to use for instantiating the scalar variables.
pytential.symbolic.primitives.make_sym_mv(name, num_components)[source]
pytential.symbolic.primitives.make_sym_surface_mv(name, ambient_dim, dim, where=None)[source]

Functions

pytential.symbolic.primitives.real
pytential.symbolic.primitives.imag
pytential.symbolic.primitives.conj
pytential.symbolic.primitives.abs
pytential.symbolic.primitives.sqrt
pytential.symbolic.primitives.sin
pytential.symbolic.primitives.cos
pytential.symbolic.primitives.tan
pytential.symbolic.primitives.asin
pytential.symbolic.primitives.acos
pytential.symbolic.primitives.atan
pytential.symbolic.primitives.atan2
pytential.symbolic.primitives.sinh
pytential.symbolic.primitives.cosh
pytential.symbolic.primitives.tanh
pytential.symbolic.primitives.asinh
pytential.symbolic.primitives.acosh
pytential.symbolic.primitives.atanh
pytential.symbolic.primitives.exp
pytential.symbolic.primitives.log

Discretization properties

class pytential.symbolic.primitives.QWeight(where=None)[source]

Bare quadrature weights (without Jacobians).

Parameters:where – A symbolic name for a pytential.discretization.Discretization
pytential.symbolic.primitives.nodes(ambient_dim, where=None)[source]

Return a pymbolic.geometric_algebra.MultiVector of node locations.

pytential.symbolic.primitives.parametrization_derivative(ambient_dim, dim, where=None)[source]

Return a pymbolic.geometric_algebra.MultiVector representing the derivative of the reference-to-global parametrization.

pytential.symbolic.primitives.parametrization_derivative_matrix(ambient_dim, dim, where=None)[source]

Return a np.array representing the derivative of the reference-to-global parametrization.

pytential.symbolic.primitives.pseudoscalar(ambient_dim, dim=None, where=None)[source]

Same as the outer product of all parametrization derivative columns.

pytential.symbolic.primitives.area_element(ambient_dim, dim=None, where=None)[source]
pytential.symbolic.primitives.sqrt_jac_q_weight(ambient_dim, dim=None, where=None)[source]
pytential.symbolic.primitives.normal(ambient_dim, dim=None, where=None)[source]

Exterior unit normals.

pytential.symbolic.primitives.mean_curvature(ambient_dim, dim=None, where=None)[source]

(Numerical) mean curvature.

pytential.symbolic.primitives.first_fundamental_form(ambient_dim, dim=None, where=None)[source]
pytential.symbolic.primitives.second_fundamental_form(ambient_dim, dim=None, where=None)[source]

Compute the second fundamental form of a surface. This is in reference to the reference-to-global mapping in use for each element.

Note

Some references assume that the second fundamental form is computed with respect to an orthonormal basis, which this is not.

pytential.symbolic.primitives.shape_operator(ambient_dim, dim=None, where=None)[source]

Elementary numerics

class pytential.symbolic.primitives.NumReferenceDerivative(ref_axes, operand, where=None)[source]

An operator that takes a derivative of operand with respect to the the element reference coordinates.

Parameters:
  • ref_axes

    a tuple of tuples indicating indices of coordinate axes of the reference element to the number of derivatives which will be taken. For example, the value ((0, 2), (1, 1)) indicates that Each axis must occur at most once. The tuple must be sorted by the axis index.

    May also be a singile integer i, which is viewed as equivalent to ((i, 1),).

  • where – A symbolic name for a pytential.discretization.Discretization
class pytential.symbolic.primitives.NodeSum(operand)[source]

Implements a global sum over all discretization nodes.

class pytential.symbolic.primitives.NodeMax(operand)[source]

Implements a global maximum over all discretization nodes.

class pytential.symbolic.primitives.ElementwiseSum(operand, where=None)[source]

Returns a vector of DOFs with all entries on each element set to the sum of DOFs on that element.

class pytential.symbolic.primitives.ElementwiseMax(operand, where=None)[source]

Returns a vector of DOFs with all entries on each element set to the maximum of DOFs on that element.

pytential.symbolic.primitives.integral(ambient_dim, dim, operand, where=None)[source]

A volume integral of operand.

class pytential.symbolic.primitives.Ones(where=None)[source]

A DOF-vector that is constant one on the whole discretization.

pytential.symbolic.primitives.ones_vec(dim, where=None)[source]
pytential.symbolic.primitives.area(ambient_dim, dim, where=None)[source]
pytential.symbolic.primitives.mean(ambient_dim, dim, operand, where=None)[source]
class pytential.symbolic.primitives.IterativeInverse(expression, rhs, variable_name, extra_vars={}, where=None)[source]

Geometric Calculus (based on Geometric/Clifford Algebra)

class pytential.symbolic.primitives.Derivative[source]

Conventional Calculus

pytential.symbolic.primitives.dd_axis(axis, ambient_dim, operand)[source]

Return the derivative along (XYZ) axis axis (in ambient_dim-dimensional space) of operand.

pytential.symbolic.primitives.d_dx()
pytential.symbolic.primitives.d_dy()
pytential.symbolic.primitives.d_dz()
pytential.symbolic.primitives.grad_mv(ambient_dim, operand)[source]

Return the ambient_dim-dimensional gradient of operand as a pymbolic.geometric_algebra.MultiVector.

pytential.symbolic.primitives.grad(ambient_dim, operand)[source]

Return the ambient_dim-dimensional gradient of operand as a numpy.ndarray.

pytential.symbolic.primitives.laplace(ambient_dim, operand)[source]

Layer potentials

class pytential.symbolic.primitives.IntG(kernel, density, qbx_forced_limit, source=None, target=None, kernel_arguments=None, **kwargs)[source]
\[\int_\Gamma g_k(x-y) \sigma(y) dS_y\]

where \(\sigma\) is density.

target_derivatives and later arguments should be considered keyword-only.

Parameters:
  • kernel – a kernel as accepted by sumpy.kernel.to_kernel_and_args(), likely a sumpy.kernel.Kernel.
  • qbx_forced_limit

    +1 if the output is required to originate from a QBX center on the “+” side of the boundary. -1 for the other side. Evaluation at a target with a value of +/- 1 in qbx_forced_limit will fail if no QBX center is found.

    +2 may be used to allow evaluation QBX center on the “+” side of the (but disallow evaluation using a center on the “-” side). Potential evaluation at the target still succeeds if no applicable QBX center is found. (-2 for the analogous behavior on the “-” side.)

    None may be used to avoid expressing a side preference for close evaluation.

    'avg' may be used as a shorthand to evaluate this potential as an average of the +1 and the -1 value.

  • kernel_arguments – A dictionary mapping named sumpy.kernel.Kernel arguments (see sumpy.kernel.Kernel.get_args() and sumpy.kernel.Kernel.get_source_args()) to expressions that determine them
  • source – The symbolic name of the source discretization. This name is bound to a concrete pytential.source.LayerPotentialSourceBase by pytential.bind().
  • target – The symbolic name of the set of targets. This name gets assigned to a concrete target set by pytential.bind().

kwargs has the same meaning as kernel_arguments can be used as a more user-friendly interface.

pytential.symbolic.primitives.int_g_dsource(ambient_dim, dsource, kernel, density, qbx_forced_limit, source=None, target=None, kernel_arguments=None, **kwargs)[source]
\[\int_\Gamma \operatorname{dsource} \dot \nabla_y \dot g(x-y) \sigma(y) dS_y\]

where \(\sigma\) is density, and dsource, a multivector. Note that the first product in the integrand is a geometric product.

pytential.symbolic.primitives.dsource

A pymbolic.geometric_algebra.MultiVector.

pytential.symbolic.primitives.S(kernel, density, qbx_forced_limit=<class 'pytential.symbolic.primitives._unspecified'>, source=None, target=None, kernel_arguments=None, **kwargs)[source]
pytential.symbolic.primitives.Sp(kernel, *args, **kwargs)[source]
pytential.symbolic.primitives.Spp(kernel, *args, **kwargs)[source]
pytential.symbolic.primitives.D(kernel, *args, **kwargs)[source]
pytential.symbolic.primitives.Dp(kernel, *args, **kwargs)[source]
pytential.symbolic.primitives.normal_derivative(ambient_dim, operand, dim=None, where=None)[source]
pytential.symbolic.primitives.tangential_derivative(ambient_dim, operand, dim=None, where=None)[source]

“Conventional” Vector Calculus

pytential.symbolic.primitives.tangential_onb(ambient_dim, dim=None, where=None)[source]

Return a matrix of shape (ambient_dim, dim) with orthogonal columns spanning the tangential space of the surface of where.

pytential.symbolic.primitives.xyz_to_tangential(xyz_vec, where=None)[source]
pytential.symbolic.primitives.tangential_to_xyz(tangential_vec, where=None)[source]
pytential.symbolic.primitives.project_to_tangential(xyz_vec, where=None)[source]
pytential.symbolic.primitives.cross(vec_a, vec_b)[source]
pytential.symbolic.primitives.n_dot(vec, where=None)[source]
pytential.symbolic.primitives.n_cross(vec, where=None)[source]
pytential.symbolic.primitives.curl(vec)[source]
pytential.symbolic.primitives.pretty(expr)[source]

Pretty-printing expressions

pytential.symbolic.primitives.pretty(expr)[source]

Binding an operator to a discretization

class pytential.GeometryCollection(places, auto_where=None)[source]

A mapping from symbolic identifiers (“place IDs”, typically strings) to ‘geometries’, where a geometry can be a pytential.source.PotentialSource or a pytential.target.TargetBase. This class is meant to hold a specific combination of sources and targets serve to host caches of information derived from them, e.g. FMM trees of subsets of them, as well as related common subexpressions such as metric terms.

__getitem__()[source]
get_discretization()[source]
get_cache()[source]
Parameters:
  • places – a scalar, tuple of or mapping of symbolic names to geometry objects. Supported objects are PotentialSource, TargetBase and Discretization.
  • auto_where – location identifier for each geometry object, used to denote specific discretizations, e.g. in the case where places is a LayerPotentialSourceBase. By default, we assume DEFAULT_SOURCE and DEFAULT_TARGET for sources and targets, respectively.
pytential.bind(places, expr, auto_where=None)[source]
Parameters:
  • places – a pytential.symbolic.execution.GeometryCollection. Alternatively, any list or mapping that is a valid argument for its constructor can also be used.
  • auto_where – for simple source-to-self or source-to-target evaluations, find ‘where’ attributes automatically.

PDE operators

Scalar PDEs

class pytential.symbolic.pde.scalar.L2WeightedPDEOperator(kernel, use_l2_weighting)[source]
class pytential.symbolic.pde.scalar.DirichletOperator(kernel, loc_sign, alpha=None, use_l2_weighting=False, kernel_arguments=None)[source]

IE operator and field representation for solving Dirichlet boundary value problems with scalar kernels (e.g. sumpy.kernel.LaplaceKernel, sumpy.kernel.HelmholtzKernel, sumpy.kernel.YukawaKernel)

Note

When testing this as a potential matcher, note that it can only access potentials that come from charge distributions having no net charge. (This is true at least in 2D.)

is_unique_only_up_to_constant()[source]
representation(u, map_potentials=None, qbx_forced_limit=None)[source]
operator(u)[source]
Parameters:
  • loc_sign – +1 for exterior, -1 for interior
  • alpha – the coefficient for the combined-field representation Set to 0 for Laplace.
class pytential.symbolic.pde.scalar.NeumannOperator(kernel, loc_sign, alpha=None, use_improved_operator=True, laplace_kernel=0, use_l2_weighting=False, kernel_arguments=None)[source]

IE operator and field representation for solving Dirichlet boundary value problems with scalar kernels (e.g. sumpy.kernel.LaplaceKernel, sumpy.kernel.HelmholtzKernel, sumpy.kernel.YukawaKernel)

Note

When testing this as a potential matcher, note that it can only access potentials that come from charge distributions having no net charge. (This is true at least in 2D.)

is_unique_only_up_to_constant()[source]
representation(u, map_potentials=None, qbx_forced_limit=None, **kwargs)[source]
operator(u)[source]
Parameters:
  • loc_sign – +1 for exterior, -1 for interior
  • alpha – the coefficient for the combined-field representation Set to 0 for Laplace.
  • use_improved_operator – Whether to use the least singular operator available

Maxwell’s equations

pytential.symbolic.pde.maxwell.get_sym_maxwell_point_source(kernel, jxyz, k)[source]

Return a symbolic expression that, when bound to a pytential.source.PointPotentialSource will yield a field satisfying Maxwell’s equations.

Uses the sign convention \(\exp(-1 \omega t)\) for the time dependency.

This will return an object of six entries, the first three of which represent the electric, and the second three of which represent the magnetic field. This satisfies the time-domain Maxwell’s equations as verified by sumpy.point_calculus.frequency_domain_maxwell().

pytential.symbolic.pde.maxwell.get_sym_maxwell_plane_wave(amplitude_vec, v, omega, epsilon=1, mu=1, where=None)[source]

Return a symbolic expression that, when bound to a pytential.source.PointPotentialSource will yield a field satisfying Maxwell’s equations.

Parameters:
  • amplitude_vec – should be orthogonal to v. If it is not, it will be orthogonalized.
  • v – a three-vector representing the phase velocity of the wave (may be an object array of variables or a vector of concrete numbers) While v may mathematically be complex-valued, this function is for now only tested for real values.
  • omega – Accepts the “Helmholtz k” to be compatible with other parts of this module.

Uses the sign convention \(\exp(-1 \omega t)\) for the time dependency.

This will return an object of six entries, the first three of which represent the electric, and the second three of which represent the magnetic field. This satisfies the time-domain Maxwell’s equations as verified by sumpy.point_calculus.frequency_domain_maxwell().

class pytential.symbolic.pde.maxwell.PECChargeCurrentMFIEOperator(k=Variable('k'))[source]

Magnetic Field Integral Equation operator with PEC boundary conditions, under the assumption of no surface charges.

See contrib/notes/mfie.tm in the repository for a derivation.

The arguments loc below decide between the exterior and the interior MFIE. The “exterior” (loc=1) MFIE enforces a zero field on the interior of the integration surface, whereas the “interior” MFIE (loc=-1) enforces a zero field on the exterior.

Uses the sign convention \(\exp(-1 \omega t)\) for the time dependency.

for the sinusoidal time dependency.

j_operator(loc, Jt)[source]
j_rhs(Hinc_xyz)[source]
rho_operator(loc, rho)[source]
rho_rhs(Jt, Einc_xyz)[source]
scattered_volume_field(Jt, rho, qbx_forced_limit=None)[source]

This will return an object array of six entries, the first three of which represent the electric, and the second three of which represent the magnetic field. This satisfies the time-domain Maxwell’s equations as verified by sumpy.point_calculus.frequency_domain_maxwell().

Internal affairs

Mappers

How a symbolic operator gets executed