CarpetX: add prolongation_type 'grmhd-auto'

EinsteinToolkit · Oct 6, 2024 · ae35b4a · ae35b4a
1 parent e20789f
commit ae35b4a
Show file tree

Hide file tree

Showing 2 changed files with 99 additions and 1 deletion.
diff --git a/CarpetX/param.ccl b/CarpetX/param.ccl
@@ -426,6 +426,7 @@ KEYWORD prolongation_type "Prolongation type"
   "hermite" :: "Hermite-interpolate in vertex centred and conserve in cell centred directions"
   "natural" :: "interpolate in vertex centred and conserve in cell centred directions, using the same order"
   "poly-cons3lfb" :: "interpolate polynomially in vertex centred directions and conserve with 3rd order accuracy and a linear fallback in cell centred directions"
+  "grmhd-auto" :: "choose different prolongation operators for different staggered variables"
 } "natural"
 
 CCTK_INT prolongation_order "Prolongation order"

diff --git a/CarpetX/src/driver.cxx b/CarpetX/src/driver.cxx
@@ -542,7 +542,7 @@ amrex::Interpolater *get_interpolator(const std::array<int, dim> indextype) {
     hermite,
     natural,
     poly_cons3lfb,
-    poly_eno3lfb,
+    grmhd_auto,
   };
   static interp_t interp = [&]() {
     if (CCTK_EQUALS(prolongation_type, "interpolate"))
@@ -559,11 +559,108 @@ amrex::Interpolater *get_interpolator(const std::array<int, dim> indextype) {
       return interp_t::natural;
     else if (CCTK_EQUALS(prolongation_type, "poly-cons3lfb"))
       return interp_t::poly_cons3lfb;
+    else if (CCTK_EQUALS(prolongation_type, "grmhd-auto"))
+      return interp_t::grmhd_auto;
     else
       assert(0);
   }();
 
   switch (interp) {
+  case interp_t::grmhd_auto:
+    switch ((indextype[0] << 2) | (indextype[1] << 1) | (indextype[2] << 0)) {
+
+    // vertex-centered
+    case 0b000:
+      switch (prolongation_order) {
+      case 1:
+        return &prolongate_poly_3d_rf2_c000_o1;
+      case 3:
+        return &prolongate_poly_3d_rf2_c000_o3;
+      case 5:
+        return &prolongate_poly_3d_rf2_c000_o5;
+      case 7:
+        return &prolongate_poly_3d_rf2_c000_o7;
+      }
+      break;
+
+    // edge-centered
+    case 0b001:
+      switch (prolongation_order) {
+      case 1:
+        return &prolongate_hermite_3d_rf2_c001_o1;
+      case 3:
+        return &prolongate_hermite_3d_rf2_c001_o3;
+      case 5:
+        return &prolongate_hermite_3d_rf2_c001_o5;
+      }
+      break;
+    case 0b010:
+      switch (prolongation_order) {
+      case 1:
+        return &prolongate_hermite_3d_rf2_c010_o1;
+      case 3:
+        return &prolongate_hermite_3d_rf2_c010_o3;
+      case 5:
+        return &prolongate_hermite_3d_rf2_c010_o5;
+      }
+      break;
+    case 0b100:
+      switch (prolongation_order) {
+      case 1:
+        return &prolongate_hermite_3d_rf2_c100_o1;
+      case 3:
+        return &prolongate_hermite_3d_rf2_c100_o3;
+      case 5:
+        return &prolongate_hermite_3d_rf2_c100_o5;
+      }
+      break;
+
+    // face-centered
+    case 0b011:
+      switch (prolongation_order) {
+      case 1:
+        return &prolongate_poly_cons3lfb_3d_rf2_c011_o1;
+      case 3:
+        return &prolongate_poly_cons3lfb_3d_rf2_c011_o3;
+      case 5:
+        return &prolongate_poly_cons3lfb_3d_rf2_c011_o5;
+      }
+      break;
+    case 0b101:
+      switch (prolongation_order) {
+      case 1:
+        return &prolongate_poly_cons3lfb_3d_rf2_c101_o1;
+      case 3:
+        return &prolongate_poly_cons3lfb_3d_rf2_c101_o3;
+      case 5:
+        return &prolongate_poly_cons3lfb_3d_rf2_c101_o5;
+      }
+      break;
+    case 0b110:
+      switch (prolongation_order) {
+      case 1:
+        return &prolongate_poly_cons3lfb_3d_rf2_c110_o1;
+      case 3:
+        return &prolongate_poly_cons3lfb_3d_rf2_c110_o3;
+      case 5:
+        return &prolongate_poly_cons3lfb_3d_rf2_c110_o5;
+      }
+      break;
+
+    // cell-centered
+    case 0b111:
+      switch (prolongation_order) {
+      case 1:
+        return &prolongate_eno_3d_rf2_c111_o1;
+      case 3:
+        return &prolongate_eno_3d_rf2_c111_o3;
+      case 5:
+        return &prolongate_eno_3d_rf2_c111_o5;
+      }
+      break;
+    }
+    break;
+
   case interp_t::interpolate:
 
     switch ((indextype[0] << 2) | (indextype[1] << 1) | (indextype[2] << 0)) {