|
SpECTRE
v2025.08.19
|
Functions that are intended to be used by external programs, e.g. to interpolate data in volume files to target points. More...
Classes | |
| struct | ObservationId |
| Identifies an observation by its ID in the volume data file. More... | |
| struct | ObservationStep |
| Identifies an observation by its index in the ordered list of observations. Negative indices are counted from the end of the list. More... | |
| struct | PointwiseInterpolator |
| Interpolates data in volume files to target points by reading the volume data into memory. More... | |
| struct | SelectObservation |
| Determines the selected observation ID in the volume data file, given multiple possible ways to specify the observation. More... | |
| struct | SpacetimeInterpolator |
| Interpolates data in volume files in both space and time. More... | |
Typedefs | |
| using | ObservationVariant = std::variant< ObservationId, ObservationStep, double > |
Functions | |
| template<size_t Dim> | |
| std::vector< std::vector< double > > | interpolate_to_points (const std::variant< std::vector< std::string >, std::string > &volume_files_or_glob, const std::string &subfile_name, const ObservationVariant &observation, const std::vector< std::string > &tensor_components, const std::array< std::vector< double >, Dim > &target_points, bool extrapolate_into_excisions=false, bool error_on_missing_points=false, std::optional< size_t > num_threads=std::nullopt) |
| Interpolate data in volume files to target points. More... | |
| template<typename Tags > | |
| auto | get_tensor_components () |
| Collect tensor component names from Tags list. | |
| template<typename Tags , typename DataType > | |
| auto | make_tagged_tuple (std::vector< DataType > interpolated_data) |
| Convert tensor components to a tagged_tuple. | |
| template<typename ResultDataType , size_t Dim, typename Frame > | |
| void | interpolate_to_points (gsl::not_null< std::vector< ResultDataType > * > result, const std::variant< std::vector< std::string >, std::string > &volume_files_or_glob, const std::string &subfile_name, const ObservationVariant &observation, const std::vector< std::string > &tensor_components, const tnsr::I< DataVector, Dim, Frame > &target_points, bool extrapolate_into_excisions=false, bool error_on_missing_points=false, std::optional< size_t > num_threads=std::nullopt) |
| Interpolates data in volume files to target points. More... | |
| template<size_t Dim, typename Frame > | |
| std::vector< DataVector > | interpolate_to_points (const std::variant< std::vector< std::string >, std::string > &volume_files_or_glob, const std::string &subfile_name, const ObservationVariant &observation, const std::vector< std::string > &tensor_components, const tnsr::I< DataVector, Dim, Frame > &target_points, bool extrapolate_into_excisions=false, bool error_on_missing_points=false, std::optional< size_t > num_threads=std::nullopt) |
| Interpolates data in volume files to target points. More... | |
| template<typename Tags , size_t Dim, typename Frame > | |
| tuples::tagged_tuple_from_typelist< Tags > | interpolate_to_points (const std::variant< std::vector< std::string >, std::string > &volume_files_or_glob, const std::string &subfile_name, const ObservationVariant &observation, const tnsr::I< DataVector, Dim, Frame > &target_points, bool extrapolate_into_excisions=false, const bool error_on_missing_points=false, std::optional< size_t > num_threads=std::nullopt) |
| Interpolates data in volume files to target points. More... | |
Functions that are intended to be used by external programs, e.g. to interpolate data in volume files to target points.
| std::vector< std::vector< double > > spectre::Exporter::interpolate_to_points | ( | const std::variant< std::vector< std::string >, std::string > & | volume_files_or_glob, |
| const std::string & | subfile_name, | ||
| const ObservationVariant & | observation, | ||
| const std::vector< std::string > & | tensor_components, | ||
| const std::array< std::vector< double >, Dim > & | target_points, | ||
| bool | extrapolate_into_excisions = false, |
||
| bool | error_on_missing_points = false, |
||
| std::optional< size_t > | num_threads = std::nullopt |
||
| ) |
Interpolate data in volume files to target points.
| Dim | Dimension of the domain |
| volume_files_or_glob | The list of H5 files, or a glob pattern |
| subfile_name | The name of the subfile in the H5 files containing the volume data |
| observation | Either the observation ID as a size_t, or the index of the observation in the volume files to interpolate as an int (a value of 0 would be the first observation, and a value of -1 would be the last observation). |
| tensor_components | The tensor components to interpolate, e.g. "Lapse", "Shift_x", "Shift_y", "Shift_z", "SpatialMetric_xx", etc. Look into the H5 file to see what components are available. |
| target_points | The points to interpolate to, in inertial coordinates. |
| extrapolate_into_excisions | Enables extrapolation into excision regions of the domain (default is false). This can be useful to fill the excision region with (constraint-violating but smooth) data so it can be imported into moving puncture codes. Specifically, we implement the strategy used in [69] adjusted for distorted excisions: we choose uniformly spaced radial anchor points spaced as \(\Delta r = 0.3 r_\mathrm{AH}\) in the grid frame (where the excision is spherical), then map the anchor points to the distorted frame (where we have the target point) and do a 7th order polynomial extrapolation into the excision region. |
| error_on_missing_points | If true, an error will be thrown if any of the target points are outside the domain. If false, the result will be filled with NaNs for points outside the domain (default is false). |
| num_threads | The number of threads to use if OpenMP is linked in. If not specified, OpenMP will determine the number of threads automatically. It's also possible to set the number of threads using the environment variable OMP_NUM_THREADS. It's an error to specify num_threads if OpenMP is not linked in. Set num_threads to 1 to disable OpenMP. |
Returns: std::vector<std::vector<double>> The interpolated data. The first dimension corresponds to the selected tensor components, and the second dimension corresponds to the target points.
| std::vector< DataVector > spectre::Exporter::interpolate_to_points | ( | const std::variant< std::vector< std::string >, std::string > & | volume_files_or_glob, |
| const std::string & | subfile_name, | ||
| const ObservationVariant & | observation, | ||
| const std::vector< std::string > & | tensor_components, | ||
| const tnsr::I< DataVector, Dim, Frame > & | target_points, | ||
| bool | extrapolate_into_excisions = false, |
||
| bool | error_on_missing_points = false, |
||
| std::optional< size_t > | num_threads = std::nullopt |
||
| ) |
Interpolates data in volume files to target points.
These are overloads of the interpolate_to_points function that work with Tensor types and tags, rather than the raw C++ types that are used in Exporter.hpp so it can be used by external programs.
The Tags template parameter is a typelist of tags that should be read from the volume files. The dataset names to read are constructed from the tag names. Here is an example of how to use this function:
| tuples::tagged_tuple_from_typelist< Tags > spectre::Exporter::interpolate_to_points | ( | const std::variant< std::vector< std::string >, std::string > & | volume_files_or_glob, |
| const std::string & | subfile_name, | ||
| const ObservationVariant & | observation, | ||
| const tnsr::I< DataVector, Dim, Frame > & | target_points, | ||
| bool | extrapolate_into_excisions = false, |
||
| const bool | error_on_missing_points = false, |
||
| std::optional< size_t > | num_threads = std::nullopt |
||
| ) |
Interpolates data in volume files to target points.
These are overloads of the interpolate_to_points function that work with Tensor types and tags, rather than the raw C++ types that are used in Exporter.hpp so it can be used by external programs.
The Tags template parameter is a typelist of tags that should be read from the volume files. The dataset names to read are constructed from the tag names. Here is an example of how to use this function:
| void spectre::Exporter::interpolate_to_points | ( | gsl::not_null< std::vector< ResultDataType > * > | result, |
| const std::variant< std::vector< std::string >, std::string > & | volume_files_or_glob, | ||
| const std::string & | subfile_name, | ||
| const ObservationVariant & | observation, | ||
| const std::vector< std::string > & | tensor_components, | ||
| const tnsr::I< DataVector, Dim, Frame > & | target_points, | ||
| bool | extrapolate_into_excisions = false, |
||
| bool | error_on_missing_points = false, |
||
| std::optional< size_t > | num_threads = std::nullopt |
||
| ) |
Interpolates data in volume files to target points.
These are overloads of the interpolate_to_points function that work with Tensor types and tags, rather than the raw C++ types that are used in Exporter.hpp so it can be used by external programs.
The Tags template parameter is a typelist of tags that should be read from the volume files. The dataset names to read are constructed from the tag names. Here is an example of how to use this function: