Fixes after rebase

Author: TaigoFr

Examples/KerrBH/Main_KerrBH.cpp

@@ -46,7 +46,8 @@ int runGRChombo(int argc, char *argv[])
     // Note: 'interpolator' needs to be in scope when gr_amr.run() is called,
     // otherwise pointer is lost
     AMRInterpolator<Lagrange<4>> interpolator(
-        gr_amr, sim_params.origin, sim_params.dx, sim_params.verbosity);
+        gr_amr, sim_params.origin, sim_params.dx, sim_params.boundary_params,
+        sim_params.verbosity);
     gr_amr.set_interpolator(&interpolator);
 
     using Clock = std::chrono::steady_clock;

Examples/KerrBH/params_cheap.txt

@@ -82,7 +82,7 @@ plot_interval = 0
 num_plot_vars = 0
 dt_multiplier = 0.25
 stop_time = 10.0
-max_steps = 2
+max_steps = 1
 
 #Lapse evolution
 lapse_power = 1.0
@@ -114,4 +114,4 @@ extraction_levels = 1 1 0
 num_points_phi = 50
 num_points_theta = 52
 extraction_extrapolation_order = 3
-extraction_extrapolation_radii = 2 0 1
+extraction_extrapolation_radii = 0 1 2

Source/AMRInterpolator/SurfaceExtraction.hpp

@@ -194,7 +194,7 @@ template <class SurfaceGeometry> class SurfaceExtraction
                         const std::vector<double> &a_integrals,
                         const std::string &a_label = "") const;
 
-protected:
+  protected:
     // returns empty vector if no extrapolation is done
     std::vector<double> richardson_extrapolation(
         const std::vector<std::vector<double>> &integrals) const;

Source/AMRInterpolator/SurfaceExtraction.impl.hpp

@@ -403,7 +403,8 @@ void SurfaceExtraction<SurfaceGeometry>::write_integrals(
 {
     if (procID() == 0)
     {
-        std::vector<double> extrapolations = richardson_extrapolation(a_integrals);
+        std::vector<double> extrapolations =
+            richardson_extrapolation(a_integrals);
         bool extrapolation_done = extrapolations.size() > 0;
 
         const int num_integrals_per_surface = a_integrals.size();
@@ -439,7 +440,8 @@ void SurfaceExtraction<SurfaceGeometry>::write_integrals(
             for (int iintegral = 0; iintegral < num_integrals_per_surface;
                  ++iintegral)
             {
-                for (int isurface = 0; isurface < m_params.num_surfaces; ++isurface)
+                for (int isurface = 0; isurface < m_params.num_surfaces;
+                     ++isurface)
                 {
                     int idx = isurface * num_integrals_per_surface + iintegral;
                     if (a_labels.empty())
@@ -451,8 +453,9 @@ void SurfaceExtraction<SurfaceGeometry>::write_integrals(
                 }
                 if (extrapolation_done)
                 {
-                    int idx = m_params.num_surfaces * num_integrals_per_surface +
-                              iintegral;
+                    int idx =
+                        m_params.num_surfaces * num_integrals_per_surface +
+                        iintegral;
                     header1_strings[idx] = a_labels[iintegral];
                     header2_strings[idx] = "Infinity";
                 }
@@ -468,17 +471,18 @@ void SurfaceExtraction<SurfaceGeometry>::write_integrals(
         // make vector of data for writing
         std::vector<double> data_for_writing(num_integrals_per_surface *
                                              total_num_surfaces);
-         for (int iintegral = 0; iintegral < num_integrals_per_surface; ++iintegral)
-         {
+        for (int iintegral = 0; iintegral < num_integrals_per_surface;
+             ++iintegral)
+        {
             for (int isurface = 0; isurface < m_params.num_surfaces; ++isurface)
             {
                 int idx = isurface * num_integrals_per_surface + iintegral;
                 data_for_writing[idx] = a_integrals[iintegral][isurface];
             }
             if (extrapolation_done)
             {
-                int idx =
-                    m_params.num_surfaces * num_integrals_per_surface + iintegral;
+                int idx = m_params.num_surfaces * num_integrals_per_surface +
+                          iintegral;
                 data_for_writing[idx] = extrapolations[iintegral];
             }
         }
@@ -535,8 +539,6 @@ SurfaceExtraction<SurfaceGeometry>::richardson_extrapolation(
     const std::vector<std::vector<double>> &integrals) const
 {
     CH_TIME("SurfaceExtraction::richardson_extrapolation");
-    int num_comps = integrals.size();
-    int num_surfaces = m_params.num_surfaces;
 
     // already validated the radii
     int extrapolation_order = m_params.radii_idxs_for_extrapolation.size();
@@ -546,6 +548,7 @@ SurfaceExtraction<SurfaceGeometry>::richardson_extrapolation(
         return std::vector<double>();
     }
 
+    int num_comps = integrals.size();
     std::vector<double> extrapolations(num_comps);
 
     if (extrapolation_order == 3)

Source/CCZ4/ADMQuantities.hpp

@@ -38,8 +38,8 @@ class ADMQuantities
 
     ADMQuantities(const std::array<double, CH_SPACEDIM> &a_center, double a_dx,
                   int a_c_Madm = -1, int a_c_Jadm = -1, double a_G_Newton = 1.0)
-        : m_deriv(a_dx), m_center(a_center), m_G_Newton(a_G_Newton), m_dir(Z),
-          m_c_Madm(a_c_Madm), m_c_Jadm(a_c_Jadm)
+        : m_deriv(a_dx), m_center(a_center), m_G_Newton(a_G_Newton),
+          m_c_Madm(a_c_Madm), m_c_Jadm(a_c_Jadm), m_dir(Z)
     {
     }
 
@@ -59,10 +59,7 @@ class ADMQuantities
 
         // Surface element for integration
         Coordinates<data_t> coords(current_cell, m_deriv.m_dx, m_center);
-        data_t r = coords.get_radius();
         Tensor<1, data_t> x = {coords.x, coords.y, coords.z};
-        // This is multiplied by r^2 as SphericalExtraction assumes it is
-        // normalised as such.
         Tensor<1, data_t> dS_U = x;
 
         data_t dS_norm = 0.;