SpECTRE
v2024.12.16
|
The building blocks used to describe the computational domain. More...
Namespaces | |
namespace | domain |
Holds entities related to the computational domain. | |
namespace | domain::BoundaryConditions |
Domain support for applying boundary conditions | |
namespace | domain::CoordinateMaps::Tags |
Tags for the coordinate maps. | |
namespace | domain::creators |
Defines classes that create Domains. | |
namespace | domain::creators::time_dependence |
Classes and functions for adding time dependence to a domain. | |
namespace | domain::FunctionsOfTime |
Contains functions of time to support the dual frame system. | |
namespace | domain::Tags |
Tags for the domain. | |
Classes | |
class | Block< VolumeDim > |
A Block<VolumeDim> is a region of a VolumeDim-dimensional computational domain that defines the root node of a tree which is used to construct the Elements that cover a region of the computational domain. More... | |
struct | domain::CoordinateMaps::Tags::CoordinateMap< VolumeDim, SourceFrame, TargetFrame > |
The coordinate map from source to target coordinates. More... | |
class | domain::creators::BinaryCompactObject< UseWorldtube > |
A general domain for two compact objects. More... | |
class | domain::creators::CylindricalBinaryCompactObject |
A general domain for two compact objects based on cylinders. More... | |
class | DomainCreator< VolumeDim > |
Base class for creating Domains from an option string. More... | |
struct | domain::OptionTags::DomainCreator< Dim > |
The input file tag for the DomainCreator to use. More... | |
struct | domain::Tags::Domain< VolumeDim > |
The Domain. More... | |
struct | domain::Tags::InitialExtents< Dim > |
The number of grid points per dimension for all elements in each block of the initial computational domain. More... | |
struct | domain::Tags::InitialRefinementLevels< Dim > |
The initial refinement level per dimension for all elements in each block of the initial computational domain. More... | |
struct | domain::Tags::ObjectCenter< Label > |
The grid frame center of the given object. More... | |
class | Domain< VolumeDim > |
A wrapper around a vector of Blocks that represent the computational domain. More... | |
struct | PairOfFaces |
Each member in PairOfFaces holds the global corner ids of a block face. PairOfFaces is used in setting up periodic boundary conditions by identifying the two faces with each other. More... | |
class | VolumeCornerIterator< VolumeDim > |
Iterates over the corners of a VolumeDim-dimensional cube. More... | |
class | FaceCornerIterator< VolumeDim > |
Iterates over the 2^(VolumeDim-1) logical corners of the face of a VolumeDim-dimensional cube in the given direction. More... | |
struct | ElementLogicalCoordHolder< Dim > |
Holds element logical coordinates of an arbitrary set of points on a single Element . The arbitrary set of points is assumed to be a subset of a larger set of points spanning multiple Element s, and this class holds offsets that index into that larger set of points. More... | |
class | ElementMap< Dim, TargetFrame > |
The CoordinateMap for the Element from the Logical frame to the TargetFrame More... | |
class | ExcisionSphere< VolumeDim > |
The excision sphere information of a computational domain. The excision sphere is assumed to be a coordinate sphere in the grid frame. More... | |
struct | domain::Tags::UnnormalizedFaceNormal< VolumeDim, Frame > |
The unnormalized face normal one form. More... | |
struct | domain::Tags::InterfaceCompute< Tags::BoundaryDirectionsExterior< VolumeDim >, UnnormalizedFaceNormalCompute< VolumeDim, Frame > > |
Specialisation of UnnormalizedFaceNormal for the external boundaries which inverts the normals. Since ExternalBoundariesDirections are meant to represent ghost elements, the normals should correspond to the normals in said element, which are inverted with respect to the current element. More... | |
class | domain::FunctionsOfTime::FunctionOfTime |
Base class for FunctionsOfTime. More... | |
class | domain::FunctionsOfTime::PiecewisePolynomial< MaxDeriv > |
A function that has a piecewise-constant MaxDeriv th derivative. More... | |
class | domain::FunctionsOfTime::QuaternionFunctionOfTime< MaxDeriv > |
A FunctionOfTime that stores quaternions for the rotation map. More... | |
struct | domain::Tags::InterfaceCompute< DirectionsTag, Tag > |
Compute tag for representing items computed on a set of interfaces. Can be retrieved using Tags::Interface<DirectionsTag, Tag> More... | |
struct | domain::Tags::Slice< DirectionsTag, Tag > |
Compute tag for representing a compute item that slices data from the volume to a set of interfaces. More... | |
struct | domain::Tags::InterfaceMesh< VolumeDim > |
Computes the VolumeDim-1 dimensional mesh on an interface from the volume mesh. Tags::InterfaceCompute<Dirs, InterfaceMesh<VolumeDim>> is retrievable as Tags::Interface<Dirs, Mesh<VolumeDim>> from the DataBox. More... | |
struct | domain::Tags::BoundaryCoordinates< VolumeDim, MovingMesh > |
Computes the coordinates in the frame Frame on the faces defined by Direction . Intended to be prefixed by a Tags::InterfaceCompute to define the directions on which to compute the coordinates. More... | |
struct | domain::Tags::JacobianDiagnostic< Dim > |
A diagnostic comparing the analytic and numerical Jacobians for a map. See domain::jacobian_diagnostic for details. More... | |
struct | domain::Tags::JacobianDiagnosticCompute< Dim, TargetFrame > |
Computes the Jacobian diagnostic, which compares the analytic Jacobian (provided by some coordinate map) to a numerical Jacobian computed using numerical partial derivatives. The coordinates must be in the target frame of the map. See domain::jacobian_diagnostic for details of the calculation. More... | |
struct | domain::Tags::MinimumGridSpacing< Dim, Frame > |
The minimum coordinate distance between grid points. More... | |
struct | domain::Tags::SizeOfElement< VolumeDim > |
The inertial-coordinate size of an element along each of its logical directions. More... | |
class | domain::BlockId |
Index a block of the computational domain. More... | |
class | BlockNeighbor< VolumeDim > |
Information about the neighbor of a host Block in a particular direction. More... | |
class | Direction< VolumeDim > |
A particular Side along a particular coordinate Axis. More... | |
struct | DirectionHash< Dim > |
Provides a perfect hash if the size of the hash table is 2 * Dim . To take advantage of this, use the FixedHashMap class. More... | |
class | DirectionMap< Dim, T > |
An optimized map with Direction keys. More... | |
class | Element< VolumeDim > |
A spectral element with knowledge of its neighbors. More... | |
class | ElementId< VolumeDim > |
An ElementId uniquely labels an Element. More... | |
class | Neighbors< VolumeDim > |
Information about the neighbors of a host Element in a particular direction. More... | |
class | OrientationMap< VolumeDim > |
A mapping of the logical coordinate axes of a host to the logical coordinate axes of a neighbor of the host. More... | |
class | SegmentId |
A SegmentId labels a segment of the interval | |
struct | domain::Tags::Element< VolumeDim > |
The Element associated with the DataBox. More... | |
struct | domain::Tags::Mesh< VolumeDim > |
The computational grid of the Element in the DataBox. More... | |
struct | domain::Tags::ElementMap< VolumeDim, TargetFrame > |
The coordinate map from the ElementLogical frame to the TargetFrame. More... | |
struct | domain::Tags::Coordinates< Dim, Frame > |
The coordinates in a given frame. More... | |
struct | domain::Tags::MappedCoordinates< MapTag, SourceCoordsTag, CoordinatesTag > |
The coordinates in the target frame of MapTag . The SourceCoordsTag 's frame must be the source frame of MapTag More... | |
struct | domain::Tags::InverseJacobian< Dim, SourceFrame, TargetFrame > |
The inverse Jacobian from the source frame to the target frame. More... | |
struct | domain::Tags::InverseJacobianCompute< MapTag, SourceCoordsTag > |
Computes the inverse Jacobian of the map held by MapTag at the coordinates held by SourceCoordsTag . The coordinates must be in the source frame of the map. More... | |
struct | domain::Tags::Jacobian< Dim, SourceFrame, TargetFrame > |
The Jacobian from the source frame to the target frame. More... | |
struct | domain::Tags::JacobianCompute< Dim, SourceFrame, TargetFrame > |
Computes the Jacobian of the map from the InverseJacobian<Dim, SourceFrame, TargetFrame> tag. More... | |
struct | domain::Tags::DetInvJacobian< SourceFrame, TargetFrame > |
The determinant of the inverse Jacobian from the source frame to the target frame. More... | |
struct | domain::Tags::DetInvJacobianCompute< Dim, SourceFrame, TargetFrame > |
Computes the determinant of the inverse Jacobian. More... | |
struct | domain::Tags::DetJacobian< SourceFrame, TargetFrame > |
The determinant of the Jacobian from the source frame to the target frame. More... | |
struct | domain::Tags::DetTimesInvJacobian< Dim, SourceFrame, TargetFrame > |
The inverse Jacobian times the determinant of the Jacobian. More... | |
struct | domain::Tags::VariablesBoundaryData |
Base tag for boundary data needed for updating the variables. More... | |
struct | domain::Tags::InternalDirections< VolumeDim > |
The set of directions to neighboring Elements. More... | |
struct | domain::Tags::BoundaryDirectionsInterior< VolumeDim > |
The set of directions which correspond to external boundaries. Used for representing data on the interior side of the external boundary faces. More... | |
struct | domain::Tags::BoundaryDirectionsExterior< VolumeDim > |
The set of directions which correspond to external boundaries. To be used to represent data which exists on the exterior side of the external boundary faces. More... | |
struct | domain::Tags::Interface< DirectionsTag, Tag > |
Tag which is either a SimpleTag for quantities on an interface, base tag to a compute item which acts on tags on an interface, or base tag to a compute item which slices a tag from the volume to an interface. More... | |
struct | domain::Tags::Direction< VolumeDim > |
Direction to an interface More... | |
struct | domain::OptionTags::ElementDistribution |
struct | domain::Tags::ElementDistribution |
Tag that holds method for how to distribute the elements on the given resources. More... | |
struct | domain::Tags::LogicalCoordinates< VolumeDim > |
The logical coordinates in the Element. More... | |
struct | domain::Actions::CheckFunctionsOfTimeAreReady< Dim > |
Check that functions of time are up-to-date. More... | |
struct | domain::CheckFunctionsOfTimeAreReadyPostprocessor< Dim > |
Dense-output postprocessor to check that functions of time are up-to-date. More... | |
Macros | |
#define | INSTANTIATE_MAPS_SIMPLE_FUNCTIONS(_, data) |
Generate instantiations of member functions of the CoordinateMap class template. More... | |
#define | INSTANTIATE_MAPS_DATA_TYPE_FUNCTIONS(_, data) |
Generate instantiations of member functions of the CoordinateMap class template. More... | |
#define | INSTANTIATE_MAPS_FUNCTIONS(MAPS_TUPLE, SOURCE_FRAME, TARGET_FRAMES_TUPLE, TYPES_TUPLE) |
Generate instantiations of member functions of the CoordinateMap class template. More... | |
Enumerations | |
enum class | ShellWedges { ShellWedges::All , ShellWedges::FourOnEquator , ShellWedges::OneAlongMinusX } |
The number of wedges to include in the Sphere domain. More... | |
enum class | Side : uint8_t { Uninitialized = 0 << detail::side_shift , Lower = 1 << detail::side_shift , Upper = 2 << detail::side_shift , Self = 3 << detail::side_shift } |
A label for the side of a manifold. More... | |
Functions | |
template<typename SourceFrame , typename TargetFrame , typename... Maps> | |
auto | domain::make_coordinate_map (Maps &&... maps) -> CoordinateMap< SourceFrame, TargetFrame, std::decay_t< Maps >... > |
Creates a CoordinateMap of maps... | |
template<typename SourceFrame , typename TargetFrame , typename... Maps> | |
auto | domain::make_coordinate_map_base (Maps &&... maps) -> std::unique_ptr< CoordinateMapBase< SourceFrame, TargetFrame, CoordinateMap< SourceFrame, TargetFrame, std::decay_t< Maps >... >::dim > > |
Creates a std::unique_ptr<CoordinateMapBase> of maps... | |
template<typename SourceFrame , typename TargetFrame , typename Arg0 , typename... Args> | |
auto | domain::make_vector_coordinate_map_base (Arg0 &&arg_0, Args &&... remaining_args) -> std::vector< std::unique_ptr< CoordinateMapBase< SourceFrame, TargetFrame, std::decay_t< Arg0 >::dim > > > |
Creates a std::vector<std::unique_ptr<CoordinateMapBase>> containing the result of make_coordinate_map_base applied to each argument passed in. | |
template<typename SourceFrame , typename TargetFrame , size_t Dim, typename Map , typename... Maps> | |
auto | domain::make_vector_coordinate_map_base (std::vector< Map > maps, const Maps &... remaining_maps) -> std::vector< std::unique_ptr< CoordinateMapBase< SourceFrame, TargetFrame, Dim > > > |
Creates a std::vector<std::unique_ptr<CoordinateMapBase>> containing the result of make_coordinate_map_base applied to each element of the vector of maps composed with the rest of the arguments passed in. | |
template<typename SourceFrame , typename TargetFrame , typename... Maps, typename NewMap > | |
CoordinateMap< SourceFrame, TargetFrame, Maps..., NewMap > | domain::push_back (CoordinateMap< SourceFrame, TargetFrame, Maps... > old_map, NewMap new_map) |
Creates a CoordinateMap by appending the new map to the end of the old maps. | |
template<typename SourceFrame , typename TargetFrame , typename... Maps, typename NewMap > | |
CoordinateMap< SourceFrame, TargetFrame, NewMap, Maps... > | domain::push_front (CoordinateMap< SourceFrame, TargetFrame, Maps... > old_map, NewMap new_map) |
Creates a CoordinateMap by prepending the new map to the beginning of the old maps. | |
template<size_t VolumeDim> | |
void | set_internal_boundaries (gsl::not_null< std::vector< DirectionMap< VolumeDim, BlockNeighbor< VolumeDim > > > * > neighbors_of_all_blocks, const std::vector< std::array< size_t, two_to_the(VolumeDim)> > &corners_of_all_blocks) |
Sets up the BlockNeighbors using the corner numbering scheme provided by the user to deduce the correct neighbors and orientations. Does not set up periodic boundary conditions. | |
template<size_t VolumeDim> | |
void | set_internal_boundaries (gsl::not_null< std::vector< DirectionMap< VolumeDim, BlockNeighbor< VolumeDim > > > * > neighbors_of_all_blocks, const std::vector< std::unique_ptr< domain::CoordinateMapBase< Frame::BlockLogical, Frame::Inertial, VolumeDim > > > &maps) |
Sets up the BlockNeighbors using the corner numbering scheme implied by the maps provided by the user to deduce the correct neighbors and orientations. More... | |
template<size_t VolumeDim> | |
void | set_identified_boundaries (const std::vector< PairOfFaces > &identifications, const std::vector< std::array< size_t, two_to_the(VolumeDim)> > &corners_of_all_blocks, gsl::not_null< std::vector< DirectionMap< VolumeDim, BlockNeighbor< VolumeDim > > > * > neighbors_of_all_blocks) |
Sets up additional BlockNeighbors corresponding to any identifications of faces provided by the user. Can be used for manually setting up periodic boundary conditions. | |
template<size_t VolumeDim> | |
auto | indices_for_rectilinear_domains (const Index< VolumeDim > &domain_extents, const std::vector< Index< VolumeDim > > &block_indices_to_exclude={}) -> std::vector< Index< VolumeDim > > |
The multi-indices that identify the individual Blocks in the lattice. | |
template<size_t VolumeDim> | |
auto | corners_for_rectilinear_domains (const Index< VolumeDim > &domain_extents, const std::vector< Index< VolumeDim > > &block_indices_to_exclude={}) -> std::vector< std::array< size_t, two_to_the(VolumeDim)> > |
The corners for a rectilinear domain made of n-cubes. More... | |
std::array< OrientationMap< 3 >, 6 > | orientations_for_sphere_wrappings () |
An array of the orientations of the six blocks that make up a Sphere. More... | |
size_t | which_wedge_index (const ShellWedges &which_wedges) |
The first index in the list "UpperZ", "LowerZ", "UpperY", "LowerY", "UpperX" "LowerX" that is included in which_wedges . It is 0 for ShellWedges::All , 2 for ShellWedges::FourOnEquator , and 5 for ShellWedges::OneAlongMinusX . | |
std::vector< domain::CoordinateMaps::Wedge< 3 > > | sph_wedge_coordinate_maps (double inner_radius, double outer_radius, double inner_sphericity, double outer_sphericity, bool use_equiangular_map, const std::optional< std::pair< double, std::array< double, 3 > > > &offset_options=std::nullopt, bool use_half_wedges=false, const std::vector< double > &radial_partitioning={}, const std::vector< domain::CoordinateMaps::Distribution > &radial_distribution={domain::CoordinateMaps::Distribution::Linear}, ShellWedges which_wedges=ShellWedges::All, double opening_angle=M_PI_2) |
std::vector< domain::CoordinateMaps::Frustum > | frustum_coordinate_maps (double length_inner_cube, double length_outer_cube, bool equiangular_map_at_outer, bool equiangular_map_at_inner, const std::array< double, 3 > &origin_preimage={{0.0, 0.0, 0.0}}, domain::CoordinateMaps::Distribution radial_distribution=domain::CoordinateMaps::Distribution::Linear, std::optional< double > distribution_value=std::nullopt, double sphericity=0.0, double opening_angle=M_PI_2) |
These are the ten Frustums used in the DomainCreators for binary compact objects. The Frustums partition the volume defined by two bounding surfaces: The inner surface is the surface of the two joined inner cubes enveloping the two compact objects, while the outer is the surface of the outer cube. More... | |
std::vector< std::array< size_t, 8 > > | corners_for_radially_layered_domains (size_t number_of_layers, bool include_central_block, const std::array< size_t, 8 > ¢ral_block_corners={{1, 2, 3, 4, 5, 6, 7, 8}}, ShellWedges which_wedges=ShellWedges::All) |
The corners for a domain with radial layers. More... | |
std::vector< std::array< size_t, 8 > > | corners_for_biradially_layered_domains (size_t number_of_radial_layers, size_t number_of_biradial_layers, bool include_central_block_lhs, bool include_central_block_rhs, const std::array< size_t, 8 > ¢ral_block_corners_lhs={ {1, 2, 3, 4, 5, 6, 7, 8}}) |
The corners for a domain with biradial layers. More... | |
template<typename TargetFrame > | |
auto | cyl_wedge_coordinate_maps (double inner_radius, double outer_radius, double lower_z_bound, double upper_z_bound, bool use_equiangular_map, const std::vector< double > &radial_partitioning={}, const std::vector< double > &partitioning_in_z={}, const std::vector< domain::CoordinateMaps::Distribution > &radial_distribution={domain::CoordinateMaps::Distribution::Linear}, const std::vector< domain::CoordinateMaps::Distribution > &distribution_in_z={domain::CoordinateMaps::Distribution::Linear}) -> std::vector< std::unique_ptr< domain::CoordinateMapBase< Frame::BlockLogical, TargetFrame, 3 > > > |
These are the CoordinateMaps used in the Cylinder DomainCreator. More... | |
auto | cyl_wedge_coord_map_center_blocks (double inner_radius, double lower_z_bound, double upper_z_bound, bool use_equiangular_map, const std::vector< double > &partitioning_in_z={}, const std::vector< domain::CoordinateMaps::Distribution > &distribution_in_z={domain::CoordinateMaps::Distribution::Linear}, CylindricalDomainParityFlip parity_flip=CylindricalDomainParityFlip::none) -> std::vector< domain::CoordinateMaps::ProductOf3Maps< domain::CoordinateMaps::Interval, domain::CoordinateMaps::Interval, domain::CoordinateMaps::Interval > > |
Same as cyl_wedge_coordinate_maps , but only the center square blocks,. More... | |
auto | cyl_wedge_coord_map_surrounding_blocks (double inner_radius, double outer_radius, double lower_z_bound, double upper_z_bound, bool use_equiangular_map, double inner_circularity, const std::vector< double > &radial_partitioning={}, const std::vector< double > &partitioning_in_z={}, const std::vector< domain::CoordinateMaps::Distribution > &radial_distribution={domain::CoordinateMaps::Distribution::Linear}, const std::vector< domain::CoordinateMaps::Distribution > &distribution_in_z={domain::CoordinateMaps::Distribution::Linear}, CylindricalDomainParityFlip parity_flip=CylindricalDomainParityFlip::none) -> std::vector< domain::CoordinateMaps::ProductOf2Maps< domain::CoordinateMaps::Wedge< 2 >, domain::CoordinateMaps::Interval > > |
Same as cyl_wedge_coordinate_maps, but only the surrounding wedge blocks. More... | |
std::vector< std::array< size_t, 8 > > | corners_for_cylindrical_layered_domains (size_t number_of_shells, size_t number_of_discs) |
The corners for a cylindrical domain split into discs with radial shells. More... | |
template<size_t VolumeDim> | |
std::array< size_t, two_to_the(VolumeDim)> | discrete_rotation (const OrientationMap< VolumeDim > &orientation, const std::array< size_t, two_to_the(VolumeDim)> &corners_of_aligned) |
Permutes the corner numbers of an n-cube. More... | |
template<typename TargetFrame , size_t VolumeDim> | |
auto | maps_for_rectilinear_domains (const Index< VolumeDim > &domain_extents, const std::array< std::vector< double >, VolumeDim > &block_demarcations, const std::vector< Index< VolumeDim > > &block_indices_to_exclude={}, const std::vector< OrientationMap< VolumeDim > > &orientations_of_all_blocks={}, bool use_equiangular_map=false) -> std::vector< std::unique_ptr< domain::CoordinateMapBase< Frame::BlockLogical, TargetFrame, VolumeDim > > > |
The CoordinateMaps for a rectilinear domain of n-cubes. More... | |
template<size_t VolumeDim> | |
Domain< VolumeDim > | rectilinear_domain (const Index< VolumeDim > &domain_extents, const std::array< std::vector< double >, VolumeDim > &block_demarcations, const std::vector< Index< VolumeDim > > &block_indices_to_exclude={}, const std::vector< OrientationMap< VolumeDim > > &orientations_of_all_blocks={}, const std::array< bool, VolumeDim > &dimension_is_periodic=make_array< VolumeDim >(false), const std::vector< PairOfFaces > &identifications={}, bool use_equiangular_map=false) |
Create a rectilinear Domain of multicubes. More... | |
template<size_t Dim> | |
auto | element_logical_coordinates (const std::vector< ElementId< Dim > > &element_ids, const std::vector< BlockLogicalCoords< Dim > > &block_coord_holders) -> std::unordered_map< ElementId< Dim >, ElementLogicalCoordHolder< Dim > > |
Given a set of points in block logical coordinates and their BlockIds , as returned from the function block_logical_coordinates , determines which Element s in a list of ElementId s contains each point, and determines the element logical coordinates of each point. More... | |
template<size_t VolumeDim> | |
tnsr::I< DataVector, VolumeDim, Frame::ElementLogical > | interface_logical_coordinates (const Mesh< VolumeDim - 1 > &mesh, const Direction< VolumeDim > &direction) |
Defines functions interface_logical_coordinates. More... | |
template<size_t Dim, typename Fr > | |
void | domain::jacobian_diagnostic (const gsl::not_null< tnsr::i< DataVector, Dim, typename Frame::ElementLogical > * > jacobian_diag, const Jacobian< DataVector, Dim, typename Frame::ElementLogical, Fr > &analytic_jacobian, const TensorMetafunctions::prepend_spatial_index< tnsr::I< DataVector, Dim, Fr >, Dim, UpLo::Lo, typename Frame::ElementLogical > &numeric_jacobian_transpose) |
A diagnostic comparing the analytic and numerical Jacobians for a map. More... | |
template<size_t Dim, typename Frame > | |
double | minimum_grid_spacing (const Index< Dim > &extents, const tnsr::I< DataVector, Dim, Frame > &coords) |
Finds the minimum coordinate distance between grid points. | |
template<size_t Dim> | |
size_t | index_to_slice_at (const Index< Dim > &extents, const Direction< Dim > &direction, const size_t offset=0) |
Finds the index in the perpendicular dimension of an element boundary. More... | |
template<size_t VolumeDim> | |
std::vector< ElementId< VolumeDim > > | initial_element_ids (size_t block_id, std::array< size_t, VolumeDim > initial_ref_levs, size_t grid_index=0) |
Create the ElementId s of the a single Block. | |
template<size_t VolumeDim> | |
std::vector< ElementId< VolumeDim > > | initial_element_ids (const std::vector< std::array< size_t, VolumeDim > > &initial_refinement_levels, size_t grid_index=0) |
Create the ElementId s of the initial computational domain. | |
constexpr size_t | maximum_number_of_neighbors (const size_t dim) |
Returns the maximum number of neighbors an element can have in dim dimensions. More... | |
constexpr size_t | maximum_number_of_neighbors_per_direction (const size_t dim) |
Returns the maximum number of neighbors in each direction an element can have in dim dimensions. More... | |
template<size_t VolumeDim, typename T > | |
std::array< tt::remove_cvref_wrap_t< T >, VolumeDim > | discrete_rotation (const OrientationMap< VolumeDim > &rotation, std::array< T, VolumeDim > source_coords) |
OrientationMap s define an active rotation of the logical axes that bring the axes of a host block into alignment with the logical axes of the neighbor block. discrete_rotation applies this active rotation on the coordinates as opposed to the axes. For a two-dimensional example, consider a host block and a neighbor block, where the OrientationMap between them is discrete_rotation interprets the OrientationMap passed to it. | |
template<size_t VolumeDim> | |
tnsr::Ij< double, VolumeDim, Frame::NoFrame > | discrete_rotation_jacobian (const OrientationMap< VolumeDim > &orientation) |
Computes the Jacobian of the transformation that is computed by discrete_rotation() More... | |
template<size_t VolumeDim> | |
tnsr::Ij< double, VolumeDim, Frame::NoFrame > | discrete_rotation_inverse_jacobian (const OrientationMap< VolumeDim > &orientation) |
Computes the inverse Jacobian of the transformation that is computed by discrete_rotation() | |
template<size_t VolumeDim> | |
Mesh< VolumeDim - 1 > | orient_mesh_on_slice (const Mesh< VolumeDim - 1 > &mesh_on_slice, size_t sliced_dim, const OrientationMap< VolumeDim > &orientation_of_neighbor) |
Orient a sliced Mesh to the logical frame of a neighbor element with the given orientation. | |
template<typename CacheTag , typename SimpleAction , typename Metavariables , typename ArrayIndex , typename Component , typename... Args> | |
bool | domain::functions_of_time_are_ready_simple_action_callback (Parallel::GlobalCache< Metavariables > &cache, const ArrayIndex &array_index, const Component *component_p, const double time, const std::optional< std::unordered_set< std::string > > &functions_to_check, Args &&... args) |
Check that functions of time are up-to-date. More... | |
template<typename CacheTag , typename ThreadedAction , typename Metavariables , typename ArrayIndex , typename Component , typename... Args> | |
bool | domain::functions_of_time_are_ready_threaded_action_callback (Parallel::GlobalCache< Metavariables > &cache, const ArrayIndex &array_index, const Component *component_p, const double time, const std::optional< std::unordered_set< std::string > > &functions_to_check, Args &&... args) |
Check that functions of time are up-to-date. More... | |
template<typename CacheTag , size_t Dim, typename Metavariables , typename ArrayIndex , typename Component > | |
bool | domain::functions_of_time_are_ready_algorithm_callback (Parallel::GlobalCache< Metavariables > &cache, const ArrayIndex &array_index, const Component *component_p, const double time, const std::optional< std::unordered_set< std::string > > &functions_to_check=std::nullopt) |
Check that functions of time are up-to-date. More... | |
template<size_t Dim, typename Fr > | |
auto | block_logical_coordinates (const Domain< Dim > &domain, const tnsr::I< DataVector, Dim, Fr > &x, double time=std::numeric_limits< double >::signaling_NaN(), const domain::FunctionsOfTimeMap &functions_of_time={}) -> std::vector< BlockLogicalCoords< Dim > > |
Computes the block logical coordinates and the containing BlockId of a set of points, given coordinates in a particular frame. More... | |
template<size_t Dim, typename Fr > | |
std::optional< tnsr::I< double, Dim, ::Frame::BlockLogical > > | block_logical_coordinates_single_point (const tnsr::I< double, Dim, Fr > &input_point, const Block< Dim > &block, double time=std::numeric_limits< double >::signaling_NaN(), const domain::FunctionsOfTimeMap &functions_of_time={}) |
Computes the block logical coordinates and the containing BlockId of a set of points, given coordinates in a particular frame. More... | |
template<size_t VolumeDim, typename TargetFrame > | |
void | unnormalized_face_normal (const gsl::not_null< tnsr::i< DataVector, VolumeDim, TargetFrame > * > result, const Mesh< VolumeDim - 1 > &interface_mesh, const InverseJacobian< DataVector, VolumeDim, Frame::ElementLogical, TargetFrame > &inv_jacobian_on_interface, const Direction< VolumeDim > &direction) |
Compute the outward grid normal on a face of an Element. More... | |
template<size_t VolumeDim, typename TargetFrame > | |
tnsr::i< DataVector, VolumeDim, TargetFrame > | unnormalized_face_normal (const Mesh< VolumeDim - 1 > &interface_mesh, const InverseJacobian< DataVector, VolumeDim, Frame::ElementLogical, TargetFrame > &inv_jacobian_on_interface, const Direction< VolumeDim > &direction) |
Compute the outward grid normal on a face of an Element. More... | |
template<size_t VolumeDim, typename TargetFrame > | |
void | unnormalized_face_normal (gsl::not_null< tnsr::i< DataVector, VolumeDim, TargetFrame > * > result, const Mesh< VolumeDim - 1 > &interface_mesh, const ElementMap< VolumeDim, TargetFrame > &map, const Direction< VolumeDim > &direction) |
Compute the outward grid normal on a face of an Element. More... | |
template<size_t VolumeDim, typename TargetFrame > | |
tnsr::i< DataVector, VolumeDim, TargetFrame > | unnormalized_face_normal (const Mesh< VolumeDim - 1 > &interface_mesh, const ElementMap< VolumeDim, TargetFrame > &map, const Direction< VolumeDim > &direction) |
Compute the outward grid normal on a face of an Element. More... | |
template<size_t VolumeDim, typename TargetFrame > | |
void | unnormalized_face_normal (gsl::not_null< tnsr::i< DataVector, VolumeDim, TargetFrame > * > result, const Mesh< VolumeDim - 1 > &interface_mesh, const domain::CoordinateMapBase< Frame::ElementLogical, TargetFrame, VolumeDim > &map, const Direction< VolumeDim > &direction) |
Compute the outward grid normal on a face of an Element. More... | |
template<size_t VolumeDim, typename TargetFrame > | |
tnsr::i< DataVector, VolumeDim, TargetFrame > | unnormalized_face_normal (const Mesh< VolumeDim - 1 > &interface_mesh, const domain::CoordinateMapBase< Frame::ElementLogical, TargetFrame, VolumeDim > &map, const Direction< VolumeDim > &direction) |
Compute the outward grid normal on a face of an Element. More... | |
template<size_t VolumeDim> | |
void | unnormalized_face_normal (gsl::not_null< tnsr::i< DataVector, VolumeDim, Frame::Inertial > * > result, const Mesh< VolumeDim - 1 > &interface_mesh, const ElementMap< VolumeDim, Frame::Grid > &logical_to_grid_map, const domain::CoordinateMapBase< Frame::Grid, Frame::Inertial, VolumeDim > &grid_to_inertial_map, double time, const std::unordered_map< std::string, std::unique_ptr< domain::FunctionsOfTime::FunctionOfTime > > &functions_of_time, const Direction< VolumeDim > &direction) |
Compute the outward grid normal on a face of an Element. More... | |
template<size_t VolumeDim> | |
tnsr::i< DataVector, VolumeDim, Frame::Inertial > | unnormalized_face_normal (const Mesh< VolumeDim - 1 > &interface_mesh, const ElementMap< VolumeDim, Frame::Grid > &logical_to_grid_map, const domain::CoordinateMapBase< Frame::Grid, Frame::Inertial, VolumeDim > &grid_to_inertial_map, double time, const std::unordered_map< std::string, std::unique_ptr< domain::FunctionsOfTime::FunctionOfTime > > &functions_of_time, const Direction< VolumeDim > &direction) |
Compute the outward grid normal on a face of an Element. More... | |
template<size_t Dim, typename Fr > | |
void | domain::jacobian_diagnostic (const gsl::not_null< tnsr::i< DataVector, Dim, Frame::ElementLogical > * > jacobian_diag, const ::Jacobian< DataVector, Dim, Frame::ElementLogical, Fr > &analytic_jacobian, const tnsr::I< DataVector, Dim, Fr > &mapped_coords, const ::Mesh< Dim > &mesh) |
A diagnostic comparing the analytic and numerical Jacobians for a map. More... | |
template<size_t Dim, typename Fr > | |
tnsr::i< DataVector, Dim, Frame::ElementLogical > | domain::jacobian_diagnostic (const ::Jacobian< DataVector, Dim, Frame::ElementLogical, Fr > &analytic_jacobian, const tnsr::I< DataVector, Dim, Fr > &mapped_coords, const ::Mesh< Dim > &mesh) |
A diagnostic comparing the analytic and numerical Jacobians for a map. More... | |
template<typename DataType , size_t Dim, typename CoordsFrame > | |
void | domain::radially_compressed_coordinates (gsl::not_null< tnsr::I< DataType, Dim, CoordsFrame > * > result, const tnsr::I< DataType, Dim, CoordsFrame > &coordinates, double inner_radius, double outer_radius, CoordinateMaps::Distribution compression) |
Coordinates suitable for visualizing large radii by compressing them logarithmically or inversely. More... | |
template<typename DataType , size_t Dim, typename CoordsFrame > | |
tnsr::I< DataType, Dim, CoordsFrame > | domain::radially_compressed_coordinates (const tnsr::I< DataType, Dim, CoordsFrame > &coordinates, double inner_radius, double outer_radius, CoordinateMaps::Distribution compression) |
Coordinates suitable for visualizing large radii by compressing them logarithmically or inversely. More... | |
template<size_t VolumeDim> | |
std::array< double, VolumeDim > | size_of_element (const ElementMap< VolumeDim, Frame::Inertial > &logical_to_inertial_map) |
Compute the inertial-coordinate size of an element along each of its logical directions. More... | |
template<size_t VolumeDim> | |
std::array< double, VolumeDim > | size_of_element (const ElementMap< VolumeDim, Frame::Grid > &logical_to_grid_map, const domain::CoordinateMapBase< Frame::Grid, Frame::Inertial, VolumeDim > &grid_to_inertial_map, double time, const std::unordered_map< std::string, std::unique_ptr< domain::FunctionsOfTime::FunctionOfTime > > &functions_of_time) |
Compute the inertial-coordinate size of an element along each of its logical directions. More... | |
template<typename VectorType , size_t VolumeDim> | |
void | orient_variables (gsl::not_null< VectorType * > result, const VectorType &variables, const Index< VolumeDim > &extents, const OrientationMap< VolumeDim > &orientation_of_neighbor) |
Orient a DataVector , ComplexDataVector , std::vector<double> , or std::vector<std::complex<double>> to the data-storage order of a neighbor element with the given orientation. More... | |
template<typename VectorType , size_t VolumeDim> | |
VectorType | orient_variables (const VectorType &variables, const Index< VolumeDim > &extents, const OrientationMap< VolumeDim > &orientation_of_neighbor) |
Orient a DataVector , ComplexDataVector , std::vector<double> , or std::vector<std::complex<double>> to the data-storage order of a neighbor element with the given orientation. More... | |
template<typename VectorType , size_t VolumeDim> | |
void | orient_variables_on_slice (gsl::not_null< VectorType * > result, const VectorType &variables_on_slice, const Index< VolumeDim - 1 > &slice_extents, size_t sliced_dim, const OrientationMap< VolumeDim > &orientation_of_neighbor) |
Orient a DataVector , ComplexDataVector , std::vector<double> , or std::vector<std::complex<double>> to the data-storage order of a neighbor element with the given orientation. More... | |
template<typename VectorType , size_t VolumeDim> | |
VectorType | orient_variables_on_slice (const VectorType &variables_on_slice, const Index< VolumeDim - 1 > &slice_extents, size_t sliced_dim, const OrientationMap< VolumeDim > &orientation_of_neighbor) |
Orient a DataVector , ComplexDataVector , std::vector<double> , or std::vector<std::complex<double>> to the data-storage order of a neighbor element with the given orientation. More... | |
template<size_t VolumeDim, typename TagsList > | |
Variables< TagsList > | orient_variables (const Variables< TagsList > &variables, const Index< VolumeDim > &extents, const OrientationMap< VolumeDim > &orientation_of_neighbor) |
Orient variables to the data-storage order of a neighbor element with the given orientation. | |
template<size_t VolumeDim, typename TagsList > | |
Variables< TagsList > | orient_variables_on_slice (const Variables< TagsList > &variables_on_slice, const Index< VolumeDim - 1 > &slice_extents, const size_t sliced_dim, const OrientationMap< VolumeDim > &orientation_of_neighbor) |
Orient variables to the data-storage order of a neighbor element with the given orientation. | |
template<size_t VolumeDim> | |
void | logical_coordinates (gsl::not_null< tnsr::I< DataVector, VolumeDim, Frame::ElementLogical > * > logical_coords, const Mesh< VolumeDim > &mesh) |
Compute the logical coordinates of a Mesh in an Element. More... | |
template<size_t VolumeDim> | |
tnsr::I< DataVector, VolumeDim, Frame::ElementLogical > | logical_coordinates (const Mesh< VolumeDim > &mesh) |
Compute the logical coordinates of a Mesh in an Element. More... | |
The building blocks used to describe the computational domain.
The VolumeDim-dimensional computational Domain is constructed from a set of non-overlapping Blocks. Each Block is a distorted VolumeDim-dimensional hypercube. Each codimension-1 boundary of a Block is either part of the external boundary of the computational domain, or is identical to a boundary of one other Block. Each Block is subdivided into one or more Elements that may be changed dynamically if AMR is enabled.
The VolumeDim-dimensional computational Domain is constructed from a set of non-overlapping Blocks. Each Block is a distorted VolumeDim-dimensional hypercube. Each codimension-1 boundary of a Block is either part of the external boundary of the computational domain, or is identical to a boundary of one other Block. Each Block is subdivided into one or more Elements that may be changed dynamically if AMR is enabled.
#define INSTANTIATE_MAPS_DATA_TYPE_FUNCTIONS | ( | _, | |
data | |||
) |
Generate instantiations of member functions of the CoordinateMap
class template.
Called as follows:
The first tuple passed to GENERATE_INSTANTIATIONS
has a bunch of tuples in it that is the list of maps being composed. The reason for defining the type aliases Affine2d
and Affine3d
is that otherwise the number of maps being composed is calculated incorrectly. The second tuple contains the source frames for the map. The third tuple passed to GENERATE_INSTANTIATIONS
contains the target frames to instantiate for, typically Frame::Grid
and Frame::Inertial
. The last tuple is the data types for which to instantiate the functions, usually double
and DataVector
.
Instantiates:
call_impl
inv_jacobian_impl
jacobian_impl
coords_frame_velocity_jacobians_impl
#define INSTANTIATE_MAPS_FUNCTIONS | ( | MAPS_TUPLE, | |
SOURCE_FRAME, | |||
TARGET_FRAMES_TUPLE, | |||
TYPES_TUPLE | |||
) |
Generate instantiations of member functions of the CoordinateMap
class template.
Called as follows:
The first tuple passed to GENERATE_INSTANTIATIONS
has a bunch of tuples in it that is the list of maps being composed. The reason for defining the type aliases Affine2d
and Affine3d
is that otherwise the number of maps being composed is calculated incorrectly. The second tuple contains the source frames for the map. The third tuple passed to GENERATE_INSTANTIATIONS
contains the frames to instantiate for, typically Frame::Grid
and Frame::Inertial
.
Instantiates:
get_to_grid_frame_impl
inverse_impl
class CoordinateMap
call_impl
inv_jacobian_impl
jacobian_impl
coords_frame_velocity_jacobians_impl
#define INSTANTIATE_MAPS_SIMPLE_FUNCTIONS | ( | _, | |
data | |||
) |
Generate instantiations of member functions of the CoordinateMap
class template.
Called as follows:
The first tuple passed to GENERATE_INSTANTIATIONS
has a bunch of tuples in it that is the list of maps being composed. The reason for defining the type aliases Affine2d
and Affine3d
is that otherwise the number of maps being composed is calculated incorrectly. The second tuple contains the source frames for the map. The third tuple passed to GENERATE_INSTANTIATIONS
contains the target frames to instantiate for, typically Frame::Grid
and Frame::Inertial
.
Instantiates:
get_to_grid_frame_impl
inverse_impl
class CoordinateMap
|
strong |
|
strong |
A label for the side of a manifold.
Lower and Upper are with respect to the logical coordinate whose axis is normal to the side, i.e. beyond the Upper (Lower) side, the logical coordinate is increasing (decreasing).
Self
is used to mark when an ElementId
does not have a direction.
Currently Direction
assumes this enum uses no more than 2 bits.
auto block_logical_coordinates | ( | const Domain< Dim > & | domain, |
const tnsr::I< DataVector, Dim, Fr > & | x, | ||
double | time = std::numeric_limits< double >::signaling_NaN() , |
||
const domain::FunctionsOfTimeMap & | functions_of_time = {} |
||
) | -> std::vector< BlockLogicalCoords< Dim > > |
Computes the block logical coordinates and the containing BlockId
of a set of points, given coordinates in a particular frame.
Returns a std::vector<std::optional<IdPair<BlockId,coords>>>, where the vector runs over the points and is indexed in the same order as the input coordinates x
. For each point, the IdPair
holds the block logical coords of that point and the BlockId
of the Block
that contains that point. The std::optional is invalid if the point is not in any Block. If a point is on a shared boundary of two or more Block
s, it is returned only once, and is considered to belong to the Block
with the smaller BlockId
.
The block_logical_coordinates_single_point
function will search the passed in block for the passed in coordinate and return the logical coordinates of that point. It will return a std::nullopt
if it can't find the point in that block.
element_logical_coordinates
.block_logical_coordinates
with x in Frame::Distorted
ignores all Block
s that lack a distorted frame, and it will return std::nullopt for points that lie outside all distorted-frame-endowed Block
s. This is what is expected for typical use cases. This means that block_logical_coordinates
does not assume that grid and distorted frames are equal in Block
s that lack a distorted frame. std::optional< tnsr::I< double, Dim, ::Frame::BlockLogical > > block_logical_coordinates_single_point | ( | const tnsr::I< double, Dim, Fr > & | input_point, |
const Block< Dim > & | block, | ||
double | time = std::numeric_limits< double >::signaling_NaN() , |
||
const domain::FunctionsOfTimeMap & | functions_of_time = {} |
||
) |
Computes the block logical coordinates and the containing BlockId
of a set of points, given coordinates in a particular frame.
Returns a std::vector<std::optional<IdPair<BlockId,coords>>>, where the vector runs over the points and is indexed in the same order as the input coordinates x
. For each point, the IdPair
holds the block logical coords of that point and the BlockId
of the Block
that contains that point. The std::optional is invalid if the point is not in any Block. If a point is on a shared boundary of two or more Block
s, it is returned only once, and is considered to belong to the Block
with the smaller BlockId
.
The block_logical_coordinates_single_point
function will search the passed in block for the passed in coordinate and return the logical coordinates of that point. It will return a std::nullopt
if it can't find the point in that block.
element_logical_coordinates
.block_logical_coordinates
with x in Frame::Distorted
ignores all Block
s that lack a distorted frame, and it will return std::nullopt for points that lie outside all distorted-frame-endowed Block
s. This is what is expected for typical use cases. This means that block_logical_coordinates
does not assume that grid and distorted frames are equal in Block
s that lack a distorted frame. std::vector< std::array< size_t, 8 > > corners_for_biradially_layered_domains | ( | size_t | number_of_radial_layers, |
size_t | number_of_biradial_layers, | ||
bool | include_central_block_lhs, | ||
bool | include_central_block_rhs, | ||
const std::array< size_t, 8 > & | central_block_corners_lhs = { {1, 2, 3, 4, 5, 6, 7, 8}} |
||
) |
The corners for a domain with biradial layers.
Generates the corners for a BBH-like Domain which is made of one or more layers of Blocks fully enveloping two interior volumes. The number_of_radial_layers
gives the number of layers that fully envelop each interior volume with six Blocks each. The number_of_biradial_layers
gives the number of layers that fully envelop both volumes at once, using ten Blocks per layer as opposed to six. The central_block_corners_lhs
are used as seed values to generate the corners for the surrounding Blocks.
std::vector< std::array< size_t, 8 > > corners_for_cylindrical_layered_domains | ( | size_t | number_of_shells, |
size_t | number_of_discs | ||
) |
The corners for a cylindrical domain split into discs with radial shells.
Generates the corners for a Domain which is made of one or more stacked discs consisting of layers of Blocks enveloping an interior square prism. The number_of_shells
specifies how many of these layers of Blocks to have in each disc.
The number_of_discs
specifies how many discs make up the domain. The very basic cylinder with one shell and one layer serves as a base to generate the corners for subsequent shells first and discs second.
std::vector< std::array< size_t, 8 > > corners_for_radially_layered_domains | ( | size_t | number_of_layers, |
bool | include_central_block, | ||
const std::array< size_t, 8 > & | central_block_corners = {{1, 2, 3, 4, 5, 6, 7, 8}} , |
||
ShellWedges | which_wedges = ShellWedges::All |
||
) |
The corners for a domain with radial layers.
Generates the corners for a Domain which is made of one or more layers of Blocks fully enveloping an interior volume, e.g. Sphere.
number_of_layers | specifies how many layers of Blocks to have in the final domain. |
include_central_block | set to true where the interior volume is filled with a central Block, and false where the interior volume is left empty. |
central_block_corners | are used as seed values to generate the corners for the surrounding Blocks. |
which_wedges | can be used to exclude a subset of the wedges. |
auto corners_for_rectilinear_domains | ( | const Index< VolumeDim > & | domain_extents, |
const std::vector< Index< VolumeDim > > & | block_indices_to_exclude = {} |
||
) | -> std::vector< std::array< size_t, two_to_the(VolumeDim)> > |
The corners for a rectilinear domain made of n-cubes.
The domain_extents
argument holds the number of blocks to have in each dimension. The blocks all have aligned orientations by construction. The block_indices_to_exclude
argument allows the user to selectively exclude blocks from the resulting domain. This allows for the creation of non-trivial shapes such as the net for a tesseract.
auto cyl_wedge_coord_map_center_blocks | ( | double | inner_radius, |
double | lower_z_bound, | ||
double | upper_z_bound, | ||
bool | use_equiangular_map, | ||
const std::vector< double > & | partitioning_in_z = {} , |
||
const std::vector< domain::CoordinateMaps::Distribution > & | distribution_in_z = {domain::CoordinateMaps::Distribution::Linear} , |
||
CylindricalDomainParityFlip | parity_flip = CylindricalDomainParityFlip::none |
||
) | -> std::vector< domain::CoordinateMaps::ProductOf3Maps< domain::CoordinateMaps::Interval, domain::CoordinateMaps::Interval, domain::CoordinateMaps::Interval > > |
Same as cyl_wedge_coordinate_maps
, but only the center square blocks,.
If CylindricalDomainParityFlip::z_direction
is specified, then the returned maps describe a cylinder with lower_z_bound
corresponding to logical coordinate upper_zeta
and upper_z_bound
corresponding to logical coordinate lower_zeta
, and thus the resulting maps are left-handed. CylindricalDomainParityFlip::z_direction
is therefore useful only when composing with another map that is also left-handed, so that the composed coordinate system is right-handed.
Returned as a vector of the coordinate maps so that they can be composed with other maps later.
auto cyl_wedge_coord_map_surrounding_blocks | ( | double | inner_radius, |
double | outer_radius, | ||
double | lower_z_bound, | ||
double | upper_z_bound, | ||
bool | use_equiangular_map, | ||
double | inner_circularity, | ||
const std::vector< double > & | radial_partitioning = {} , |
||
const std::vector< double > & | partitioning_in_z = {} , |
||
const std::vector< domain::CoordinateMaps::Distribution > & | radial_distribution = {domain::CoordinateMaps::Distribution::Linear} , |
||
const std::vector< domain::CoordinateMaps::Distribution > & | distribution_in_z = {domain::CoordinateMaps::Distribution::Linear} , |
||
CylindricalDomainParityFlip | parity_flip = CylindricalDomainParityFlip::none |
||
) | -> std::vector< domain::CoordinateMaps::ProductOf2Maps< domain::CoordinateMaps::Wedge< 2 >, domain::CoordinateMaps::Interval > > |
Same as cyl_wedge_coordinate_maps, but only the surrounding wedge blocks.
If CylindricalDomainParityFlip::z_direction
is specified, then the returned maps describe a cylinder with lower_z_bound
corresponding to logical coordinate upper_zeta
and upper_z_bound
corresponding to logical coordinate lower_zeta
, and thus the resulting maps are left-handed. CylindricalDomainParityFlip::z_direction
is therefore useful only when composing with another map that is also left-handed, so that the composed coordinate system is right-handed.
Returned as a vector of the coordinate maps so that they can be composed with other maps later.
auto cyl_wedge_coordinate_maps | ( | double | inner_radius, |
double | outer_radius, | ||
double | lower_z_bound, | ||
double | upper_z_bound, | ||
bool | use_equiangular_map, | ||
const std::vector< double > & | radial_partitioning = {} , |
||
const std::vector< double > & | partitioning_in_z = {} , |
||
const std::vector< domain::CoordinateMaps::Distribution > & | radial_distribution = {domain::CoordinateMaps::Distribution::Linear} , |
||
const std::vector< domain::CoordinateMaps::Distribution > & | distribution_in_z = {domain::CoordinateMaps::Distribution::Linear} |
||
) | -> std::vector< std::unique_ptr< domain::CoordinateMapBase< Frame::BlockLogical, TargetFrame, 3 > > > |
These are the CoordinateMaps used in the Cylinder DomainCreator.
The radial_partitioning
specifies the radial boundaries of sub-shells between inner_radius
and outer_radius
, while partitioning_in_z
specifies the z-boundaries, splitting the cylinder into stacked 3-dimensional disks. The circularity of the shell wedges changes from 0 to 1 within the innermost sub-shell.
Set the radial_distribution
to select the radial distribution of grid points in the cylindrical shells. The innermost shell must have domain::CoordinateMaps::Distribution::Linear
because it changes the circularity. The distribution along the z-axis for each circular disc is specified through distribution_in_z
.
std::array< size_t, two_to_the(VolumeDim)> discrete_rotation | ( | const OrientationMap< VolumeDim > & | orientation, |
const std::array< size_t, two_to_the(VolumeDim)> & | corners_of_aligned | ||
) |
Permutes the corner numbers of an n-cube.
Returns the correct ordering of global corner numbers for a rotated block in an otherwise aligned edifice of blocks, given the OrientationMap a block aligned with the edifice has relative to this one, and given the corner numbering the rotated block would have if it were aligned. This is useful in creating domains for testing purposes, e.g. RotatedIntervals, RotatedRectangles, and RotatedBricks.
tnsr::Ij< double, VolumeDim, Frame::NoFrame > discrete_rotation_jacobian | ( | const OrientationMap< VolumeDim > & | orientation | ) |
Computes the Jacobian of the transformation that is computed by discrete_rotation()
double
because the Jacobian is spatially constant. auto element_logical_coordinates | ( | const std::vector< ElementId< Dim > > & | element_ids, |
const std::vector< BlockLogicalCoords< Dim > > & | block_coord_holders | ||
) | -> std::unordered_map< ElementId< Dim >, ElementLogicalCoordHolder< Dim > > |
Given a set of points in block logical coordinates and their BlockIds
, as returned from the function block_logical_coordinates
, determines which Element
s in a list of ElementId
s contains each point, and determines the element logical coordinates of each point.
Returns a std::unordered_map from ElementId
s to ElementLogicalCoordHolder
s. It is expected that only a subset of the points will be found in the given Element
s. Boundary points: If a point is on the boundary of an Element, it is considered contained in that Element only if it is on the lower bound of the Element, or if it is on the upper bound of the element and that upper bound coincides with the upper bound of the containing Block. This means that each boundary point is contained in one and only one Element. We assume that the input block_coord_holders associates a point on a Block boundary with only a single Block, the one with the smaller BlockId, which is always the lower-bounding Block.
std::vector< domain::CoordinateMaps::Frustum > frustum_coordinate_maps | ( | double | length_inner_cube, |
double | length_outer_cube, | ||
bool | equiangular_map_at_outer, | ||
bool | equiangular_map_at_inner, | ||
const std::array< double, 3 > & | origin_preimage = {{0.0, 0.0, 0.0}} , |
||
domain::CoordinateMaps::Distribution | radial_distribution = domain::CoordinateMaps::Distribution::Linear , |
||
std::optional< double > | distribution_value = std::nullopt , |
||
double | sphericity = 0.0 , |
||
double | opening_angle = M_PI_2 |
||
) |
These are the ten Frustums used in the DomainCreators for binary compact objects. The Frustums partition the volume defined by two bounding surfaces: The inner surface is the surface of the two joined inner cubes enveloping the two compact objects, while the outer is the surface of the outer cube.
length_inner_cube | The side length of the cubes enveloping the two shells. |
length_outer_cube | The side length of the outer cube. |
equiangular_map_at_outer | Whether to apply a tangent map in the angular directions at the outer boundary. |
equiangular_map_at_inner | Whether to apply a tangent map in the angular directions at the inner boundary. |
origin_preimage | The center of the two joined inner cubes is moved away from the origin and to this point, origin_preimage. |
radial_distribution | The gridpoint distribution in the radial direction, possibly dependent on the value passed to distribution_value . |
distribution_value | Used by radial_distribution . |
sphericity | Determines whether the outer surface is a cube (value of 0), a sphere (value of 1) or somewhere in between. |
opening_angle | determines the gridpoint distribution used in the Frustums such that they conform to the outer sphere of Wedges with the same value for opening_angle . |
bool domain::functions_of_time_are_ready_algorithm_callback | ( | Parallel::GlobalCache< Metavariables > & | cache, |
const ArrayIndex & | array_index, | ||
const Component * | component_p, | ||
const double | time, | ||
const std::optional< std::unordered_set< std::string > > & | functions_to_check = std::nullopt |
||
) |
Check that functions of time are up-to-date.
Check that functions of time in CacheTag
with names in functions_to_check
are ready at time time
. If functions_to_check
is a std::nullopt
, checks all functions in CacheTag
. If any function is not ready, schedules a Parallel::PerformAlgorithmCallback
or Parallel::Actions::PerformAlgorithmOnElement<false>
callback with the GlobalCache.
bool domain::functions_of_time_are_ready_simple_action_callback | ( | Parallel::GlobalCache< Metavariables > & | cache, |
const ArrayIndex & | array_index, | ||
const Component * | component_p, | ||
const double | time, | ||
const std::optional< std::unordered_set< std::string > > & | functions_to_check, | ||
Args &&... | args | ||
) |
Check that functions of time are up-to-date.
Check that functions of time in CacheTag
with names in functions_to_check
are ready at time time
. If functions_to_check
is a std::nullopt
, checks all functions in CacheTag
. If any function is not ready, schedules a Parallel::SimpleActionCallback
with the GlobalCache which calls the simple action passed in as a template parameter. The Args
are forwareded to the callback.
bool domain::functions_of_time_are_ready_threaded_action_callback | ( | Parallel::GlobalCache< Metavariables > & | cache, |
const ArrayIndex & | array_index, | ||
const Component * | component_p, | ||
const double | time, | ||
const std::optional< std::unordered_set< std::string > > & | functions_to_check, | ||
Args &&... | args | ||
) |
Check that functions of time are up-to-date.
Check that functions of time in CacheTag
with names in functions_to_check
are ready at time time
. If functions_to_check
is a std::nullopt
, checks all functions in CacheTag
. If any function is not ready, schedules a Parallel::ThreadedActionCallback
with the GlobalCache which calls the threaded action passed in as a template parameter. The Args
are forwareded to the callback.
size_t index_to_slice_at | ( | const Index< Dim > & | extents, |
const Direction< Dim > & | direction, | ||
const size_t | offset = 0 |
||
) |
Finds the index in the perpendicular dimension of an element boundary.
Optionally provide an offset
to find an index offset from the element boundary.
tnsr::I< DataVector, VolumeDim, Frame::ElementLogical > interface_logical_coordinates | ( | const Mesh< VolumeDim - 1 > & | mesh, |
const Direction< VolumeDim > & | direction | ||
) |
Defines functions interface_logical_coordinates.
Compute the logical coordinates on a face of an Element.
Returns: element logical-frame vector holding coordinates
tnsr::i< DataVector, Dim, Frame::ElementLogical > domain::jacobian_diagnostic | ( | const ::Jacobian< DataVector, Dim, Frame::ElementLogical, Fr > & | analytic_jacobian, |
const tnsr::I< DataVector, Dim, Fr > & | mapped_coords, | ||
const ::Mesh< Dim > & | mesh | ||
) |
A diagnostic comparing the analytic and numerical Jacobians for a map.
Specifically, returns
, where
void domain::jacobian_diagnostic | ( | const gsl::not_null< tnsr::i< DataVector, Dim, Frame::ElementLogical > * > | jacobian_diag, |
const ::Jacobian< DataVector, Dim, Frame::ElementLogical, Fr > & | analytic_jacobian, | ||
const tnsr::I< DataVector, Dim, Fr > & | mapped_coords, | ||
const ::Mesh< Dim > & | mesh | ||
) |
A diagnostic comparing the analytic and numerical Jacobians for a map.
Specifically, returns
, where
void domain::jacobian_diagnostic | ( | const gsl::not_null< tnsr::i< DataVector, Dim, typename Frame::ElementLogical > * > | jacobian_diag, |
const Jacobian< DataVector, Dim, typename Frame::ElementLogical, Fr > & | analytic_jacobian, | ||
const TensorMetafunctions::prepend_spatial_index< tnsr::I< DataVector, Dim, Fr >, Dim, UpLo::Lo, typename Frame::ElementLogical > & | numeric_jacobian_transpose | ||
) |
A diagnostic comparing the analytic and numerical Jacobians for a map.
Specifically, returns
, where
tnsr::I< DataVector, VolumeDim, Frame::ElementLogical > logical_coordinates | ( | const Mesh< VolumeDim > & | mesh | ) |
Compute the logical coordinates of a Mesh in an Element.
The logical coordinates are the collocation points associated with the Spectral::Basis and Spectral::Quadrature of the mesh
. The Spectral::Basis determines the domain of the logical coordinates, while the Spectral::Quadrature determines their distribution. For Legendre or Chebyshev bases, the logical coordinates are in the interval
void logical_coordinates | ( | gsl::not_null< tnsr::I< DataVector, VolumeDim, Frame::ElementLogical > * > | logical_coords, |
const Mesh< VolumeDim > & | mesh | ||
) |
Compute the logical coordinates of a Mesh in an Element.
The logical coordinates are the collocation points associated with the Spectral::Basis and Spectral::Quadrature of the mesh
. The Spectral::Basis determines the domain of the logical coordinates, while the Spectral::Quadrature determines their distribution. For Legendre or Chebyshev bases, the logical coordinates are in the interval
auto maps_for_rectilinear_domains | ( | const Index< VolumeDim > & | domain_extents, |
const std::array< std::vector< double >, VolumeDim > & | block_demarcations, | ||
const std::vector< Index< VolumeDim > > & | block_indices_to_exclude = {} , |
||
const std::vector< OrientationMap< VolumeDim > > & | orientations_of_all_blocks = {} , |
||
bool | use_equiangular_map = false |
||
) | -> std::vector< std::unique_ptr< domain::CoordinateMapBase< Frame::BlockLogical, TargetFrame, VolumeDim > > > |
The CoordinateMaps for a rectilinear domain of n-cubes.
Allows for both Affine and Equiangular maps.
|
constexpr |
Returns the maximum number of neighbors an element can have in dim
dimensions.
|
constexpr |
Returns the maximum number of neighbors in each direction an element can have in dim
dimensions.
VectorType orient_variables | ( | const VectorType & | variables, |
const Index< VolumeDim > & | extents, | ||
const OrientationMap< VolumeDim > & | orientation_of_neighbor | ||
) |
Orient a DataVector
, ComplexDataVector
, std::vector<double>
, or std::vector<std::complex<double>>
to the data-storage order of a neighbor element with the given orientation.
The vector may represent more than one tensor component over the grid represented by extents
.
variables.size()
).In most cases the Variables
version of orient_variables
should be called. However, in some cases the tags and thus the type of the data being sent is determined at runtime. In these cases the std::vector
version of orient_variables
is useful. A concrete example of this is when hybridizing DG with finite difference methods, where sometimes the data sent is both the variables for reconstruction and the fluxes for either the DG or finite difference scheme, while at other points only one of these three is sent.
void orient_variables | ( | gsl::not_null< VectorType * > | result, |
const VectorType & | variables, | ||
const Index< VolumeDim > & | extents, | ||
const OrientationMap< VolumeDim > & | orientation_of_neighbor | ||
) |
Orient a DataVector
, ComplexDataVector
, std::vector<double>
, or std::vector<std::complex<double>>
to the data-storage order of a neighbor element with the given orientation.
The vector may represent more than one tensor component over the grid represented by extents
.
variables.size()
).In most cases the Variables
version of orient_variables
should be called. However, in some cases the tags and thus the type of the data being sent is determined at runtime. In these cases the std::vector
version of orient_variables
is useful. A concrete example of this is when hybridizing DG with finite difference methods, where sometimes the data sent is both the variables for reconstruction and the fluxes for either the DG or finite difference scheme, while at other points only one of these three is sent.
VectorType orient_variables_on_slice | ( | const VectorType & | variables_on_slice, |
const Index< VolumeDim - 1 > & | slice_extents, | ||
size_t | sliced_dim, | ||
const OrientationMap< VolumeDim > & | orientation_of_neighbor | ||
) |
Orient a DataVector
, ComplexDataVector
, std::vector<double>
, or std::vector<std::complex<double>>
to the data-storage order of a neighbor element with the given orientation.
The vector may represent more than one tensor component over the grid represented by extents
.
variables.size()
).In most cases the Variables
version of orient_variables
should be called. However, in some cases the tags and thus the type of the data being sent is determined at runtime. In these cases the std::vector
version of orient_variables
is useful. A concrete example of this is when hybridizing DG with finite difference methods, where sometimes the data sent is both the variables for reconstruction and the fluxes for either the DG or finite difference scheme, while at other points only one of these three is sent.
void orient_variables_on_slice | ( | gsl::not_null< VectorType * > | result, |
const VectorType & | variables_on_slice, | ||
const Index< VolumeDim - 1 > & | slice_extents, | ||
size_t | sliced_dim, | ||
const OrientationMap< VolumeDim > & | orientation_of_neighbor | ||
) |
Orient a DataVector
, ComplexDataVector
, std::vector<double>
, or std::vector<std::complex<double>>
to the data-storage order of a neighbor element with the given orientation.
The vector may represent more than one tensor component over the grid represented by extents
.
variables.size()
).In most cases the Variables
version of orient_variables
should be called. However, in some cases the tags and thus the type of the data being sent is determined at runtime. In these cases the std::vector
version of orient_variables
is useful. A concrete example of this is when hybridizing DG with finite difference methods, where sometimes the data sent is both the variables for reconstruction and the fluxes for either the DG or finite difference scheme, while at other points only one of these three is sent.
std::array< OrientationMap< 3 >, 6 > orientations_for_sphere_wrappings | ( | ) |
An array of the orientations of the six blocks that make up a Sphere.
A Block or Blocks can be wrapped in an outer layer of Blocks surrounding the original Block(s). In the BBH Domain, this occurs several times, using both Wedges and Frustums. This standardizes the ordering of the orientations for both.
tnsr::I< DataType, Dim, CoordsFrame > domain::radially_compressed_coordinates | ( | const tnsr::I< DataType, Dim, CoordsFrame > & | coordinates, |
double | inner_radius, | ||
double | outer_radius, | ||
CoordinateMaps::Distribution | compression | ||
) |
Coordinates suitable for visualizing large radii by compressing them logarithmically or inversely.
Rescales the coordinates
for inner_radius
. The coordinates are compressed from outer_radius
so the compressed outer radius is a multiple of the inner radius and increases with the outer radius as well, but exponentials are tamed.
The radial compression map domain::CoordinateMaps::Interval
map, which is also used to distribute grid points radially. Therefore, radial grid points will be distributed linearly in the radially compressed coordinates if you use the same compression
distribution that you used to distribute radial grid points in the CoordsFrame
.
void domain::radially_compressed_coordinates | ( | gsl::not_null< tnsr::I< DataType, Dim, CoordsFrame > * > | result, |
const tnsr::I< DataType, Dim, CoordsFrame > & | coordinates, | ||
double | inner_radius, | ||
double | outer_radius, | ||
CoordinateMaps::Distribution | compression | ||
) |
Coordinates suitable for visualizing large radii by compressing them logarithmically or inversely.
Rescales the coordinates
for inner_radius
. The coordinates are compressed from outer_radius
so the compressed outer radius is a multiple of the inner radius and increases with the outer radius as well, but exponentials are tamed.
The radial compression map domain::CoordinateMaps::Interval
map, which is also used to distribute grid points radially. Therefore, radial grid points will be distributed linearly in the radially compressed coordinates if you use the same compression
distribution that you used to distribute radial grid points in the CoordsFrame
.
Domain< VolumeDim > rectilinear_domain | ( | const Index< VolumeDim > & | domain_extents, |
const std::array< std::vector< double >, VolumeDim > & | block_demarcations, | ||
const std::vector< Index< VolumeDim > > & | block_indices_to_exclude = {} , |
||
const std::vector< OrientationMap< VolumeDim > > & | orientations_of_all_blocks = {} , |
||
const std::array< bool, VolumeDim > & | dimension_is_periodic = make_array< VolumeDim >(false) , |
||
const std::vector< PairOfFaces > & | identifications = {} , |
||
bool | use_equiangular_map = false |
||
) |
Create a rectilinear Domain of multicubes.
Useful for constructing domains for testing non-trivially connected rectilinear domains made up of cubes. We refer to a domain of this type as an edifice. The domain_extents
provides the size (in the number of blocks) of the initial aligned edifice to construct. The block_indices_to_exclude
parameter is used in refining the shape of the edifice from a cube to sometime more non-trivial, such as an L-shape or the net of a tesseract. The block_demarcations
and use_equiangular_map
parameters determine the CoordinateMaps to be used. orientations_of_all_blocks
contains the OrientationMap of the edifice relative to each block.
The identifications
parameter is used when identifying the faces of blocks in an edifice. This is used to identify the 1D boundaries in the 2D net for a 3D cube to construct a domain with topology S2. Note: If the user wishes to rotate the blocks as well as manually identify their faces, the user must provide the PairOfFaces corresponding to the rotated corners.
void set_internal_boundaries | ( | gsl::not_null< std::vector< DirectionMap< VolumeDim, BlockNeighbor< VolumeDim > > > * > | neighbors_of_all_blocks, |
const std::vector< std::unique_ptr< domain::CoordinateMapBase< Frame::BlockLogical, Frame::Inertial, VolumeDim > > > & | maps | ||
) |
Sets up the BlockNeighbors using the corner numbering scheme implied by the maps provided by the user to deduce the correct neighbors and orientations.
std::array< double, VolumeDim > size_of_element | ( | const ElementMap< VolumeDim, Frame::Grid > & | logical_to_grid_map, |
const domain::CoordinateMapBase< Frame::Grid, Frame::Inertial, VolumeDim > & | grid_to_inertial_map, | ||
double | time, | ||
const std::unordered_map< std::string, std::unique_ptr< domain::FunctionsOfTime::FunctionOfTime > > & | functions_of_time | ||
) |
Compute the inertial-coordinate size of an element along each of its logical directions.
For each logical direction, compute the distance (in inertial coordinates) between the element's lower and upper faces in that logical direction. The distance is measured between centers of the faces, with the centers defined in the logical coordinates. Note that for curved elements, this is an approximate measurement of size.
Because this quantity is defined in terms of specific coordinates, it is not well represented by a Tensor
, so we use a std::array
.
std::array< double, VolumeDim > size_of_element | ( | const ElementMap< VolumeDim, Frame::Inertial > & | logical_to_inertial_map | ) |
Compute the inertial-coordinate size of an element along each of its logical directions.
For each logical direction, compute the distance (in inertial coordinates) between the element's lower and upper faces in that logical direction. The distance is measured between centers of the faces, with the centers defined in the logical coordinates. Note that for curved elements, this is an approximate measurement of size.
Because this quantity is defined in terms of specific coordinates, it is not well represented by a Tensor
, so we use a std::array
.
std::vector< domain::CoordinateMaps::Wedge< 3 > > sph_wedge_coordinate_maps | ( | double | inner_radius, |
double | outer_radius, | ||
double | inner_sphericity, | ||
double | outer_sphericity, | ||
bool | use_equiangular_map, | ||
const std::optional< std::pair< double, std::array< double, 3 > > > & | offset_options = std::nullopt , |
||
bool | use_half_wedges = false , |
||
const std::vector< double > & | radial_partitioning = {} , |
||
const std::vector< domain::CoordinateMaps::Distribution > & | radial_distribution = {domain::CoordinateMaps::Distribution::Linear} , |
||
ShellWedges | which_wedges = ShellWedges::All , |
||
double | opening_angle = M_PI_2 |
||
) |
These are the CoordinateMaps of the Wedge<3>s used in the Sphere and binary compact object DomainCreators. This function can also be used to wrap the Sphere in a cube made of six Wedge<3>s.
inner_radius | Radius of the inner boundary of the shell, or the radius circumscribing the inner cube of a sphere. |
outer_radius | Outer radius of the shell or sphere. |
inner_sphericity | Specifies if the wedges form a spherical inner boundary (1.0) or a cubical inner boundary (0.0). |
outer_sphericity | Specifies if the wedges form a spherical outer boundary (1.0) or a cubical outer boundary (0.0). |
offset_options | A pair of values with the first being half the length of the cube that would form the outer boundary and the second being the offset to apply to the wedges. |
use_equiangular_map | Toggles the equiangular map of the Wedge map. |
use_half_wedges | When true , the wedges in the +z,-z,+y,-y directions are cut in half along their xi-axes. The resulting ten CoordinateMaps are used for the outermost Blocks of the BBH Domain. |
radial_partitioning | Specifies the radial boundaries of sub-shells between inner_radius and outer_radius . If the inner and outer sphericities are different, the innermost shell does the transition. |
radial_distribution | Select the radial distribution of grid points in the spherical shells. |
which_wedges | Select a subset of wedges. |
opening_angle | sets the combined opening angle of the two half wedges that open up along the y-z plane. The endcap wedges are then given an angle of pi minus this opening angle. This parameter only has an effect if use_half_wedges is set to true . |
void unnormalized_face_normal | ( | const gsl::not_null< tnsr::i< DataVector, VolumeDim, TargetFrame > * > | result, |
const Mesh< VolumeDim - 1 > & | interface_mesh, | ||
const InverseJacobian< DataVector, VolumeDim, Frame::ElementLogical, TargetFrame > & | inv_jacobian_on_interface, | ||
const Direction< VolumeDim > & | direction | ||
) |
Compute the outward grid normal on a face of an Element.
Computes the grid-frame normal by taking the logical-frame unit one-form in the given Direction and mapping it to the grid frame with the given map, or the given inverse Jacobian.
tnsr::i< DataVector, VolumeDim, TargetFrame > unnormalized_face_normal | ( | const Mesh< VolumeDim - 1 > & | interface_mesh, |
const domain::CoordinateMapBase< Frame::ElementLogical, TargetFrame, VolumeDim > & | map, | ||
const Direction< VolumeDim > & | direction | ||
) |
Compute the outward grid normal on a face of an Element.
Computes the grid-frame normal by taking the logical-frame unit one-form in the given Direction and mapping it to the grid frame with the given map, or the given inverse Jacobian.
tnsr::i< DataVector, VolumeDim, Frame::Inertial > unnormalized_face_normal | ( | const Mesh< VolumeDim - 1 > & | interface_mesh, |
const ElementMap< VolumeDim, Frame::Grid > & | logical_to_grid_map, | ||
const domain::CoordinateMapBase< Frame::Grid, Frame::Inertial, VolumeDim > & | grid_to_inertial_map, | ||
double | time, | ||
const std::unordered_map< std::string, std::unique_ptr< domain::FunctionsOfTime::FunctionOfTime > > & | functions_of_time, | ||
const Direction< VolumeDim > & | direction | ||
) |
Compute the outward grid normal on a face of an Element.
Computes the grid-frame normal by taking the logical-frame unit one-form in the given Direction and mapping it to the grid frame with the given map, or the given inverse Jacobian.
tnsr::i< DataVector, VolumeDim, TargetFrame > unnormalized_face_normal | ( | const Mesh< VolumeDim - 1 > & | interface_mesh, |
const ElementMap< VolumeDim, TargetFrame > & | map, | ||
const Direction< VolumeDim > & | direction | ||
) |
Compute the outward grid normal on a face of an Element.
Computes the grid-frame normal by taking the logical-frame unit one-form in the given Direction and mapping it to the grid frame with the given map, or the given inverse Jacobian.
tnsr::i< DataVector, VolumeDim, TargetFrame > unnormalized_face_normal | ( | const Mesh< VolumeDim - 1 > & | interface_mesh, |
const InverseJacobian< DataVector, VolumeDim, Frame::ElementLogical, TargetFrame > & | inv_jacobian_on_interface, | ||
const Direction< VolumeDim > & | direction | ||
) |
Compute the outward grid normal on a face of an Element.
Computes the grid-frame normal by taking the logical-frame unit one-form in the given Direction and mapping it to the grid frame with the given map, or the given inverse Jacobian.
void unnormalized_face_normal | ( | gsl::not_null< tnsr::i< DataVector, VolumeDim, Frame::Inertial > * > | result, |
const Mesh< VolumeDim - 1 > & | interface_mesh, | ||
const ElementMap< VolumeDim, Frame::Grid > & | logical_to_grid_map, | ||
const domain::CoordinateMapBase< Frame::Grid, Frame::Inertial, VolumeDim > & | grid_to_inertial_map, | ||
double | time, | ||
const std::unordered_map< std::string, std::unique_ptr< domain::FunctionsOfTime::FunctionOfTime > > & | functions_of_time, | ||
const Direction< VolumeDim > & | direction | ||
) |
Compute the outward grid normal on a face of an Element.
Computes the grid-frame normal by taking the logical-frame unit one-form in the given Direction and mapping it to the grid frame with the given map, or the given inverse Jacobian.
void unnormalized_face_normal | ( | gsl::not_null< tnsr::i< DataVector, VolumeDim, TargetFrame > * > | result, |
const Mesh< VolumeDim - 1 > & | interface_mesh, | ||
const domain::CoordinateMapBase< Frame::ElementLogical, TargetFrame, VolumeDim > & | map, | ||
const Direction< VolumeDim > & | direction | ||
) |
Compute the outward grid normal on a face of an Element.
Computes the grid-frame normal by taking the logical-frame unit one-form in the given Direction and mapping it to the grid frame with the given map, or the given inverse Jacobian.
void unnormalized_face_normal | ( | gsl::not_null< tnsr::i< DataVector, VolumeDim, TargetFrame > * > | result, |
const Mesh< VolumeDim - 1 > & | interface_mesh, | ||
const ElementMap< VolumeDim, TargetFrame > & | map, | ||
const Direction< VolumeDim > & | direction | ||
) |
Compute the outward grid normal on a face of an Element.
Computes the grid-frame normal by taking the logical-frame unit one-form in the given Direction and mapping it to the grid frame with the given map, or the given inverse Jacobian.