Revamped data reader test. Redoing KdS tests

tests/KdS_test/ADMFixedBGVars.hpp

@@ -0,0 +1,42 @@
+/* GRChombo
+ * Copyright 2012 The GRChombo collaboration.
+ * Please refer to LICENSE in GRChombo's root directory.
+ *
+ * Modified by GRBoondi 2024
+ * Please refer to LICENSE in GRBoondi's root directory.
+ */
+
+#ifndef ADMFIXEDBGVARS_HPP_
+#define ADMFIXEDBGVARS_HPP_
+
+#include "Tensor.hpp"
+#include "UserVariables.hpp"
+#include "VarsTools.hpp"
+
+/// Namespace for ADM vars for fixed BG evolution
+namespace ADMFixedBGVars
+{
+/// Vars object for ADM vars used in FixedBG evolution
+template <class data_t> struct Vars
+{
+    // ADM vars needed in matter only rhs (ok for Proca and SF)
+    Tensor<2, data_t> K_tensor;
+    data_t K;
+
+    data_t lapse;
+    Tensor<1, data_t> shift;
+    Tensor<2, data_t> gamma;
+
+    Tensor<1, data_t> d1_lapse;
+    Tensor<1, Tensor<1, data_t>> d1_shift;
+    Tensor<2, Tensor<1, data_t>> d1_gamma;
+
+    // Optional second derivatives of the vars
+    Tensor<2, data_t> d2_lapse;
+    Tensor<1, Tensor<2, data_t>> d2_shift;
+    Tensor<2, Tensor<2, data_t>> d2_gamma;
+};
+
+}; // namespace ADMFixedBGVars
+
+#endif /* ADMFIXEDBGVARS_HPP_ */

tests/KdS_test/AssignFixedBGtoBSSNVars.hpp

@@ -0,0 +1,64 @@
+/* GRChombo
+ * Copyright 2012 The GRChombo collaboration.
+ * Please refer to LICENSE in GRChombo's root directory.
+ */
+
+#ifndef ASSIGNFIXEDBGTOBSSNVARS_HPP_
+#define ASSIGNFIXEDBGTOBSSNVARS_HPP_
+
+#include "ADMFixedBGVars.hpp"
+#include "BoxLoops.hpp"
+#include "CCZ4RHS.hpp"
+#include "Cell.hpp"
+#include "Coordinates.hpp"
+#include "MetricVariablesInterface.hpp"
+#include "Tensor.hpp"
+#include "UserVariables.hpp" //This files needs NUM_VARS - total number of components
+#include "simd.hpp"
+
+template <class background_t> class AssignFixedBGtoBSSNVars
+{
+  public:
+    AssignFixedBGtoBSSNVars(const background_t a_background, const double a_dx,
+                            const std::array<double, CH_SPACEDIM> a_center)
+        : m_dx(a_dx), m_center(a_center), m_background(a_background)
+    {
+    }
+
+    template <class data_t> void compute(Cell<data_t> current_cell) const
+    {
+        // compute the ADM vars
+        const Coordinates<data_t> coords(current_cell, m_dx, m_center);
+        ADMFixedBGVars::Vars<data_t> adm_vars;
+        m_background.compute_metric_background(adm_vars, coords);
+
+        // assign values to the BSSN vars
+        using namespace TensorAlgebra;
+        CCZ4RHS<>::Vars<data_t> bssn_vars;
+
+        bssn_vars.K = adm_vars.K;
+        bssn_vars.lapse = adm_vars.lapse;
+        FOR1(i) { bssn_vars.shift[i] = adm_vars.shift[i]; }
+        bssn_vars.chi = compute_determinant_sym(adm_vars.gamma);
+        bssn_vars.chi = pow(bssn_vars.chi, -1.0 / 3.0);
+        FOR2(i, j)
+        {
+            bssn_vars.h[i][j] = adm_vars.gamma[i][j] * bssn_vars.chi;
+            bssn_vars.A[i][j] = adm_vars.K_tensor[i][j] * bssn_vars.chi -
+                                1.0 / 3.0 * bssn_vars.K * bssn_vars.h[i][j];
+        }
+
+        // set the field to something arbitrary
+        data_t r = coords.get_radius();
+        current_cell.store_vars(0.1 * sin(coords.x) * exp(-r) * coords.z,
+                                c_Avec1);
+        current_cell.store_vars(bssn_vars);
+    }
+
+  protected:
+    const double m_dx;
+    const std::array<double, CH_SPACEDIM> m_center;
+    const background_t m_background;
+};
+
+#endif /* ASSIGNFIXEDBGTOBSSNVARS_HPP_ */

tests/KdS_test/BaseProcaFieldTest.hpp

