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, dcoll: DiscretizationCollection, dd: DOFDesc | None = None) 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 (seedt_geometric_factors()
).- Returns:
a
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: DiscretizationCollection, dd: DOFDesc | None = None) Sequence[float] [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.
Mesh size utilities#
- grudge.dt_utils.dt_geometric_factors(dcoll: DiscretizationCollection, dd: DOFDesc | None = None) 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.
- grudge.dt_utils.h_max_from_volume(dcoll: DiscretizationCollection, dim=None, dd: DOFDesc | None = None) int | float | complex | generic [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: DiscretizationCollection, dim=None, dd: DOFDesc | None = None) int | float | complex | generic [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: pyopencl.tools.AllocatorBase | None = None, wait_event_queue_length: int | None = None, force_device_scalars: bool = False)[source]#
Inherits from
meshmode.array_context.PyOpenCLArrayContext
. Extends it to understandgrudge
-specific transform metadata. (Of which there isnβt any, for now.)
- class grudge.array_context.PytatoPyOpenCLArrayContext(queue, allocator=None, *, compile_trace_callback: Callable[[Any, str, Any], None] | None = None)[source]#
Inherits from
meshmode.array_context.PytatoPyOpenCLArrayContext
. Extends it to understandgrudge
-specific transform metadata. (Of which there isnβt any, for now.)
- class grudge.array_context.MPIPyOpenCLArrayContext(mpi_communicator, queue: pyopencl.CommandQueue, *, allocator: pyopencl.tools.AllocatorBase | None = None, wait_event_queue_length: int | None = None, force_device_scalars: bool = False)[source]#
An array context for using distributed computation with
pyopencl
eager evaluation.
- grudge.array_context.get_reasonable_array_context_class(lazy: bool = True, distributed: bool = True, fusion: bool | None = None) Type[ArrayContext] [source]#
Returns a reasonable
PyOpenCLArrayContext
currently supported given the constraints of lazy and distributed.
Miscellaneous tools#
- grudge.tools.build_jacobian(actx: ArrayContext, f: Callable[[ArrayOrContainerT], ArrayOrContainerT], base_state: ArrayOrContainerT, stepsize: float) ndarray [source]#
Returns a Jacobian matrix of f determined by a one-sided finite difference approximation with stepsize.
- Parameters:
f β a callable with a single argument, any array or
arraycontext.ArrayContainer
supported by actx.base_state β The point at which the Jacobian is to be calculated. May be any array or
arraycontext.ArrayContainer
supported by actx.
- Returns:
a two-dimensional
numpy.ndarray
.
- grudge.tools.map_subarrays(f: Callable[[Any], Any], in_shape: Tuple[int, ...], out_shape: Tuple[int, ...], ary: Any, *, return_nested: bool = False) Any [source]#
Apply a function f to subarrays of shape in_shape of an
numpy.ndarray
, typically (but not necessarily) of dtypeobject
. Return annumpy.ndarray
with the corresponding subarrays replaced by the return values of f, and with the shape adapted to reflect out_shape.Similar to
numpy.vectorize
.Example 1: given a function f that maps arrays of shape
(2, 2)
to scalars and an input array ary of shape(3, 2, 2)
, the call:map_subarrays(f, (2, 2), (), ary)
will produce an array of shape
(3,)
containing the results of calling f on the 3 subarrays of shape(2, 2)
in ary.Example 2: given a function f that maps arrays of shape
(2,)
to arrays of shape(2, 2)
and an input array ary of shape(3, 2)
, the call:map_subarrays(f, (2,), (2, 2), ary)
will produce an array of shape
(3, 2, 2)
containing the results of calling f on the 3 subarrays of shape(2,)
in ary. The call:map_subarrays(f, (2,), (2, 2), ary, return_nested=True)
will instead produce an array of shape
(3,)
with each entry containing an array of shape(2, 2)
.- Parameters:
f β the function to be called.
in_shape β the shape of any inputs to f.
out_shape β the shape of the result of calling f on an array of shape in_shape.
ary β a
numpy.ndarray
instance.return_nested β if out_shape is nontrivial, this flag indicates whether to return a nested array (containing one entry for each application of f), or to return a single, higher-dimensional array.
- Returns:
an array with the subarrays of shape in_shape replaced with subarrays of shape out_shape resulting from the application of f.
- grudge.tools.rec_map_subarrays(f: Callable[[Any], Any], in_shape: Tuple[int, ...], out_shape: Tuple[int, ...], ary: ArrayOrContainer, *, scalar_cls: type | Tuple[type] | None = None, return_nested: bool = False) ArrayOrContainer [source]#
Like
map_subarrays()
, but with support forarraycontext.ArrayContainer
s.- Parameters:
scalar_cls β An array container of this type will be considered a scalar and arrays of it will be passed to f without further destructuring.