Fixed *_RADIATIVE_BOUNDARY_3D macros which accidentally applied rho_x

on all the faces. Bug was spotted by Mihaela.


git-svn-id: http://svn.cactuscode.org/arrangements/CactusBase/Boundary/trunk@245 6a38eb6e-646e-4a02-a296-d141613ad6c4

1 files changed, 48 insertions, 31 deletions

diff --git a/src/RadiationBoundary.c b/src/RadiationBoundary.c
index 98ea3b9..7ed9d15 100644
--- a/src/RadiationBoundary.c
+++ b/src/RadiationBoundary.c
@@ -270,8 +270,7 @@ CCTK_INT BndRadiative(const cGH *GH, CCTK_INT num_vars, CCTK_INT *vars,
     gdim = CCTK_GroupDimI(gi);
     if (gdim > max_gdim)
     {
-      width_alldirs = 
-	(CCTK_INT *) realloc(width_alldirs, 2*gdim*sizeof(CCTK_INT));
+      width_alldirs = realloc(width_alldirs, 2*gdim*sizeof(CCTK_INT));
       max_gdim = gdim;
     }
 
@@ -304,7 +303,7 @@ CCTK_INT BndRadiative(const cGH *GH, CCTK_INT num_vars, CCTK_INT *vars,
 
     /* Apply the boundary condition */
     if ((retval = ApplyBndRadiative(GH, 0, width_alldirs, dir, limit, speed,
-				    vars[i], prev_time_level, j)) < 0)
+                                    vars[i], prev_time_level, j)) < 0)
     {
       CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING,
                  "ApplyBndRadiative() returned %d", retval);
