sofia package

Subpackages

Module contents

sofia.requires(*core_deps: str)

Build a decorator raising ImportError if given core dependencies are missing.

class sofia.Domain(transformation: CoordSys = global)

Bases: Basic

The position space where a field lives.

dim(coordinate: Symbol) → int

Return the dimension index of a coordinate.

vec(*components: Field)
property identity

Return a sum over the outer products of the unit vectors with themselves.

default_assumptions = {}
class sofia.Field(name, domain, *args, **kwargs)

Bases: Expr

Symbolic field living on Domain.

is_commutative = True
copy(new_name: str, append_name: bool = False)
as_dummy(name=None)

Return the expression with any objects having structurally bound symbols replaced with unique, canonical symbols within the object in which they appear and having only the default assumption for commutativity being True. When applied to a symbol a new symbol having only the same commutativity will be returned.

Examples

>>> from sympy import Integral, Symbol
>>> from sympy.abc import x
>>> r = Symbol('r', real=True)
>>> Integral(r, (r, x)).as_dummy()
Integral(_0, (_0, x))
>>> _.variables[0].is_real is None
True
>>> r.as_dummy()
_r

Notes

Any object that has structurally bound variables should have a property, bound_symbols that returns those symbols appearing in the object.

default_assumptions = {'commutative': True}
class sofia.VectorField(name: str | Symbol, domain: Domain, staggered: str | tuple[Staggered, ...] | None = None, **kwargs)
class sofia.VectorField(*terms: Vector, **options)

Bases: VectorAdd

A vector field composed of one scalar Field per axis.

It is a sympy vector so it can be used in any vector expression.

Parameters:
  • name – Each component field is named f”{name}_{axis}” (where axis is “x”, “y” or “z”).

  • domain – Domain the component fields live on.

  • staggered – Where the components sit within a cell on a Grid. This can be given either as a Staggered per component or as the string “faces”, which means each component lives on the face normal to its axis.

  • **kwargs – Forwarded to every component Field.

Examples

>>> from sofia import Grid
>>> G = Grid(ndim=2, name="G")
>>> u = VectorField("u", G)
>>> u
u_x*global.ihat + u_y*global.jhat
>>> u.fields
(u_x, u_y)
>>> [f.staggered for f in VectorField("m", G, staggered="faces").fields]
[XFACE, YFACE]
>>> u << 0
CodeBlock(
EvalAssignment(u_x, 0),
EvalAssignment(u_y, 0)
)
is_commutative = True
property fields: tuple[Field, ...]

The scalar field of each component.

assignments(other) → dict[Field, Expr]

Return the value each component field gets when other is assigned.

default_assumptions = {'commutative': True}
class sofia.TensorField(name: str | Symbol, domain: Domain, symmetric: bool = False, staggered: str | tuple[tuple[Staggered | None, ...], ...] | None = None, **kwargs)
class sofia.TensorField(*terms: Dyadic, **options)

Bases: DyadicAdd

An order 2 tensor field (dyadic) composed of one scalar Field per pair of axes.

It is a sympy dyadic, so it can be used in any dyadic expression.

Parameters:
  • name – Each component field is named f”{name}_{a}{b}” (e.g. T_xy).

  • domain – Domain the component fields live on.

  • symmetric – Use one field for T_ab and T_ba, so the tensor is symmetric by construction.

  • staggered – Where the components sit within a cell on a Grid. This can be given either as one Staggered per component, as rows (staggered[a][b] for component T_ab), or as the string “mac”, which puts the diagonal components at cell centers and the off-diagonal ones on the edges.

  • **kwargs – Forwarded to every component Field.

Examples

>>> from sofia import Grid
>>> G = Grid(ndim=2, name="G")
>>> T = TensorField("T", G)
>>> T.fields
(T_xx, T_xy, T_yx, T_yy)
>>> T.dot(G.unit_vectors[0])
T_xx*global.ihat + T_yx*global.jhat
>>> S = TensorField("S", G, symmetric=True)
>>> S.fields
(S_xx, S_xy, S_yy)
>>> S << 2
CodeBlock(
EvalAssignment(S_xx, 2),
EvalAssignment(S_xy, 2),
EvalAssignment(S_yy, 2)
)
is_commutative = True
property fields: tuple[Field, ...]

The distinct component fields.

property symmetric: bool

Whether T_ab and T_ba are the same field for every pair of axes.

assignments(other) → dict[Field, Expr]

Return the value each component field gets when other is assigned.

