Update CSW SetInitialData

Update CSW SetInitialData Fix test
sxs-collaboration · Dec 14, 2024 · ac800d7 · ac800d7
1 parent fd14167
commit ac800d7
Show file tree

Hide file tree

Showing 4 changed files with 103 additions and 22 deletions.
diff --git a/src/Evolution/Systems/CurvedScalarWave/Actions/CMakeLists.txt b/src/Evolution/Systems/CurvedScalarWave/Actions/CMakeLists.txt
@@ -4,12 +4,12 @@
 spectre_target_sources(
   ${LIBRARY}
   PRIVATE
-  NumericInitialData.cpp
+  SetInitialData.cpp
   )
 
 spectre_target_headers(
   ${LIBRARY}
   INCLUDE_DIRECTORY ${CMAKE_SOURCE_DIR}/src
   HEADERS
-  NumericInitialData.hpp
+  SetInitialData.hpp
   )
diff --git a/...ScalarWave/Actions/NumericInitialData.cpp → ...rvedScalarWave/Actions/SetInitialData.cpp b/...ScalarWave/Actions/NumericInitialData.cpp → ...rvedScalarWave/Actions/SetInitialData.cpp
@@ -1,9 +1,10 @@
 // Distributed under the MIT License.
 // See LICENSE.txt for details.
 
-#include "Evolution/Systems/CurvedScalarWave/Actions/NumericInitialData.hpp"
+#include "Evolution/Systems/CurvedScalarWave/Actions/SetInitialData.hpp"
 
 #include <boost/functional/hash.hpp>
+#include <optional>
 #include <string>
 #include <utility>
 #include <variant>
@@ -15,9 +16,12 @@ namespace CurvedScalarWave {
 NumericInitialData::NumericInitialData(
     std::string file_glob, std::string subfile_name,
     std::variant<double, importers::ObservationSelector> observation_value,
-    const bool enable_interpolation, ScalarVars selected_variables)
-    : importer_options_(std::move(file_glob), std::move(subfile_name),
-                        observation_value, enable_interpolation),
+    std::optional<double> observation_value_epsilon, bool enable_interpolation,
+    ScalarVars selected_variables
+    )
+    : importer_options_(
+          std::move(file_glob), std::move(subfile_name), observation_value,
+          observation_value_epsilon.value_or(1.0e-12), enable_interpolation),
       selected_variables_(std::move(selected_variables)) {}
 
 NumericInitialData::NumericInitialData(CkMigrateMessage* msg)

diff --git a/...ScalarWave/Actions/NumericInitialData.hpp → ...rvedScalarWave/Actions/SetInitialData.hpp b/...ScalarWave/Actions/NumericInitialData.hpp → ...rvedScalarWave/Actions/SetInitialData.hpp
@@ -4,6 +4,7 @@
 #pragma once
 
 #include <cstddef>
+#include <optional>
 #include <pup.h>
 #include <string>
 #include <variant>
@@ -36,8 +37,13 @@
 #include "Utilities/TMPL.hpp"
 #include "Utilities/TaggedTuple.hpp"
 
