Fix load_vars in InitialFluidData

dinatraykova · Apr 4, 2024 · b1b9f9b · b1b9f9b
1 parent a5bb17e
commit b1b9f9b
Show file tree

Hide file tree

Showing 3 changed files with 52 additions and 49 deletions.
diff --git a/Examples/Fluid_Kerr/InitialFluidData.hpp b/Examples/Fluid_Kerr/InitialFluidData.hpp
@@ -35,12 +35,46 @@ class InitialFluidData : public EoS
     };
 
     //! The constructor
-    InitialFluidData(params_t a_params, double a_dx);
-
+    InitialFluidData(params_t a_params, double a_dx)
+        : m_dx(a_dx), m_params(a_params)
+    {
+    }
     //! Function to compute the value of all the initial vars on the grid
     //! The constructor
     template <class data_t> void compute(Cell<data_t> current_cell) const;
 
+    //! Structure containing the rhs variables for the matter fields
+    template <class data_t> struct Vars
+    {
+        data_t D;
+        Tensor<1, data_t> Sj;
+        data_t tau;
+        //        data_t Jt;
+        Tensor<1, data_t> vi;
+        data_t rho;
+        data_t eps;
+        //        data_t nn;
+
+        /// Defines the mapping between members of Vars and Chombo grid
+        /// variables (enum in User_Variables)
+        template <typename mapping_function_t>
+        void enum_mapping(mapping_function_t mapping_function)
+        {
+            using namespace VarsTools; // define_enum_mapping is part of
+            // VarsTools
+            define_enum_mapping(mapping_function, c_D, D);
+            define_enum_mapping(mapping_function, GRInterval<c_Sj1, c_Sj3>(),
+                                Sj);
+            define_enum_mapping(mapping_function, c_tau, tau);
+            //            define_enum_mapping(mapping_function, c_Jt, Jt);
+            define_enum_mapping(mapping_function, GRInterval<c_vi1, c_vi3>(),
+                                vi);
+            define_enum_mapping(mapping_function, c_rho, rho);
+            define_enum_mapping(mapping_function, c_eps, eps);
+            //           define_enum_mapping(mapping_function, c_nn, nn);
+        }
+    };
+
   protected:
     double m_dx;
     const params_t m_params;

diff --git a/Examples/Fluid_Kerr/InitialFluidData.impl.hpp b/Examples/Fluid_Kerr/InitialFluidData.impl.hpp
@@ -10,18 +10,13 @@
 #ifndef INITIALFLUIDDATA_IMPL_HPP_
 #define INITIALFLUIDDATA_IMPL_HPP_
 
-inline InitialFluidData::InitialFluidData(params_t a_params, double a_dx)
-    : m_dx(a_dx), m_params(a_params)
-{
-}
-
 template <class data_t>
 void InitialFluidData::compute(Cell<data_t> current_cell) const
 {
     // load vars
-    // FluidCCZ4RHS<PerfectFluid<>>::Vars<data_t> vars;
-    // VarsTools::assign(vars, 0.);
     const auto metric_vars = current_cell.template load_vars<MetricVars>();
+    auto matter_vars = current_cell.template load_vars<Vars>();
+    VarsTools::assign(matter_vars, 0.);
 
     // where am i?
     Coordinates<data_t> coords(current_cell, m_dx, m_params.center);
@@ -32,44 +27,35 @@ void InitialFluidData::compute(Cell<data_t> current_cell) const
     double y = coords.y;
     double z = coords.z;
 
-    Tensor<1, data_t> vi, Sj;
+    //    Tensor<1, data_t> vi, Sj;
     data_t chi_regularised = simd_max(metric_vars.chi, 1e-6);
 
-    FOR(i) { vi[i] = 0.; }
-    data_t D = 0.;
-    data_t tau = 0.;
-    data_t rho = 0.;
-    data_t eps = 0.;
-
     // calculate the field value
-    rho =
+    matter_vars.rho =
         m_params.rho0 * (exp(-pow(rr / m_params.awidth, 2.0))) + m_params.delta;
     data_t v2 = 0.;
-    FOR(i, j) v2 += metric_vars.h[i][j] * vi[i] * vi[j] / chi_regularised;
+    FOR(i, j)
+    v2 += metric_vars.h[i][j] * matter_vars.vi[i] * matter_vars.vi[j] /
+          chi_regularised;
 
-    data_t P_over_rho = (1. + eps) / 3.;
-    // data_t P_over_rho = 0.;
-    // EoS::compute_eos(P_over_rho, vars);
+    data_t P_over_rho = 0.;
+    EoS::compute_eos(P_over_rho, matter_vars);
 
     data_t WW = 1. / (1. - v2);
-    data_t hh = 1. + eps + P_over_rho;
+    data_t hh = 1. + matter_vars.eps + P_over_rho;
 
-    D = rho * sqrt(WW);
-    tau = rho * (hh * WW - P_over_rho) - D;
+    matter_vars.D = matter_vars.rho * sqrt(WW);
+    matter_vars.tau = matter_vars.rho * (hh * WW - P_over_rho) - matter_vars.D;
     FOR(i)
     {
-        Sj[i] = 0.;
+        matter_vars.Sj[i] = 0.;
         FOR(j)
-        Sj[i] += rho * hh * WW * metric_vars.h[i][j] * vi[j] / chi_regularised;
+        matter_vars.Sj[i] += matter_vars.rho * hh * WW * metric_vars.h[i][j] *
+                             matter_vars.vi[j] / chi_regularised;
     }
 
     // store the vars
-    // current_cell.store_vars(vars);
-    current_cell.store_vars(rho, c_rho);
-    current_cell.store_vars(vi, GRInterval<c_vi1, c_vi3>());
-    current_cell.store_vars(D, c_D);
-    current_cell.store_vars(Sj, GRInterval<c_Sj1, c_Sj3>());
-    current_cell.store_vars(tau, c_tau);
+    current_cell.store_vars(matter_vars);
 }
 
 #endif /* INITIALFLUIDDATA_IMPL_HPP_ */
diff --git a/Examples/Fluid_Kerr/PerfectFluid.hpp b/Examples/Fluid_Kerr/PerfectFluid.hpp
@@ -77,23 +77,6 @@ template <class eos_t = DefaultEoS> class PerfectFluid
         }
     };
 
-    //! Structure containing the rhs variables for the matter fields requiring
-    //!  2nd derivs
-    //    template <class data_t> struct Diff2Vars
-    //{
-    //    data_t rho;
-
-    /// Defines the mapping between members of Vars and Chombo grid
-    ///  variables (enum in User_Variables)
-    //  template <typename mapping_function_t>
-    //  void enum_mapping(mapping_function_t mapping_function)
-    //  {
-    //      VarsTools::define_enum_mapping(mapping_function, c_rho, rho);
-    //        VarsTools::define_enum_mapping(
-    //              mapping_function, GRInterval<c_Avec1, c_Avec3>(), Avec);
-    //      }
-    //};
-
     //! The function which calculates the EM Tensor, given the vars and
     //! derivatives
     template <class data_t, template <typename> class vars_t>