handle vector groups of vectors and higher order tensorsHEADsvnmaster

this patch adds support for groups of the form

CCTK_REAL vvel[3] "group of vectors"
{
  vvelx, vvely, vvelz
}

it changes how the vector/tensor component is constructed from the variable
index. It assumes that Cactus order variables as
vvelx[0], vvelx[1], vvelx[2], vvely[0], ...


git-svn-id: http://svn.cactuscode.org/arrangements/CactusNumerical/ReflectionSymmetry/trunk@61 082bdb00-0f4f-0410-b49e-b1835e5f2039

2 files changed, 30 insertions, 22 deletions

diff --git a/src/apply.c b/src/apply.c
index 0c09b8e..9c39379 100644
--- a/src/apply.c
+++ b/src/apply.c
@@ -302,19 +302,23 @@ BndReflectVI (cGH const * restrict const cctkGH,
     assert (numvars % numcomps == 0 && 
             (numvars == numcomps * vectorlength || numvars == vectorlength));
     ttype = VECTOR;
-    tcomponent = (vi - firstvar) % numcomps;
+    if(numvars == vectorlength) {
+      tcomponent = (vi - firstvar);
+    } else {
+      tcomponent = (vi - firstvar) / vectorlength;
+    }
   } else if (CCTK_EQUALS (tensortypealias, "4u")
              || CCTK_EQUALS (tensortypealias, "4d"))
   {
     /* 4-vector */
     const int numcomps = 4;
     assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
-    if ((vi - firstvar) % numcomps == 0) {
+    if ((vi - firstvar) / vectorlength == 0) {
       ttype = SCALAR;
       tcomponent = 0;
     } else {
       ttype = VECTOR;
-      tcomponent = (vi - firstvar) % numcomps - 1;
+      tcomponent = (vi - firstvar) / vectorlength - 1;
     }
   } else if (CCTK_EQUALS (tensortypealias, "uu_sym")
              || CCTK_EQUALS (tensortypealias, "dd_sym"))
@@ -323,7 +327,7 @@ BndReflectVI (cGH const * restrict const cctkGH,
     const int numcomps = 6;
     assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
     ttype = SYMTENSOR;
-    tcomponent = (vi - firstvar) % numcomps;
+    tcomponent = (vi - firstvar) / vectorlength;
   } else if (CCTK_EQUALS (tensortypealias, "uu")
              || CCTK_EQUALS (tensortypealias, "ud")
              || CCTK_EQUALS (tensortypealias, "du")
@@ -333,29 +337,29 @@ BndReflectVI (cGH const * restrict const cctkGH,
     const int numcomps = 9;
     assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
     ttype = TENSOR;
-    tcomponent = (vi - firstvar) % numcomps;
+    tcomponent = (vi - firstvar) / vectorlength;
   } else if (CCTK_EQUALS (tensortypealias, "4uu_sym")
              || CCTK_EQUALS (tensortypealias, "4dd_sym"))
   {
     /* symmetric 4-tensor */
     const int numcomps = 10;
     assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
-    if ((vi - firstvar) % numcomps == 0) {
+    if ((vi - firstvar) / vectorlength == 0) {
       ttype = SCALAR;
       tcomponent = 0;
-    } else if ((vi - firstvar) % numcomps <= 3) {
+    } else if ((vi - firstvar) / vectorlength <= 3) {
       ttype = VECTOR;
-      tcomponent = (vi - firstvar) % numcomps - 1;
+      tcomponent = (vi - firstvar) / vectorlength - 1;
     } else {
       ttype = SYMTENSOR;
-      tcomponent = (vi - firstvar) % numcomps - 4;
+      tcomponent = (vi - firstvar) / vectorlength - 4;
     }
   } else if (CCTK_EQUALS (tensortypealias, "ddd_sym")) {
     /* 3rd rank tensor, symmetric in last 2 indices */
     const int numcomps = 18;
     assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
     ttype = SYMTENSOR3;
-    tcomponent = (vi - firstvar) % numcomps;
+    tcomponent = (vi - firstvar) / vectorlength;
   } else if (CCTK_EQUALS (tensortypealias, "weylscalars_real")) {
     /* Weyl scalars, stored as 10 real values.  NOTE: This assumes
        that Psi_0 comes first, which is NOT the default with
@@ -363,7 +367,7 @@ BndReflectVI (cGH const * restrict const cctkGH,
     const int numcomps = 10;
     assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
     ttype = WEYLSCALARS_REAL;
-    tcomponent = (vi - firstvar) % numcomps;
+    tcomponent = (vi - firstvar) / vectorlength;
   } else if (CCTK_EQUALS (tensortypealias, "ManualCartesian")) {
       /* Reflection symmetries specified by hand */
       ttype = MANUALCARTESIAN;
diff --git a/src/interpolate.c b/src/interpolate.c
index 571c048..bed4cae 100644
--- a/src/interpolate.c
+++ b/src/interpolate.c
@@ -244,19 +244,23 @@ ReflectionSymmetry_Interpolate (CCTK_POINTER_TO_CONST restrict const cctkGH,
         assert (numvars % numcomps == 0 && 
                 (numvars == numcomps * vectorlength || numvars == vectorlength));
         ttype = VECTOR;
-        tcomponent = (vi - firstvar) % numcomps;
+        if(numvars == vectorlength) {
+          tcomponent = (vi - firstvar);
+        } else {
+          tcomponent = (vi - firstvar) / vectorlength;
+        }
       } else if (CCTK_EQUALS (tensortypealias, "4u")
                  || CCTK_EQUALS (tensortypealias, "4d"))
       {
         /* 4-vector */
         const int numcomps = 4;
         assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
-        if ((vi - firstvar) % numcomps == 0) {
+        if ((vi - firstvar) / vectorlength == 0) {
           ttype = SCALAR;
           tcomponent = 0;
         } else {
           ttype = VECTOR;
-          tcomponent = (vi - firstvar) % numcomps - 1;
+          tcomponent = (vi - firstvar) / vectorlength - 1;
         }
       } else if (CCTK_EQUALS (tensortypealias, "uu_sym")
                  || CCTK_EQUALS (tensortypealias, "dd_sym"))
@@ -265,7 +269,7 @@ ReflectionSymmetry_Interpolate (CCTK_POINTER_TO_CONST restrict const cctkGH,
         const int numcomps = 6;
         assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
         ttype = SYMTENSOR;
-        tcomponent = (vi - firstvar) % numcomps;
+        tcomponent = (vi - firstvar) / vectorlength;
       } else if (CCTK_EQUALS (tensortypealias, "uu")
                  || CCTK_EQUALS (tensortypealias, "ud")
                  || CCTK_EQUALS (tensortypealias, "du")
@@ -275,35 +279,35 @@ ReflectionSymmetry_Interpolate (CCTK_POINTER_TO_CONST restrict const cctkGH,
         const int numcomps = 9;
         assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
         ttype = TENSOR;
-        tcomponent = (vi - firstvar) % numcomps;
+        tcomponent = (vi - firstvar) / vectorlength;
       } else if (CCTK_EQUALS (tensortypealias, "4uu_sym")
                  || CCTK_EQUALS (tensortypealias, "4dd_sym"))
       {
         /* symmetric 4-tensor */
         const int numcomps = 10;
         assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
-        if ((vi - firstvar) % numcomps == 0) {
+        if ((vi - firstvar) / vectorlength == 0) {
           ttype = SCALAR;
           tcomponent = 0;
-        } else if ((vi - firstvar) % numcomps <= 3) {
+        } else if ((vi - firstvar) / vectorlength <= 3) {
           ttype = VECTOR;
-          tcomponent = (vi - firstvar) % numcomps - 1;
+          tcomponent = (vi - firstvar) / vectorlength - 1;
         } else {
           ttype = SYMTENSOR;
-          tcomponent = (vi - firstvar) % numcomps - 4;
+          tcomponent = (vi - firstvar) / vectorlength - 4;
         }
       } else if (CCTK_EQUALS (tensortypealias, "ddd_sym")) {
         /* 3rd rank tensor, symmetric in last 2 indices */
         const int numcomps = 18;
         assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
         ttype = SYMTENSOR3;
-        tcomponent = (vi - firstvar) % numcomps;
+        tcomponent = (vi - firstvar) / vectorlength;
       } else if (CCTK_EQUALS (tensortypealias, "weylscalars_real")) {
         /* Weyl scalars, stored as 10 real numbers */
         const int numcomps = 10;
         assert (numvars % numcomps == 0 && numvars == numcomps * vectorlength);
         ttype = WEYLSCALARS_REAL;
-        tcomponent = (vi - firstvar) % numcomps;
+        tcomponent = (vi - firstvar) / vectorlength;
       } else if (CCTK_EQUALS (tensortypealias, "ManualCartesian")) {
         /* Reflection symmetries specified by hand */
         ttype = MANUALCARTESIAN;