@@ -545,9 +544,9 @@ int BndRadiativeGI (const cGH *GH,
   first_vi_from = CCTK_FirstVarIndexI (gi_from);
   if (first_vi_to >= 0 && first_vi_from >= 0)
   {
-    retval = OldApplyBndRadiative (GH, -1, (const CCTK_INT *) stencil, 0, var0, 
-				   speed, first_vi_to, first_vi_from,
-				   CCTK_NumVarsInGroupI (gi_to));
+    retval = OldApplyBndRadiative (GH, -1, (const CCTK_INT *) stencil, 0, var0,
+                                   speed, first_vi_to, first_vi_from,
+                                   CCTK_NumVarsInGroupI (gi_to));
   }
   else
   {
@@ -1172,6 +1171,11 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
    @vtype      int
    @vio        in
    @endvar
+   @var        rho
+   @vdesc      wave speed multiplied by courant factor
+   @vtype      CCTK_REAL
+   @vio        in
+   @endvar
    @var        var_to
    @vdesc      target variable array
    @vtype      cctk_type []
@@ -1192,6 +1196,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
                                     xyz,                                      \
                                     radius,                                   \
                                     offset,                                   \
+                                    rho,                                      \
                                     var_to,                                   \
                                     var_from,                                 \
                                     cctk_type)                                \
@@ -1228,7 +1233,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
           H *= dtv * (0.25*(_to[_1] + _to[_2] + _from[_1] + _from[_2]) - var0)\
                + 0.5*(  _radius[_1] * (_to[_1] - _from[_1])                   \
                       + _radius[_2] * (_to[_2] - _from[_2]))                  \
-               + 0.25*(_to[_2] - _to[_1] + _from[_2] - _from[_1]) * rho[0]    \
+               + 0.25*(_to[_2] - _to[_1] + _from[_2] - _from[_1]) * rho       \
                      *(  SQR (_radius[_1]) / _xyz[_1]                         \
                        + SQR (_radius[_2]) / _xyz[_2]);                       \
           dtvvar0H = dtvvar0 + H;                                             \
@@ -1237,14 +1242,14 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
         _to[_0] = (cctk_type) (                                               \
                    (dtvvar0H * (  _xyz[_0] * SQR (_radius0_inv)               \
                                 + _xyz[_1] * SQR (_radius1_inv))              \
-                    - _to[_1]   * ( rho[0] + _xyz[_1]*_radius1_inv            \
+                    - _to[_1]   * ( rho + _xyz[_1]*_radius1_inv               \
                                                      *(1 + dtvh*_radius1_inv))\
-                    + _from[_0] * ( rho[0] + _xyz[_0]*_radius0_inv            \
+                    + _from[_0] * ( rho + _xyz[_0]*_radius0_inv               \
                                                      *(1 - dtvh*_radius0_inv))\
-                    - _from[_1] * ( rho[0] - _xyz[_1]*_radius1_inv            \
+                    - _from[_1] * ( rho - _xyz[_1]*_radius1_inv               \
                                                      *(1 - dtvh*_radius1_inv))\
                    )                                                          \
-                   /              (-rho[0] + _xyz[_0]*_radius0_inv            \
+                   /              (-rho + _xyz[_0]*_radius0_inv               \
                                                      *(1 + dtvh*_radius0_inv))\
                   );                                                          \
         _radius--; _xyz--; _to--; _from--;                                    \
@@ -1282,6 +1287,11 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
    @vtype      int
    @vio        in
    @endvar
+   @var        rho
+   @vdesc      wave speed multiplied by courant factor
+   @vtype      CCTK_REAL
+   @vio        in
+   @endvar
    @var        var_to
    @vdesc      target variable array
    @vtype      cctk_type []
@@ -1302,6 +1312,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
                                     xyz,                                      \
                                     radius,                                   \
                                     offset,                                   \
+                                    rho,                                      \
                                     var_to,                                   \
                                     var_from,                                 \
                                     cctk_type)                                \
@@ -1338,7 +1349,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
           H *= dtv * (0.25*(_to[_1] + _to[_2] + _from[_1] + _from[_2]) - var0)\
                + 0.5 *(  _radius[_1] * (_to[_1] - _from[_1])                  \
                        + _radius[_2] * (_to[_2] - _from[_2]))                 \
-               + 0.25*(_to[_1] - _to[_2] + _from[_1] - _from[_2]) * rho[0]    \
+               + 0.25*(_to[_1] - _to[_2] + _from[_1] - _from[_2]) * rho       \
                      * (  SQR (_radius[_1]) / _xyz[_1]                        \
                         + SQR (_radius[_2]) / _xyz[_2]);                      \
           dtvvar0H = dtvvar0 + H;                                             \
@@ -1347,14 +1358,14 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
         _to[_0] = (cctk_type) (                                               \
                    (dtvvar0H * (  _xyz[_0] * (SQR (_radius0_inv))             \
                                 + _xyz[_1] * (SQR (_radius1_inv)))            \
-                    + _to[_1]   * ( rho[0] - _xyz[_1]*_radius1_inv            \
+                    + _to[_1]   * ( rho - _xyz[_1]*_radius1_inv               \
                                                      *(1 + dtvh*_radius1_inv))\
-                    + _from[_0] * (-rho[0] + _xyz[_0]*_radius0_inv            \
+                    + _from[_0] * (-rho + _xyz[_0]*_radius0_inv               \
                                                      *(1 - dtvh*_radius0_inv))\
-                    + _from[_1] * ( rho[0] + _xyz[_1]*_radius1_inv            \
+                    + _from[_1] * ( rho + _xyz[_1]*_radius1_inv               \
                                                      *(1 - dtvh*_radius1_inv))\
                    )                                                          \
-                   /              ( rho[0] + _xyz[_0]*_radius0_inv            \
+                   /              ( rho + _xyz[_0]*_radius0_inv               \
                                                      *(1 + dtvh*_radius0_inv))\
                   );                                                          \
         _radius++; _xyz++; _to++; _from++;                                    \
@@ -1400,6 +1411,11 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
    @vtype      int [ dim ]
    @vio        in
    @endvar
+   @var        rho
+   @vdesc      wave speed multiplied by courant factor
+   @vtype      CCTK_REAL []
+   @vio        in
+   @endvar
    @var        var_to
    @vdesc      target variable array
    @vtype      cctk_type []
@@ -1426,6 +1442,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
                            stencil,                                           \
                            coords,                                            \
                            offset,                                            \
+                           rho,                                               \
                            var_to,                                            \
                            var_from,                                          \
                            dim,                                               \
@@ -1444,7 +1461,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
   if (doBC[0])                                                                \
   {                                                                           \
     LOWER_RADIATIVE_BOUNDARY_3D (stencil[0], lsh[1], lsh[2],                  \
-                                 coords[0], coords[MAXDIM], offset[0],        \
+                                 coords[0], coords[MAXDIM], offset[0], rho[0],\
                                  var_to, var_from, cctk_type);                \
   }                                                                           \
                                                                               \
@@ -1452,7 +1469,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
   if (doBC[1])                                                                \
   {                                                                           \
     UPPER_RADIATIVE_BOUNDARY_3D (lsh[0] - stencil[1], 0, 0,                   \
-                                 coords[0], coords[MAXDIM], offset[0],        \
+                                 coords[0], coords[MAXDIM], offset[0], rho[0],\
                                  var_to, var_from, cctk_type);                \
   }                                                                           \
                                                                               \
@@ -1460,7 +1477,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
   if (doBC[2])                                                                \
   {                                                                           \
     LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], stencil[2], lsh[2],                  \
-                                 coords[1], coords[MAXDIM], offset[1],        \
+                                 coords[1], coords[MAXDIM], offset[1], rho[1],\
                                  var_to, var_from, cctk_type);                \
   }                                                                           \
                                                                               \
@@ -1468,7 +1485,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
   if (doBC[3])                                                                \
   {                                                                           \
     UPPER_RADIATIVE_BOUNDARY_3D (0, lsh[1] - stencil[3], 0,                   \
-                                 coords[1], coords[MAXDIM], offset[1],        \
+                                 coords[1], coords[MAXDIM], offset[1], rho[1],\
                                  var_to, var_from, cctk_type);                \
   }                                                                           \
                                                                               \
@@ -1476,7 +1493,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
   if (doBC[4])                                                                \
   {                                                                           \
     LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], lsh[1], stencil[4],                  \
-                                 coords[2], coords[MAXDIM], offset[2],        \
+                                 coords[2], coords[MAXDIM], offset[2], rho[2],\
                                  var_to, var_from, cctk_type);                \
   }                                                                           \
                                                                               \
@@ -1484,7 +1501,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN)
   if (doBC[5])                                                                \
   {                                                                           \
     UPPER_RADIATIVE_BOUNDARY_3D (0, 0, lsh[2] - stencil[5],                   \
-                                 coords[2], coords[MAXDIM], offset[2],        \
+                                 coords[2], coords[MAXDIM], offset[2], rho[2],\
                                  var_to, var_from, cctk_type);                \
   }                                                                           \
 }
@@ -1717,26 +1734,26 @@ static int ApplyBndRadiative (const cGH *GH,
     switch (CCTK_VarTypeI (var_to))
     {
       case CCTK_VARIABLE_CHAR:
-        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset,
+        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, rho,
                             GH->data[var_to][timelvl_to],
                             GH->data[var_from][timelvl_from], gdim, CCTK_CHAR);
         break;
 
       case CCTK_VARIABLE_INT:
-        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset,
+        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, rho,
                             GH->data[var_to][timelvl_to],
                             GH->data[var_from][timelvl_from], gdim, CCTK_INT);
         break;
 
       case CCTK_VARIABLE_REAL:
-        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset,
+        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, rho,
                             GH->data[var_to][timelvl_to],
                             GH->data[var_from][timelvl_from], gdim, CCTK_REAL);
         break;
 
 #ifdef CCTK_INT2
       case CCTK_VARIABLE_INT2:
-        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset,
+        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, rho,
                             GH->data[var_to][timelvl_to],
                             GH->data[var_from][timelvl_from], gdim, CCTK_INT2);
         break;