@@ -0,0 +1,78 @@
+/* GRBoondi 2024
+ * Please refer to LICENSE in GRBoondi's root directory.
+ */
+
+/*
+   Base class for the Proca field. Must be inherited by the Proca field class
+defined by the user
+
+NOTE:
+       We use the method of 'Curiously Recurring Template Pattern' to allow
+arbitrary templated modifications to the theory
+       https://en.wikipedia.org/wiki/Curiously_recurring_template_pattern
+
+       This allows us to use virtual functions that can be specified by derived
+classes. e.g. see ProcaKerrBH Example
+*/
+#ifndef BASEPROCAFIELDTEST_H_INCLUDED
+#define BASEPROCAFIELDTEST_H_INCLUDED
+
+#include "ADMFixedBGVars.hpp" //For metric variables
+#include "ADMProcaVars.hpp"   //For matter variables
+#include "CCZ4RHS.hpp"
+#include "DefaultBackground.hpp"      //Minkowski background as default
+#include "FourthOrderDerivatives.hpp" //For calculating derivatives
+#include "Tensor.hpp"                 //For performing tensorial operations
+#include "TensorAlgebra.hpp"
+#include "UserVariables.hpp" //For user-defined variables (e.g. see EMKerrBH)
+#include "VarsTools.hpp"
+#include "simd.hpp" //for SIMD operations
+
+#include "MetricVariablesInterface.hpp"
+#include "ProcaField.hpp"
+
+template <class background_t, class modification_t> class BaseProcaFieldTest
+{
+  protected:
+    const background_t m_background;
+
+    struct params_t
+    {
+    };
+
+  public:
+    // constructor, inputs are matter params
+    BaseProcaFieldTest(background_t a_background)
+        : m_background{a_background} {};
+
+    template <class data_t> using Vars = ADMProcaVars::MatterVars<data_t>;
+
+    template <class data_t>
+    using Diff2Vars = ADMProcaVars::Diff2MatterVars<data_t>;
+
+    // we set G=0, so this doenst really matter, need it just for using
+    // MatterCCZ4
+    template <class data_t, template <typename> class vars_t>
+    emtensor_t<data_t> compute_emtensor(
+        const vars_t<data_t> &vars,          //!< the value of the variables
+        const vars_t<Tensor<1, data_t>> &d1, //!< the value of the 1st derivs
+        const Tensor<2, data_t> &h_UU, //!< the inverse metric (raised indices)
+        const Tensor<3, data_t> &chris_ULL) const
+    {
+        return emtensor_t<data_t>{};
+    }; //!< the conformal christoffel symbol
+
+    template <class data_t, template <typename> class vars_t,
+              template <typename> class diff2_vars_t,
+              template <typename> class rhs_vars_t>
+    void add_matter_rhs(
+        rhs_vars_t<data_t> &total_rhs,             // RHS terms for all vars
+        const vars_t<data_t> &matter_vars,         // the value fo the variables
+        const vars_t<Tensor<1, data_t>> &d1,       // value of 1st derivs
+        const diff2_vars_t<Tensor<2, data_t>> &d2, // 2nd derivs
+        const vars_t<data_t> &advec // value of the beta^i d_i(var) terms
+    ) const;
+};
+
+#include "BaseProcaFieldTest.impl.hpp"
+#endif // BASEPROCAFIELDTEST_H_INCLUDED

tests/KdS_test/BaseProcaFieldTest.impl.hpp

