SpECTRE
v2024.12.16
|
Functions and classes specific to the Discontinuous Galerkin algorithm. More...
Enumerations | |
enum class | dg::Formulation { StrongInertial , WeakInertial , StrongLogical } |
The DG formulation to use. More... | |
Functions | |
template<typename InboxTag , size_t Dim, typename TemporalIdType , typename... InboxTags> | |
bool | dg::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> | |
void | dg::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 | dg::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 > | dg::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 > | dg::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... BoundaryCorrectionTags> | |
void | dg::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 | dg::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 | dg::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 > | dg::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 | dg::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 > | dg::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. | |
Functions and classes specific to the Discontinuous Galerkin algorithm.
|
strong |
The DG formulation to use.
StrongInertial
formulation is also known as the integrate then transform formulation. The "Inertial" part of the name refers to the fact that the integration is done over the physical/inertial coordinates, while the "strong" part refers to the fact that the boundary correction terms are zero if the solution is continuous at the interfaces. See [188] for an overview.WeakInertial
formulation is also known as the integrate then transform formulation. The "Inertial" part of the name refers to the fact that the integration is done over the physical/inertial coordinates, while the "weak" part refers to the fact that the boundary correction terms are non-zero even if the solution is continuous at the interfaces. See [188] for an overview.StrongLogical
formulation is also known as the transform then integrate formulation. The "logical" part of the name refers to the fact that the integration is done over the logical coordinates, while the "strong" part refers to the fact that the boundary correction terms are zero if the solution is continuous at the interfaces. This formulation arises from the StrongInertial
formulation by moving the Jacobians that appear when computing the divergence of fluxes into the divergence so the divergence is computed in logical coordinates:
dg::metric_identity_det_jac_times_inv_jac
for details and for functions that compute the Jacobians so they satisfy the metric identities numerically (which may or may not be useful or necessary). bool dg::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
.
The InboxTag
must hold a container indexed by the TemporalIdType
, representing received neighbor data at a given time. Each element in the container, i.e. the received neighbor data at a given time, must be another container indexed by the dg::MortarId<Dim>
. The value it holds is not relevant for this function. Here's an example for such a type:
void dg::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.
The lifting operation takes the (d-1)-dimensional flux term at the interface and computes the corresponding d-dimensional term in the volume. SpECTRE implements an efficient DG method in which each interface grid point contributes only to that same grid point of the volume.
SpECTRE implements a DG method with a diagonalized mass matrix (also known as a mass-lumping scheme). This choice gives a large reduction in the computational cost of the lifting operation, however, the scheme is slightly less accurate, especially when the grid is deformed by non-trivial Jacobians. For more details on the diagonalization of the mass matrix and its implications, [188], especially Section 3.
auto dg::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.
The lifting operation takes the (d-1)-dimensional flux term at the interface and computes the corresponding d-dimensional term in the volume. SpECTRE implements an efficient DG method in which each interface grid point contributes only to that same grid point of the volume.
SpECTRE implements a DG method with a diagonalized mass matrix (also known as a mass-lumping scheme). This choice gives a large reduction in the computational cost of the lifting operation, however, the scheme is slightly less accurate, especially when the grid is deformed by non-trivial Jacobians. For more details on the diagonalization of the mass matrix and its implications, [188], especially Section 3.
void dg::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.
The metric identities are given by
We want to compute
The discretized form, with the Jacobian determinant
where indices
In 1d we have:
In 2d we have:
In 3d we have:
Again, this ensures that the metric identities are satisfied discretely. That is,
numerically.
The reason for calculating
which should be identically zero if
The subtraction technique is most commonly used in finite difference codes.
void dg::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.
Uses dg::metric_identity_jacobian()
to compute the determinant of the Jacobian times the inverse Jacobian. By taking the determinant of this product, we can isolate
We assume the determinant of the Jacobian is positive, which means logical and inertial coordinates have the same handedness. With this assumption, we have
We can now compute the inverse Jacobian using:
This guarantees that multiplying the determinant of the Jacobian by the inverse Jacobian gives a result that satisfies the metric identities. We also compute the Jacobian by inverting the inverse Jacobian, which guarantees they are (numerical) inverses of each other.
jacobian
must be the Jacobian to use for computing the determinant of the Jacobian times the inverse Jacobian so that it satisfies the metric identities. The jacobian
can be computed analytically or numerically, either is fine. On output the jacobian
is the inverse of the inverse Jacobian that satisfies the metric identities. Mesh< Dim > dg::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.
orientation
passed to domain::Initialization::create_initial_mesh
, for example.