default_assumptions = {'commutative': True}
class sofia.GridField(name: str, domain: Grid, staggered: Staggered = Center, discretisation_order: int = 0, axes: Tuple | tuple[Symbol | int, ...] = None, _shifts: Tuple | tuple[Symbol | float | int, ...] | None = None)

Bases: Field

Symbolic representation of a field living on a structured grid.

A GridField is a symbol that carries the information needed to discretise expressions it appears in: the grid it is defined on, which of the grid axes it varies along, where within a cell its values sit, and how far it is shifted from its natural location.

Parameters:
  • name – Name of the field, used for its symbolic and generated-code names.

  • domain – The Grid the field is defined on.

  • staggered – Location of the field values within a cell, as a Staggered member. Defaults to Staggered.CENTER.

  • discretisation_order – Order of accuracy for derivatives and interpolation of this field. Must be even. The default of 0 means the order of domain is used instead.

  • axes – Grid axes the field varies along, given as indices or coordinates. Defaults to all of them. Passing a subset gives a field of lower dimensionality than the grid; passing () gives a scalar.

Examples

>>> import sympy as sp
>>> from sofia import Grid, Field, Staggered
>>> G = Grid(ndim=2)
>>> u = Field("u", domain=G)
>>> u.shape
(G.Nx, G.Ny)
>>> u.axes
(0, 1)
>>> Field("v", domain=G, axes=(G.x,)).shape
(G.Nx,)
>>> Field("dt", domain=G, axes=()).shape
()

Shifts:

>>> u.s(1, 0).shifts
(1, 0)
>>> u.s(1, 0).s(2, 0).shifts
(3, 0)
>>> u.s(j=-1).shifts
(0, -1)
is_commutative = True
copy(new_name: str, append_name: bool = False)
with_order(order)
to_staggered_origin()
reset_shifts()
s(*shifts, **kwargs)
indexify(bare=False)
fractional_offset(wrt) → Expr

Combine shift and staggered offset of a field.

interpolate()
default_assumptions = {'commutative': True}
class sofia.Grid(ndim: Literal[1, 2, 3], name: str = 'G', discretisation_order: int = 2, boundaries: Numeric | tuple[Numeric | tuple[Numeric, Numeric], ...] | None = None, _C: CoordSys | None = None, _shape: tuple[Symbol | int, ...] | None = None, *, transformations: __annotationlib_name_1__ | Callable | None = None, transformations_inverse: __annotationlib_name_2__ | Callable | None = None)

Bases: Domain

Symbolic representation of a structured grid.

This class defines a discretized grid in arbitrary dimensional space using symbolic quantities. A Grid instance serves multiple purposes when building a simulation code, but most directly it can be used to conveniently get symbolic coordinates and indices, as well as symbolic grid shape and size specifications.

Parameters:
  • ndim – Number of dimensions of the grid.

  • name – Name of the grid. This name is used to generate symbolic names for the grid coordinates, indices, and other symbolic quantities.

  • discretisation_order – Order of the discretisation scheme. Must be an even positive integer.

  • boundaries – Boundary specifications for each dimension in index space. Each dimension can be specified as either a single value (for periodic boundaries) or a pair of values (for non-periodic boundaries). The values represent the insets from first and last grid nodes. For example, a boundary specification of (0.5, -2.5) means that the boundary is at i = 0.5 and i = N - 1 - 2.5, where N is the number of grid points in that dimension.

  • transformation – A callable function that defines a transformation from the base coordinate system to a new coordinate system.

property uniform_axes

Axes whose scale factor is constant – these stay circulant/spectral.

property boundaries_absolute

Physical domain edges in index space.

property exterior_nodes

Number of nodes per dim that lie outside the physcal domain (ghost cells).

dim(coordinate_or_index: Symbol) → int

Return the dimension index of a coordinate.

coord_to_ind(coordinate: Symbol) → Symbol
dds_eval(L, boundary_left, boundary_right)
maps_coord_to_ind()
maps_ind_to_coord()
concrete(shape: tuple[int, ...] | int, extent: tuple[Numeric | tuple[Numeric, Numeric] | list[Numeric], ...] | list[Numeric | tuple[Numeric, Numeric] | list[Numeric]] | Numeric | None = None, subs: dict | None = None) → dict[Symbol, int | float]

Return the substitutions corresponding to the concrete numbers of the grid.

Parameters:
  • shape – Grid shape per dimension.

  • extent – Physical domain extents per dimension, given as concrete floats.

  • subs – Any additional substitutions.

Examples

