aboutsummaryrefslogtreecommitdiff
path: root/src/ScalarBoundary.c
diff options
context:
space:
mode:
Diffstat (limited to 'src/ScalarBoundary.c')
-rw-r--r--src/ScalarBoundary.c238
1 files changed, 187 insertions, 51 deletions
diff --git a/src/ScalarBoundary.c b/src/ScalarBoundary.c
index cebd85e..1a7cb03 100644
--- a/src/ScalarBoundary.c
+++ b/src/ScalarBoundary.c
@@ -37,6 +37,13 @@ static int ApplyBndScalar (const cGH *GH,
CCTK_REAL scalar,
int first_var,
int num_vars);
+static int OldApplyBndScalar (const cGH *GH,
+ int stencil_dir,
+ const int *stencil_alldirs,
+ int dir,
+ CCTK_REAL scalar,
+ int first_var,
+ int num_vars);
/********************************************************************
******************** External Routines ************************
@@ -51,6 +58,7 @@ static int ApplyBndScalar (const cGH *GH,
the Scalar boundary condition
@enddesc
@calls ApplyBndScalar
+
@var GH
@vdesc Pointer to CCTK grid hierarchy
@vtype const cGH *
@@ -72,6 +80,11 @@ static int ApplyBndScalar (const cGH *GH,
@vtype int
@vio in
@endvar
+ @var widths
+ @vdesc array of boundary widths for each variable
+ @vtype int
+ @vio in
+ @endvar
@var table_handles
@vdesc array of table handles which hold extra arguments
@vtype int
@@ -80,28 +93,35 @@ static int ApplyBndScalar (const cGH *GH,
@returntype int
@returndesc
return code of @seeroutine ApplyBndScalar
+ -21 error reading boundary width array from table
+ -22 wrong size boundary width array in table
@endreturndesc
@@*/
-int BndScalar(const cGH *GH, int num_vars, int *vars, int *faces, int *tables)
+int BndScalar(const cGH *GH, int num_vars, int *vars, int *faces, int *widths,
+ int *tables)
{
- int i, j, k, gdim, max_gdim, err, retval;
+ int i, j, k, gi, gdim, max_gdim, err, retval;
/* variables to pass to ApplyBndScalar */
- /*int stencil_dir;*/ /* width of stencil in direction dir, not needed */
- int *stencil_alldirs; /* width of stencil in all directions */
- int dir; /* direction in which to apply bc */
+ int *width_alldirs; /* width of stencil in all directions */
+ int dir; /* direction in which to apply bc */
CCTK_REAL scalar;
- retval = 0; stencil_alldirs = NULL; max_gdim = 0;
+ retval = 0; width_alldirs = NULL; max_gdim = 0;
/* loop through variables, j at a time */
for (i=0; i<num_vars; i+=j) {
-
/* find other adjacent vars which are selected for identical bcs */
j=1;
- while (i+j<num_vars && vars[i+j]==vars[i]+j && tables[i+j]==tables[i] &&
- faces[i+j]==faces[i])
+ /* Since GFs are allowed to have different staggering, the best we
+ can do is find variables of the same group which are selected
+ for identical bcs. If all GFs had the same staggering then we
+ could groups many GFs together. */
+ gi = CCTK_GroupIndexFromVarI (vars[i]);
+ while (i+j<num_vars && CCTK_GroupIndexFromVarI(vars[i+j])==gi &&
+ tables[i+j]==tables[i] && faces[i+j]==faces[i] &&
+ widths[i+j]==widths[i])
{
++j;
}
@@ -112,28 +132,12 @@ int BndScalar(const cGH *GH, int num_vars, int *vars, int *faces, int *tables)
CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING,
"Faces specification %d for Scalar boundary conditions on "
"%s is not implemented yet. "
- "Applying radiative bcs to all (external) faces.", faces[i],
+ "Applying Scalar bcs to all (external) faces.", faces[i],
CCTK_VarName(vars[i]));
}
dir = 0; /* apply bc to all faces */
/* Set up default arguments for ApplyBndScalar */
- /* Allocate memory and populate stencil width array */
- gdim = CCTK_GroupDimFromVarI(vars[i]);
- if (!stencil_alldirs)
- {
- stencil_alldirs = (int *) malloc(gdim*sizeof(int));
- max_gdim = gdim;
- } else if (gdim > max_gdim)
- {
- realloc(stencil_alldirs, gdim*sizeof(int));
- max_gdim = gdim;
- }
- for (k=0; k<gdim; ++k)
- {
- stencil_alldirs[k] = 1;
- }
- /* Defaults for remainder of arguments */
scalar = 0.;
/* Look on table for possible non-default arguments
@@ -148,20 +152,58 @@ int BndScalar(const cGH *GH, int num_vars, int *vars, int *faces, int *tables)
"Invalid table handle passed for Scalar boundary "
"conditions for %s. Using all default values.",
CCTK_VarName(vars[i]));
+ }
+
+ /* Determine boundary width on all faces */
+ /* allocate memory for buffer */
+ gdim = CCTK_GroupDimI(gi);
+ if (!width_alldirs)
+ {
+ width_alldirs = (int *) malloc(2*gdim*sizeof(int));
+ max_gdim = gdim;
+ } else if (gdim > max_gdim)
+ {
+ width_alldirs = realloc(width_alldirs, gdim*sizeof(int));
+ max_gdim = gdim;
+ }
+
+ /* fill it with values, either from table or the boundary_width
+ parameter */
+ if (widths[i]<0)
+ {
+ err = Util_TableGetIntArray(tables[i], gdim, width_alldirs,
+ "BOUNDARY WIDTH");
+ if (err<0)
+ {
+ CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING,
+ "Error %d when reading boundary width array from table "
+ "for %s",
+ err, CCTK_VarName(vars[i]));
+ return -21;
+ } else if (err!=2*gdim)
+ {
+ CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING,
+ "Boundary width array for %s has %d elements, but %d "
+ "expected", CCTK_VarName(vars[i]), err, 2*gdim);
+ return -22;
+ }
} else
{
- /* Stencil width array */
- Util_TableGetIntArray(tables[i], gdim, stencil_alldirs, "STENCIL WIDTH");
+ for (k=0; k<2*gdim; ++k)
+ {
+ width_alldirs[k] = widths[i];
+ }
}
- if ((retval = ApplyBndScalar(GH, 0, stencil_alldirs, dir, scalar, vars[i],
+ /* Apply the boundary condition */
+ if ((retval = ApplyBndScalar(GH, 0, width_alldirs, dir, scalar, vars[i],
j)) < 0)
{
CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING,
"ApplyBndScalar() returned %d", retval);
}
}
- free(stencil_alldirs);
+ free(width_alldirs);
return retval;
}
@@ -347,7 +389,7 @@ int BndScalarVI (const cGH *GH,
if (vi >= 0 && vi < CCTK_NumVars ())
{
- retval = ApplyBndScalar (GH, -1, stencil, 0, scalar, vi, 1);
+ retval = OldApplyBndScalar (GH, -1, stencil, 0, scalar, vi, 1);
}
else
{
@@ -498,7 +540,7 @@ int BndScalarGI (const cGH *GH,
first_vi = CCTK_FirstVarIndexI (gi);
if (first_vi >= 0)
{
- retval = ApplyBndScalar (GH, -1, stencil, 0, scalar, first_vi,
+ retval = OldApplyBndScalar (GH, -1, stencil, 0, scalar, first_vi,
CCTK_NumVarsInGroupI (gi));
}
else
@@ -918,28 +960,28 @@ void CCTK_FCALL CCTK_FNAME (BndScalarVN)
if (gdim > 0) \
{ \
/* lower x */ \
- SCALAR_BOUNDARY_TYPED (doBC[0], stencil[0], lssh[1], lssh[2], \
+ SCALAR_BOUNDARY_TYPED (doBC[0], widths[0], lssh[1], lssh[2], \
i, j, k, left_cctk_type, right_cctk_type); \
/* upper x */ \
- SCALAR_BOUNDARY_TYPED (doBC[1], stencil[0], lssh[1], lssh[2], \
+ SCALAR_BOUNDARY_TYPED (doBC[1], widths[1], lssh[1], lssh[2], \
lssh[0]-i-1, j, k, left_cctk_type, right_cctk_type);\
} \
if (gdim > 1) \
{ \
/* lower y */ \
- SCALAR_BOUNDARY_TYPED (doBC[2], lssh[0], stencil[1], lssh[2], \
+ SCALAR_BOUNDARY_TYPED (doBC[2], lssh[0], widths[2], lssh[2], \
i, j, k, left_cctk_type, right_cctk_type); \
/* upper y */ \
- SCALAR_BOUNDARY_TYPED (doBC[3], lssh[0], stencil[1], lssh[2], \
+ SCALAR_BOUNDARY_TYPED (doBC[3], lssh[0], widths[3], lssh[2], \
i, lssh[1]-j-1, k, left_cctk_type, right_cctk_type);\
} \
if (gdim > 2) \
{ \
/* lower z */ \
- SCALAR_BOUNDARY_TYPED (doBC[4], lssh[0], lssh[1], stencil[2], \
+ SCALAR_BOUNDARY_TYPED (doBC[4], lssh[0], lssh[1], widths[4], \
i, j, k, left_cctk_type, right_cctk_type); \
/* upper z */ \
- SCALAR_BOUNDARY_TYPED (doBC[5], lssh[0], lssh[1], stencil[2], \
+ SCALAR_BOUNDARY_TYPED (doBC[5], lssh[0], lssh[1], widths[5], \
i, j, lssh[2]-k-1, left_cctk_type, right_cctk_type);\
} \
}
@@ -967,13 +1009,13 @@ void CCTK_FCALL CCTK_FNAME (BndScalarVN)
@vtype const cGH *
@vio in
@endvar
- @var stencil_dir
- @vdesc stencil width in direction dir
+ @var width_dir
+ @vdesc boundary width in direction dir
@vtype int
@vio in
@endvar
- @var stencil_alldirs
- @vdesc stencil widths for all directions
+ @var in_widths
+ @vdesc boundary widths for all directions
@vtype int [ dimension of variable(s) ]
@vio in
@endvar
@@ -1010,13 +1052,13 @@ void CCTK_FCALL CCTK_FNAME (BndScalarVN)
0 for success
-1 if abs(direction) is greater than variables' dimension
-2 if variable dimension is not supported
- -3 if NULL pointer passed as stencil array
+ -3 if NULL pointer passed as boundary width array
-4 if variable type is not supported
@endreturndesc
@@*/
static int ApplyBndScalar (const cGH *GH,
- int stencil_dir,
- const int *stencil_alldirs,
+ int width_dir,
+ const int *in_widths,
int dir,
CCTK_REAL scalar,
int first_var,
@@ -1025,7 +1067,8 @@ static int ApplyBndScalar (const cGH *GH,
int i, j, k;
int gindex, gdim;
int var, timelvl;
- int doBC[2*MAXDIM], dstag[MAXDIM], lsh[MAXDIM], lssh[MAXDIM], stencil[MAXDIM];
+ int doBC[2*MAXDIM], dstag[MAXDIM], lsh[MAXDIM], lssh[MAXDIM],
+ widths[2*MAXDIM];
SymmetryGHex *sGHex;
/* check the direction parameter */
@@ -1050,19 +1093,20 @@ static int ApplyBndScalar (const cGH *GH,
return (-2);
}
- /* set up stencil width array */
+ /* set up boundary width array */
if (dir)
{
- stencil[abs (dir) - 1] = stencil_dir;
+ widths[2*(abs(dir)-1)] = width_dir;
+ widths[2*(abs(dir)-1)+1] = width_dir;
}
- else if (stencil_alldirs)
+ else if (in_widths)
{
- memcpy (stencil, stencil_alldirs, gdim * sizeof (int));
+ memcpy (widths, in_widths, 2* gdim * sizeof (int));
}
else
{
CCTK_WARN (1, "ApplyBndScalar: NULL pointer passed "
- "for stencil width array");
+ "for boundary width array");
return (-3);
}
@@ -1176,3 +1220,95 @@ static int ApplyBndScalar (const cGH *GH,
return(0);
}
+
+
+/*@@
+ @routine OldApplyBndScalar
+ @date 5 May 2003
+ @author David Rideout
+ @desc
+ The new boundary API expects a 2d-element array for the
+ boundary_widths (d=dimension of grid variable), while
+ the old API expects a d-element array. This function
+ converts the old array to the new format.
+ @enddesc
+
+ @var GH
+ @vdesc Pointer to CCTK grid hierarchy
+ @vtype const cGH *
+ @vio in
+ @endvar
+ @var stencil_dir
+ @vdesc boundary width in direction dir
+ @vtype int
+ @vio in
+ @endvar
+ @var stencil_alldirs
+ @vdesc boundary widths for all directions
+ @vtype int [ dimension of variable(s) ]
+ @vio in
+ @endvar
+ @var dir
+ @vdesc direction to set boundaries (0 for setting all directions)
+ @vtype int
+ @vio in
+ @endvar
+ @var scalar
+ @vdesc scalar value to set the boundaries to
+ @vtype CCTK_REAL
+ @vio in
+ @endvar
+ @var first_var
+ @vdesc index of first variable to apply boundary conditions to
+ @vtype int
+ @vio in
+ @endvar
+ @var num_vars
+ @vdesc number of variables
+ @vtype int
+ @vio in
+ @endvar
+
+ @calls CCTK_GroupIndexFromVarI
+ ApplyBndScalar
+ @returntype int
+ @returndesc
+ returncode from @seeroutine ApplyBndScalar
+ @endreturndesc
+@@*/
+
+static int OldApplyBndScalar(const cGH *GH,
+ int stencil_dir,
+ const int *stencil_alldirs,
+ int dir,
+ CCTK_REAL scalar,
+ int first_var,
+ int num_vars)
+{
+ int retval, *boundary_widths, dim, i;
+ static int warned;
+
+ /* Convert stencil_alldirs to new format */
+ dim = CCTK_GroupDimFromVarI (first_var);
+ boundary_widths = (int *) malloc(2*dim*sizeof(int));
+ for (i=0; i<2*dim; ++i)
+ {
+ boundary_widths[i] = stencil_alldirs[i/2];
+ }
+
+ /* Bug people for using the old interface */
+ if (!warned)
+ {
+ CCTK_VWarn (2, __LINE__, __FILE__, CCTK_THORNSTRING,
+ "Copied older d-element array of boundary widths into the "
+ "newer 2d-element format. Please use the new boundary "
+ "interface to avoid this.");
+ warned = 1;
+ }
+
+ /* Call ApplyBnd... with new boundary width array */
+ retval = ApplyBndScalar(GH, stencil_dir, boundary_widths, dir, scalar,
+ first_var, num_vars);
+
+ return retval;
+}