Add a new group "Bvec" to HydroBase, containing the magnetic field.

Same as the velocity "vel", this is a vector group with three
components.

There are also two new parameters "initial_Bvec" with default value
"none", specifying that Bvec does not have storage.  This keeps things
backward compatible.  Another parameter "Bvec_evolution_method"
specifies its method; the default value is also "none".  These
parameters need to be extended by initial data and evolution thorns
that handle magnetic fields.


git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinBase/HydroBase/trunk@18 57fe0bb3-ccba-405f-9b23-de0201f165b7

4 files changed, 89 insertions, 1 deletions

diff --git a/interface.ccl b/interface.ccl
index 164a8ff..2806eac 100755
--- a/interface.ccl
+++ b/interface.ccl
@@ -16,3 +16,4 @@ CCTK_REAL vel[3]  type = GF Timelevels = 3 tags='ProlongationParameter="HydroBas
 
 CCTK_REAL Y_e     type = GF Timelevels = 3 tags='ProlongationParameter="HydroBase::prolongation_type" tensortypealias="Scalar" interpolator="matter"' "electron fraction"
 
+CCTK_REAL Bvec[3] type = GF Timelevels = 3 tags='ProlongationParameter="HydroBase::prolongation_type" tensortypealias="U" interpolator="matter"' "Magnetic field components B^i"
diff --git a/param.ccl b/param.ccl
index 5aca620..50e4cc7 100755
--- a/param.ccl
+++ b/param.ccl
@@ -41,3 +41,15 @@ KEYWORD Y_e_evolution_method "Evolution method for Y_e"
   "none" :: "Evolution for Y_e is disabled"
 } "none"
 
+# Settings for magnetic field B^i
+
+KEYWORD initial_Bvec "Initial value for Bvec"
+{
+  "none" :: "inactive"
+  "zero" :: "initially set to zero"
+} "none"
+
+KEYWORD Bvec_evolution_method "Evolution method for Bvec"
+{
+  "none" :: "Evolution for Bvec is disabled"
+} "none"
diff --git a/schedule.ccl b/schedule.ccl
index 57de555..9d9a956 100755
--- a/schedule.ccl
+++ b/schedule.ccl
@@ -10,6 +10,10 @@ if (timelevels == 1)
   {
     STORAGE: Y_e[1]
   }
+  if (!CCTK_EQUALS(initial_Bvec, "none"))
+  {
+    STORAGE: Bvec[1]
+  }
 }
 else if (timelevels == 2)
 {
@@ -21,6 +25,10 @@ else if (timelevels == 2)
   {
     STORAGE: Y_e[2]
   }
+  if (!CCTK_EQUALS(initial_Bvec, "none"))
+  {
+    STORAGE: Bvec[2]
+  }
 }
 else if (timelevels == 3)
 {
@@ -32,6 +40,10 @@ else if (timelevels == 3)
   {
     STORAGE: Y_e[3]
   }
+  if (!CCTK_EQUALS(initial_Bvec, "none"))
+  {
+    STORAGE: Bvec[3]
+  }
 }
 
 schedule group HydroBase_Initial at CCTK_INITIAL    \
@@ -107,7 +119,7 @@ if (CCTK_EQUALS(initial_hydro, "zero"))
   schedule HydroBase_Zero in HydroBase_Initial
   {
     LANG: C
-  } "Set up zero (vacuum without atmosphere) hydro initial data"
+  } "Set up vacuum hydro initial data"
 }
 
 if (CCTK_Equals(initial_Y_e, "one"))
@@ -117,3 +129,11 @@ if (CCTK_Equals(initial_Y_e, "one"))
     LANG: C
   } "Set electron fraction to 1"
 }
+
+if (CCTK_Equals(initial_Bvec, "zero"))
+{
+  SCHEDULE HydroBase_Bvec_zero in HydroBase_Initial
+  {
+    LANG: C
+  } "Set magnetic field to 0"
+}
diff --git a/src/Initialisation.c b/src/Initialisation.c
index a739091..fcd71f5 100644
--- a/src/Initialisation.c
+++ b/src/Initialisation.c
@@ -111,3 +111,58 @@ void HydroBase_Y_e_one (CCTK_ARGUMENTS)
                "Unsupported parameter value for InitBase::initial_data_setup_method");
   }
 }
+
+
+
+void HydroBase_Bvec_zero (CCTK_ARGUMENTS)
+{
+  DECLARE_CCTK_ARGUMENTS;
+  DECLARE_CCTK_PARAMETERS;
+  
+  int const np = cctk_lsh[0] * cctk_lsh[1] * cctk_lsh[2];
+  
+#pragma omp parallel for
+  for (int i=0; i<np; ++i) {
+    Bvec[i]      = 0.0;
+    Bvec[i+  np] = 0.0;
+    Bvec[i+2*np] = 0.0;
+  }
+  
+  if (CCTK_EQUALS (initial_data_setup_method, "init_some_levels") ||
+      CCTK_EQUALS (initial_data_setup_method, "init_single_levels"))
+  {
+    /* do nothing */
+  }
+  else if (CCTK_EQUALS (initial_data_setup_method, "init_all_levels"))
+  {
+    
+    if (CCTK_ActiveTimeLevels(cctkGH, "HydroBase::Bvec") >= 2) {
+#pragma omp parallel for
+      for (int i=0; i<np; ++i) {
+        Bvec_p[i]      = 0.0;
+        Bvec_p[i+  np] = 0.0;
+        Bvec_p[i+2*np] = 0.0;
+      }
+    }
+    
+    if (CCTK_ActiveTimeLevels(cctkGH, "HydroBase::Bvec") >= 3) {
+#pragma omp parallel for
+      for (int i=0; i<np; ++i) {
+        Bvec_p_p[i]      = 0.0;
+        Bvec_p_p[i+  np] = 0.0;
+        Bvec_p_p[i+2*np] = 0.0;
+      }
+    }
+    
+    if (CCTK_ActiveTimeLevels(cctkGH, "HydroBase::Bvec") >= 4) {
+      CCTK_WARN (CCTK_WARN_ABORT,
+                 "Too many active time levels for HydroBase::Bvec");
+    }
+    
+  }
+  else
+  {
+    CCTK_WARN (CCTK_WARN_ABORT,
+               "Unsupported parameter value for InitBase::initial_data_setup_method");
+  }
+}