SpECTRE  v2023.09.07
evolution::dg Namespace Reference

Functionality for evolving hyperbolic partial differential equations using the discontinuous Galerkin method. More...

Namespaces

namespace  Actions
 Actions for using the discontinuous Galerkin to evolve hyperbolic partial differential equations.
 
namespace  Initialization
 Functionality for initializing the discontinuous Galerkin to evolve hyperbolic partial differential equations.
 
namespace  subcell
 Implementation of a generic finite volume/conservative finite difference subcell limiter.
 
namespace  Tags
 Tags used for DG evolution scheme.
 

Classes

struct  ApplyBoundaryCorrections
 Apply corrections from boundary communication. More...
 
struct  BackgroundGrVars
 Allocate or assign background general relativity quantities needed for evolution systems run on a curved spacetime without solving Einstein equations (e.g. ValenciaDivclean, ForceFree). More...
 
struct  BoundaryMessage
 [Charm++ Message] (https://charm.readthedocs.io/en/latest/charm%2B%2B/manual.html#messages) intended to be used in receive_data calls on the elements to send boundary data from one element on one node, to a different element on a (potentially) different node. More...
 
class  MortarData
 Data on the mortar used to compute the boundary correction for the DG scheme. More...
 
struct  using_subcell
 If Metavars has a SubcellOptions member struct and SubcellOptions::subcell_enabled is true then inherits from std::true_type, otherwise inherits from std::false_type. More...
 

Functions

template<typename Metavariables , typename DbTagsList , typename... InboxTags>
bool receive_boundary_data_global_time_stepping (const gsl::not_null< db::DataBox< DbTagsList > * > box, const gsl::not_null< tuples::TaggedTuple< InboxTags... > * > inboxes)
 Receive boundary data for global time-stepping. Returns true if all necessary data has been received.
 
template<typename Metavariables , bool DenseOutput, typename DbTagsList , typename... InboxTags>
bool receive_boundary_data_local_time_stepping (const gsl::not_null< db::DataBox< DbTagsList > * > box, const gsl::not_null< tuples::TaggedTuple< InboxTags... > * > inboxes)
 Receive boundary data for local time-stepping. Returns true if all necessary data has been received. More...
 
template<size_t Dim, typename DtTagsList >
void interpolate_dt_terms_gauss_points (const gsl::not_null< Variables< DtTagsList > * > dt_vars, const Mesh< Dim > &volume_mesh, const Direction< Dim > &direction, const Variables< DtTagsList > &dt_corrections)
 Interpolate the Bjorhus/time derivative corrections to the volume time derivatives in the specified direction. More...
 
template<size_t Dim, typename DtTagsList , typename BoundaryCorrectionTagsList >
void lift_boundary_terms_gauss_points (const gsl::not_null< Variables< DtTagsList > * > dt_vars, const Scalar< DataVector > &volume_det_inv_jacobian, const Mesh< Dim > &volume_mesh, const Direction< Dim > &direction, const Variables< BoundaryCorrectionTagsList > &boundary_corrections, const Scalar< DataVector > &magnitude_of_face_normal, const Scalar< DataVector > &face_det_jacobian)
 Lift the boundary corrections to the volume time derivatives in the specified direction. More...
 
template<size_t Dim, typename DtTagsList , typename BoundaryCorrectionTagsList >
void lift_boundary_terms_gauss_points (const gsl::not_null< Variables< DtTagsList > * > dt_vars, const Scalar< DataVector > &volume_det_inv_jacobian, const Mesh< Dim > &volume_mesh, const size_t dimension, const Variables< BoundaryCorrectionTagsList > &upper_boundary_corrections, const Scalar< DataVector > &upper_magnitude_of_face_normal, const Scalar< DataVector > &upper_face_det_jacobian, const Variables< BoundaryCorrectionTagsList > &lower_boundary_corrections, const Scalar< DataVector > &lower_magnitude_of_face_normal, const Scalar< DataVector > &lower_face_det_jacobian)
 Lift both the upper and lower (in logical coordinates) boundary corrections to the volume time derivatives in the specified logical dimension. More...
 
template<size_t Dim>
bool operator== (const BoundaryMessage< Dim > &lhs, const BoundaryMessage< Dim > &rhs)
 
template<size_t Dim>
bool operator!= (const BoundaryMessage< Dim > &lhs, const BoundaryMessage< Dim > &rhs)
 
template<size_t Dim>
std::ostreamoperator<< (std::ostream &os, const BoundaryMessage< Dim > &message)
 
template<size_t Dim>
bool operator!= (const MortarData< Dim > &lhs, const MortarData< Dim > &rhs)
 
template<size_t Dim>
std::ostreamoperator<< (std::ostream &os, const MortarData< Dim > &mortar_data)
 
template<typename VolumeVarsTagsList , typename FaceVarsTagsList , size_t Dim>
void project_contiguous_data_to_boundary (const gsl::not_null< Variables< FaceVarsTagsList > * > face_fields, const Variables< VolumeVarsTagsList > &volume_fields, const Mesh< Dim > &volume_mesh, const Direction< Dim > &direction)
 Projects a Variables of volume data to a contiguous subset of a boundary Variables More...
 
template<typename TagsToProjectList , typename VolumeVarsTagsList , typename FaceVarsTagsList , size_t Dim>
void project_tensors_to_boundary (const gsl::not_null< Variables< FaceVarsTagsList > * > face_fields, const Variables< VolumeVarsTagsList > &volume_fields, const Mesh< Dim > &volume_mesh, const Direction< Dim > &direction)
 Projects a subset of the tensors in the volume_fields onto the face. More...
 
template<typename Symm , typename IndexList , size_t Dim>
void project_tensor_to_boundary (const gsl::not_null< Tensor< DataVector, Symm, IndexList > * > face_field, const Tensor< DataVector, Symm, IndexList > &volume_field, const Mesh< Dim > &volume_mesh, const Direction< Dim > &direction)
 Projects a tensor to the face. More...
 

Variables

template<typename Metavars >
constexpr bool using_subcell_v = using_subcell<Metavars>::value
 If Metavars has a SubcellOptions member struct and SubcellOptions::subcell_enabled is true then is true, otherwise false. More...
 

Detailed Description

Functionality for evolving hyperbolic partial differential equations using the discontinuous Galerkin method.

Function Documentation

◆ interpolate_dt_terms_gauss_points()

template<size_t Dim, typename DtTagsList >
void evolution::dg::interpolate_dt_terms_gauss_points ( const gsl::not_null< Variables< DtTagsList > * >  dt_vars,
const Mesh< Dim > &  volume_mesh,
const Direction< Dim > &  direction,
const Variables< DtTagsList > &  dt_corrections 
)

Interpolate the Bjorhus/time derivative corrections to the volume time derivatives in the specified direction.

The general interpolation term (for the \(+\xi\)-dimension) is:

\begin{align*} \partial_t u_{\alpha\breve{\imath}\breve{\jmath}\breve{k}}=\cdots +\ell^{\mathrm{Gauss-Lobatto}}_{N} \left(\xi_{\breve{\imath}}^{\mathrm{Gauss}}\right) \partial_t u^{\mathrm{BC}}_{\alpha\breve{\jmath}\breve{k}}, \end{align*}

where \(\breve{\imath}\), \(\breve{\jmath}\), and \(\breve{k}\) are indices in the logical \(\xi\), \(\eta\), and \(\zeta\) dimensions. \(\partial_t u^{\mathrm{BC}}\) is the time derivative correction, and the function Spectral::boundary_interpolation_term() is used to compute and cache the terms from the lifting terms.

◆ lift_boundary_terms_gauss_points() [1/2]

template<size_t Dim, typename DtTagsList , typename BoundaryCorrectionTagsList >
void evolution::dg::lift_boundary_terms_gauss_points ( const gsl::not_null< Variables< DtTagsList > * >  dt_vars,
const Scalar< DataVector > &  volume_det_inv_jacobian,
const Mesh< Dim > &  volume_mesh,
const Direction< Dim > &  direction,
const Variables< BoundaryCorrectionTagsList > &  boundary_corrections,
const Scalar< DataVector > &  magnitude_of_face_normal,
const Scalar< DataVector > &  face_det_jacobian 
)

Lift the boundary corrections to the volume time derivatives in the specified direction.

The general lifting term (for the \(\xi\)-dimension) is:

\begin{align*} \partial_t u_{\alpha\breve{\imath}\breve{\jmath}\breve{k}}=\cdots -\frac{\ell_{\breve{\imath}}\left(\xi=1\right)} {w_{\breve{\imath}}J_{\breve{\imath}\breve{\jmath}\breve{k}}} \left[J\sqrt{ \frac{\partial\xi}{\partial x^i} \gamma^{ij} \frac{\partial\xi}{\partial x^j}} \left(G_{\alpha} + D_{\alpha}\right) \right]_{\breve{\jmath}\breve{k}}\left(\xi=1\right), \end{align*}

where \(\breve{\imath}\), \(\breve{\jmath}\), and \(\breve{k}\) are indices in the logical \(\xi\), \(\eta\), and \(\zeta\) dimensions. The \(G+D\) terms correspond to the boundary_corrections function argument, and the function Spectral::boundary_lifting_term() is used to compute and cache the terms from the lifting terms \(\ell_{\breve{\imath}}(\xi=\pm1)/w_{\breve{\imath}}\).

Note
that normal vectors are pointing out of the element.

◆ lift_boundary_terms_gauss_points() [2/2]

template<size_t Dim, typename DtTagsList , typename BoundaryCorrectionTagsList >
void evolution::dg::lift_boundary_terms_gauss_points ( const gsl::not_null< Variables< DtTagsList > * >  dt_vars,
const Scalar< DataVector > &  volume_det_inv_jacobian,
const Mesh< Dim > &  volume_mesh,
const size_t  dimension,
const Variables< BoundaryCorrectionTagsList > &  upper_boundary_corrections,
const Scalar< DataVector > &  upper_magnitude_of_face_normal,
const Scalar< DataVector > &  upper_face_det_jacobian,
const Variables< BoundaryCorrectionTagsList > &  lower_boundary_corrections,
const Scalar< DataVector > &  lower_magnitude_of_face_normal,
const Scalar< DataVector > &  lower_face_det_jacobian 
)

Lift both the upper and lower (in logical coordinates) boundary corrections to the volume time derivatives in the specified logical dimension.

The upper and lower boundary corrections in the logical dimension are lifted together in order to reduce the amount of striding through data that is needed and to improve cache-friendliness.

The general lifting term (for the \(\xi\)-dimension) is:

\begin{align*} \partial_t u_{\alpha\breve{\imath}\breve{\jmath}\breve{k}}=\cdots -\frac{\ell_{\breve{\imath}}\left(\xi=1\right)} {w_{\breve{\imath}}J_{\breve{\imath}\breve{\jmath}\breve{k}}} \left[J\sqrt{ \frac{\partial\xi}{\partial x^i} \gamma^{ij} \frac{\partial\xi}{\partial x^j}} \left(G_{\alpha} + D_{\alpha}\right) \right]_{\breve{\jmath}\breve{k}}\left(\xi=1\right) - \frac{\ell_{\breve{\imath}}\left(\xi=-1\right)} {w_{\breve{\imath}}J_{\breve{\imath}\breve{\jmath}\breve{k}}} \left[J\sqrt{ \frac{\partial\xi}{\partial x^i} \gamma^{ij} \frac{\partial\xi}{\partial x^j}} \left(G_{\alpha} + D_{\alpha}\right) \right]_{\breve{\jmath}\breve{k}}\left(\xi=-1\right), \end{align*}

where \(\breve{\imath}\), \(\breve{\jmath}\), and \(\breve{k}\) are indices in the logical \(\xi\), \(\eta\), and \(\zeta\) dimensions. The \(G+D\) terms correspond to the upper_boundary_corrections and lower_boundary_corrections function arguments, and the function Spectral::boundary_lifting_term() is used to compute and cache the terms from the lifting terms \(\ell_{\breve{\imath}}(\xi=\pm1)/w_{\breve{\imath}}\).

Note
that normal vectors are pointing out of the element and therefore both terms have the same sign.

◆ project_contiguous_data_to_boundary()

template<typename VolumeVarsTagsList , typename FaceVarsTagsList , size_t Dim>
void evolution::dg::project_contiguous_data_to_boundary ( const gsl::not_null< Variables< FaceVarsTagsList > * >  face_fields,
const Variables< VolumeVarsTagsList > &  volume_fields,
const Mesh< Dim > &  volume_mesh,
const Direction< Dim > &  direction 
)

Projects a Variables of volume data to a contiguous subset of a boundary Variables

The volume_fields are all projected into the face_fields in the direction direction. The tags in VolumeVarsTagsList must be a contiguous subset of the tags in FaceVarsTagsList. That is, FaceVarsTagsList must be equivalent to tmpl::append<Before, VolumeVarsTagsList, After> where Before and After are tmpl::lists of arbitrary size. This is because the projection is applied on all of the tensor components of the volume variables and is written into contiguous memory on the boundary.

In general, this function will be used for projecting all the evolved variables or all the volume fluxes to the faces. The function evolution::dg::project_tensors_to_boundary() should be used for projecting individual tensors to the face.

Note
This function works for both Gauss and Gauss-Lobatto uniform meshes.

◆ project_tensor_to_boundary()

template<typename Symm , typename IndexList , size_t Dim>
void evolution::dg::project_tensor_to_boundary ( const gsl::not_null< Tensor< DataVector, Symm, IndexList > * >  face_field,
const Tensor< DataVector, Symm, IndexList > &  volume_field,
const Mesh< Dim > &  volume_mesh,
const Direction< Dim > &  direction 
)

Projects a tensor to the face.

Note
This function works for both Gauss and Gauss-Lobatto uniform meshes.

◆ project_tensors_to_boundary()

template<typename TagsToProjectList , typename VolumeVarsTagsList , typename FaceVarsTagsList , size_t Dim>
void evolution::dg::project_tensors_to_boundary ( const gsl::not_null< Variables< FaceVarsTagsList > * >  face_fields,
const Variables< VolumeVarsTagsList > &  volume_fields,
const Mesh< Dim > &  volume_mesh,
const Direction< Dim > &  direction 
)

Projects a subset of the tensors in the volume_fields onto the face.

The tensors to project are listed in the TagsToProjectList.

Note
This function works for both Gauss and Gauss-Lobatto uniform meshes.

◆ receive_boundary_data_local_time_stepping()

template<typename Metavariables , bool DenseOutput, typename DbTagsList , typename... InboxTags>
bool evolution::dg::receive_boundary_data_local_time_stepping ( const gsl::not_null< db::DataBox< DbTagsList > * >  box,
const gsl::not_null< tuples::TaggedTuple< InboxTags... > * >  inboxes 
)

Receive boundary data for local time-stepping. Returns true if all necessary data has been received.

Setting DenseOutput to true receives data required for output at Tags::Time instead of Tags::Next<::Tags::TimeStepId>.

Variable Documentation

◆ using_subcell_v

template<typename Metavars >
constexpr bool evolution::dg::using_subcell_v = using_subcell<Metavars>::value
constexpr

If Metavars has a SubcellOptions member struct and SubcellOptions::subcell_enabled is true then is true, otherwise false.

Note
This check is intentionally not inside the DgSubcell library so that executables that do not use subcell do not need to link against it.