Add noise to boundary condition.

git-svn-id: http://svn.cactuscode.org/arrangements/CactusNumerical/Noise/trunk@15 936e0995-0e4f-0410-aade-aa24bf7baf3d

4 files changed, 66 insertions, 67 deletions

diff --git a/schedule.ccl b/schedule.ccl
index c29205d..13bb5be 100644
--- a/schedule.ccl
+++ b/schedule.ccl
@@ -14,7 +14,7 @@ if (apply_bc_noise) {
      LANG: C
   } "Add noise to boundary condition"
 
-  SCHEDULE bc_noise AT POSTRESTRICT
+  SCHEDULE bc_noise AT CCTK_POSTRESTRICT
   {
      LANG: C
   } "Add noise to boundary condition"
diff --git a/src/bc_noise.c b/src/bc_noise.c
index e6c25ed..1b34721 100644
--- a/src/bc_noise.c
+++ b/src/bc_noise.c
@@ -1,8 +1,15 @@
+/* $Header$ */
+
+#include <assert.h>
+#include <stdio.h>
+#include <stdlib.h>
+#include <string.h>
 
 
 
 #include "noise.h"
 
+#include "cctk.h"
 #include "cctk_Parameters.h"
 #include "cctk_Arguments.h"
 #include "cctk_FortranString.h"
@@ -11,9 +18,11 @@
 
 
 
+/* #define DEBUG_BOUNDARY 1 */
+
 static int ApplyBndNoise (const cGH *GH,
 			  int stencil_dir,
-			  const int *stencil_alldirs,
+			  const int *stencil,
 			  int dir,
 			  int first_var,
 			  int num_vars);
@@ -29,10 +38,10 @@ BndNoiseVI (const cGH *GH, const int *stencil, int vi)
 }
 
 void
-CCTK_FCALL CCTK_FNAME (BndNoiseVI) (int *ierr, const cGH *GH,
+CCTK_FCALL CCTK_FNAME (BndNoiseVI) (int *ierr, const cGH **GH,
 				    const int *stencil, const int *vi)
 {
-  *ierr = BndNoiseVI (GH, stencil, *vi);
+  *ierr = BndNoiseVI (*GH, stencil, *vi);
 }
 
 int
@@ -45,11 +54,11 @@ BndNoiseVN (const cGH *GH, const int *stencil, const char *vn)
 }
 
 void
