SpECTRE
v2024.09.29
|
Functions and classes for finite difference methods. More...
Namespaces | |
namespace | reconstruction |
Variable and flux vector splitting reconstruction schemes for finite difference methods. | |
Enumerations | |
enum class | DerivativeOrder : int { OneHigherThanRecons = -1 , OneHigherThanReconsButFiveToFour = -2 , Two = 2 , Four = 4 , Six = 6 , Eight = 8 , Ten = 10 } |
Controls which FD derivative order is used. More... | |
Functions | |
std::ostream & | operator<< (std::ostream &os, DerivativeOrder der_order) |
template<size_t Dim> | |
void | low_pass_filter (gsl::not_null< gsl::span< double > * > filtered_data, const gsl::span< const double > &volume_vars, const DirectionMap< Dim, gsl::span< const double > > &ghost_cell_vars, const Mesh< Dim > &volume_mesh, size_t number_of_variables, size_t fd_order, double epsilon) |
Apply a low-pass filter to the data. More... | |
template<size_t Dim> | |
void | kreiss_oliger_filter (gsl::not_null< gsl::span< double > * > filtered_data, const gsl::span< const double > &volume_vars, const DirectionMap< Dim, gsl::span< const double > > &ghost_cell_vars, const Mesh< Dim > &volume_mesh, size_t number_of_variables, size_t fd_order, double epsilon) |
Apply Kreiss-Oliger dissipation [115] of order fd_order to the variables. More... | |
template<size_t Dim, typename FluxesTags > | |
void | set_cartesian_neighbor_cell_centered_fluxes (const gsl::not_null< DirectionMap< Dim, Variables< FluxesTags > > * > flux_neighbor_data, const DirectionalIdMap< Dim, evolution::dg::subcell::GhostData > &all_ghost_data, const Mesh< Dim > &subcell_mesh, const size_t ghost_zone_size, const size_t number_of_rdmp_values_in_ghost_data) |
Fill the flux_neighbor_data with pointers into the all_ghost_data . More... | |
template<size_t Dim, typename... EvolvedVarsTags, typename FluxesTags = tmpl::list<::Tags::Flux< EvolvedVarsTags, tmpl::size_t<Dim>, Frame::Inertial>...>> | |
void | cartesian_high_order_flux_corrections (const gsl::not_null< std::optional< std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > > * > high_order_corrections, const std::optional< Variables< FluxesTags > > &cell_centered_fluxes, const std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > &second_order_boundary_corrections, const fd::DerivativeOrder &fd_derivative_order, const DirectionalIdMap< Dim, evolution::dg::subcell::GhostData > &all_ghost_data, const Mesh< Dim > &subcell_mesh, const size_t ghost_zone_size, const std::array< gsl::span< std::uint8_t >, Dim > &reconstruction_order={}, const size_t number_of_rdmp_values_in_ghost_data=0) |
Computes the high-order Cartesian flux corrections if necessary. More... | |
template<size_t Dim, typename ReconstructionTags > | |
void | neighbor_data_as_variables (const gsl::not_null< DirectionalIdMap< Dim, Variables< ReconstructionTags > > * > vars_neighbor_data, const DirectionalIdMap< Dim, evolution::dg::subcell::GhostData > &all_ghost_data, const size_t ghost_zone_size, const Mesh< Dim > &subcell_mesh) |
Given the type-erased neighbor data for reconstruction stored in a DataVector , have Variables point into them. More... | |
template<size_t StencilSize> | |
std::array< std::array< double, StencilSize >, StencilSize > | non_uniform_1d_weights (const std::deque< double > ×) |
Returns the weights for a 1D non-uniform finite difference stencil. More... | |
template<size_t Dim> | |
void | logical_partial_derivatives (const gsl::not_null< std::array< gsl::span< double >, Dim > * > logical_derivatives, const gsl::span< const double > &volume_vars, const DirectionMap< Dim, gsl::span< const double > > &ghost_cell_vars, const Mesh< Dim > &volume_mesh, size_t number_of_variables, size_t fd_order) |
Compute the logical partial derivatives using cell-centered finite difference derivatives. More... | |
template<typename DerivativeTags , size_t Dim, typename DerivativeFrame > | |
void | partial_derivatives (gsl::not_null< Variables< db::wrap_tags_in< Tags::deriv, DerivativeTags, tmpl::size_t< Dim >, DerivativeFrame > > * > partial_derivatives, const gsl::span< const double > &volume_vars, const DirectionMap< Dim, gsl::span< const double > > &ghost_cell_vars, const Mesh< Dim > &volume_mesh, size_t number_of_variables, size_t fd_order, const InverseJacobian< DataVector, Dim, Frame::ElementLogical, DerivativeFrame > &inverse_jacobian) |
Compute the partial derivative on the DerivativeFrame using the inverse_jacobian . More... | |
template<DerivativeOrder DerivOrder, size_t Dim, typename... EvolvedVarsTags> | |
void | cartesian_high_order_fluxes_using_nodes (const gsl::not_null< std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > * > high_order_boundary_corrections_in_logical_direction, const std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > &second_order_boundary_corrections_in_logical_direction, const Variables< tmpl::list< ::Tags::Flux< EvolvedVarsTags, tmpl::size_t< Dim >, Frame::Inertial >... > > &cell_centered_inertial_flux, const DirectionMap< Dim, Variables< tmpl::list<::Tags::Flux< EvolvedVarsTags, tmpl::size_t< Dim >, Frame::Inertial >... > > > &ghost_cell_inertial_flux, const Mesh< Dim > &subcell_mesh, const size_t number_of_ghost_cells, const std::array< gsl::span< std::uint8_t >, Dim > &reconstruction_order={}) |
Computes a high-order boundary correction | |
template<size_t Dim, typename... EvolvedVarsTags> | |
void | cartesian_high_order_fluxes_using_nodes (const gsl::not_null< std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > * > high_order_boundary_corrections_in_logical_direction, const std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > &second_order_boundary_corrections_in_logical_direction, const Variables< tmpl::list< ::Tags::Flux< EvolvedVarsTags, tmpl::size_t< Dim >, Frame::Inertial >... > > &cell_centered_inertial_flux, const DirectionMap< Dim, Variables< tmpl::list<::Tags::Flux< EvolvedVarsTags, tmpl::size_t< Dim >, Frame::Inertial >... > > > &ghost_cell_inertial_flux, const Mesh< Dim > &subcell_mesh, const size_t number_of_ghost_cells, const DerivativeOrder derivative_order, const std::array< gsl::span< std::uint8_t >, Dim > &reconstruction_order={}) |
Computes a high-order boundary correction | |
Functions and classes for finite difference methods.
|
strong |
Controls which FD derivative order is used.
void fd::cartesian_high_order_flux_corrections | ( | const gsl::not_null< std::optional< std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > > * > | high_order_corrections, |
const std::optional< Variables< FluxesTags > > & | cell_centered_fluxes, | ||
const std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > & | second_order_boundary_corrections, | ||
const fd::DerivativeOrder & | fd_derivative_order, | ||
const DirectionalIdMap< Dim, evolution::dg::subcell::GhostData > & | all_ghost_data, | ||
const Mesh< Dim > & | subcell_mesh, | ||
const size_t | ghost_zone_size, | ||
const std::array< gsl::span< std::uint8_t >, Dim > & | reconstruction_order = {} , |
||
const size_t | number_of_rdmp_values_in_ghost_data = 0 |
||
) |
Computes the high-order Cartesian flux corrections if necessary.
The cell_centered_fluxes
is stored in the tag evolution::dg::subcell::Tags::CellCenteredFlux
, fd_derivative_order
is from evolution::dg::subcell::Tags::SubcellOptions
(.finite_difference_derivative_order()
), the all_ghost_data
is stored in the tag evolution::dg::subcell::Tags::GhostDataForReconstruction
, the ghost_zone_size
should come from the FD reconstructor.
By default we assume no RDMP data is in the ghost_data
buffer. In the future we will want to update how we store the data in order to eliminate more memory allocations and copies, in which case that value will be non-zero.
high_order_corrections
must either not have a value or have all elements be of the same size as second_order_boundary_corrections[0].number_of_grid_points()
, where we've assumed second_order_boundary_corrections
is the same in all directions. void fd::cartesian_high_order_fluxes_using_nodes | ( | const gsl::not_null< std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > * > | high_order_boundary_corrections_in_logical_direction, |
const std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > & | second_order_boundary_corrections_in_logical_direction, | ||
const Variables< tmpl::list< ::Tags::Flux< EvolvedVarsTags, tmpl::size_t< Dim >, Frame::Inertial >... > > & | cell_centered_inertial_flux, | ||
const DirectionMap< Dim, Variables< tmpl::list<::Tags::Flux< EvolvedVarsTags, tmpl::size_t< Dim >, Frame::Inertial >... > > > & | ghost_cell_inertial_flux, | ||
const Mesh< Dim > & | subcell_mesh, | ||
const size_t | number_of_ghost_cells, | ||
const DerivativeOrder | derivative_order, | ||
const std::array< gsl::span< std::uint8_t >, Dim > & | reconstruction_order = {} |
||
) |
Computes a high-order boundary correction
The correction to the second-order boundary correction is given by
where
where
This is a generalization of the correction presented in [37].
This high-order flux can be fed into a flux limiter, e.g. to guarantee positivity.
G^{(4)}
, which differs by a minus sign. This amounts to a minus sign change in front of the void fd::cartesian_high_order_fluxes_using_nodes | ( | const gsl::not_null< std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > * > | high_order_boundary_corrections_in_logical_direction, |
const std::array< Variables< tmpl::list< EvolvedVarsTags... > >, Dim > & | second_order_boundary_corrections_in_logical_direction, | ||
const Variables< tmpl::list< ::Tags::Flux< EvolvedVarsTags, tmpl::size_t< Dim >, Frame::Inertial >... > > & | cell_centered_inertial_flux, | ||
const DirectionMap< Dim, Variables< tmpl::list<::Tags::Flux< EvolvedVarsTags, tmpl::size_t< Dim >, Frame::Inertial >... > > > & | ghost_cell_inertial_flux, | ||
const Mesh< Dim > & | subcell_mesh, | ||
const size_t | number_of_ghost_cells, | ||
const std::array< gsl::span< std::uint8_t >, Dim > & | reconstruction_order = {} |
||
) |
Computes a high-order boundary correction
The correction to the second-order boundary correction is given by
where
where
This is a generalization of the correction presented in [37].
This high-order flux can be fed into a flux limiter, e.g. to guarantee positivity.
G^{(4)}
, which differs by a minus sign. This amounts to a minus sign change in front of the void fd::kreiss_oliger_filter | ( | gsl::not_null< gsl::span< double > * > | filtered_data, |
const gsl::span< const double > & | volume_vars, | ||
const DirectionMap< Dim, gsl::span< const double > > & | ghost_cell_vars, | ||
const Mesh< Dim > & | volume_mesh, | ||
size_t | number_of_variables, | ||
size_t | fd_order, | ||
double | epsilon | ||
) |
Apply Kreiss-Oliger dissipation [115] of order fd_order
to the variables.
Define the operators
where the subscript
where
We arrive at the following operators:
void fd::logical_partial_derivatives | ( | const gsl::not_null< std::array< gsl::span< double >, Dim > * > | logical_derivatives, |
const gsl::span< const double > & | volume_vars, | ||
const DirectionMap< Dim, gsl::span< const double > > & | ghost_cell_vars, | ||
const Mesh< Dim > & | volume_mesh, | ||
size_t | number_of_variables, | ||
size_t | fd_order | ||
) |
Compute the logical partial derivatives using cell-centered finite difference derivatives.
Up to 8th order stencils are supported.
stride
parameter in the derivative stencils to make testing performance easier in the future.void fd::low_pass_filter | ( | gsl::not_null< gsl::span< double > * > | filtered_data, |
const gsl::span< const double > & | volume_vars, | ||
const DirectionMap< Dim, gsl::span< const double > > & | ghost_cell_vars, | ||
const Mesh< Dim > & | volume_mesh, | ||
size_t | number_of_variables, | ||
size_t | fd_order, | ||
double | epsilon | ||
) |
Apply a low-pass filter to the data.
The filter fits a Legendre polynomial of degree equal to fd_order
to each grid point and its neighboring points, then subtracts out the highest mode contribution at the grid point. This is inspired by a Heaviside filter for Legendre polynomial-based spectral methods.
The filter at different orders is given by:
void fd::neighbor_data_as_variables | ( | const gsl::not_null< DirectionalIdMap< Dim, Variables< ReconstructionTags > > * > | vars_neighbor_data, |
const DirectionalIdMap< Dim, evolution::dg::subcell::GhostData > & | all_ghost_data, | ||
const size_t | ghost_zone_size, | ||
const Mesh< Dim > & | subcell_mesh | ||
) |
Given the type-erased neighbor data for reconstruction stored in a DataVector
, have Variables
point into them.
This function is helpful for reconstruction, especially when wanting to apply different reconstruction methods to different tags. This can happen, for example, when doing positivity-preserving reconstruction. The density should remain positive, but negative velocities are fine.
void fd::partial_derivatives | ( | gsl::not_null< Variables< db::wrap_tags_in< Tags::deriv, DerivativeTags, tmpl::size_t< Dim >, DerivativeFrame > > * > | partial_derivatives, |
const gsl::span< const double > & | volume_vars, | ||
const DirectionMap< Dim, gsl::span< const double > > & | ghost_cell_vars, | ||
const Mesh< Dim > & | volume_mesh, | ||
size_t | number_of_variables, | ||
size_t | fd_order, | ||
const InverseJacobian< DataVector, Dim, Frame::ElementLogical, DerivativeFrame > & | inverse_jacobian | ||
) |
Compute the partial derivative on the DerivativeFrame
using the inverse_jacobian
.
Logical partial derivatives are first computed using the fd::logical_partial_derivatives()
function.
void fd::set_cartesian_neighbor_cell_centered_fluxes | ( | const gsl::not_null< DirectionMap< Dim, Variables< FluxesTags > > * > | flux_neighbor_data, |
const DirectionalIdMap< Dim, evolution::dg::subcell::GhostData > & | all_ghost_data, | ||
const Mesh< Dim > & | subcell_mesh, | ||
const size_t | ghost_zone_size, | ||
const size_t | number_of_rdmp_values_in_ghost_data | ||
) |
Fill the flux_neighbor_data
with pointers into the all_ghost_data
.
The all_ghost_data
is stored in the tag evolution::dg::subcell::Tags::GhostDataForReconstruction
, and the ghost_zone_size
should come from the FD reconstructor.