revamped data reader test

tests/SimpleDataReader_test/GNUmakefile

@@ -0,0 +1,29 @@
+# -*- Mode: Makefile -*-
+
+### This makefile produces an executable for each name in the `ebase'
+###  variable using the libraries named in the `LibNames' variable.
+
+# Included makefiles need an absolute path to the Chombo installation
+# CHOMBO_HOME := Please set the CHOMBO_HOME locally (e.g. export CHOMBO_HOME=... in bash)
+# GRCHOMBO_SOURCE := Set locally (e.g. export GRCHOMBO_SOURCE=path/to/GRChombo/Source in bash)
+
+GRBoondi_SOURCE = $(shell pwd)/../../Source
+
+ebase := Main_datareader
+
+LibNames := AMRTimeDependent AMRTools BoxTools
+
+src_dirs := $(GRCHOMBO_SOURCE)/utils \
+            $(GRCHOMBO_SOURCE)/simd  \
+            $(GRCHOMBO_SOURCE)/BoxUtils  \
+            $(GRCHOMBO_SOURCE)/TaggingCriteria  \
+            $(GRCHOMBO_SOURCE)/GRChomboCore  \
+            $(GRCHOMBO_SOURCE)/AMRInterpolator  \
+            $(GRBoondi_SOURCE)/Background  \
+			$(GRBoondi_SOURCE)/Diagnostics \
+			$(GRBoondi_SOURCE)/Core \
+            $(GRBoondi_SOURCE)/Matter \
+            $(GRBoondi_SOURCE)/Lagrangians \
+            $(GRBoondi_SOURCE)/utils
+
+include $(CHOMBO_HOME)/mk/Make.test

tests/SimpleDataReader_test/Main_datareader.cpp

@@ -0,0 +1,83 @@
+
+// GRBoondi
+#include "DataContainer.hpp"
+#include "DataManipulation.hpp"
+#include "SetupFunctions.hpp" //For setting up MPI processes
+#include "SimpleDataReader.hpp"
+
+#define SMALL_ERR 1e-4
+
+int main(int argc, char *argv[])
+{
+    mainSetup(argc, argv); // setup MPI processes
+
+    int reader_failed { 0 };
+    int interp_failed { 0 };
+
+    SimpleDataReader<double> reader{"data_test.dat"};
+    DataContainer<double> data = reader.get_data();
+
+    // Test the nearest neighbor algorithm
+    std::cout << "###########################################" << std::endl;
+    std::cout << "Testing nearest neighbor algorithm" << std::endl;
+    std::cout << "###########################################" << std::endl;
+
+    //get the coordinates and data values
+    std::vector<std::vector<double>> coords{data.get_coords()};
+    std::vector<double> vals{data.get_data()};
+
+    //set the query point
+    std::vector<double> query_point{2.5, 2.5};
+
+    //find the nearest neighbors
+    std::pair<std::vector<double>, std::vector<std::vector<double>>>
+        nearest_neighbors{
+            DataManipulation::find_nearest_neighbors(coords, query_point, 4)};
+
+    //extract nearest neighbor coordinates
+    std::vector<std::vector<double>> nearest_neighbor_coords{nearest_neighbors.second};
+
+    //compare against known neighbors
+    std::vector<std::vector<double>> ref_data { {2.,2.},{2.,3.},{3.,2.},{3.,3.} };
+    if ( nearest_neighbor_coords == ref_data )
+    {
+        std::cout << "Data reader test ..... PASSED" << std::endl;
+    } else
+    {
+        std::cout << "Data reader test ..... FAILED" << std::endl;
+        reader_failed = 1;
+    }
+
+    // Test the interpolation algorithm
+    std::cout << "\n\n###########################################" << std::endl;
+    std::cout << "Testing interpolation algorithm" << std::endl;
+    std::cout << "###########################################" << std::endl;
+
+    // create 3D point
+    std::vector<double> point1{nearest_neighbors.second[0][0],
+                               nearest_neighbors.second[0][1],
+                               vals[nearest_neighbors.first[0]]};
+    std::vector<double> point2{nearest_neighbors.second[1][0],
+                               nearest_neighbors.second[1][1],
+                               vals[nearest_neighbors.first[1]]};
+    std::vector<double> point3{nearest_neighbors.second[2][0],
+                               nearest_neighbors.second[2][1],
+                               vals[nearest_neighbors.first[2]]};
+
+    //get interpolated value
+    double interp_data { DataManipulation::lin_interp_2d(point1, point2, point3, query_point) };
+
+    //compare against known value
+    double ref_interp_data { 4. };
+    if (interp_data == ref_interp_data)
+    {
+        std::cout << "Interpolation test ..... PASSED" << std::endl;
+    } else 
+    {
+        std::cout << "Interpolation test ..... FAILED" << std::endl;
+        interp_failed = 1;
+    }
+
+    mainFinalize(); // cleanup MPI processes
+    return interp_failed && reader_failed;
+}

