Fix for mixed slicing, registering static and geodesic for mixed slicing

and formatting


git-svn-id: http://svn.einsteintoolkit.org/cactus/EinsteinBase/CoordGauge/trunk@33 edbb7e70-9571-45d5-a481-0a560a9b4751

1 files changed, 186 insertions, 121 deletions

diff --git a/src/Slicing.c b/src/Slicing.c
index 720f4a0..dd1584c 100644
--- a/src/Slicing.c
+++ b/src/Slicing.c
@@ -7,15 +7,17 @@ static char *rcsid = "$Id$";
    @desc 
      A flexible calling structure to calling slicing functions:
      * a slicing is registered with its NAME, with a NUMBER derived from
-       the position in the parameter and indicating its priority, with a FUNCTION, 
-       which tells if the slicing should be used in the next iteration. 
-     * The slicing is REGISTERED in a startup routine of the thorn providing the 
-       slicing.
+       the position in the parameter and indicating its priority, with a 
+       FUNCTION, which tells if the slicing should be used in the next 
+       iteration. 
+     * The slicing is REGISTERED in a startup routine of the thorn 
+       providing the slicing.
      * This structure (HandledData) is created during STARTUP
      * This strucure is initialized during CCTK_INITIAL
      * This structure is evaluted during a CCTK_PRESTEP and an integer
-       active_slicing_handle defines which slicing is next. All slicings have to check
-       if they match this handle and only then take their turn.
+       active_slicing_handle defines which slicing is next. 
+       All slicings have to check if they match this handle and only 
+       then take their turn.
    @enddesc 
  @@*/
 
@@ -29,7 +31,7 @@ static char *rcsid = "$Id$";
 #include "Slicing.h"
 #include "ErrorCodes.h"
 
-/*$#define ESLIC_DEBUG$*/
+/*#define ESLIC_DEBUG*/
 
 
 /* Local data holding info on slicings..*/
@@ -66,27 +68,29 @@ int Einstein_RegisterSlicing(const char *slice_name)
   /*Check that this Slicing has not been registered, yet */
   handle = Util_GetHandle(Eslicings,slice_name,NULL);
   