@@ -1744,7 +1761,7 @@ static int ApplyBndRadiative (const cGH *GH,
 
 #ifdef CCTK_INT4
       case CCTK_VARIABLE_INT4:
-        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset,
+        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, rho,
                             GH->data[var_to][timelvl_to],
                             GH->data[var_from][timelvl_from], gdim, CCTK_INT4);
         break;
@@ -1752,7 +1769,7 @@ static int ApplyBndRadiative (const cGH *GH,
 
 #ifdef CCTK_INT8
       case CCTK_VARIABLE_INT8:
-        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset,
+        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, rho,
                             GH->data[var_to][timelvl_to],
                             GH->data[var_from][timelvl_from], gdim, CCTK_INT8);
         break;
@@ -1760,7 +1777,7 @@ static int ApplyBndRadiative (const cGH *GH,
 
 #ifdef CCTK_REAL4
       case CCTK_VARIABLE_REAL4:
-        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset,
+        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, rho,
                             GH->data[var_to][timelvl_to],
                             GH->data[var_from][timelvl_from], gdim, CCTK_REAL4);
         break;
@@ -1768,7 +1785,7 @@ static int ApplyBndRadiative (const cGH *GH,
 
 #ifdef CCTK_REAL8
       case CCTK_VARIABLE_REAL8:
-        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset,
+        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, rho,
                             GH->data[var_to][timelvl_to],
                             GH->data[var_from][timelvl_from], gdim, CCTK_REAL8);
         break;
@@ -1776,7 +1793,7 @@ static int ApplyBndRadiative (const cGH *GH,
 
 #ifdef CCTK_REAL16
       case CCTK_VARIABLE_REAL16:
-        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset,
+        RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, rho,
                             GH->data[var_to][timelvl_to],
                             GH->data[var_from][timelvl_from], gdim,CCTK_REAL16);
         break;