tests/SimpleDataReader_test/ProcaField.hpp

@@ -0,0 +1,67 @@
+#ifndef PROCAFIELD_H_INCLUDED
+#define PROCAFIELD_H_INCLUDED
+
+/*
+This class adds the simplest L2 lagrangian to the base equations of motion
+*/
+#include "ADMFixedBGVars.hpp"
+#include "ADMProcaVars.hpp"
+#include "BaseProcaField.hpp"
+#include "DefaultG.hpp"
+#include "KerrSchild.hpp"
+#include "L2_simp.hpp"
+
+// Note: base class BaseProcaField uses CRTP, so pass ProcaField itself as
+// template argument
+class ProcaField : public BaseProcaField<KerrSchild, ProcaField>
+{
+
+  protected:
+    template <class data_t>
+    using MatterVars = typename ADMProcaVars::MatterVars<data_t>;
+
+    template <class data_t>
+    using MatterVarsD2 = typename ADMProcaVars::Diff2MatterVars<data_t>;
+
+    template <class data_t>
+    using MetricVars = typename ADMFixedBGVars::template Vars<data_t>;
+
+    using L2_t = L2<DefaultG>;
+
+  public:
+    struct params_t
+    {
+        double mass;
+        double alpha2;
+        double vector_damping;
+    };
+
+    KerrSchild m_background;
+    params_t m_params;
+
+    ProcaField(KerrSchild a_background, params_t a_params)
+        : BaseProcaField<KerrSchild, ProcaField>(a_background),
+          m_background(a_background), m_params(a_params){};
+
+    template <class data_t>
+    void compute_emtensor_modification(
+        emtensor_t<data_t>
+            &base_emtensor, // pass by reference to allow modifications
+        const MatterVars<data_t> &vars, const MetricVars<data_t> &metric_vars,
+        const MatterVars<Tensor<1, data_t>> &d1,
+        const MatterVarsD2<Tensor<2, data_t>> &d2, // 2nd derivs
+        const MatterVars<data_t> &advec // value of the beta^i d_i(var) terms
+    ) const {};
+
+    template <class data_t, template <typename> class rhs_vars_t>
+    void matter_rhs_modification(
+        rhs_vars_t<data_t> &total_rhs,         // RHS terms for all vars
+        const MatterVars<data_t> &matter_vars, // the value fo the variables
+        const MetricVars<data_t> &metric_vars,
+        const MatterVars<Tensor<1, data_t>> &d1,   // value of 1st derivs
+        const MatterVarsD2<Tensor<2, data_t>> &d2, // 2nd derivs
+        const MatterVars<data_t> &advec // value of the beta^i d_i(var) terms
+    ) const {};
+};
+
+#endif // PROCAFIELD_H_INCLUDED

tests/SimpleDataReader_test/SimulationParameters.hpp

@@ -0,0 +1,24 @@
+#ifndef SIMULATIONPARAMETERS_HPP_
+#define SIMULATIONPARAMETERS_HPP_
+
+#include "KerrSchild.hpp"
+#include "L2_simp.hpp"
+#include "ProcaField.hpp"
+#include "ProcaSimulationParameters.hpp"
+
+class SimulationParameters : public ProcaSimulationParameters
+{
+
+  public:
+    SimulationParameters(GRParmParse &pp) : ProcaSimulationParameters(pp)
+    {
+        read_params(pp);
+        check_params();
+    }
+
+    void read_params(GRParmParse &pp) {}
+
+    void check_params() {}
+};
+
+#endif /* SIMULATIONPARAMETERS_HPP_ */

tests/SimpleDataReader_test/UserVariables.hpp

@@ -0,0 +1,44 @@
+/* GRChombo
+ * Copyright 2012 The GRChombo collaboration.
+ * Please refer to LICENSE in GRChombo's root directory.
+ */
+
+#ifndef USERVARIABLES_HPP
+#define USERVARIABLES_HPP
+
+#include "DiagnosticVariables.hpp"
+
+/// This enum gives the index of every variable stored in the grid
+enum
+{
+    c_phi, // scalar part of proca field
+
+    c_Avec1, // spatial part of proca field
+    c_Avec2,
+    c_Avec3,
+
+    c_Evec1, // conjugate momentum of proca field
+    c_Evec2,
+    c_Evec3,
+
+    c_Z, // auxiliary scalar field
+
+    NUM_VARS
+};
+
+namespace UserVariables
+{
+static const std::array<std::string, NUM_VARS> variable_names = {
+    "phi",
+
+    "Avec1", "Avec2", "Avec3",
+
+    "Evec1", "Evec2", "Evec3",
+
+    "Z"};
+
+} // namespace UserVariables
+
+#include "UserVariables.inc.hpp"
+
+#endif /* USERVARIABLES_HPP */