-  if (handle<0) {
-  
+  if (handle<0) 
+  {
     /* Allocate new slicing structure */
-    new_slicing = (struct Einstein_slicing*) malloc(sizeof(struct Einstein_slicing));
+    new_slicing = (struct Einstein_slicing*) 
+      malloc(sizeof(struct Einstein_slicing));
 
-    if (new_slicing) {
+    if (new_slicing) 
+    {
       
-      /* store new_slicing in the StoreHandledData-array and returns the handle for it */
+      /* store new_slicing */
       handle = Util_NewHandle(&Eslicings, slice_name, new_slicing);
     
 #ifdef ESLIC_DEBUG
-      printf("REGISTER_SLICING: name: >>%s<<  new handle %d \n",slice_name,handle);
+      printf("REGISTER_SLICING: >>%s<<  new handle %d \n",slice_name,handle);
 #endif
 
-      /*Initialize the Einstein_slicing structure */
+      /* Initialize the Einstein_slicing structure */
 
       /* we know the slicing name */
       new_slicing->name = (char *)slice_name; 
 
       /* Function pointer , evaluates if to perform this slicings, allowing 
-	 for iteration checks, gridfucntion behavior, etc.  Will be set by
+	 for iteration checks, gridfunction behavior, etc.  Will be set by
          Einstein_RegisterTimeToSlice */
       new_slicing->timetoslice = NULL;
       
@@ -100,20 +104,17 @@ int Einstein_RegisterSlicing(const char *slice_name)
       {
 	/* Memory failure */
 	CCTK_WARN(0,"cannot allocate memory for new slicing");
-	handle = -2;  /* FIXME: some decent MACROS here */
+	handle = -2; 
       }
   }
   else
-    { 
-      /* Extension handle already exists */
-      err =(char*)malloc(256*sizeof(char));
-      sprintf(err,"REGISTER_SLICING: new handle %s (# %d) allready in use \n",
+  { 
+    /* Extension handle already exists */
+    CCTK_VWarn(0,__LINE__,__FILE__,CCTK_THORNSTRING,
+	      "New handle %s (%d) already in use",
 	      slice_name,handle);
-      CCTK_WARN(0,err);
-      if (err) free(err);
-
-      handle = -1;
-    }
+    handle = -1;
+  }
 
   return handle;
 }
@@ -155,7 +156,7 @@ int CCTK_FCALL CCTK_FNAME(Einstein_RegisterSlicing)(ONE_FORTSTRING_ARG)
       iteration. The user can register any function to check for iteration 
       conditions, grid function behavior, etc. 
 
-      This has to be in C currently. We cannot wrapp this easily (TOM).
+      This has to be in C currently. We cannot wrap this easily (TOM).
    @enddesc 
    @calls     
    @calledby   The user
@@ -173,17 +174,17 @@ int Einstein_RegisterTimeToSlice(int handle, int (*func)(cGH *))
 
   slicing = Util_GetHandledData(Eslicings, handle);
   
-  if (slicing)  {
+  if (slicing)  
+  {
     slicing->timetoslice=func;
     return_code = 1;
   }
-  else {
+  else 
+  {
     /* handle didn't return the data pointer */ 
-    err=(char*)malloc(128*sizeof(char));
-    sprintf(err,"could not obtain slicing-structure for handle %d",handle);
-    CCTK_WARN(1,err);
-    if (err) free (err);
-
+    CCTK_VWarn(1,__LINE__,__FILE__,CCTK_THORNSTRING,
+	      "Could not obtain slicing structure for handle %d",
+	      handle);
     return_code = 0;
   }   
 
@@ -230,7 +231,10 @@ void Einstein_ActivateSlicing(cGH *GH)
 
   /* Get our grid scalar pointer: initialize with -1*/
   i = CCTK_VarIndex("Einstein::active_slicing_handle");
-  if (i<0) CCTK_WARN(0,"INTERNAL ERROR: cannot find grid scalar: active_slicing_handle");
+  if (i<0) 
+  {
+    CCTK_WARN(0,"Cannot find grid scalar: active_slicing_handle");
+  }
   active_slicing_handle = (CCTK_INT *)CCTK_VarDataPtrI(GH,0,i);
   *active_slicing_handle = -1;
 
@@ -240,39 +244,43 @@ void Einstein_ActivateSlicing(cGH *GH)
     { 
       handle= Util_GetHandle(Eslicings,slicing,NULL);
 
-      if (handle<0) {
+      if (handle<0) 
+      {
 
 	/* If slicing is set to "none", only warn level 2 (could be initial data only) */
 	if (CCTK_Equals(slicing,"none")) 
 	{
-	  CCTK_WARN(2,"No slicing set. Prob. fatal when you enter evolution loop (itlast>0).\n");
+	  CCTK_WARN(2,"No slicing set. Probably fatal for evolution");
 	}
         else if (CCTK_Equals(slicing,"static"))
 	{
           handle2=Einstein_RegisterSlicing("static");
-          if (handle2<0) CCTK_WARN(1,"Cannot register static slicing");
+          if (handle2<0) 
+	  {
+	    CCTK_WARN(1,"Cannot register static slicing");
+	  }
 	}
         else if (CCTK_Equals(slicing,"geodesic"))
 	{
           handle2=Einstein_RegisterSlicing("geodesic");
-          if (handle2<0) CCTK_WARN(1,"Cannot register geodesic slicing");
+          if (handle2<0) 
+	  {
+	    CCTK_WARN(1,"Cannot register geodesic slicing");
+	  }
 	} 
 	else 
 	{
-	  err = (char*)malloc(256*sizeof(char));
-	  sprintf(err,"%s%s%s%s",
-		  "WARNING: cannot get handle for slicing: >>",
-		  slicing,
-		  "<< \nTypo? Slicing registered by its provider ",
-		  "with Einstein_RegisterSlicing? \nProb. fatal when you continue to evolve.\n");
-	  CCTK_WARN(1,err);
-	  if (err) free(err);
+	  CCTK_VWarn(1,__LINE__,__FILE__,CCTK_THORNSTRING,
+		     "Cannot get handle for slicing %s",slicing);
 	}
 	return;
       }
 
       slic  = Util_GetHandledData(Eslicings,handle);
-      if (!slic) CCTK_WARN(0,"INTERNAL ERROR: cannot access handle structure!");
+      if (!slic)
+      {
+	CCTK_WARN(0,"Cannot access handle structure!");
+      }
       slic->param_active = 1;
 
       /*Some error checking: */
@@ -307,13 +315,30 @@ void Einstein_ActivateSlicing(cGH *GH)
       priority    = 1;      
       while (item)
       {
+        if (CCTK_Equals(item,"static"))
+	{
+	  printf("Registering Static\n");
+          handle2=Einstein_RegisterSlicing("static");
+          if (handle2<0) 
+	  {
+	    CCTK_WARN(1,"Cannot register static slicing");
+	  }
+	}
+        else if (CCTK_Equals(item,"geodesic"))
+	{
+          handle2=Einstein_RegisterSlicing("geodesic");
+          if (handle2<0) 
+	  {
+	    CCTK_WARN(1,"Cannot register geodesic slicing");
+	  }
+	} 
+
 	slic = Util_GetHandledData(Eslicings,Util_GetHandle(Eslicings,item,NULL));
 
-	if (!slic) { 
-	  err=(char*)malloc(128*sizeof(char));
-	  sprintf(err,"ERROR: no registered slicing found for item >%s<\n",item);
-	  CCTK_WARN(0,err);
-	  if (err) free(err);
+	if (!slic) 
+        { 
+	  CCTK_VWarn(0,__LINE__,__FILE__,CCTK_THORNSTRING,
+		     "No registered slicing for <%s>",item);
 	}
 
 	slic->param_active = priority;
@@ -342,24 +367,28 @@ void Einstein_ActivateSlicing(cGH *GH)
 
 #ifdef ESLIC_DEBUG   
   if (CCTK_Equals(slicing_verbose,"yes")) 
-    { 
-      /* Some diagnostic output */
-      for (h=0;h<num_slicings;h++) {
-	slic = Util_GetHandledData(Eslicings,h);
-	if (!slic) CCTK_WARN(0,"INTERNAL ERROR: cannot get slicing for handle.");
-	
-	printf(" Einstein_ActivateSlicing:  >>%s<<   priority: %d   %s \n",
-	       slic->name,
-	       slic->param_active,
-	       (char*)((slic->timetoslice)?"FUNC":"nofunc "));
+  { 
+    int h;
+    /* Some diagnostic output */
+    for (h=0;h<num_slicings;h++) 
+    {
+      slic = Util_GetHandledData(Eslicings,h);
+      if (!slic) 
+      {
+	CCTK_WARN(0,"Cannot get slicing for handle.");
       }
+      printf(" Einstein_ActivateSlicing:  >>%s<<   priority: %d   %s \n",
+	     slic->name,
+	     slic->param_active,
+	     (char*)((slic->timetoslice)?"FUNC":"nofunc "));
     }
+  }
 #endif
 
-  
-
 }
 
+
+
  /*@@
    @routine    Einstein_GetSlicingHandle
    @date       Thu Jul 22 11:30:47 1999
@@ -437,7 +466,8 @@ int CCTK_FCALL CCTK_FNAME(Einstein_GetSlicingHandle)(ONE_FORTSTRING_ARG)
    @endhistory 
 @@*/
 
-void Einstein_SetNextSlicing(cGH *GH) {
+void Einstein_SetNextSlicing(cGH *GH) 
+{
   DECLARE_CCTK_PARAMETERS
 
   struct    Einstein_slicing *slic;
@@ -453,86 +483,121 @@ void Einstein_SetNextSlicing(cGH *GH) {
 
   /* Get our grid scalar pointer */
   i = CCTK_VarIndex("Einstein::active_slicing_handle");
-  if (i<0) CCTK_WARN(0,"INTERNAL ERROR: cannot find grid scalar: active_slicing_handle");
+  if (i<0) 
+  {
+    CCTK_WARN(0,"Cannot find grid scalar: active_slicing_handle");
+  }
   active_slicing_handle = (CCTK_INT *)CCTK_VarDataPtrI(GH,0,i);
 
-
-  /* First check if we run mixed NOT: a single slicing specified,
-     that's easy: */
   if (!CCTK_Equals(slicing,"mixed")) 
+  {
+    h = Util_GetHandle(Eslicings,slicing,NULL);
+    if (h<0) 
     {
-
-      h = Util_GetHandle(Eslicings,slicing,NULL);
-      if (h<0) CCTK_WARN(0,"ERROR: slicing parameter specifies a non-activated/registered slicing! \n");
-      
-      slic = (struct Einstein_slicing *)Util_GetHandledData(Eslicings,h);
-      if (!slic) 
-	CCTK_WARN(0,"INTERNAL ERROR: Slicing registry out of sync! No slicing found!");
+      CCTK_WARN(0,"Slicing parameter specifies non-activated/registered slicing!");
+    }
       
-      /* Set the active handle to the only slicing and return*/
-      /*$*active_slicing_handle = h;$*/
-      h_yes = h;
-
+    slic = (struct Einstein_slicing *)Util_GetHandledData(Eslicings,h);
+    if (!slic) 
+    {
+      CCTK_WARN(0,"Slicing registry out of sync! No slicing found!");
     }
+      
+    /* Set the active handle to the only slicing and return*/
+    /*$*active_slicing_handle = h;$*/
+    h_yes = h;
+    
+  }
   /* We run mixed: Loop over all registered handles, 
      check if active, get the timetoslice function and decide which slicing
      to use  */
   else
+  {
+    for (h=0;h<num_slicings;h++) 
     {
-      for (h=0;h<num_slicings;h++) 
-	{
 #ifdef ESLIC_DEBUG	  
-	  printf("HANDLE LOOP: %d of %d\n",h,num_slicings);
+      printf("HANDLE LOOP: %d of %d\n",h,num_slicings);
 #endif	  
-	  slic = (struct Einstein_slicing *)Util_GetHandledData(Eslicings,h);
-	  if (!slic) 
-	    CCTK_WARN(0,"INTERNAL ERROR: Slicing registry out of sync! No slicing found!");
-	  
-	  /* we have two competing flags: the priority and the 
-	     timetoslice function. We pick two handles for
-	     1) highest priority (lowest param_active number) 
-	        for timetoslice = yes  (1)
-	     2) highest priority (lowest param_active number) 
-                for timetoslice = don't care (0)
-	     set the active_slicing_handle to one of the two handles in this order
-	  */
-#ifdef ESLIC_DEBUG
-	  printf(" handle: %d pactive: %d yes: %d egal: %d \n",h,slic->param_active,h_yes, h_egal);
-#endif
-	  if ((slic->param_active>0)) {
-	    if (slic->timetoslice!=NULL) {
-	      doicare = slic->timetoslice(GH);
-	    }
-	    if ((doicare== SLICING_YES)     &&(slic->param_active<h_yes))  h_yes = h; 
-	    if ((doicare== SLICING_DONTCARE)&&(slic->param_active<h_egal)) h_egal= h;
-	    if ((doicare== SLICING_NO)      &&(slic->param_active<h_no ))  h_no  = h;
-	    
+      slic = (struct Einstein_slicing *)Util_GetHandledData(Eslicings,h);
+      if (!slic)
+      {
+	CCTK_WARN(0,"Slicing registry out of sync! No slicing found!");
+      }
+      else
+      {
+	/* we have two competing flags: the priority and the 
+	   timetoslice function. We pick two handles for
+	   1) highest priority (lowest param_active number) 
+	   for timetoslice = yes  (1)
+	   2) highest priority (lowest param_active number) 
+	   for timetoslice = don't care (0)
+	   set the active_slicing_handle to one of the two handles 
+	   in this order
+	*/
+	if ((slic->param_active>0)) 
+	{
+	  doicare = 0;
+	  if (slic->timetoslice!=NULL) 
+	  {
+	    doicare = slic->timetoslice(GH);
+	  }
+	  if ((doicare== SLICING_YES) && (slic->param_active<h_yes))  
+	  {
+	    h_yes = h; 
+	  }
+	  if ((doicare== SLICING_DONTCARE) && (slic->param_active<h_egal)) 
+	  {
+	    h_egal= h;
+	  }
+	  if ((doicare== SLICING_NO) && (slic->param_active<h_no ))  
+	  {
+	    h_no  = h;
 	  }
 	}
+#ifdef ESLIC_DEBUG
+	printf(" handle:%d name:%s param_active:%d h_yes:%d h_egal:%d \n",
+	       h,slic->name,slic->param_active,h_yes, h_egal);
+#endif
+      }
     }
+  }
   
   /* First try to set the scalar to the highest priority that 
-     cares (h_yes), else to the highest prio, that doesn;t 
+     cares (h_yes), else to the highest priority that doesn;t 
      care (h_egal). If that fails, we are in trouble.  */
 
-  if (h_yes!=num_slicings) *active_slicing_handle = h_yes;
+  if (h_yes!=num_slicings) 
+  {
+    *active_slicing_handle = h_yes;
+  }
   else 
-    if (h_egal!=num_slicings) *active_slicing_handle = h_egal;
-    else printf("ERROR: no active slicing found: hyes/hegal/hno: %d %d %d \n",
-		h_yes,h_no,h_egal);
+  {
+    if (h_egal != num_slicings)
+    {
+      *active_slicing_handle = h_egal;
+    }
+    else 
+    {
+      printf("ERROR: no active slicing found: hyes/hegal/hno: %d %d %d \n",
+	     h_yes,h_egal,h_no);
+    }
+  }
   
   if (CCTK_Equals(slicing_verbose,"yes")) 
-    { 
-      info = (char*)malloc(256*sizeof(char));
-      slic = (struct  Einstein_slicing *)
-	Util_GetHandledData(Eslicings,*active_slicing_handle);
-      if (!slic) CCTK_WARN(0,"INTERNAL ERROR: Cannot find slicing handle");
-
-      sprintf(info,"next active_slicing_handle: %s  (handle# %d)",
-	      slic->name,*active_slicing_handle);
-      CCTK_INFO(info);
-      if (info) free(info);
+  { 
+    info = (char*)malloc(256*sizeof(char));
+    slic = (struct  Einstein_slicing *)
+      Util_GetHandledData(Eslicings,*active_slicing_handle);
+    if (!slic) 
+    {
+      CCTK_WARN(0,"Cannot find slicing handle");
     }
+    
+    sprintf(info,"Next slicing: %s (%d)",
+	    slic->name,*active_slicing_handle);
+    CCTK_INFO(info);
+    if (info) free(info);
+  }