Add average atomic mass: HydroBase::Abar. Ticket #1649.svn

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

5 files changed, 48 insertions, 1 deletions

diff --git a/interface.ccl b/interface.ccl
index fa4fcd3..9e176f1 100755
--- a/interface.ccl
+++ b/interface.ccl
@@ -38,6 +38,7 @@ CCTK_REAL vel[3]  type = GF Timelevels = 3 tags='ProlongationParameter="HydroBas
 CCTK_REAL w_lorentz  type = GF Timelevels = 3 tags='ProlongationParameter="HydroBase::prolongation_type" tensortypealias="Scalar" interpolator="matter"' "Lorentz Factor"
 
 CCTK_REAL Y_e         type = GF Timelevels = 3 tags='ProlongationParameter="HydroBase::prolongation_type" tensortypealias="Scalar" interpolator="matter"' "Electron Fraction"
+CCTK_REAL Abar        type = GF Timelevels = 3 tags='ProlongationParameter="HydroBase::prolongation_type" tensortypealias="Scalar" interpolator="matter"' "Average atomic mass [atomic mass unit]"
 CCTK_REAL temperature type = GF Timelevels = 3 tags='ProlongationParameter="HydroBase::prolongation_type" tensortypealias="Scalar" interpolator="matter"' "Temperature [MeV]"
 CCTK_REAL entropy     type = GF Timelevels = 3 tags='ProlongationParameter="HydroBase::prolongation_type" tensortypealias="Scalar" interpolator="matter"' "Specific Entropy [k_b/baryon]"
 
diff --git a/param.ccl b/param.ccl
index e1001c5..dc1c4ce 100755
--- a/param.ccl
+++ b/param.ccl
@@ -42,11 +42,23 @@ KEYWORD initial_Y_e "Initial value for Y_e"
   "read from file" :: "Read the initial data using the IOUtil file reader.  Note that this only allows you to read the data from a file, it does not actually do it.  You still have to programme the IOUtil file reader accordingly."
 } "none"
 
+KEYWORD initial_Abar "Initial value for Abar"
+{
+  "none" :: "inactive"
+  "zero" :: "initially set to zero"
+  "read from file" :: "Read the initial data using the IOUtil file reader.  Note that this only allows you to read the data from a file, it does not actually do it.  You still have to programme the IOUtil file reader accordingly."
+} "none"
+
 KEYWORD Y_e_evolution_method "Evolution method for Y_e"
 {
   "none" :: "Evolution for Y_e is disabled"
 } "none"
 
+KEYWORD Abar_evolution_method "Evolution method for Abar"
+{
+  "none" :: "Evolution for Abar is disabled"
+} "none"
+
 KEYWORD temperature_evolution_method "Evolution method for temperature"
 {
   "none" :: "Evolution for temperature is disabled"
diff --git a/schedule.ccl b/schedule.ccl
index 7cf56a3..3a296fe 100755
--- a/schedule.ccl
+++ b/schedule.ccl
@@ -11,6 +11,10 @@ if (timelevels == 1)
   {
     STORAGE: Y_e[1]
   }
+  if (!CCTK_EQUALS(initial_Abar, "none"))
+  {
+    STORAGE: Abar[1]
+  }
   if (!CCTK_EQUALS(initial_Bvec, "none"))
   {
     STORAGE: Bvec[1]
@@ -43,6 +47,10 @@ else if (timelevels == 2)
   {
     STORAGE: Y_e[2]
   }
+  if (!CCTK_EQUALS(initial_Abar, "none"))
+  {
+    STORAGE: Abar[2]
+  }
   if (!CCTK_EQUALS(initial_Bvec, "none"))
   {
     STORAGE: Bvec[2]
@@ -75,6 +83,10 @@ else if (timelevels == 3)
   {
     STORAGE: Y_e[3]
   }
+  if (!CCTK_EQUALS(initial_Abar, "none"))
+  {
+    STORAGE: Abar[3]
+  }
   if (!CCTK_EQUALS(initial_Bvec, "none"))
   {
     STORAGE: Bvec[3]
diff --git a/src/Initialisation.c b/src/Initialisation.c
index 9aa31ec..762b402 100644
--- a/src/Initialisation.c
+++ b/src/Initialisation.c
@@ -20,7 +20,14 @@ void HydroBase_Zero (CCTK_ARGUMENTS)
     w_lorentz[i] = 1.0;
     eps[i     ]  = 0.0;
   }
-  
+
+  if (CCTK_ActiveTimeLevels(cctkGH, "HydroBase::Abar") >= 1) {
+#pragma omp parallel for
+    for (int i=0; i<np; ++i) {
+      Abar[i     ] = 0.0;
+    }
+  }
+
   if (CCTK_ActiveTimeLevels(cctkGH, "HydroBase::temperature") >= 1) {
 #pragma omp parallel for
     for (int i=0; i<np; ++i) {
@@ -55,6 +62,13 @@ void HydroBase_Zero (CCTK_ARGUMENTS)
       }
     }
 
+    if (CCTK_ActiveTimeLevels(cctkGH, "HydroBase::Abar") >= 2) {
+#pragma omp parallel for
+      for (int i=0; i<np; ++i) {
+        Abar_p[i     ] = 0.0;
+      }
+    }
+
     if (CCTK_ActiveTimeLevels(cctkGH, "HydroBase::temperature") >= 2) {
 #pragma omp parallel for
       for (int i=0; i<np; ++i) {
@@ -86,6 +100,13 @@ void HydroBase_Zero (CCTK_ARGUMENTS)
                  "Too many active time levels for HydroBase variables");
     }
 
+    if (CCTK_ActiveTimeLevels(cctkGH, "HydroBase::Abar") >= 3) {
+#pragma omp parallel for
+      for (int i=0; i<np; ++i) {
+        Abar_p_p[i     ] = 0.0;
+      }
+    }
+
     if (CCTK_ActiveTimeLevels(cctkGH, "HydroBase::temperature") >= 3) {
 #pragma omp parallel for
       for (int i=0; i<np; ++i) {
diff --git a/src/ParamCheck.c b/src/ParamCheck.c
index 7e3f77e..6a83dac 100644
--- a/src/ParamCheck.c
+++ b/src/ParamCheck.c
@@ -35,6 +35,7 @@ void HydroBase_ParamCheck (CCTK_ARGUMENTS)
     {"initial_Bvec", initial_Bvec, "HydroBase::Bvec[1]"},
     {"initial_Bvec", initial_Bvec, "HydroBase::Bvec[2]"},
     {"initial_Y_e", initial_Y_e, "HydroBase::Y_e"},
+    {"initial_Abar", initial_Abar, "HydroBase::Abar"},
     {"initial_temperature", initial_temperature, "HydroBase::temperature"},
     {"initial_entropy", initial_entropy, "HydroBase::entropy"},
   };