-CCTK_FCALL CCTK_FNAME (BndNoiseVN) (int *ierr, const cGH *GH,
+CCTK_FCALL CCTK_FNAME (BndNoiseVN) (int *ierr, const cGH **GH,
 				    const int *stencil, ONE_FORTSTRING_ARG)
 {
   ONE_FORTSTRING_CREATE (vn);
-  *ierr = BndNoiseVN (GH, stencil, vn);
+  *ierr = BndNoiseVN (*GH, stencil, vn);
   free (vn);
 }
 
@@ -79,10 +88,10 @@ BndNoiseGI (const cGH *GH, const int *stencil, int gi)
 
 
 void
-CCTK_FCALL CCTK_FNAME (BndNoiseGI) (int *ierr, const cGH *GH,
+CCTK_FCALL CCTK_FNAME (BndNoiseGI) (int *ierr, const cGH **GH,
 				    const int *stencil, const int *gi)
 {
-  *ierr = BndNoiseGI (GH, stencil, *gi);
+  *ierr = BndNoiseGI (*GH, stencil, *gi);
 }
 
 
@@ -108,12 +117,12 @@ BndNoiseGN (const cGH *GH, const int *stencil, const char *gn)
 
 void CCTK_FCALL CCTK_FNAME (BndNoiseGN)
      (int *ierr,
-      const cGH *GH,
+      const cGH **GH,
       const int *stencil,
       ONE_FORTSTRING_ARG)
 {
   ONE_FORTSTRING_CREATE (gn)
-    *ierr = BndNoiseGN (GH, stencil, gn);
+    *ierr = BndNoiseGN (*GH, stencil, gn);
   free (gn);
 }
 
@@ -123,29 +132,25 @@ void CCTK_FCALL CCTK_FNAME (BndNoiseGN)
 
 static int ApplyBndNoise (const cGH *GH,
 			  int stencil_dir,
-			  const int *stencil_alldirs,
+			  const int *stencil,
 			  int dir,
 			  int first_var,
 			  int num_vars)
 {
+  DECLARE_CCTK_PARAMETERS;
   int i, j, k;
   int var, vtypesize, gindex, gdim, timelvl;
-  int doBC[2*MAXDIM], dstag[MAXDIM], lsh[MAXDIM], lssh[MAXDIM],
-    stencil[MAXDIM];
+  int doBC[2*MAXDIM], dstag[MAXDIM], lsh[MAXDIM], lssh[MAXDIM];
   SymmetryGHex *sGHex;
-  CCTK_REAL amplitude;
-  CCTK_REAL offset;
   int type;
 
-  stencil[0] = 1; stencil[1] = 1; stencil[2] = 1;
+  /* This argument is unused an undocumented; better make sure people
+     don't try to use it for something.  */
+  assert (stencil_dir == -1);
 
   /* get the group index of the variables */
   gindex = CCTK_GroupIndexFromVarI (first_var);
 
-  amplitude = *((CCTK_REAL*) CCTK_ParameterGet("amplitude", "Noise", &type))
-    / RAND_MAX;
-  offset = 0.5 * amplitude;
-
   /* get the number of dimensions and the size of the variables' type */
   gdim      = CCTK_GroupDimI (gindex);
   vtypesize = CCTK_VarTypeSize (CCTK_VarTypeI (first_var));
@@ -171,7 +176,7 @@ static int ApplyBndNoise (const cGH *GH,
      so that we can use the INDEX_3D macro later on */
   for (i = gdim; i < MAXDIM; i++)
     {
-      lsh[i]  = 0;
+      lsh[i]  = 1;
       lssh[i] = 1;
     }
 
@@ -210,8 +215,8 @@ static int ApplyBndNoise (const cGH *GH,
 	  doBC[i*2+1] &= GH->cctk_lsh[i] > 1 && GH->cctk_bbox[i*2+1];
 	  if (dir != 0)
 	    {
-	      doBC[i*2]   &= (dir < 0 && (i + 1 == abs (dir)));
-	      doBC[i*2+1] &= (dir > 0 && (i + 1 == abs (dir)));
+	      doBC[i*2]   &= i+1 == -dir;
+	      doBC[i*2+1] &= i+1 ==  dir;
 	    }
 	}
 
@@ -284,16 +289,21 @@ static int ApplyBndNoise (const cGH *GH,
 static void
 add_bc_noise_to_var (int idx, const char* optstring, void* cctkGH)
 {
-  int i, j, k, ijk;
-  CCTK_REAL* data;
+  DECLARE_CCTK_PARAMETERS;
   cGH* GH = cctkGH;
   int sw[3];
 
-  sw[0] = 1;
-  sw[1] = 1;
-  sw[2] = 1;
-
-  BndNoiseVI(GH, sw, idx);
+  /* Change type from CCTK_INT to int */
+  sw[0] = noise_stencil[0];
+  sw[1] = noise_stencil[1];
+  sw[2] = noise_stencil[2];
+
+  if (noise_boundaries[0]) ApplyBndNoise (GH, -1, sw, -1, idx, 1);
+  if (noise_boundaries[1]) ApplyBndNoise (GH, -1, sw, +1, idx, 1);
+  if (noise_boundaries[2]) ApplyBndNoise (GH, -1, sw, -2, idx, 1);
+  if (noise_boundaries[3]) ApplyBndNoise (GH, -1, sw, +2, idx, 1);
+  if (noise_boundaries[4]) ApplyBndNoise (GH, -1, sw, -3, idx, 1);
+  if (noise_boundaries[5]) ApplyBndNoise (GH, -1, sw, +3, idx, 1);
 }
 
 
@@ -303,6 +313,7 @@ bc_noise(CCTK_ARGUMENTS)
   DECLARE_CCTK_ARGUMENTS
   DECLARE_CCTK_PARAMETERS
 
+/*   Boundary_MakeSureThatTheSelectionIsEmpty(); */
   if (CCTK_TraverseString(bc_vars, add_bc_noise_to_var, cctkGH,
       		    CCTK_GROUP_OR_VAR) < 0)
     {
diff --git a/src/id_noise.c b/src/id_noise.c
index 9951435..4d4e8e1 100644
--- a/src/id_noise.c
+++ b/src/id_noise.c
@@ -18,21 +18,16 @@
 #include "cctk_Parameters.h"
 
 void
-add_noise_to_group (int idx, const char* optstring, void* cctkGH)
+add_noise_to_var (int idx, const char* optstring, void* cctkGH)
 {
+  DECLARE_CCTK_PARAMETERS;
   int i, j, k, ijk;
   CCTK_REAL* data;
-  CCTK_REAL amplitude;
-  CCTK_REAL offset;
   cGH* GH = cctkGH;
   int type;
 
   data = (CCTK_REAL*) CCTK_VarDataPtrI(GH, 0, idx);
 
-  amplitude = *((CCTK_REAL*) CCTK_ParameterGet("amplitude", "Noise", &type))
-    / RAND_MAX;
-  offset = 0.5 * amplitude;
-
   for (k=1; k< GH->cctk_lsh[2]-1; ++k)
     {
       for (j=1; j< GH->cctk_lsh[1]-1; ++j)
@@ -54,7 +49,7 @@ id_noise(CCTK_ARGUMENTS)
   DECLARE_CCTK_ARGUMENTS
   DECLARE_CCTK_PARAMETERS
 
-  if (CCTK_TraverseString(id_vars, add_noise_to_group, cctkGH,
+  if (CCTK_TraverseString(id_vars, add_noise_to_var, cctkGH,
 			  CCTK_GROUP_OR_VAR) < 0)
     {
       CCTK_WARN (1, "Failed to parse 'IDRandom::noisy_id_vars' parameter");
diff --git a/src/noise.h b/src/noise.h
index 7cb839a..749b411 100644
--- a/src/noise.h
+++ b/src/noise.h
@@ -3,50 +3,43 @@
 #ifndef NOISE_H
 #define NOISE_H
 
+#include <assert.h>
 #include <stdlib.h>
+#include <string.h>
 #include <time.h>
 #include "cctk.h"
 
 
 /* constants */
 
-#define TRUE 1
-#define FALSE 0
-
 #define MAXDIM  3
-#
+
 
 /* macros */
 
-#define RAND_VAL random()*amplitude - offset
-
-#define INDEX_3D(lsh, i, j, k)      ((i) + (lsh)[0]*((j) + (lsh)[1]*(k)))
-
-
-#define BOUNDARY_NOISE(doBC,                                                 \
-                       iend, jend, kend,                                     \
-                       ii, jj, kk)            \
-{                                                                             \
-if (doBC)                                        \
-  {                                              \
-    for (k = 0; k < kend; k++)                   \
-    {                                                                         \
-      for (j = 0; j < jend; j++)                                              \
-      {                                                                       \
-        for (i = 0; i < iend; i++)                                            \
-        {                                                                     \
-          int _index = CCTK_GFINDEX3D(GH, ii, jj, kk);          \
-          CCTK_REAL* v= (CCTK_REAL *) GH->data[var][timelvl] + _index; \
-          v[0] += RAND_VAL;  \
-        }                                                  \
+#define RAND_VAL ((random()*(1.0/RAND_MAX)-0.5)*amplitude)
+
+#define BOUNDARY_NOISE(doBC,                                            \
+                       iend, jend, kend,                                \
+                       ii, jj, kk)                                      \
+{                                                                       \
+  if (doBC)                                                             \
+  {                                                                     \
+    CCTK_REAL* v= CCTK_VarDataPtrI(GH, timelvl, var);                   \
+    for (k = 0; k < kend; k++)                                          \
+    {                                                                   \
+      for (j = 0; j < jend; j++)                                        \
+      {                                                                 \
+        for (i = 0; i < iend; i++)                                      \
+        {                                                               \
+          const int _index = CCTK_GFINDEX3D(GH, (ii), (jj), (kk));      \
+          v[_index] += RAND_VAL;                                        \
+        }                                                               \
       }                                                                 \
     }                                                                   \
-  }                                                      \
+  }                                                                     \
 }
 
 
 
 #endif /* !define(NOISE_H) */
-
-
-