symmetry conditions for quadrant/bitant and 1D/2D support

git-svn-id: http://svn.cactuscode.org/arrangements/CactusBase/CartGrid3D/trunk@98 c78560ca-4b45-4335-b268-5f3340f3cb52

1 files changed, 92 insertions, 361 deletions

diff --git a/src/Symmetry.c b/src/Symmetry.c
index e833eaf..0b0209c 100644
--- a/src/Symmetry.c
+++ b/src/Symmetry.c
@@ -1,395 +1,126 @@
-/*@@
-   @file      Symmetry.c
-   @date      Mon Mar 15 15:09:00 1999
-   @author    Gerd Lanfermann
-   @desc 
-     This file contains the routines for registering and applying symmetry 
-     boundary conditions
-   @enddesc 
- @@*/
 
+
+#include <stdio.h>
 #include <assert.h>
 #include <stdlib.h>
-#include <stdio.h>
+#include <string.h>
 
 #include "cctk.h"
-#include "cctk_Parameters.h"
-#include "Symmetry.h" 
-#include "cctk_FortranString.h"
-
-/*#define DEBUG_BOUND*/
-
- /*@@
-   @routine    SetCartSymmetry
-   @date       Mon Mar 15 15:10:58 1999
-   @author     Gerd Lanfermann
-   @desc 
-   This routine sets the GH extension (EinsteinBoundGHex *bGHex), 
-   which describes the symmetry  boundary type of each GF. Takes 
-   the name of the GF ("implementation::gfname") and the symmetry operators 
-   sx,sy,sz and inserts them in the array bGHex. 
-   These values will looked up by ApplySym
-   @enddesc 
-   @calls     
-   @calledby   
-   @history 
- 
-   @endhistory 
-
-@@*/
-
-void SetCartSymmetry(cGH *GH, int *sym, const char *imp_gf) {
-  
-  DECLARE_CCTK_PARAMETERS
-
-  SymmetryGHex *sGHex; 
-  int index;
 
-  /* Pointer to the SymmetryGHextension */
-  sGHex  = (SymmetryGHex *)GH->extensions[CCTK_GHExtensionHandle("Symmetry")];
-
-  /* now that we have the same, get the index to the GFs */
-  index = CCTK_VarIndex(imp_gf);
-  if (index<0) 
-  {
-    char *message;
-    message = (char *)malloc((100*strlen(imp_gf))*sizeof(char));
-    sprintf(message,"Grid function %s has no index",imp_gf);
-    CCTK_WARN(0,message);
-    free(message);
-  }
-  
-  /* Reference the hash table in the GHex and tell it what kind of
-     symmetry is being applied 
-     (depending on sym and the grid layout) 
-     If there is no symmetry necessary,set ESYM_NOSYM
-     When we apply a symmetry and find ESYM_UNSET, something went wrong! 
-   */
-
-  if (CCTK_Equals(domain,"full")) 
-  {
-
-#ifdef DEBUG_BOUND
-    printf(" Registered full grid symmetries for -%s- in SetCartSymmetry\n",imp_gf);
-#endif
+/*#define SYM_DEBUG*/
 
-    sGHex->GFSym[index][0] = GFSYM_NOSYM;
-    sGHex->GFSym[index][2] = GFSYM_NOSYM;
-    sGHex->GFSym[index][4] = GFSYM_NOSYM;
+/* index convention:
+   i ~ x ~ 0
+   j ~ y ~ 1
+   k ~ z ~ 2 */
+   
+int CartApplySym3Di(cGH *GH, int gdim, int *doSym, 
+		   int *lssh, int *ghostz, int *sym, CCTK_REAL *var) 
+{
 
-  } 
-  else if (CCTK_Equals(domain,"octant")) 
-  {
+  int i,j,k;
 
-#ifdef DEBUG_BOUND
-    printf(" Registered octant symmetries for -%s- in SetCartSymmetry\n",imp_gf);
+#ifdef SYM_DEBUG
+  printf(" doSym: %d %d /  %d %d /  %d %d \n",
+	 doSym[0],doSym[1],
+	 doSym[2],doSym[3],
+	 doSym[4],doSym[5]);
+  printf(" lssh: %d %d %d sym: %d %d %d \n", 
+	 lssh[0],lssh[1],lssh[2], sym[0], sym[2], sym[4] );
+  printf(" ghostz %d %d %d \n",ghostz[0],ghostz[1],ghostz[2]);
 #endif
 
-    sGHex->GFSym[index][0] = sym[0];
-    sGHex->GFSym[index][2] = sym[1];
-    sGHex->GFSym[index][4] = sym[2];
-
-  } 
-  else if (CCTK_Equals(domain,"quadrant")) 
-  {
-
-#ifdef DEBUG_BOUND
-    printf("Registered quadrant symmetries for -%s- in SetCartSymmetry\n",imp_gf);
-#endif    
-
-    sGHex->GFSym[index][0] = sym[0];
-    sGHex->GFSym[index][2] = sym[1];
-    sGHex->GFSym[index][4] = GFSYM_NOSYM;
-
-  } 
-  else if (CCTK_Equals(domain,"bitant")) 
+  if (doSym[0] == 1)
   {
-
-#ifdef DEBUG_BOUND
-    printf("Registered bitant symmetries for -%s- in SetCartSymmetry\n",imp_gf);
-#endif    
-
-    sGHex->GFSym[index][4] = sym[2];
-    sGHex->GFSym[index][0] = GFSYM_NOSYM;
-    sGHex->GFSym[index][2] = GFSYM_NOSYM;
+    for(k=0; k < lssh[2]; k++)
+    {
+      for(j=0; j < lssh[1]; j++)
+      {
+        for(i=0; i < ghostz[0]; i++)
+        {
+          var[CCTK_GFINDEX3D(GH,i,j,k)] = 
+	    sym[0]*var[CCTK_GFINDEX3D(GH,2*ghostz[0]-1-i,j,k)];
+        }
+      }
+    }
   }
-
-  /* All untouched GFSym[][] will hold GFSYM_UNSET */
-
-  
-
-}
-
-void FMODIFIER FORTRAN_NAME(SetCartSymmetry)(cGH *GH,  int *sym, ONE_FORTSTRING_ARG) {
-
-  ONE_FORTSTRING_CREATE(imp_gf) 
-  SetCartSymmetry(GH, sym, imp_gf);
-  free(imp_gf);
-} 
-
-
-
-
- /*@@
-   @routine    ApplySymmetry
-   @date       Mon Mar 15 15:16:28 1999
-   @author     Gerd Lanfermann
-   @desc 
-      Takes the GH and the name of the group.
-
-      Routine applies the symmetry BC to the GFs of a  group as set by SetSym 
-      in GHExtension.  It is called by C or Fortran (via a wrapper) by 
-      passing the group name. The assignment of the symmetries is carried out
-      by the F routine FortranSym.
-   @enddesc 
-   @calls   FortranSym  
-   @calledby   
-   @history 
- 
-   @endhistory 
-
-@@*/
-
-
-
-
-int CartSymBCGroupI(cGH *GH, int groupnum) {
-
-  void FORTRAN_NAME(SymmetryCondition)(int *, CCTK_REAL *, int *, int *, int *);
-  
-  DECLARE_CCTK_PARAMETERS
-
-  SymmetryGHex *sGHex; 
-  int first,last,index,j;
-  int time;
-  int doSym[6];
-
-  /* Get out if we are sure no symmetries should be applied */
-  if (CCTK_Equals(domain,"full")) return 0;
-
-  /* Get the pointer to the Symmetry GH extension */
-  sGHex  = (SymmetryGHex*)GH->extensions[CCTK_GHExtensionHandle("Symmetry")];
-
-
-  if (groupnum < 0) 
+ if (doSym[2] == 1)
   {
-    char *message=NULL;
-    message = (char *)malloc(300*sizeof(char));
-    sprintf(message,"Invalid group number %d",groupnum);
-    CCTK_WARN(0,message);
-    free(message);
+    for(i=0; i < lssh[0]; i++)
+      {
+	for(k=0; k < lssh[2]; k++)
+	  {
+	    for(j=0; j < ghostz[1]; j++)
+	      {
+		var[CCTK_GFINDEX3D(GH,i,j,k)] = 
+		  sym[2]*var[CCTK_GFINDEX3D(GH,i,2*ghostz[1]-1-j,k)];
+	      }
+	  }
+      }
   }
-
-  /*get the index of the first GF in the group and how many Vars there are*/
-  first = CCTK_FirstVarIndexI(groupnum);
-  last  = first+CCTK_NumVarsInGroupI(groupnum)-1;
-  assert (first>=0 && last>=first);
-
-  /* Use next time level, if present */
-  time = CCTK_NumTimeLevelsFromVarI(first) - 1;
-  if (time < 0)
+ if (doSym[4] == 1)
   {
-    time = 0;
+    for(i=0; i < lssh[0]; i++)
+      {
+	for(j=0; j < lssh[1]; j++)
+	  {
+	    for(k=0; k < ghostz[2]; k++)
+	      {
+		var[CCTK_GFINDEX3D(GH,i,j,k)] = 
+		  sym[4]*var[CCTK_GFINDEX3D(GH,i,j,2*ghostz[2]-1-k)];
+	      }
+	  }
+      }
   }
+ return(0);
+}
+   
+int CartApplySym2Di(cGH *GH, int gdim, int *doSym, 
+		   int *lssh, int *ghostz, int *sym, CCTK_REAL *var) 
+{
+  int i,j;
 
-  /* loop over the variables in the group */
-  for (index=first; index<=last; index++) {
-
-    /* and check that we actually have a grid function (and not a scalar)*/
-    if (CCTK_GroupTypeFromVarI(index) == GROUP_GF) {
-
-      /*at this point, there should be NO ESYM_UNSET anymore 
-        if there is, we forgot to register the symmetries for a GF */
-      if ((sGHex->GFSym[index][0]==GFSYM_UNSET)||
-	  (sGHex->GFSym[index][2]==GFSYM_UNSET)||
-	  (sGHex->GFSym[index][4]==GFSYM_UNSET)) {
-	char *message=NULL;
-	message = (char *)malloc(300*sizeof(char));
-	sprintf(message,"Cannot apply symmetry to >%s< (# %d) without registered symmetries",
-		CCTK_GroupNameFromVarI(index),index);
-	CCTK_WARN(1,message);
-	free(message);
-      }
-      
-      /* whether we want to apply sym depends on several things: */
-      for (j=0;j<3;j++) 
+  if (doSym[0] == 1)
+  {
+    for(j=0; j < lssh[1]; j++)
+    {
+      for(i=0; i < ghostz[0]; i++)
       {
-	doSym[2*j+1]=0;
-	if ( (GH->cctk_lsh[j]>1) && (GH->cctk_bbox[2*j]==1) &&
-	    (sGHex->GFSym[index][2*j] != GFSYM_UNSET) &&
-	    (sGHex->GFSym[index][2*j] != GFSYM_NOSYM) 
-	     )
-	  doSym[2*j] = 1;
-	else 
-	  doSym[2*j] = 0;
+	var[CCTK_GFINDEX2D(GH,i,j)] = 
+	  sym[0]*var[CCTK_GFINDEX2D(GH,2*ghostz[0]-1-i,j)];
       }
-
-      /* Call the Fortran Symmetry Routine */
-      FORTRAN_NAME(SymmetryCondition)(
-	     GH->cctk_lsh,            
-	     /* xyz-size of PE local grid */
-	     GH->data[index][time],         
-	     /* pointer to start of data array for GF[index]*/
-	     GH->cctk_nghostzones,            
-	     /* number of ghost zones */
-	     sGHex->GFSym[index],  
-	     /* the symmetries for this GF */
-	     doSym
-	     /* flags whether to apply syms */
-	     );                       
     }
   }
-
-  return 0;
-
-  
-
-}
-
-
-
-
-
-int CartSymBCVarI(cGH *GH, int varnum) {
-
-  void FORTRAN_NAME(SymmetryCondition)(int *, CCTK_REAL *, int *, int *, int *);
   
-  DECLARE_CCTK_PARAMETERS
-
-  SymmetryGHex *sGHex; 
-  int j;
-  int time;
-  int doSym[6];
-
-  /* Get out if we are sure no symmetries should be applied */
-  if (CCTK_Equals(domain,"full")) return 0;
-
-  assert (varnum>=0);
-
-  /* Use next time level, if present */
-  time = CCTK_NumTimeLevelsFromVarI(varnum) - 1;
-  if (time < 0)
+  if (doSym[2] == 1)
+  {
+    for(i=0; i < lssh[0]; i++)
     {
-      time = 0;
-    }
-
-  /* Get the pointer to the Symmetry GH extension */
-  sGHex  = (SymmetryGHex*)GH->extensions[CCTK_GHExtensionHandle("Symmetry")];
-
-    /* and check that we actually have a grid function (and not a scalar)*/
-    if (CCTK_GroupTypeFromVarI(varnum) == GROUP_GF) {
-
-
-
-      /*at this point, there should be NO ESYM_UNSET anymore 
-        if there is, we forgot to register the symmetries for a GF */
-      if ((sGHex->GFSym[varnum][0]==GFSYM_UNSET)||
-	  (sGHex->GFSym[varnum][2]==GFSYM_UNSET)||
-	  (sGHex->GFSym[varnum][4]==GFSYM_UNSET)) {
-	char *message=NULL;
-	message = (char *)malloc(300*sizeof(char));
-	sprintf(message,"Cannot apply symmetry to >%s< (# %d) without registered symmetries",
-		CCTK_GroupNameFromVarI(varnum),varnum);
-	CCTK_WARN(1,message);
-	free(message);
-      }
-      
-      /* whether we want to apply sym depends on several things: */
-      for (j=0;j<3;j++) 
+      for(j=0; j < ghostz[1]; j++)
       {
-	doSym[2*j+1]=0;
-	if ( (GH->cctk_lsh[j]>1) && (GH->cctk_bbox[2*j]==1) &&
-	    (sGHex->GFSym[varnum][2*j] != GFSYM_UNSET) &&
-	    (sGHex->GFSym[varnum][2*j] != GFSYM_NOSYM) 
-	     )
-	  doSym[2*j] = 1;
-	else 
-	  doSym[2*j] = 0;
+	var[CCTK_GFINDEX2D(GH,i,j)] = 
+	  sym[2]*var[CCTK_GFINDEX2D(GH,i,2*ghostz[1]-1-j)];
       }
-
-      /* Call the Fortran Symmetry Routine */
-      FORTRAN_NAME(SymmetryCondition)(
-	     GH->cctk_lsh,            
-	     /* xyz-size of PE local grid */
-	     GH->data[varnum][time],         
-	     /* pointer to start of data array for GF[varnum]*/
-	     GH->cctk_nghostzones,            
-	     /* number of ghost zones */
-	     sGHex->GFSym[varnum],  
-	     /* the symmetries for this GF */
-	     doSym
-	     /* flags whether to apply syms */
-	     );                       
     }
-
-    return 0;
-
-    
-
-}
-
-
-int CartSymBCVar(cGH *GH, char *name) 
-{
-
-  int index;
-  index = CCTK_VarIndex(name);
-  CartSymBCVarI(GH,index);
-  return 0;
-
+  }
+  
+  return(0);
 }
 
-int CartSymBCGroup(cGH *GH, char *name) 
+int CartApplySym1Di(cGH *GH, int gdim, int *doSym, 
+		   int *lssh, int *ghostz, int *sym, CCTK_REAL *var) 
 {
+  int i;
 
-  int index;
-  int retval;
-  index = CCTK_GroupIndex(name);
-  retval = CartSymBCGroupI(GH,index);
-  return retval;
-
-}
-
-void FMODIFIER FORTRAN_NAME(CartSymBCVar)(int *ierr, cGH *GH, ONE_FORTSTRING_ARG) 
-{ 
-
-  ONE_FORTSTRING_CREATE(name)  
-  *ierr = CartSymBCVar(GH,name);
-  free(name);
-}
-
-void FMODIFIER FORTRAN_NAME(CartSymBCGroup)(int *ierr, cGH *GH, ONE_FORTSTRING_ARG) 
-{ 
-
-  ONE_FORTSTRING_CREATE(name)  
-  *ierr = CartSymBCGroup(GH,name);
-  free(name);
-}
-
-void FMODIFIER FORTRAN_NAME(ApplySymmetry)(cGH *GH, ONE_FORTSTRING_ARG) 
-{ 
-
-  ONE_FORTSTRING_CREATE(name)  
-  CartSymBCGroup(GH,name);
-  free(name);
-}
-
-void FMODIFIER FORTRAN_NAME(CartSymBCVarI)(int *ierr, cGH *GH, int *index) 
-{ 
-  *ierr = CartSymBCVarI(GH,*index);
-}
+  if (doSym[0] == 1)
+  {
+    for(i=0; i < ghostz[0]; i++)
+    {
+      var[CCTK_GFINDEX1D(GH,i)] = 
+	sym[0]*var[CCTK_GFINDEX1D(GH,2*ghostz[0]-1-i)];
+    }
+  }
 
-void FMODIFIER FORTRAN_NAME(CartSymBCGRoupI)(int *ierr, cGH *GH, int *index) 
-{ 
-  *ierr = CartSymBCGroupI(GH,*index);
+  return(0);
 }
 
-
-
-
-
-
-
-
-
-