>>> from sofia import Grid
>>> G = Grid(ndim=2)
>>> G.concrete(shape=(100, 200), extent=((0, 1), (0, 2)))
{G.Nx: 100, G.Ny: 200, G.Lx: 1, G.Ly: 2, G.x0: 0, G.y0: 0, global.x: i/100, global.y: j/100}
default_assumptions = {}
class sofia.NeighbourSum(*args)

Bases: Function

Sum a per-pair expression over the neighbours of each particle.

Parameters:
  • expr – Per-pair expression. field is the value at the particle itself, field.neighbour the value at the neighbour. Can be scalar or vector.

  • pos – Particle position, determines which particles count as neighbours. Only pairs closer than the domain’s smoothing_length contribute to the sum.

  • exclude_self – If True, skip the particle itself.

property cutoff: Basic

Condition for a neighbour to contribute to the sum.

A neighbour contributes if it is closer than smoothing_length to the particle.

default_assumptions = {}
class sofia.ParticleDomain(ndim: Literal[1, 2, 3] = 3, name: str = 'P', _n_particles: Expr | int | None = None, _extent: Tuple | None = None, _smoothing_length: Expr | float | None = None)

Bases: Domain

Symbolic representation of a particle domain.

A domain where fields are indexed by particle. Particles live in ndim-dimensional Cartesian space.

Parameters:
  • ndim – Spatial dimensionality.

  • name – Name used to identify the particle domain.

property smoothing_length: Expr

Cutoff radius of every NeighbourSum on this domain.

concrete(n_particles: int | None = None, extent: tuple[Numeric | tuple[Numeric, Numeric] | list[Numeric], ...] | list[Numeric | tuple[Numeric, Numeric] | list[Numeric]] | Numeric | None = None, smoothing_length: float | None = None, subs: dict | None = None) → dict[Symbol, int | float]

Return substitutions corresponding to the concrete numbers of the domain.

default_assumptions = {}
class sofia.ParticleField(name: str, domain: ParticleDomain)

Bases: Field

Symbolic field living on a ParticleDomain.

indexify()
default_assumptions = {'commutative': True}
is_commutative = True
class sofia.Scalar(name: str | Symbol, *, dtype: DTypeLike | None = None, **assumptions: bool | None)

Bases: AtomicExpr

Symbolic scalar.

Identified by name and assumptions, compared after sympy has drawn their consequences. If dtype is boolean, a BooleanScalar is returned instead.

Parameters:
  • name – Identifier the generated code declares this scalar under.

  • dtype – Data type. By default it is the backend’s float dtype.

  • assumptions – Sympy assumptions, restricted to SUPPORTED_ASSUMPTIONS.

is_commutative = True
is_symbol = True
is_Symbol = True
kind = NumberKind
default_assumptions = {'commutative': True}
class sofia.Staggered(offsets: dict[int, Expr | float] | None = None, name=None)

Bases: Basic

Fractional per-axis offsets of a field relative to cell centers.

CENTER: ClassVar[Staggered] = Center
XFACE: ClassVar[Staggered] = XFACE
YFACE: ClassVar[Staggered] = YFACE
ZFACE: ClassVar[Staggered] = ZFACE
staggered_offsets(ndim: int) → tuple[Expr, ...]
XY_EDGE = XY_EDGE
XZ_EDGE = XZ_EDGE
YZ_EDGE = YZ_EDGE
default_assumptions = {}
sofia.fields(names: str, domain: Domain, **kwargs: object) → tuple[Field, ...]

Create one or more Fields on a shared domain.

Equivalent to sympy.symbols(), but for Fields. Always returns a tuple for better type support.

Parameters:
  • names – Field names, comma- or whitespace-separated; see sympy.symbols().

  • domain – Domain the fields live on.

  • **kwargs – Forwarded to the Field subclass registered for domain.

Returns:

A tuple of the created Field(s).

Return type:

tuple[Field, …]

Examples

>>> (u,) = fields("u", domain=Grid(ndim=1))
>>> v, w = fields("v w", domain=Grid(ndim=2))
sofia.print_environment_info() → None

Print the environment info for reporting bugs, mimicking the jax equivalent.

sofia.scalars(names: str, **kwargs: object) → tuple[Scalar, ...]

Create one or more Scalars.

Equivalent to sympy.symbols(), but for Scalars. Always returns a tuple for better type support.

Parameters:
  • names – Scalar names, comma- or whitespace-separated; see sympy.symbols().

  • **kwargs – Keyword arguments passed to Scalar.

Returns:

A tuple of the created Scalar(s).

Return type:

tuple[Scalar, …]

Examples

>>> (a,) = scalars("a")
>>> b, c = scalars("b c", positive=True)