-namespace CurvedScalarWave {
+/// \cond
+namespace Tags {
+struct Time;
+}  // namespace Tags
+/// \endcond
 
+namespace CurvedScalarWave {
 /*!
  * \brief Numeric initial data loaded from volume data files
  */
@@ -75,10 +81,7 @@ class NumericInitialData : public evolution::initial_data::InitialData {
   };
 
   using options =
-      tmpl::list<importers::OptionTags::FileGlob,
-                 importers::OptionTags::Subgroup,
-                 importers::OptionTags::ObservationValue,
-                 importers::OptionTags::EnableInterpolation, Variables>;
+      tmpl::push_back<importers::ImporterOptions::tags_list, Variables>;
 
   static constexpr Options::String help =
       "Numeric initial data loaded from volume data files";
@@ -104,6 +107,7 @@ class NumericInitialData : public evolution::initial_data::InitialData {
   NumericInitialData(
       std::string file_glob, std::string subfile_name,
       std::variant<double, importers::ObservationSelector> observation_value,
+      std::optional<double> observation_value_epsilon,
       bool enable_interpolation, ScalarVars selected_variables);
 
   const importers::ImporterOptions& importer_options() const {
@@ -159,11 +163,11 @@ namespace Actions {
  * \brief Dispatch loading numeric initial data from files.
  *
  * Place this action before
- * CurvedScalarWave::Actions::SetNumericInitialData in the action list.
+ * ScalarTensor::Actions::SetNumericInitialData in the action list.
  * See importers::Actions::ReadAllVolumeDataAndDistribute for details, which is
  * invoked by this action.
  */
-struct ReadNumericInitialData {
+struct SetInitialData {
   using const_global_cache_tags =
       tmpl::list<evolution::initial_data::Tags::InitialData>;
 
@@ -175,15 +179,40 @@ struct ReadNumericInitialData {
       const tuples::TaggedTuple<InboxTags...>& /*inboxes*/,
       Parallel::GlobalCache<Metavariables>& cache,
       const ArrayIndex& /*array_index*/, const ActionList /*meta*/,
+      const ParallelComponent* const parallel_component) {
+    // Dispatch to the correct `apply` overload based on type of initial data
+    using initial_data_classes =
+        tmpl::at<typename Metavariables::factory_creation::factory_classes,
+                 evolution::initial_data::InitialData>;
+    return call_with_dynamic_type<Parallel::iterable_action_return_t,
+                                  initial_data_classes>(
+        &db::get<evolution::initial_data::Tags::InitialData>(box),
+        [&box, &cache, &parallel_component](const auto* const initial_data) {
+          return apply(make_not_null(&box), *initial_data, cache,
+                       parallel_component);
+        });
+  }
+
+ private:
+  static constexpr size_t Dim = 3;
+
+  // Numeric initial data
+  template <typename DbTagsList, typename Metavariables,
+            typename ParallelComponent>
+  static Parallel::iterable_action_return_t apply(
+      const gsl::not_null<db::DataBox<DbTagsList>*> /*box*/,
+      const NumericInitialData& initial_data,
+      Parallel::GlobalCache<Metavariables>& cache,
       const ParallelComponent* const /*meta*/) {
     // Select the subset of the available variables that we want to read from
     // the volume data file
-    const auto& initial_data = dynamic_cast<const NumericInitialData&>(
-        db::get<evolution::initial_data::Tags::InitialData>(box));
     tuples::tagged_tuple_from_typelist<db::wrap_tags_in<
         importers::Tags::Selected, NumericInitialData::all_vars>>
         selected_fields{};
     initial_data.select_for_import(make_not_null(&selected_fields));
+    // Dispatch loading the variables from the volume data file
+    // - Not using `ckLocalBranch` here to make sure the simple action
+    //   invocation is asynchronous.
     auto& reader_component = Parallel::get_parallel_component<
         importers::ElementDataReader<Metavariables>>(cache);
     Parallel::simple_action<importers::Actions::ReadAllVolumeDataAndDistribute<
@@ -192,13 +221,58 @@ struct ReadNumericInitialData {
         initial_data.volume_data_id(), std::move(selected_fields));
     return {Parallel::AlgorithmExecution::Continue, std::nullopt};
   }
+
+  // "AnalyticData"-type initial data
+  template <typename DbTagsList, typename InitialData, typename Metavariables,
+            typename ParallelComponent>
+  static Parallel::iterable_action_return_t apply(
+      const gsl::not_null<db::DataBox<DbTagsList>*> box,
+      const InitialData& initial_data,
+      Parallel::GlobalCache<Metavariables>& /*cache*/,
+      const ParallelComponent* const /*meta*/) {
+    // Get scalar variables from analytic data / solution
+    const auto& [coords, mesh, inv_jacobian] = [&box]() {
+      return std::forward_as_tuple(
+          db::get<domain::Tags::Coordinates<Dim, Frame::Inertial>>(*box),
+          db::get<domain::Tags::Mesh<Dim>>(*box),
+          db::get<domain::Tags::InverseJacobian<Dim, Frame::ElementLogical,
+                                                Frame::Inertial>>(*box));
+    }();
+    auto vars = evolution::Initialization::initial_data(
+        initial_data, coords, db::get<::Tags::Time>(*box),
+        tmpl::append<tmpl::list<>,
+                     // Don't use the scalar gradient
+                     tmpl::list<CurvedScalarWave::Tags::Psi,
+                                CurvedScalarWave::Tags::Pi>>{});
+
+    // Move scalar variables and compute gradient
+    db::mutate<CurvedScalarWave::Tags::Psi, CurvedScalarWave::Tags::Pi,
+               CurvedScalarWave::Tags::Phi<3>>(
+        [&vars](const gsl::not_null<Scalar<DataVector>*> psi_scalar,
+                const gsl::not_null<Scalar<DataVector>*> pi_scalar,
+                const gsl::not_null<tnsr::i<DataVector, 3>*> phi_scalar,
+                const Mesh<3>& local_mesh,
+                const InverseJacobian<DataVector, 3_st, Frame::ElementLogical,
+                                      Frame::Inertial>& local_inv_jacobian) {
+          *psi_scalar = std::move(get<CurvedScalarWave::Tags::Psi>(vars));
+          *pi_scalar = std::move(get<CurvedScalarWave::Tags::Pi>(vars));
+          // Set Phi to the numerical spatial derivative of the scalar
+          partial_derivative(phi_scalar, *psi_scalar, local_mesh,
+                             local_inv_jacobian);
+        },
+        box, mesh, inv_jacobian);
+
+    // No need to import numeric initial data, so we terminate the phase by
+    // pausing the algorithm on this element
+    return {Parallel::AlgorithmExecution::Pause, std::nullopt};
+  }
 };
 
 /*!
  * \brief Receive numeric initial data loaded by
- * CurvedScalarWave::Actions::ReadNumericInitialData.
+ * CurvedScalarWave::Actions::SetInitialData.
  */
-struct SetNumericInitialData {
+struct ReceiveNumericInitialData {
   static constexpr size_t Dim = 3;
   using inbox_tags =
       tmpl::list<importers::Tags::VolumeData<NumericInitialData::all_vars>>;

diff --git a/tests/Unit/Evolution/Systems/CurvedScalarWave/Actions/Test_NumericInitialData.cpp b/tests/Unit/Evolution/Systems/CurvedScalarWave/Actions/Test_NumericInitialData.cpp
@@ -21,7 +21,7 @@
 #include "Domain/CoordinateMaps/Wedge.hpp"
 #include "Domain/Structure/ElementId.hpp"
 #include "Domain/Tags.hpp"
-#include "Evolution/Systems/CurvedScalarWave/Actions/NumericInitialData.hpp"
+#include "Evolution/Systems/CurvedScalarWave/Actions/SetInitialData.hpp"
 #include "Framework/ActionTesting.hpp"
 #include "Framework/TestCreation.hpp"
 #include "IO/Importers/Actions/ReadVolumeData.hpp"
@@ -61,8 +61,8 @@ struct MockElementArray {
                                             Frame::Inertial>>>>>,
       Parallel::PhaseActions<
           Parallel::Phase::Testing,
-          tmpl::list<CurvedScalarWave::Actions::ReadNumericInitialData,
-                     CurvedScalarWave::Actions::SetNumericInitialData>>>;
+          tmpl::list<CurvedScalarWave::Actions::SetInitialData,
+                     CurvedScalarWave::Actions::ReceiveNumericInitialData>>>;
 };
 
 struct MockReadVolumeData {
@@ -147,7 +147,8 @@ void test_numeric_initial_data(const NumericInitialData& initial_data,
   const auto map =
       domain::make_coordinate_map<Frame::ElementLogical, Frame::Inertial>(
           domain::CoordinateMaps::Wedge<3>{
-              wave_radius / 2., wave_radius * 2., 1., 1., {}, true});
+              wave_radius / 2., wave_radius * 2., 1., 1.,
+              OrientationMap<3>::create_aligned(), true});
   const auto logical_coords = logical_coordinates(mesh);
   const auto coords = map(logical_coords);
   const auto inv_jacobian = map.inv_jacobian(logical_coords);
@@ -204,13 +205,15 @@ SPECTRE_TEST_CASE("Unit.Evolution.Systems.CurvedScalarWave.SetInitialData",
   test_numeric_initial_data(NumericInitialData{"TestInitialData.h5",
                                                "VolumeData",
                                                0.,
+                                               {1.0e-9},
                                                false,
                                                {"CustomPsi", "CustomPi"}},
                             "NumericInitialData:\n"
                             "  FileGlob: TestInitialData.h5\n"
                             "  Subgroup: VolumeData\n"
                             "  ObservationValue: 0.\n"
-                            "  Interpolate: False\n"
+                            "  ObservationValueEpsilon: 1e-9\n"
+                            "  ElementsAreIdentical: False\n"
                             "  Variables:\n"
                             "    Psi: CustomPsi\n"
                             "    Pi: CustomPi");