Skip to content

Commit

Permalink
Merge pull request #5577 from isaaclegred/FixDirichletFreeOutflow
Browse files Browse the repository at this point in the history
Fix boundary condition bugs for temperature reconstruction
  • Loading branch information
nilsdeppe authored Oct 24, 2023
2 parents f17af07 + 0b7f8f4 commit c587a22
Show file tree
Hide file tree
Showing 5 changed files with 98 additions and 75 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -221,6 +221,7 @@ class DirichletFreeOutflow final : public BoundaryCondition {
using fd_interior_primitive_variables_tags =
tmpl::list<hydro::Tags::RestMassDensity<DataVector>,
hydro::Tags::ElectronFraction<DataVector>,
hydro::Tags::Temperature<DataVector>,
hydro::Tags::Pressure<DataVector>,
hydro::Tags::SpecificInternalEnergy<DataVector>,
hydro::Tags::LorentzFactor<DataVector>,
Expand All @@ -241,7 +242,7 @@ class DirichletFreeOutflow final : public BoundaryCondition {
const gsl::not_null<tnsr::iaa<DataVector, 3, Frame::Inertial>*> phi,
const gsl::not_null<Scalar<DataVector>*> rest_mass_density,
const gsl::not_null<Scalar<DataVector>*> electron_fraction,
const gsl::not_null<Scalar<DataVector>*> pressure,
const gsl::not_null<Scalar<DataVector>*> temperature,
const gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
lorentz_factor_times_spatial_velocity,
const gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
Expand All @@ -255,6 +256,7 @@ class DirichletFreeOutflow final : public BoundaryCondition {
// interior prim vars tags
const Scalar<DataVector>& interior_rest_mass_density,
const Scalar<DataVector>& interior_electron_fraction,
const Scalar<DataVector>& interior_temperature,
const Scalar<DataVector>& interior_pressure,
const Scalar<DataVector>& interior_specific_internal_energy,
const Scalar<DataVector>& interior_lorentz_factor,
Expand Down Expand Up @@ -312,6 +314,7 @@ class DirichletFreeOutflow final : public BoundaryCondition {
Flux, typename grmhd::ValenciaDivClean::System::flux_variables>>>
cell_centered_ghost_fluxes{std::nullopt};
// Set to zero since it shouldn't be used
Scalar<DataVector> pressure{};
Scalar<DataVector> specific_internal_energy{};
tnsr::I<DataVector, 3> spatial_velocity{};
Scalar<DataVector> lorentz_factor{};
Expand All @@ -326,8 +329,8 @@ class DirichletFreeOutflow final : public BoundaryCondition {

grmhd::ValenciaDivClean::BoundaryConditions::HydroFreeOutflow::
fd_ghost_impl(
rest_mass_density, electron_fraction, pressure,
make_not_null(&specific_internal_energy),
rest_mass_density, electron_fraction, temperature,
make_not_null(&pressure), make_not_null(&specific_internal_energy),
lorentz_factor_times_spatial_velocity,
make_not_null(&spatial_velocity), make_not_null(&lorentz_factor),
magnetic_field, divergence_cleaning_field,
Expand All @@ -343,9 +346,9 @@ class DirichletFreeOutflow final : public BoundaryCondition {

// fd_interior_primitive_variables_tags
interior_rest_mass_density, interior_electron_fraction,
interior_pressure, interior_specific_internal_energy,
interior_lorentz_factor, interior_spatial_velocity,
interior_magnetic_field,
interior_temperature, interior_pressure,
interior_specific_internal_energy, interior_lorentz_factor,
interior_spatial_velocity, interior_magnetic_field,
// Note: metric vars are empty because they shouldn't be used
interior_spatial_metric, interior_lapse, interior_shift,

Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -171,7 +171,7 @@ std::optional<std::string> HydroFreeOutflow::dg_ghost(
void HydroFreeOutflow::fd_ghost(
const gsl::not_null<Scalar<DataVector>*> rest_mass_density,
const gsl::not_null<Scalar<DataVector>*> electron_fraction,
const gsl::not_null<Scalar<DataVector>*> pressure,
const gsl::not_null<Scalar<DataVector>*> temperature,
const gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
lorentz_factor_times_spatial_velocity,
const gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
Expand All @@ -193,6 +193,7 @@ void HydroFreeOutflow::fd_ghost(
// fd_interior_primitive_variables_tags
const Scalar<DataVector>& interior_rest_mass_density,
const Scalar<DataVector>& interior_electron_fraction,
const Scalar<DataVector>& interior_temperature,
const Scalar<DataVector>& interior_pressure,
const Scalar<DataVector>& interior_specific_internal_energy,
const Scalar<DataVector>& interior_lorentz_factor,
Expand All @@ -202,36 +203,38 @@ void HydroFreeOutflow::fd_ghost(
// fd_gridless_tags
const fd::Reconstructor& reconstructor) {
Variables<tmpl::push_back<typename System::variables_tag::tags_list,
SpatialVelocity, LorentzFactor,
SpatialVelocity, LorentzFactor, Pressure,
SpecificInternalEnergy, SqrtDetSpatialMetric,
SpatialMetric, InvSpatialMetric, Lapse, Shift>>
temp_vars{get(*rest_mass_density).size()};
fd_ghost_impl(
rest_mass_density, electron_fraction, pressure,
make_not_null(&get<SpecificInternalEnergy>(temp_vars)),
lorentz_factor_times_spatial_velocity,
make_not_null(&get<SpatialVelocity>(temp_vars)),
make_not_null(&get<LorentzFactor>(temp_vars)), magnetic_field,
divergence_cleaning_field, make_not_null(&get<SpatialMetric>(temp_vars)),
make_not_null(&get<InvSpatialMetric>(temp_vars)),
make_not_null(&get<SqrtDetSpatialMetric>(temp_vars)),
make_not_null(&get<Lapse>(temp_vars)),
make_not_null(&get<Shift>(temp_vars)),

direction,

// fd_interior_temporary_tags
subcell_mesh,

// fd_interior_primitive_variables_tags
interior_rest_mass_density, interior_electron_fraction, interior_pressure,
interior_specific_internal_energy, interior_lorentz_factor,
interior_spatial_velocity, interior_magnetic_field,
interior_spatial_metric, interior_lapse, interior_shift,

reconstructor.ghost_zone_size(),

cell_centered_ghost_fluxes->has_value());
fd_ghost_impl(rest_mass_density, electron_fraction, temperature,
make_not_null(&get<Pressure>(temp_vars)),
make_not_null(&get<SpecificInternalEnergy>(temp_vars)),
lorentz_factor_times_spatial_velocity,
make_not_null(&get<SpatialVelocity>(temp_vars)),
make_not_null(&get<LorentzFactor>(temp_vars)), magnetic_field,
divergence_cleaning_field,
make_not_null(&get<SpatialMetric>(temp_vars)),
make_not_null(&get<InvSpatialMetric>(temp_vars)),
make_not_null(&get<SqrtDetSpatialMetric>(temp_vars)),
make_not_null(&get<Lapse>(temp_vars)),
make_not_null(&get<Shift>(temp_vars)),

direction,

// fd_interior_temporary_tags
subcell_mesh,

// fd_interior_primitive_variables_tags
interior_rest_mass_density, interior_electron_fraction,
interior_temperature, interior_pressure,
interior_specific_internal_energy, interior_lorentz_factor,
interior_spatial_velocity, interior_magnetic_field,
interior_spatial_metric, interior_lapse, interior_shift,

reconstructor.ghost_zone_size(),

cell_centered_ghost_fluxes->has_value());
if (cell_centered_ghost_fluxes->has_value()) {
ConservativeFromPrimitive::apply(
make_not_null(&get<Tags::TildeD>(temp_vars)),
Expand All @@ -243,7 +246,7 @@ void HydroFreeOutflow::fd_ghost(

// Note: Only the spatial velocity changes.
*rest_mass_density, *electron_fraction,
get<SpecificInternalEnergy>(temp_vars), *pressure,
get<SpecificInternalEnergy>(temp_vars), get<Pressure>(temp_vars),
get<SpatialVelocity>(temp_vars), get<LorentzFactor>(temp_vars),
*magnetic_field,

Expand All @@ -270,7 +273,7 @@ void HydroFreeOutflow::fd_ghost(

get<Lapse>(temp_vars), get<Shift>(temp_vars),
get<SqrtDetSpatialMetric>(temp_vars), get<SpatialMetric>(temp_vars),
get<InvSpatialMetric>(temp_vars), *pressure,
get<InvSpatialMetric>(temp_vars), get<Pressure>(temp_vars),
get<SpatialVelocity>(temp_vars), get<LorentzFactor>(temp_vars),
*magnetic_field);
}
Expand All @@ -279,6 +282,7 @@ void HydroFreeOutflow::fd_ghost(
void HydroFreeOutflow::fd_ghost_impl(
const gsl::not_null<Scalar<DataVector>*> rest_mass_density,
const gsl::not_null<Scalar<DataVector>*> electron_fraction,
const gsl::not_null<Scalar<DataVector>*> temperature,
const gsl::not_null<Scalar<DataVector>*> pressure,
const gsl::not_null<Scalar<DataVector>*> specific_internal_energy,
const gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
Expand All @@ -305,6 +309,7 @@ void HydroFreeOutflow::fd_ghost_impl(
// fd_interior_primitive_variables_tags
const Scalar<DataVector>& interior_rest_mass_density,
const Scalar<DataVector>& interior_electron_fraction,
const Scalar<DataVector>& interior_temperature,
const Scalar<DataVector>& interior_pressure,
const Scalar<DataVector>& interior_specific_internal_energy,
const Scalar<DataVector>& interior_lorentz_factor,
Expand All @@ -322,12 +327,12 @@ void HydroFreeOutflow::fd_ghost_impl(
subcell_extents.slice_away(dim_direction).product()};

using prim_tags_without_cleaning_field =
tmpl::list<RestMassDensity, ElectronFraction, Pressure,
tmpl::list<RestMassDensity, ElectronFraction, Temperature,
LorentzFactorTimesSpatialVelocity, MagneticField>;

// Create a single large DV to reduce the number of Variables allocations
using fluxes_tags =
tmpl::list<SpecificInternalEnergy, SpatialMetric, Lapse, Shift>;
tmpl::list<Pressure, SpecificInternalEnergy, SpatialMetric, Lapse, Shift>;
const size_t buffer_size_for_fluxes =
need_tags_for_fluxes
? Variables<fluxes_tags>::number_of_independent_components
Expand Down Expand Up @@ -357,7 +362,8 @@ void HydroFreeOutflow::fd_ghost_impl(
get_boundary_val(interior_rest_mass_density);
get<ElectronFraction>(outermost_prim_vars) =
get_boundary_val(interior_electron_fraction);
get<Pressure>(outermost_prim_vars) = get_boundary_val(interior_pressure);
get<Temperature>(outermost_prim_vars) =
get_boundary_val(interior_temperature);

{
// Kill ingoing components of spatial velocity and compute Wv^i
Expand Down Expand Up @@ -401,7 +407,7 @@ void HydroFreeOutflow::fd_ghost_impl(

*rest_mass_density = get<RestMassDensity>(ghost_prim_vars);
*electron_fraction = get<ElectronFraction>(ghost_prim_vars);
*pressure = get<Pressure>(ghost_prim_vars);
*temperature = get<Temperature>(ghost_prim_vars);
*lorentz_factor_times_spatial_velocity =
get<LorentzFactorTimesSpatialVelocity>(ghost_prim_vars);
*magnetic_field = get<MagneticField>(ghost_prim_vars);
Expand All @@ -419,6 +425,7 @@ void HydroFreeOutflow::fd_ghost_impl(
static_cast<std::ptrdiff_t>(outermost_fluxes_vars.size())),
num_face_pts * buffer_size_for_fluxes * ghost_zone_size};

get<Pressure>(outermost_fluxes_vars) = get_boundary_val(interior_pressure);
get<SpecificInternalEnergy>(outermost_fluxes_vars) =
get_boundary_val(interior_specific_internal_energy);
get<SpatialMetric>(outermost_fluxes_vars) =
Expand All @@ -431,7 +438,8 @@ void HydroFreeOutflow::fd_ghost_impl(
outermost_fluxes_vars, ghost_data_extents,
dim_direction, i_ghost);
}

// Need pressure for high-order finite difference
*pressure = get<Pressure>(ghost_fluxes_vars);
*specific_internal_energy = get<SpecificInternalEnergy>(ghost_fluxes_vars);
*spatial_metric = get<SpatialMetric>(ghost_fluxes_vars);
*lapse = get<Lapse>(ghost_fluxes_vars);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -20,6 +20,7 @@
#include "Options/String.hpp"
#include "PointwiseFunctions/GeneralRelativity/Tags.hpp"
#include "PointwiseFunctions/Hydro/Tags.hpp"
#include "PointwiseFunctions/Hydro/Temperature.hpp"
#include "Utilities/TMPL.hpp"

/// \cond
Expand Down Expand Up @@ -62,6 +63,7 @@ class HydroFreeOutflow final : public BoundaryCondition {
private:
using RestMassDensity = hydro::Tags::RestMassDensity<DataVector>;
using ElectronFraction = hydro::Tags::ElectronFraction<DataVector>;
using Temperature = hydro::Tags::Temperature<DataVector>;
using Pressure = hydro::Tags::Pressure<DataVector>;
using LorentzFactorTimesSpatialVelocity =
hydro::Tags::LorentzFactorTimesSpatialVelocity<DataVector, 3>;
Expand All @@ -79,7 +81,7 @@ class HydroFreeOutflow final : public BoundaryCondition {
using Shift = gr::Tags::Shift<DataVector, 3>;

using prim_tags_for_reconstruction =
tmpl::list<RestMassDensity, ElectronFraction, Pressure,
tmpl::list<RestMassDensity, ElectronFraction, Temperature,
LorentzFactorTimesSpatialVelocity, MagneticField,
DivergenceCleaningField>;

Expand Down Expand Up @@ -175,16 +177,18 @@ class HydroFreeOutflow final : public BoundaryCondition {
tmpl::list<evolution::dg::subcell::Tags::Mesh<3>, Shift, Lapse,
SpatialMetric>;
using fd_interior_primitive_variables_tags =
tmpl::list<RestMassDensity, ElectronFraction, Pressure,
tmpl::list<RestMassDensity, ElectronFraction, Temperature,
hydro::Tags::Pressure<DataVector>,
hydro::Tags::SpecificInternalEnergy<DataVector>,
hydro::Tags::LorentzFactor<DataVector>,
hydro::Tags::SpatialVelocity<DataVector, 3>, MagneticField>;

using fd_gridless_tags = tmpl::list<fd::Tags::Reconstructor>;

static void fd_ghost(
gsl::not_null<Scalar<DataVector>*> rest_mass_density,
gsl::not_null<Scalar<DataVector>*> electron_fraction,
gsl::not_null<Scalar<DataVector>*> pressure,
gsl::not_null<Scalar<DataVector>*> temperature,
gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
lorentz_factor_times_spatial_velocity,
gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*> magnetic_field,
Expand All @@ -205,6 +209,7 @@ class HydroFreeOutflow final : public BoundaryCondition {
// interior prim vars tags
const Scalar<DataVector>& interior_rest_mass_density,
const Scalar<DataVector>& interior_electron_fraction,
const Scalar<DataVector>& interior_temperature,
const Scalar<DataVector>& interior_pressure,
const Scalar<DataVector>& interior_specific_internal_energy,
const Scalar<DataVector>& interior_lorentz_factor,
Expand All @@ -218,6 +223,7 @@ class HydroFreeOutflow final : public BoundaryCondition {
static void fd_ghost_impl(
gsl::not_null<Scalar<DataVector>*> rest_mass_density,
gsl::not_null<Scalar<DataVector>*> electron_fraction,
gsl::not_null<Scalar<DataVector>*> temperature,
gsl::not_null<Scalar<DataVector>*> pressure,
gsl::not_null<Scalar<DataVector>*> specific_internal_energy,
gsl::not_null<tnsr::I<DataVector, 3, Frame::Inertial>*>
Expand All @@ -241,6 +247,7 @@ class HydroFreeOutflow final : public BoundaryCondition {
// fd_interior_primitive_variables_tags
const Scalar<DataVector>& interior_rest_mass_density,
const Scalar<DataVector>& interior_electron_fraction,
const Scalar<DataVector>& interior_temperature,
const Scalar<DataVector>& interior_pressure,
const Scalar<DataVector>& interior_specific_internal_energy,
const Scalar<DataVector>& interior_lorentz_factor,
Expand Down
Loading

0 comments on commit c587a22

Please sign in to comment.