SpECTRE
v2024.09.29
|
Functionality related to discontinuous Galerkin schemes. More...
Classes | |
struct | SimpleBoundaryData |
Distinguishes between field data, which can be projected to a mortar, and extra data, which will not be projected. More... | |
class | SimpleMortarData |
Storage of boundary data on two sides of a mortar. More... | |
Typedefs | |
template<size_t VolumeDim> | |
using | MortarId = DirectionalId< VolumeDim > |
template<size_t MortarDim> | |
using | MortarSize = std::array< Spectral::MortarSize, MortarDim > |
template<size_t VolumeDim, typename ValueType > | |
using | MortarMap = DirectionalIdMap< VolumeDim, ValueType > |
Enumerations | |
enum class | Formulation { StrongInertial , WeakInertial , StrongLogical } |
The DG formulation to use. More... | |
Functions | |
std::ostream & | operator<< (std::ostream &os, Formulation t) |
template<typename InboxTag , size_t Dim, typename TemporalIdType , typename... InboxTags> | |
bool | has_received_from_all_mortars (const TemporalIdType &temporal_id, const Element< Dim > &element, const tuples::TaggedTuple< InboxTags... > &inboxes) |
Determines if data on all mortars has been received for the InboxTag at time temporal_id . 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, typename... VolumeTags, typename... BoundaryCorrectionTags> | |
void | lift_boundary_terms_gauss_points (const gsl::not_null< Variables< tmpl::list< VolumeTags... > > * > volume_data, const Variables< tmpl::list< BoundaryCorrectionTags... > > &boundary_corrections, const Mesh< Dim > &volume_mesh, const Direction< Dim > &direction) |
Lift the boundary corrections from the face to the volume for Gauss points in the direction perpendicular to the face. More... | |
template<size_t Dim> | |
void | metric_identity_det_jac_times_inv_jac (gsl::not_null< InverseJacobian< DataVector, Dim, Frame::ElementLogical, Frame::Inertial > * > det_jac_times_inverse_jacobian, const Mesh< Dim > &mesh, const tnsr::I< DataVector, Dim, Frame::Inertial > &inertial_coords, const Jacobian< DataVector, Dim, Frame::ElementLogical, Frame::Inertial > &jacobian) |
Compute the Jacobian determinant times the inverse Jacobian so that the result is divergence-free. More... | |
template<size_t Dim> | |
void | metric_identity_jacobian_quantities (gsl::not_null< InverseJacobian< DataVector, Dim, Frame::ElementLogical, Frame::Inertial > * > det_jac_times_inverse_jacobian, gsl::not_null< InverseJacobian< DataVector, Dim, Frame::ElementLogical, Frame::Inertial > * > inverse_jacobian, gsl::not_null< Jacobian< DataVector, Dim, Frame::ElementLogical, Frame::Inertial > * > jacobian, gsl::not_null< Scalar< DataVector > * > det_jacobian, const Mesh< Dim > &mesh, const tnsr::I< DataVector, Dim, Frame::Inertial > &inertial_coords) |
Compute the Jacobian, inverse Jacobian, and determinant of the Jacobian so that they satisfy the metric identities. More... | |
template<size_t Dim> | |
Mesh< Dim > | mortar_mesh (const Mesh< Dim > &face_mesh1, const Mesh< Dim > &face_mesh2) |
Find a mesh for a mortar capable of representing data from either of two faces. More... | |
template<size_t Dim> | |
MortarSize< Dim - 1 > | mortar_size (const ElementId< Dim > &self, const ElementId< Dim > &neighbor, size_t dimension, const OrientationMap< Dim > &orientation) |
Determine the size of the mortar (i.e., the part of the face it covers) for communicating with a neighbor. This is the size relative to the size of self , and will not generally agree with that determined by neighbor . | |
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 DataType , size_t Dim> | |
void | apply_mass_matrix (const gsl::not_null< DataType * > data, const Mesh< Dim > &mesh) |
Apply the DG mass matrix to the data, in the diagonal mass-matrix approximation ("mass-lumping") More... | |
template<size_t Dim, typename TagsList > | |
void | apply_mass_matrix (const gsl::not_null< Variables< TagsList > * > data, const Mesh< Dim > &mesh) |
Apply the DG mass matrix to the data, in the diagonal mass-matrix approximation ("mass-lumping") More... | |
template<typename DataType , size_t Dim> | |
void | apply_inverse_mass_matrix (const gsl::not_null< DataType * > data, const Mesh< Dim > &mesh) |
Apply the inverse DG mass matrix to the data. More... | |
template<size_t Dim, typename TagsList > | |
void | apply_inverse_mass_matrix (const gsl::not_null< Variables< TagsList > * > data, const Mesh< Dim > &mesh) |
Apply the inverse DG mass matrix to the data. More... | |
template<typename... BoundaryCorrectionTags> | |
void | lift_flux (const gsl::not_null< Variables< tmpl::list< BoundaryCorrectionTags... > > * > boundary_correction_terms, const size_t extent_perpendicular_to_boundary, const Scalar< DataVector > &magnitude_of_face_normal) |
Lifts the flux contribution from an interface to the volume. More... | |
template<typename... FluxTags> | |
auto | lift_flux (Variables< tmpl::list< FluxTags... > > flux, const size_t extent_perpendicular_to_boundary, const Scalar< DataVector > &magnitude_of_face_normal) -> Variables< tmpl::list< db::remove_tag_prefix< FluxTags >... > > |
Lifts the flux contribution from an interface to the volume. More... | |
template<typename Tags , size_t Dim> | |
void | project_to_mortar (const gsl::not_null< Variables< Tags > * > result, const Variables< Tags > &vars, const Mesh< Dim > &face_mesh, const Mesh< Dim > &mortar_mesh, const MortarSize< Dim > &mortar_size) |
Project variables from a face to a mortar. | |
template<typename Tags , size_t Dim> | |
Variables< Tags > | project_to_mortar (const Variables< Tags > &vars, const Mesh< Dim > &face_mesh, const Mesh< Dim > &mortar_mesh, const MortarSize< Dim > &mortar_size) |
Project variables from a face to a mortar. | |
template<typename Tags , size_t Dim> | |
void | project_from_mortar (const gsl::not_null< Variables< Tags > * > result, const Variables< Tags > &vars, const Mesh< Dim > &face_mesh, const Mesh< Dim > &mortar_mesh, const MortarSize< Dim > &mortar_size) |
Project variables from a mortar to a face. | |
template<typename Tags , size_t Dim> | |
Variables< Tags > | project_from_mortar (const Variables< Tags > &vars, const Mesh< Dim > &face_mesh, const Mesh< Dim > &mortar_mesh, const MortarSize< Dim > &mortar_size) |
Project variables from a mortar to a face. | |
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... | |
template<typename Symm , typename IndexList , size_t Dim> | |
Tensor< DataVector, Symm, IndexList > | project_tensor_to_boundary (const Tensor< DataVector, Symm, IndexList > &volume_field, const Mesh< Dim > &volume_mesh, const Direction< Dim > &direction) |
Projects a tensor to the face. More... | |
Functionality related to discontinuous Galerkin schemes.
void dg::apply_inverse_mass_matrix | ( | const gsl::not_null< DataType * > | data, |
const Mesh< Dim > & | mesh | ||
) |
Apply the inverse DG mass matrix to the data.
void dg::apply_inverse_mass_matrix | ( | const gsl::not_null< Variables< TagsList > * > | data, |
const Mesh< Dim > & | mesh | ||
) |
Apply the inverse DG mass matrix to the data.
void dg::apply_mass_matrix | ( | const gsl::not_null< DataType * > | data, |
const Mesh< Dim > & | mesh | ||
) |
Apply the DG mass matrix to the data, in the diagonal mass-matrix approximation ("mass-lumping")
The DG mass matrix is:
where
where
void dg::apply_mass_matrix | ( | const gsl::not_null< Variables< TagsList > * > | data, |
const Mesh< Dim > & | mesh | ||
) |
Apply the DG mass matrix to the data, in the diagonal mass-matrix approximation ("mass-lumping")
The DG mass matrix is:
where
where
void 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
where
void 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
where boundary_corrections
function argument, and the function Spectral::boundary_lifting_term() is used to compute and cache the terms from the lifting terms
dt_vars | The volume time derivatives |
volume_det_inv_jacobian | The inverse Jacobian determinant in the volume |
volume_mesh | The mesh of the volume |
direction | The direction of the face on which the boundary corrections are defined |
boundary_corrections | The boundary corrections |
magnitude_of_face_normal | The term |
face_det_jacobian | The volume Jacobian determinant |
void 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
where 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
void dg::lift_boundary_terms_gauss_points | ( | const gsl::not_null< Variables< tmpl::list< VolumeTags... > > * > | volume_data, |
const Variables< tmpl::list< BoundaryCorrectionTags... > > & | boundary_corrections, | ||
const Mesh< Dim > & | volume_mesh, | ||
const Direction< Dim > & | direction | ||
) |
Lift the boundary corrections from the face to the volume for Gauss points in the direction perpendicular to the face.
This function doesn't apply any Jacobians or quadrature weights. It only applies the lifting matrix boundary_corrections
and adds the result to the volume_data
.
void 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::list
s 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 dg::project_tensors_to_boundary()
should be used for projecting individual tensors to the face.
void 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.
Tensor< DataVector, Symm, IndexList > dg::project_tensor_to_boundary | ( | const Tensor< DataVector, Symm, IndexList > & | volume_field, |
const Mesh< Dim > & | volume_mesh, | ||
const Direction< Dim > & | direction | ||
) |
Projects a tensor to the face.
void 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
.