Erik Schnetter's implementation of ICN with averaging, so the intermediate steps are always at t+dt.

git-svn-id: http://svn.cactuscode.org/arrangements/CactusNumerical/MoL/trunk@24 578cdeb0-5ea1-4b81-8215-5a3b8777ee0b

1 files changed, 127 insertions, 0 deletions

diff --git a/src/ICN.c b/src/ICN.c
index 0a1cced..a5231a7 100644
--- a/src/ICN.c
+++ b/src/ICN.c
@@ -34,6 +34,7 @@ CCTK_FILEVERSION(AlphaThorns_MoL_ICN_c);
  ********************************************************************/
 
 void MoL_ICNAdd(CCTK_ARGUMENTS);
+void MoL_ICNAverage(CCTK_ARGUMENTS);
 
 /********************************************************************
  ********************* Other Routine Prototypes *********************
@@ -176,6 +177,132 @@ void MoL_ICNAdd(CCTK_ARGUMENTS)
 
 }
 
+ /*@@
+   @routine    MoL_ICNAverage
+   @date       Fri Jul 18 14:02:00 2003
+   @author     Ian Hawke, Erik Schnetter
+   @desc 
+   Averages between the current and the previous time level.
+   @enddesc 
+   @calls     
+   @calledby   
+   @history 
+ 
+   @endhistory 
+
+@@*/
+
+void MoL_ICNAverage(CCTK_ARGUMENTS)
+{
+
+  DECLARE_CCTK_ARGUMENTS;
+  DECLARE_CCTK_PARAMETERS;
+  
+  cGroupDynamicData arraydata;
+  CCTK_INT groupindex, ierr;
+  CCTK_INT arraytotalsize, arraydim;
+
+  CCTK_INT index, var;
+  CCTK_INT totalsize;
+  CCTK_REAL *OldVar;
+  CCTK_REAL *UpdateVar;
+  CCTK_REAL *RHSVar;
+  
+  /* FIXME */
+
+#ifdef MOLDOESCOMPLEX
+
+  CCTK_COMPLEX *OldComplexVar;
+  CCTK_COMPLEX *UpdateComplexVar;
+  CCTK_COMPLEX *RHSComplexVar;
+  CCTK_COMPLEX Complex_Half = CCTK_Cmplx(0.5, 0);
+
+#endif
+
+#ifdef MOLDEBUG
+  printf("Inside ICN.\nProcessor %d.\nStep %d.\nRefinement %d.\nTimestep %g.\nSpacestep %g.\nTime %g\n",
+         CCTK_MyProc(cctkGH),
+         MoL_Intermediate_Steps - *MoL_Intermediate_Step + 1,
+         *cctk_levfac,
+         CCTK_DELTA_TIME,
+         CCTK_DELTA_SPACE(0),
+         cctk_time);
+#endif  
+
+  totalsize = cctk_lsh[0] * cctk_lsh[1] * cctk_lsh[2];
+#ifdef MOLDEBUG
+  printf("MoL: the ICN routine says dt = %f.\n", CCTK_DELTA_TIME);
+#endif
+  for (var = 0; var < MoLNumEvolvedVariables; var++)
+  {
+    OldVar = (CCTK_REAL*)CCTK_VarDataPtrI(cctkGH, 1,
+                                          EvolvedVariableIndex[var]);
+    UpdateVar = (CCTK_REAL*)CCTK_VarDataPtrI(cctkGH, 0,
+                                             EvolvedVariableIndex[var]);
+    RHSVar = (CCTK_REAL*)CCTK_VarDataPtrI(cctkGH, 0, 
+                                          RHSVariableIndex[var]);
+    
+    for (index = 0; index < totalsize; index++)
+    {
+      UpdateVar[index] = 0.5 * (UpdateVar[index] + OldVar[index]);
+    }
+  }
+  
+  for (var = 0; var < MoLNumEvolvedArrayVariables; var++)
+  {
+    OldVar = (CCTK_REAL*)CCTK_VarDataPtrI(cctkGH, 1,
+                                             EvolvedArrayVariableIndex[var]);
+    UpdateVar = (CCTK_REAL*)CCTK_VarDataPtrI(cctkGH, 0,
+                                             EvolvedArrayVariableIndex[var]);
+    RHSVar = (CCTK_REAL*)CCTK_VarDataPtrI(cctkGH, 0, 
+                                          RHSArrayVariableIndex[var]);
+    
+    groupindex = CCTK_GroupIndexFromVarI(EvolvedArrayVariableIndex[var]);
+    ierr = CCTK_GroupDynamicData(cctkGH, groupindex,
+                                 &arraydata);
+    if (ierr)
+    {
+      CCTK_VWarn(0, __LINE__, __FILE__, "MoL", 
+                 "The driver does not return group information for group '%s'.", 
+                 CCTK_GroupName(groupindex));
+    }
+    arraytotalsize = 1;
+    for (arraydim = 0; arraydim < arraydata.dim; arraydim++)
+    {
+      arraytotalsize *= arraydata.lsh[arraydim];
+    }
+
+    for (index = 0; index < arraytotalsize; index++)
+    {
+      UpdateVar[index] = 0.5 * (UpdateVar[index] + OldVar[index]);
+    }
+  }
+
+  /* FIXME */
+
+#ifdef MOLDOESCOMPLEX
+  
+  for (var = 0; var < MoLNumEvolvedComplexVariables; var++)
+  {
+    OldComplexVar = (CCTK_COMPLEX*)CCTK_VarDataPtrI(cctkGH, 1,
+                                                    EvolvedComplexVariableIndex[var]);
+    UpdateComplexVar = (CCTK_COMPLEX*)CCTK_VarDataPtrI(cctkGH, 0,
+                                                       EvolvedComplexVariableIndex[var]);
+    RHSComplexVar = (CCTK_COMPLEX*)CCTK_VarDataPtrI(cctkGH, 0, 
+                                                    RHSComplexVariableIndex[var]);
+    
+    for (index = 0; index < totalsize; index++)
+    {
+      UpdateComplexVar[index] = CCTK_CmplxMul(Complex_Half, CCTK_CmplxAdd(UpdateComplexVar[index], OldComplexVar[index]));
+    }
+  }
+
+#endif
+  
+  return;
+
+}
+
 /********************************************************************
  *********************     Local Routines   *************************
  ********************************************************************/