|
| 1 | +// Distributed under the MIT License. |
| 2 | +// See LICENSE.txt for details. |
| 3 | + |
| 4 | +#pragma once |
| 5 | + |
| 6 | +#include <array> |
| 7 | +#include <cstddef> |
| 8 | +#include <memory> |
| 9 | +#include <utility> |
| 10 | +#include <vector> |
| 11 | + |
| 12 | +#include "DataStructures/DataBox/PrefixHelpers.hpp" |
| 13 | +#include "DataStructures/DataBox/Prefixes.hpp" |
| 14 | +#include "DataStructures/Tensor/TypeAliases.hpp" |
| 15 | +#include "DataStructures/Variables.hpp" |
| 16 | +#include "DataStructures/VariablesTag.hpp" |
| 17 | +#include "Domain/Structure/DirectionalIdMap.hpp" |
| 18 | +#include "Domain/Tags.hpp" |
| 19 | +#include "Evolution/DgSubcell/Tags/GhostDataForReconstruction.hpp" |
| 20 | +#include "Evolution/DgSubcell/Tags/Mesh.hpp" |
| 21 | +#include "Evolution/Systems/GeneralizedHarmonic/Tags.hpp" |
| 22 | +#include "Evolution/Systems/GrMhd/GhValenciaDivClean/FiniteDifference/Reconstructor.hpp" |
| 23 | +#include "Evolution/Systems/GrMhd/GhValenciaDivClean/System.hpp" |
| 24 | +#include "Evolution/Systems/GrMhd/GhValenciaDivClean/Tags.hpp" |
| 25 | +#include "Evolution/Systems/GrMhd/ValenciaDivClean/Tags.hpp" |
| 26 | +#include "Evolution/VariableFixing/FixToAtmosphere.hpp" |
| 27 | +#include "Evolution/VariableFixing/Tags.hpp" |
| 28 | +#include "Options/String.hpp" |
| 29 | +#include "PointwiseFunctions/GeneralRelativity/Tags.hpp" |
| 30 | +#include "PointwiseFunctions/Hydro/Tags.hpp" |
| 31 | +#include "Utilities/Serialization/CharmPupable.hpp" |
| 32 | +#include "Utilities/TMPL.hpp" |
| 33 | + |
| 34 | +/// \cond |
| 35 | +class DataVector; |
| 36 | +template <size_t Dim> |
| 37 | +class Direction; |
| 38 | +template <size_t Dim> |
| 39 | +class Element; |
| 40 | +template <size_t Dim> |
| 41 | +class ElementId; |
| 42 | +namespace EquationsOfState { |
| 43 | +template <bool IsRelativistic, size_t ThermodynamicDim> |
| 44 | +class EquationOfState; |
| 45 | +} // namespace EquationsOfState |
| 46 | +template <size_t Dim> |
| 47 | +class Mesh; |
| 48 | +namespace gsl { |
| 49 | +template <typename T> |
| 50 | +class not_null; |
| 51 | +} // namespace gsl |
| 52 | +namespace PUP { |
| 53 | +class er; |
| 54 | +} // namespace PUP |
| 55 | +namespace evolution::dg::subcell { |
| 56 | +class GhostData; |
| 57 | +} // namespace evolution::dg::subcell |
| 58 | +/// \endcond |
| 59 | + |
| 60 | +namespace grmhd::GhValenciaDivClean::fd { |
| 61 | +/*! |
| 62 | + * \brief PPM (Piecewise Parabolic Method) reconstruction on the GRMHD |
| 63 | + * primitive variables (see ::fd::reconstruction::ppm() for details) and |
| 64 | + * unlimited 3rd order (degree 2 polynomial) reconstruction on the metric |
| 65 | + * variables. |
| 66 | + * |
| 67 | + * Only the spacetime metric is reconstructed when we and the neighboring |
| 68 | + * element in the direction are doing FD. If we are doing DG and a neighboring |
| 69 | + * element is doing FD, then the spacetime metric, \f$\Phi_{iab}\f$, and |
| 70 | + * \f$\Pi_{ab}\f$ are all reconstructed since the Riemann solver on the DG |
| 71 | + * element also needs to solve for the metric variables. |
| 72 | + */ |
| 73 | +template <typename System> |
| 74 | +class PpmPrim : public Reconstructor<System> { |
| 75 | + public: |
| 76 | + static constexpr size_t dim = 3; |
| 77 | + |
| 78 | + struct AtmosphereTreatment { |
| 79 | + using type = ::VariableFixing::FixReconstructedStateToAtmosphere; |
| 80 | + static constexpr Options::String help = { |
| 81 | + "What reconstructed states to fix to their atmosphere values."}; |
| 82 | + }; |
| 83 | + struct ReconstructRhoTimesTemperature { |
| 84 | + using type = bool; |
| 85 | + static constexpr Options::String help = { |
| 86 | + "If 'true' then we reconstruct the rho*T, if 'false' we reconstruct " |
| 87 | + "T."}; |
| 88 | + }; |
| 89 | + |
| 90 | + using options = |
| 91 | + tmpl::list<AtmosphereTreatment, ReconstructRhoTimesTemperature>; |
| 92 | + static constexpr Options::String help{ |
| 93 | + "PPM (Piecewise Parabolic Method) reconstruction scheme using primitive " |
| 94 | + "variables and the metric variables."}; |
| 95 | + |
| 96 | + PpmPrim(::VariableFixing::FixReconstructedStateToAtmosphere |
| 97 | + fix_reconstructed_state_to_atmosphere, |
| 98 | + bool reconstruct_rho_times_temperature); |
| 99 | + |
| 100 | + PpmPrim() = default; |
| 101 | + PpmPrim(PpmPrim&&) = default; |
| 102 | + PpmPrim& operator=(PpmPrim&&) = default; |
| 103 | + PpmPrim(const PpmPrim&) = default; |
| 104 | + PpmPrim& operator=(const PpmPrim&) = default; |
| 105 | + ~PpmPrim() override = default; |
| 106 | + |
| 107 | + explicit PpmPrim(CkMigrateMessage* msg); |
| 108 | + |
| 109 | + WRAPPED_PUPable_decl_base_template(Reconstructor<System>, PpmPrim); |
| 110 | + |
| 111 | + auto get_clone() const -> std::unique_ptr<Reconstructor<System>> override; |
| 112 | + |
| 113 | + static constexpr bool use_adaptive_order = false; |
| 114 | + |
| 115 | + void pup(PUP::er& p) override; |
| 116 | + |
| 117 | + size_t ghost_zone_size() const override { return 2; } |
| 118 | + |
| 119 | + using reconstruction_argument_tags = |
| 120 | + tmpl::list<::Tags::Variables<hydro::grmhd_tags<DataVector>>, |
| 121 | + typename System::variables_tag, |
| 122 | + hydro::Tags::GrmhdEquationOfState, domain::Tags::Element<dim>, |
| 123 | + evolution::dg::subcell::Tags::GhostDataForReconstruction<dim>, |
| 124 | + evolution::dg::subcell::Tags::Mesh<dim>, |
| 125 | + ::Tags::VariableFixer<VariableFixing::FixToAtmosphere<dim>>>; |
| 126 | + |
| 127 | + template <size_t ThermodynamicDim, typename TagsList> |
| 128 | + void reconstruct( |
| 129 | + gsl::not_null<std::array<Variables<TagsList>, dim>*> vars_on_lower_face, |
| 130 | + gsl::not_null<std::array<Variables<TagsList>, dim>*> vars_on_upper_face, |
| 131 | + const Variables<hydro::grmhd_tags<DataVector>>& volume_prims, |
| 132 | + const Variables<typename System::variables_tag::type::tags_list>& |
| 133 | + volume_spacetime_and_cons_vars, |
| 134 | + const EquationsOfState::EquationOfState<true, ThermodynamicDim>& eos, |
| 135 | + const Element<dim>& element, |
| 136 | + const DirectionalIdMap<dim, evolution::dg::subcell::GhostData>& |
| 137 | + ghost_data, |
| 138 | + const Mesh<dim>& subcell_mesh, |
| 139 | + const VariableFixing::FixToAtmosphere<dim>& fix_to_atmosphere) const; |
| 140 | + |
| 141 | + /// Called by an element doing DG when the neighbor is doing subcell. |
| 142 | + template <size_t ThermodynamicDim, typename TagsList> |
| 143 | + void reconstruct_fd_neighbor( |
| 144 | + gsl::not_null<Variables<TagsList>*> vars_on_face, |
| 145 | + const Variables<hydro::grmhd_tags<DataVector>>& subcell_volume_prims, |
| 146 | + const Variables< |
| 147 | + grmhd::GhValenciaDivClean::Tags::spacetime_reconstruction_tags>& |
| 148 | + subcell_volume_spacetime_metric, |
| 149 | + const EquationsOfState::EquationOfState<true, ThermodynamicDim>& eos, |
| 150 | + const Element<dim>& element, |
| 151 | + const DirectionalIdMap<dim, evolution::dg::subcell::GhostData>& |
| 152 | + ghost_data, |
| 153 | + const Mesh<dim>& subcell_mesh, |
| 154 | + const VariableFixing::FixToAtmosphere<dim>& fix_to_atmosphere, |
| 155 | + Direction<dim> direction_to_reconstruct) const; |
| 156 | + |
| 157 | + bool reconstruct_rho_times_temperature() const override; |
| 158 | + |
| 159 | + private: |
| 160 | + template <typename LocalSystem> |
| 161 | + friend bool operator==(const PpmPrim<LocalSystem>& lhs, |
| 162 | + const PpmPrim<LocalSystem>& rhs); |
| 163 | + template <typename LocalSystem> |
| 164 | + friend bool operator!=(const PpmPrim<LocalSystem>& lhs, |
| 165 | + const PpmPrim<LocalSystem>& rhs); |
| 166 | + |
| 167 | + ::VariableFixing::FixReconstructedStateToAtmosphere |
| 168 | + fix_reconstructed_state_to_atmosphere_{ |
| 169 | + ::VariableFixing::FixReconstructedStateToAtmosphere::Never}; |
| 170 | + bool reconstruct_rho_times_temperature_{false}; |
| 171 | +}; |
| 172 | +} // namespace grmhd::GhValenciaDivClean::fd |
0 commit comments