@@ -0,0 +1,141 @@
+/* GRBoondi 2024
+ * Please refer to LICENSE in GRBoondi's root directory.
+ */
+
+/*
+implementation file for ProcaField.hpp
+*/
+
+#if !defined(BASEPROCAFIELDTEST_H_INCLUDED)
+#error "This file should only be included through GeneralizedProcaField.hpp"
+#endif
+
+#ifndef BASEPROCAFIELDTEST_IMPL_H_INCLUDED
+#define BASEPROCAFIELDTEST_IMPL_H_INCLUDED
+
+template <class background_t, class modification_t>
+template <class data_t, template <typename> class vars_t,
+          template <typename> class diff2_vars_t,
+          template <typename> class rhs_vars_t>
+void BaseProcaFieldTest<background_t, modification_t>::add_matter_rhs(
+    rhs_vars_t<data_t> &total_rhs,             // RHS terms for all vars
+    const vars_t<data_t> &matter_vars,         // the value fo the variables
+    const vars_t<Tensor<1, data_t>> &d1,       // value of 1st derivs
+    const diff2_vars_t<Tensor<2, data_t>> &d2, // 2nd derivs
+    const vars_t<data_t> &advec // value of the beta^i d_i(var) terms
+) const
+{
+
+    // compute non-conformal metric
+    CCZ4RHS<>::Vars<data_t> ccz4_vars;
+    ccz4_vars.lapse = matter_vars.lapse;
+    ccz4_vars.chi = matter_vars.chi;
+    ccz4_vars.K = matter_vars.K;
+    FOR1(i) { ccz4_vars.shift[i] = matter_vars.shift[i]; };
+    FOR2(i, j)
+    {
+        ccz4_vars.h[i][j] = matter_vars.h[i][j];
+        ccz4_vars.A[i][j] = matter_vars.A[i][j];
+    }
+
+    /*     ADMFixedBGVars::Vars<data_t> fixed_bg_vars {
+     * MetricVarsInterface<data_t>(ccz4_vars) };
+     */
+
+    ADMFixedBGVars::Vars<data_t> fixed_bg_vars;
+    fixed_bg_vars.lapse = matter_vars.lapse;
+    fixed_bg_vars.K = matter_vars.K;
+    FOR1(i)
+    {
+        fixed_bg_vars.shift[i] = matter_vars.shift[i];
+        fixed_bg_vars.d1_lapse[i] = d1.lapse[i];
+    };
+    FOR2(i, j)
+    {
+        fixed_bg_vars.gamma[i][j] = matter_vars.h[i][j] / matter_vars.chi;
+        fixed_bg_vars.K_tensor[i][j] =
+            (matter_vars.A[i][j] +
+             1.0 / 3.0 * matter_vars.K * matter_vars.h[i][j]) /
+            matter_vars.chi;
+        fixed_bg_vars.d1_shift[i][j] = d1.shift[i][j];
+    }
+    FOR3(i, j, k)
+    {
+        fixed_bg_vars.d1_gamma[i][j][k] =
+            d1.h[i][j][k] / matter_vars.chi -
+            matter_vars.h[i][j] / matter_vars.chi / matter_vars.chi * d1.chi[k];
+    }
+
+    auto h_UU = TensorAlgebra::compute_inverse_sym(ccz4_vars.h);
+    const auto non_phys_chris =
+        TensorAlgebra::compute_christoffel(d1.h, h_UU).ULL;
+    const Tensor<3, data_t> chris_phys = TensorAlgebra::compute_phys_chris(
+        d1.chi, ccz4_vars.chi, ccz4_vars.h, h_UU, non_phys_chris);
+
+    // calculate conformal contravariant metric and conformal christoffel
+    // symbols
+    const Tensor<2, data_t> gamma_UU =
+        TensorAlgebra::compute_inverse(fixed_bg_vars.gamma);
+
+    // covariant derivative of spatial part of Proca field
+    Tensor<2, data_t> DA;
+    FOR2(i, j)
+    {
+        DA[i][j] = d1.Avec[j][i];
+        FOR1(k) { DA[i][j] -= chris_phys[k][i][j] * matter_vars.Avec[k]; }
+    }
+
+    // evolution equations for spatial part of vector field (index down)
+    FOR1(i)
+    {
+        total_rhs.Avec[i] = -fixed_bg_vars.lapse * d1.phi[i] -
+                            matter_vars.phi * d1.lapse[i] + advec.Avec[i];
+
+        FOR1(j)
+        {
+            total_rhs.Avec[i] += -fixed_bg_vars.lapse *
+                                     fixed_bg_vars.gamma[i][j] *
+                                     matter_vars.Evec[j] +
+                                 matter_vars.Avec[j] * d1.shift[j][i];
+        };
+    };
+
+    // evolution equations for Electric vector field (index up)
+    FOR1(i)
+    {
+        total_rhs.Evec[i] =
+            fixed_bg_vars.lapse * fixed_bg_vars.K * matter_vars.Evec[i] +
+            advec.Evec[i];
+
+        FOR1(j) { total_rhs.Evec[i] += -matter_vars.Evec[j] * d1.shift[i][j]; }
+
+        FOR3(j, k, l)
+        {
+            total_rhs.Evec[i] +=
+                gamma_UU[j][k] * gamma_UU[i][l] *
+                (d1.lapse[j] * (d1.Avec[k][l] - d1.Avec[l][k]) +
+                 fixed_bg_vars.lapse * (d2.Avec[k][l][j] - d2.Avec[l][k][j]));
+
+            FOR1(m)
+            {
+                total_rhs.Evec[i] +=
+                    -fixed_bg_vars.lapse * gamma_UU[j][k] * gamma_UU[i][l] *
+                    (chris_phys[m][j][l] * (d1.Avec[k][m] - d1.Avec[m][k]) +
+                     chris_phys[m][j][k] * (d1.Avec[m][l] - d1.Avec[l][m]));
+            };
+        };
+    };
+
+    // Evolution equations for phi field totally depend on the theory, so we
+    // leave is up to the user to specify them for their model
+    total_rhs.phi = 0.;
+
+    // Evolution for auxiliary Z field is also left up to the user in how they
+    // want to add damping terms
+    total_rhs.Z = 0.;
+
+    static_cast<const modification_t *>(this)->matter_rhs_modification(
+        total_rhs, matter_vars, fixed_bg_vars, d1, d2, advec);
+};
+
+#endif // BASEPROCAFIELDTEST_IMPL_H_INCLUDED

tests/KdS_test/EmptyDiagnosticVariables.hpp

@@ -0,0 +1,20 @@
+/* GRChombo
+ * Copyright 2012 The GRChombo collaboration.
+ * Please refer to LICENSE in GRChombo's root directory.
+ */
+
+#ifndef DIAGNOSTICVARIABLES_HPP
+#define DIAGNOSTICVARIABLES_HPP
+
+// assign an enum to each variable
+enum
+{
+    NUM_DIAGNOSTIC_VARS
+};
+
+namespace DiagnosticVariables
+{
+static const std::array<std::string, NUM_DIAGNOSTIC_VARS> variable_names = {};
+}
+
+#endif /* DIAGNOSTICVARIABLES_HPP */