From 5d2004134817f2669f7e9b513b435f8ee0010dd5 Mon Sep 17 00:00:00 2001 From: rhaas Date: Fri, 14 Jun 2013 05:56:55 +0000 Subject: handle vector groups of vectors and higher order tensors 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 --- src/apply.c | 26 +++++++++++++++----------- src/interpolate.c | 26 +++++++++++++++----------- 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; -- cgit v1.2.3