Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Move NeighborPackagedData into cpp files & allow volume tags in dg boundary terms #5888

Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
2 changes: 2 additions & 0 deletions src/Evolution/DgSubcell/CMakeLists.txt
Original file line number Diff line number Diff line change
Expand Up @@ -25,6 +25,7 @@ spectre_target_headers(
Mesh.hpp
NeighborRdmpAndVolumeData.hpp
NeighborReconstructedFaceSolution.hpp
NeighborReconstructedFaceSolution.tpp
NeighborTciDecision.hpp
PerssonTci.hpp
PrepareNeighborData.hpp
Expand Down Expand Up @@ -53,6 +54,7 @@ spectre_target_sources(
Matrices.cpp
Mesh.cpp
NeighborRdmpAndVolumeData.cpp
NeighborTciDecision.cpp
PerssonTci.cpp
Projection.cpp
RdmpTci.cpp
Expand Down
109 changes: 5 additions & 104 deletions src/Evolution/DgSubcell/NeighborReconstructedFaceSolution.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -4,26 +4,15 @@
#pragma once

#include <cstddef>
#include <iterator>
#include <optional>
#include <tuple>
#include <utility>
#include <vector>

#include "DataStructures/DataBox/DataBox.hpp"
#include "DataStructures/DataBox/Access.hpp"
#include "DataStructures/DataVector.hpp"
#include "Domain/Structure/Direction.hpp"
#include "Domain/Structure/DirectionalId.hpp"
#include "Domain/Structure/DirectionalIdMap.hpp"
#include "Domain/Structure/ElementId.hpp"
#include "Evolution/DgSubcell/GhostData.hpp"
#include "Evolution/DgSubcell/RdmpTciData.hpp"
#include "Evolution/DgSubcell/Tags/DataForRdmpTci.hpp"
#include "Evolution/DgSubcell/Tags/GhostDataForReconstruction.hpp"
#include "Evolution/DgSubcell/Tags/Mesh.hpp"
#include "NumericalAlgorithms/Spectral/Mesh.hpp"
#include "Time/TimeStepId.hpp"
#include "Utilities/ErrorHandling/Assert.hpp"
#include "Utilities/Gsl.hpp"

namespace evolution::dg::subcell {
Expand Down Expand Up @@ -58,103 +47,15 @@ namespace evolution::dg::subcell {
* the last argument to
* `Metavariables::SubcellOptions::DgComputeSubcellNeighborPackagedData::apply`.
*/
template <size_t VolumeDim, typename DgComputeSubcellNeighborPackagedData,
typename DbTagsList>
template <size_t VolumeDim, typename DgComputeSubcellNeighborPackagedData>
void neighbor_reconstructed_face_solution(
const gsl::not_null<db::DataBox<DbTagsList>*> box,
const gsl::not_null<std::pair<
gsl::not_null<db::Access*> box,
gsl::not_null<std::pair<
const TimeStepId,
DirectionalIdMap<
VolumeDim,
std::tuple<Mesh<VolumeDim>, Mesh<VolumeDim - 1>,
std::optional<DataVector>, std::optional<DataVector>,
::TimeStepId, int>>>*>
received_temporal_id_and_data) {
db::mutate<subcell::Tags::GhostDataForReconstruction<VolumeDim>,
subcell::Tags::DataForRdmpTci>(
[&received_temporal_id_and_data](const auto subcell_ghost_data_ptr,
const auto rdmp_tci_data_ptr) {
const size_t number_of_evolved_vars =
rdmp_tci_data_ptr->max_variables_values.size();
for (auto& received_mortar_data :
received_temporal_id_and_data->second) {
const auto& mortar_id = received_mortar_data.first;
ASSERT(std::get<2>(received_mortar_data.second).has_value(),
"The subcell mortar data was not sent at TimeStepId "
<< received_temporal_id_and_data->first
<< " with mortar id " << mortar_id);
const DataVector& neighbor_ghost_and_subcell_data =
*std::get<2>(received_mortar_data.second);
// Compute min and max over neighbors
const size_t offset_for_min =
neighbor_ghost_and_subcell_data.size() - number_of_evolved_vars;
const size_t offset_for_max = offset_for_min - number_of_evolved_vars;
for (size_t var_index = 0; var_index < number_of_evolved_vars;
++var_index) {
rdmp_tci_data_ptr->max_variables_values[var_index] = std::max(
rdmp_tci_data_ptr->max_variables_values[var_index],
neighbor_ghost_and_subcell_data[offset_for_max + var_index]);
rdmp_tci_data_ptr->min_variables_values[var_index] = std::min(
rdmp_tci_data_ptr->min_variables_values[var_index],
neighbor_ghost_and_subcell_data[offset_for_min + var_index]);
}

ASSERT(subcell_ghost_data_ptr->find(mortar_id) ==
subcell_ghost_data_ptr->end(),
"The subcell neighbor data is already inserted. Direction: "
<< mortar_id.direction
<< " with ElementId: " << mortar_id.id);

(*subcell_ghost_data_ptr)[mortar_id] = GhostData{1};
GhostData& all_ghost_data = subcell_ghost_data_ptr->at(mortar_id);
DataVector& neighbor_data =
all_ghost_data.neighbor_ghost_data_for_reconstruction();
neighbor_data.destructive_resize(
neighbor_ghost_and_subcell_data.size() -
2 * number_of_evolved_vars);

// Copy over the neighbor data for reconstruction. We need this
// since we might be doing a step unwind and the DG algorithm deletes
// the inbox data after lifting the fluxes to the volume.
// The std::prev avoids copying over the data for the RDMP TCI, which
// is both the maximum and minimum of each evolved variable, so
// `2*number_of_evolved_vars` components.
std::copy(neighbor_ghost_and_subcell_data.begin(),
std::prev(neighbor_ghost_and_subcell_data.end(),
2 * static_cast<std::ptrdiff_t>(
number_of_evolved_vars)),
neighbor_data.begin());
}
},
box);
std::vector<DirectionalId<VolumeDim>> mortars_to_reconstruct_to{};
for (auto& received_mortar_data : received_temporal_id_and_data->second) {
const auto& mortar_id = received_mortar_data.first;
if (not std::get<3>(received_mortar_data.second).has_value()) {
mortars_to_reconstruct_to.push_back(mortar_id);
}
}
DirectionalIdMap<VolumeDim, DataVector> neighbor_reconstructed_evolved_vars =
DgComputeSubcellNeighborPackagedData::apply(*box,
mortars_to_reconstruct_to);
ASSERT(neighbor_reconstructed_evolved_vars.size() ==
mortars_to_reconstruct_to.size(),
"Should have reconstructed "
<< mortars_to_reconstruct_to.size() << " sides but reconstructed "
<< neighbor_reconstructed_evolved_vars.size() << " sides.");
// Now copy over the packaged data _into_ the inbox in order to avoid having
// to make other changes to the DG algorithm (code in
// src/Evolution/DiscontinuousGalerkin)
for (auto& received_mortar_data : received_temporal_id_and_data->second) {
const auto& mortar_id = received_mortar_data.first;
if (not std::get<3>(received_mortar_data.second).has_value()) {
ASSERT(neighbor_reconstructed_evolved_vars.find(mortar_id) !=
neighbor_reconstructed_evolved_vars.end(),
"Could not find mortar id " << mortar_id
<< " in reconstructed data map.");
std::get<3>(received_mortar_data.second) =
std::move(neighbor_reconstructed_evolved_vars.at(mortar_id));
}
}
}
received_temporal_id_and_data);
} // namespace evolution::dg::subcell
128 changes: 128 additions & 0 deletions src/Evolution/DgSubcell/NeighborReconstructedFaceSolution.tpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,128 @@
// Distributed under the MIT License.
// See LICENSE.txt for details.

#pragma once

#include <cstddef>
#include <iterator>
#include <optional>
#include <tuple>
#include <utility>
#include <vector>

#include "DataStructures/DataBox/Access.hpp"
#include "DataStructures/DataVector.hpp"
#include "Domain/Structure/Direction.hpp"
#include "Domain/Structure/DirectionalId.hpp"
#include "Domain/Structure/DirectionalIdMap.hpp"
#include "Domain/Structure/ElementId.hpp"
#include "Evolution/DgSubcell/GhostData.hpp"
#include "Evolution/DgSubcell/RdmpTciData.hpp"
#include "Evolution/DgSubcell/Tags/DataForRdmpTci.hpp"
#include "Evolution/DgSubcell/Tags/GhostDataForReconstruction.hpp"
#include "Evolution/DgSubcell/Tags/Mesh.hpp"
#include "NumericalAlgorithms/Spectral/Mesh.hpp"
#include "Time/TimeStepId.hpp"
#include "Utilities/ErrorHandling/Assert.hpp"
#include "Utilities/Gsl.hpp"

namespace evolution::dg::subcell {
template <size_t VolumeDim, typename DgComputeSubcellNeighborPackagedData>
void neighbor_reconstructed_face_solution(
const gsl::not_null<db::Access*> box,
const gsl::not_null<std::pair<
const TimeStepId,
DirectionalIdMap<
VolumeDim,
std::tuple<Mesh<VolumeDim>, Mesh<VolumeDim - 1>,
std::optional<DataVector>, std::optional<DataVector>,
::TimeStepId, int>>>*>
received_temporal_id_and_data) {
db::mutate<subcell::Tags::GhostDataForReconstruction<VolumeDim>,
subcell::Tags::DataForRdmpTci>(
[&received_temporal_id_and_data](const auto subcell_ghost_data_ptr,
const auto rdmp_tci_data_ptr) {
const size_t number_of_evolved_vars =
rdmp_tci_data_ptr->max_variables_values.size();
for (auto& received_mortar_data :
received_temporal_id_and_data->second) {
const auto& mortar_id = received_mortar_data.first;
ASSERT(std::get<2>(received_mortar_data.second).has_value(),
"The subcell mortar data was not sent at TimeStepId "
<< received_temporal_id_and_data->first
<< " with mortar id " << mortar_id);
const DataVector& neighbor_ghost_and_subcell_data =
*std::get<2>(received_mortar_data.second);
// Compute min and max over neighbors
const size_t offset_for_min =
neighbor_ghost_and_subcell_data.size() - number_of_evolved_vars;
const size_t offset_for_max = offset_for_min - number_of_evolved_vars;
for (size_t var_index = 0; var_index < number_of_evolved_vars;
++var_index) {
rdmp_tci_data_ptr->max_variables_values[var_index] = std::max(
rdmp_tci_data_ptr->max_variables_values[var_index],
neighbor_ghost_and_subcell_data[offset_for_max + var_index]);
rdmp_tci_data_ptr->min_variables_values[var_index] = std::min(
rdmp_tci_data_ptr->min_variables_values[var_index],
neighbor_ghost_and_subcell_data[offset_for_min + var_index]);
}

ASSERT(subcell_ghost_data_ptr->find(mortar_id) ==
subcell_ghost_data_ptr->end(),
"The subcell neighbor data is already inserted. Direction: "
<< mortar_id.direction
<< " with ElementId: " << mortar_id.id);

(*subcell_ghost_data_ptr)[mortar_id] = GhostData{1};
GhostData& all_ghost_data = subcell_ghost_data_ptr->at(mortar_id);
DataVector& neighbor_data =
all_ghost_data.neighbor_ghost_data_for_reconstruction();
neighbor_data.destructive_resize(
neighbor_ghost_and_subcell_data.size() -
2 * number_of_evolved_vars);

// Copy over the neighbor data for reconstruction. We need this
// since we might be doing a step unwind and the DG algorithm deletes
// the inbox data after lifting the fluxes to the volume.
// The std::prev avoids copying over the data for the RDMP TCI, which
// is both the maximum and minimum of each evolved variable, so
// `2*number_of_evolved_vars` components.
std::copy(neighbor_ghost_and_subcell_data.begin(),
std::prev(neighbor_ghost_and_subcell_data.end(),
2 * static_cast<std::ptrdiff_t>(
number_of_evolved_vars)),
neighbor_data.begin());
}
},
box);
std::vector<DirectionalId<VolumeDim>> mortars_to_reconstruct_to{};
for (auto& received_mortar_data : received_temporal_id_and_data->second) {
const auto& mortar_id = received_mortar_data.first;
if (not std::get<3>(received_mortar_data.second).has_value()) {
mortars_to_reconstruct_to.push_back(mortar_id);
}
}
DirectionalIdMap<VolumeDim, DataVector> neighbor_reconstructed_evolved_vars =
DgComputeSubcellNeighborPackagedData::apply(*box,
mortars_to_reconstruct_to);
ASSERT(neighbor_reconstructed_evolved_vars.size() ==
mortars_to_reconstruct_to.size(),
"Should have reconstructed "
<< mortars_to_reconstruct_to.size() << " sides but reconstructed "
<< neighbor_reconstructed_evolved_vars.size() << " sides.");
// Now copy over the packaged data _into_ the inbox in order to avoid having
// to make other changes to the DG algorithm (code in
// src/Evolution/DiscontinuousGalerkin)
for (auto& received_mortar_data : received_temporal_id_and_data->second) {
const auto& mortar_id = received_mortar_data.first;
if (not std::get<3>(received_mortar_data.second).has_value()) {
ASSERT(neighbor_reconstructed_evolved_vars.find(mortar_id) !=
neighbor_reconstructed_evolved_vars.end(),
"Could not find mortar id " << mortar_id
<< " in reconstructed data map.");
std::get<3>(received_mortar_data.second) =
std::move(neighbor_reconstructed_evolved_vars.at(mortar_id));
}
}
}
} // namespace evolution::dg::subcell
65 changes: 65 additions & 0 deletions src/Evolution/DgSubcell/NeighborTciDecision.cpp
Original file line number Diff line number Diff line change
@@ -0,0 +1,65 @@
// Distributed under the MIT License.
// See LICENSE.txt for details.

#include "Evolution/DgSubcell/NeighborTciDecision.hpp"

#include <cstddef>
#include <optional>
#include <tuple>
#include <utility>

#include "DataStructures/DataBox/Access.hpp"
#include "DataStructures/DataVector.hpp"
#include "Domain/Structure/Direction.hpp"
#include "Domain/Structure/DirectionalId.hpp"
#include "Domain/Structure/DirectionalIdMap.hpp"
#include "Domain/Structure/ElementId.hpp"
#include "Evolution/DgSubcell/Tags/TciStatus.hpp"
#include "NumericalAlgorithms/Spectral/Mesh.hpp"
#include "Time/TimeStepId.hpp"
#include "Utilities/ErrorHandling/Assert.hpp"
#include "Utilities/GenerateInstantiations.hpp"
#include "Utilities/Gsl.hpp"

namespace evolution::dg::subcell {
template <size_t Dim>
void neighbor_tci_decision(
const gsl::not_null<db::Access*> box,
const std::pair<
const TimeStepId,
DirectionalIdMap<
Dim, std::tuple<Mesh<Dim>, Mesh<Dim - 1>, std::optional<DataVector>,
std::optional<DataVector>, ::TimeStepId, int>>>&
received_temporal_id_and_data) {
db::mutate<subcell::Tags::NeighborTciDecisions<Dim>>(
[&received_temporal_id_and_data](const auto neighbor_tci_decisions_ptr) {
for (const auto& [directional_element_id, neighbor_data] :
received_temporal_id_and_data.second) {
ASSERT(neighbor_tci_decisions_ptr->contains(directional_element_id),
"The NeighborTciDecisions tag does not contain the neighbor "
<< directional_element_id);
neighbor_tci_decisions_ptr->at(directional_element_id) =
std::get<5>(neighbor_data);
}
},
box);
}

#define DIM(data) BOOST_PP_TUPLE_ELEM(0, data)

#define INSTANTIATION(r, data) \
template void neighbor_tci_decision( \
gsl::not_null<db::Access*> box, \
const std::pair< \
const TimeStepId, \
DirectionalIdMap< \
DIM(data), \
std::tuple<Mesh<DIM(data)>, Mesh<DIM(data) - 1>, \
std::optional<DataVector>, std::optional<DataVector>, \
::TimeStepId, int>>>& received_temporal_id_and_data);

GENERATE_INSTANTIATIONS(INSTANTIATION, (1, 2, 3))

#undef INSTANTIATION
#undef DIM
} // namespace evolution::dg::subcell
26 changes: 4 additions & 22 deletions src/Evolution/DgSubcell/NeighborTciDecision.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -8,43 +8,25 @@
#include <tuple>
#include <utility>

#include "DataStructures/DataBox/DataBox.hpp"
#include "DataStructures/DataBox/Access.hpp"
#include "DataStructures/DataVector.hpp"
#include "Domain/Structure/Direction.hpp"
#include "Domain/Structure/DirectionalId.hpp"
#include "Domain/Structure/DirectionalIdMap.hpp"
#include "Domain/Structure/ElementId.hpp"
#include "Evolution/DgSubcell/Tags/TciStatus.hpp"
#include "NumericalAlgorithms/Spectral/Mesh.hpp"
#include "Time/TimeStepId.hpp"
#include "Utilities/ErrorHandling/Assert.hpp"
#include "Utilities/Gsl.hpp"

namespace evolution::dg::subcell {
/*!
* \brief Copies the neighbors' TCI decisions into
* `subcell::Tags::NeighborTciDecisions<Dim>`
*/
template <size_t Dim, typename DbTagsList>
template <size_t Dim>
void neighbor_tci_decision(
const gsl::not_null<db::DataBox<DbTagsList>*> box,
gsl::not_null<db::Access*> box,
const std::pair<
const TimeStepId,
DirectionalIdMap<
Dim, std::tuple<Mesh<Dim>, Mesh<Dim - 1>, std::optional<DataVector>,
std::optional<DataVector>, ::TimeStepId, int>>>&
received_temporal_id_and_data) {
db::mutate<subcell::Tags::NeighborTciDecisions<Dim>>(
[&received_temporal_id_and_data](const auto neighbor_tci_decisions_ptr) {
for (const auto& [directional_element_id, neighbor_data] :
received_temporal_id_and_data.second) {
ASSERT(neighbor_tci_decisions_ptr->contains(directional_element_id),
"The NeighborTciDecisions tag does not contain the neighbor "
<< directional_element_id);
neighbor_tci_decisions_ptr->at(directional_element_id) =
std::get<5>(neighbor_data);
}
},
box);
}
received_temporal_id_and_data);
} // namespace evolution::dg::subcell
Loading
Loading