Helper functions

Estimating stable time-steps

Helper functions for estimating stable time steps for RKDG methods.

Characteristic lengthscales

grudge.dt_utils.characteristic_lengthscales(actx: arraycontext.context.ArrayContext, dcoll: grudge.discretization.DiscretizationCollection) meshmode.dof_array.DOFArray[source]

Computes the characteristic length scale \(h_{\text{loc}}\) at each node. The characteristic length scale is mainly useful for estimating the stable time step size. E.g. for a hyperbolic system, an estimate of the stable time step can be estimated as \(h_{\text{loc}} / c\), where \(c\) is the characteristic wave speed. The estimate is obtained using the following formula:

\[h_{\text{loc}} = \operatorname{min}\left(\Delta r_i\right) r_D\]

where \(\operatorname{min}\left(\Delta r_i\right)\) is the minimum node distance on the reference cell (see dt_non_geometric_factors()), and \(r_D\) is the inradius of the cell (see dt_geometric_factors()).

Returns

a frozen DOFArray containing a characteristic lengthscale for each element, at each nodal location.

Note

While a prediction of stability is only meaningful in relation to a given time integrator with a known stability region, the lengthscale provided here is not intended to be specific to any one time integrator, though the stability region of standard four-stage, fourth-order Runge-Kutta methods has been used as a guide. Any concrete time integrator will likely require scaling of the values returned by this routine.

Non-geometric quantities

grudge.dt_utils.dt_non_geometric_factors(dcoll: grudge.discretization.DiscretizationCollection, dd=None) list[source]

Computes the non-geometric scale factors following [Hesthaven_2008], section 6.4, for each element group in the dd discretization:

\[c_{ng} = \operatorname{min}\left( \Delta r_i \right),\]

where \(\Delta r_i\) denotes the distance between two distinct nodal points on the reference element.

Parameters

dd – a DOFDesc, or a value convertible to one. Defaults to the base volume discretization if not provided.

Returns

a list of float values denoting the minimum node distance on the reference element for each group.

Mesh size utilities

grudge.dt_utils.dt_geometric_factors(dcoll: grudge.discretization.DiscretizationCollection, dd=None) meshmode.dof_array.DOFArray[source]

Computes a geometric scaling factor for each cell following [Hesthaven_2008], section 6.4, defined as the inradius (radius of an inscribed circle/sphere).

Specifically, the inradius for each element is computed using the following formula from [Shewchuk_2002], Table 1, for simplicial cells (triangles/tetrahedra):

\[r_D = \frac{d V}{\sum_{i=1}^{N_{faces}} F_i},\]

where \(d\) is the topological dimension, \(V\) is the cell volume, and \(F_i\) are the areas of each face of the cell.

Parameters

dd – a DOFDesc, or a value convertible to one. Defaults to the base volume discretization if not provided.

Returns

a frozen DOFArray containing the geometric factors for each cell and at each nodal location.

grudge.dt_utils.h_max_from_volume(dcoll: grudge.discretization.DiscretizationCollection, dim=None, dd=None) float[source]

Returns a (maximum) characteristic length based on the volume of the elements. This length may not be representative if the elements have very high aspect ratios.

Parameters
  • dim – an integer denoting topological dimension. If None, the spatial dimension specified by grudge.DiscretizationCollection.dim is used.

  • dd – a DOFDesc, or a value convertible to one. Defaults to the base volume discretization if not provided.

Returns

a scalar denoting the maximum characteristic length.

grudge.dt_utils.h_min_from_volume(dcoll: grudge.discretization.DiscretizationCollection, dim=None, dd=None) float[source]

Returns a (minimum) characteristic length based on the volume of the elements. This length may not be representative if the elements have very high aspect ratios.

Parameters
  • dim – an integer denoting topological dimension. If None, the spatial dimension specified by grudge.DiscretizationCollection.dim is used.

  • dd – a DOFDesc, or a value convertible to one. Defaults to the base volume discretization if not provided.

Returns

a scalar denoting the minimum characteristic length.

Array contexts

class grudge.array_context.PyOpenCLArrayContext(queue: pyopencl.CommandQueue, allocator: Optional[pyopencl.tools.AllocatorInterface] = None, wait_event_queue_length: Optional[int] = None, force_device_scalars: bool = False)[source]

Inherits from meshmode.array_context.PyOpenCLArrayContext. Extends it to understand grudge-specific transform metadata. (Of which there isn’t any, for now.)

class grudge.array_context.PytatoPyOpenCLArrayContext(queue, allocator=None)[source]

Inherits from meshmode.array_context.PytatoPyOpenCLArrayContext. Extends it to understand grudge-specific transform metadata. (Of which there isn’t any, for now.)