From 01051497d58f949a8d1f857fee4c814fe039f888 Mon Sep 17 00:00:00 2001 From: tradke Date: Fri, 12 Sep 2003 17:32:23 +0000 Subject: Another fix for *_RADIATIVE_BOUNDARY_3D: dx was applied on all the faces. git-svn-id: http://svn.cactuscode.org/arrangements/CactusBase/Boundary/trunk@246 6a38eb6e-646e-4a02-a296-d141613ad6c4 --- src/RadiationBoundary.c | 296 +++++++++++++----------------------------------- 1 file changed, 78 insertions(+), 218 deletions(-) (limited to 'src') diff --git a/src/RadiationBoundary.c b/src/RadiationBoundary.c index 7ed9d15..9033d68 100644 --- a/src/RadiationBoundary.c +++ b/src/RadiationBoundary.c @@ -1156,50 +1156,18 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) @vtype int @vio in @endvar - @var xyz - @vdesc coordinate arrays for x, y, and z - @vtype CCTK_REAL [ MAXDIM ] - @vio in - @endvar - @var radius - @vdesc radius coordinate array - @vtype CCTK_REAL [] - @vio in - @endvar - @var offset - @vdesc offset to the next element in this direction + @var dim + @vdesc dimension to apply BC @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 [] - @vio in - @endvar - @var var_from - @vdesc source variable array - @vtype cctk_type [] - @vio in - @endvar @var cctk_type @vdesc CCTK datatypes of the source and target variable @vtype @vio in @endvar @@*/ -#define LOWER_RADIATIVE_BOUNDARY_3D(istart, jstart, kstart, \ - xyz, \ - radius, \ - offset, \ - rho, \ - var_to, \ - var_from, \ - cctk_type) \ +#define LOWER_RADIATIVE_BOUNDARY_3D(istart, jstart, kstart, dim, cctk_type) \ { \ int _i, _j, _k; \ \ @@ -1208,51 +1176,50 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) { \ for (_j = jstart-1; _j >= 0; _j--) \ { \ - int _0 = 0*offset, \ - _1 = 1*offset, \ - _2 = 2*offset; \ - int _index = CCTK_GFINDEX3D (GH, istart-1, _j, _k); \ - CCTK_REAL *_radius = radius + _index, \ - *_xyz = xyz + _index; \ - cctk_type *_to = (cctk_type *) (var_to) + _index, \ - *_from = (cctk_type *) (var_from) + _index; \ + int _0 = 0*offset[dim], \ + _1 = 1*offset[dim], \ + _2 = 2*offset[dim]; \ + int _idx = CCTK_GFINDEX3D (GH, istart-1, _j, _k); \ + CCTK_REAL *_r = coords[MAXDIM] + _idx, \ + *_xyz = coords[dim] + _idx; \ + cctk_type *_to = (cctk_type *) GH->data[var_to][timelvl_to] + _idx, \ + *_from = (cctk_type *) GH->data[var_from][timelvl_from] + _idx;\ \ \ for (_i = istart-1; _i >= 0; _i--) \ { \ - CCTK_REAL _radius0_inv = 1/_radius[_0], \ - _radius1_inv = 1/_radius[_1]; \ + CCTK_REAL _r0_inv = 1/_r[_0], _r1_inv = 1/_r[_1]; \ \ if (radpower >= 0) \ { \ CCTK_REAL H; \ \ - H = 0.25 * radpower * dxyz[0] \ - * (_xyz[_0]*SQR (_radius0_inv) + _xyz[_1]*SQR (_radius1_inv)); \ + H = 0.25 * radpower * dxyz[dim] \ + * (_xyz[_0]*SQR (_r0_inv) + _xyz[_1]*SQR (_r1_inv)); \ H = (1 + H) / (1 - H); \ 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 \ - *( SQR (_radius[_1]) / _xyz[_1] \ - + SQR (_radius[_2]) / _xyz[_2]); \ + + 0.5*( _r[_1] * (_to[_1] - _from[_1]) \ + + _r[_2] * (_to[_2] - _from[_2])) \ + + 0.25*(_to[_2] - _to[_1] + _from[_2] - _from[_1]) * rho[dim] \ + *( SQR (_r[_1]) / _xyz[_1] \ + + SQR (_r[_2]) / _xyz[_2]); \ dtvvar0H = dtvvar0 + H; \ } \ \ _to[_0] = (cctk_type) ( \ - (dtvvar0H * ( _xyz[_0] * SQR (_radius0_inv) \ - + _xyz[_1] * SQR (_radius1_inv)) \ - - _to[_1] * ( rho + _xyz[_1]*_radius1_inv \ - *(1 + dtvh*_radius1_inv))\ - + _from[_0] * ( rho + _xyz[_0]*_radius0_inv \ - *(1 - dtvh*_radius0_inv))\ - - _from[_1] * ( rho - _xyz[_1]*_radius1_inv \ - *(1 - dtvh*_radius1_inv))\ + (dtvvar0H * ( _xyz[_0] * SQR (_r0_inv) \ + + _xyz[_1] * SQR (_r1_inv)) \ + - _to[_1] * ( rho[dim] + _xyz[_1]*_r1_inv \ + *(1 + dtvh*_r1_inv)) \ + + _from[_0] * ( rho[dim] + _xyz[_0]*_r0_inv \ + *(1 - dtvh*_r0_inv)) \ + - _from[_1] * ( rho[dim] - _xyz[_1]*_r1_inv \ + *(1 - dtvh*_r1_inv)) \ ) \ - / (-rho + _xyz[_0]*_radius0_inv \ - *(1 + dtvh*_radius0_inv))\ + / (-rho[dim] + _xyz[_0]*_r0_inv \ + *(1 + dtvh*_r0_inv)) \ ); \ - _radius--; _xyz--; _to--; _from--; \ + _r--; _xyz--; _to--; _from--; \ } \ } \ } \ @@ -1272,50 +1239,18 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) @vtype int @vio in @endvar - @var xyz - @vdesc coordinate arrays for x, y, and z - @vtype CCTK_REAL [ MAXDIM ] - @vio in - @endvar - @var radius - @vdesc radius coordinate array - @vtype CCTK_REAL [] - @vio in - @endvar - @var offset - @vdesc offset to the next element in this direction + @var dim + @vdesc dimension to apply BC @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 [] - @vio in - @endvar - @var var_from - @vdesc source variable array - @vtype cctk_type [] - @vio in - @endvar @var cctk_type @vdesc CCTK datatypes of the source and target variable @vtype @vio in @endvar @@*/ -#define UPPER_RADIATIVE_BOUNDARY_3D(istart, jstart, kstart, \ - xyz, \ - radius, \ - offset, \ - rho, \ - var_to, \ - var_from, \ - cctk_type) \ +#define UPPER_RADIATIVE_BOUNDARY_3D(istart, jstart, kstart, dim, cctk_type) \ { \ int _i, _j, _k; \ \ @@ -1324,51 +1259,50 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) { \ for (_j = jstart; _j < GH->cctk_lsh[1]; _j++) \ { \ - int _0 = -0*offset, \ - _1 = -1*offset, \ - _2 = -2*offset; \ - int _index = CCTK_GFINDEX3D (GH, istart, _j, _k); \ - CCTK_REAL *_radius = radius + _index, \ - *_xyz = xyz + _index; \ - cctk_type *_to = (cctk_type *) (var_to) + _index, \ - *_from = (cctk_type *) (var_from) + _index; \ + int _0 = -0*offset[dim], \ + _1 = -1*offset[dim], \ + _2 = -2*offset[dim]; \ + int _idx = CCTK_GFINDEX3D (GH, istart, _j, _k); \ + CCTK_REAL *_r = coords[MAXDIM] + _idx, \ + *_xyz = coords[dim] + _idx; \ + cctk_type *_to = (cctk_type *) GH->data[var_to][timelvl_to] + _idx, \ + *_from = (cctk_type *) GH->data[var_from][timelvl_from] + _idx;\ \ \ for (_i = istart; _i < GH->cctk_lsh[0]; _i++) \ { \ - CCTK_REAL _radius0_inv = 1/_radius[_0], \ - _radius1_inv = 1/_radius[_1]; \ + CCTK_REAL _r0_inv = 1/_r[_0], _r1_inv = 1/_r[_1]; \ \ if (radpower >= 0) \ { \ CCTK_REAL H; \ \ - H = 0.25 * radpower * dxyz[0] \ - * (_xyz[_0]*SQR (_radius0_inv) + _xyz[_1]*SQR (_radius1_inv)); \ + H = 0.25 * radpower * dxyz[dim] \ + * (_xyz[_0]*SQR (_r0_inv) + _xyz[_1]*SQR (_r1_inv)); \ H = (1 - H) / (1 + H); \ 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 \ - * ( SQR (_radius[_1]) / _xyz[_1] \ - + SQR (_radius[_2]) / _xyz[_2]); \ + + 0.5 *( _r[_1] * (_to[_1] - _from[_1]) \ + + _r[_2] * (_to[_2] - _from[_2])) \ + + 0.25*(_to[_1] - _to[_2] + _from[_1] - _from[_2]) * rho[dim] \ + * ( SQR (_r[_1]) / _xyz[_1] \ + + SQR (_r[_2]) / _xyz[_2]); \ dtvvar0H = dtvvar0 + H; \ } \ \ _to[_0] = (cctk_type) ( \ - (dtvvar0H * ( _xyz[_0] * (SQR (_radius0_inv)) \ - + _xyz[_1] * (SQR (_radius1_inv))) \ - + _to[_1] * ( rho - _xyz[_1]*_radius1_inv \ - *(1 + dtvh*_radius1_inv))\ - + _from[_0] * (-rho + _xyz[_0]*_radius0_inv \ - *(1 - dtvh*_radius0_inv))\ - + _from[_1] * ( rho + _xyz[_1]*_radius1_inv \ - *(1 - dtvh*_radius1_inv))\ + (dtvvar0H * ( _xyz[_0] * (SQR (_r0_inv)) \ + + _xyz[_1] * (SQR (_r1_inv))) \ + + _to[_1] * ( rho[dim] - _xyz[_1]*_r1_inv \ + *(1 + dtvh*_r1_inv)) \ + + _from[_0] * (-rho[dim] + _xyz[_0]*_r0_inv \ + *(1 - dtvh*_r0_inv)) \ + + _from[_1] * ( rho[dim] + _xyz[_1]*_r1_inv \ + *(1 - dtvh*_r1_inv)) \ ) \ - / ( rho + _xyz[_0]*_radius0_inv \ - *(1 + dtvh*_radius0_inv))\ + / ( rho[dim] + _xyz[_0]*_r0_inv \ + *(1 + dtvh*_r0_inv)) \ ); \ - _radius++; _xyz++; _to++; _from++; \ + _r++; _xyz++; _to++; _from++; \ } \ } \ } \ @@ -1391,118 +1325,62 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) @vtype int [ dim ] @vio in @endvar - @var doBC - @vdesc flags telling whether to apply BC in a given direction or not - @vtype int [ 2*dim ] - @vio in - @endvar @var stencil @vdesc stencils in every direction @vtype int [ 2*dim ] @vio in @endvar - @var coords - @vdesc coordinate arrays for x, y, z, and the radius - @vtype CCTK_REAL [ MAXDIM+1 ] - @vio in - @endvar - @var offset - @vdesc offsets to the next element in each dimension - @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 [] - @vio in - @endvar - @var var_from - @vdesc source variable array - @vtype cctk_type [] - @vio in - @endvar - @var dim - @vdesc dimension of the source and target variable - @vtype int - @vio in - @endvar @var cctk_type @vdesc CCTK datatypes of the source and target variable @vtype @vio in @endvar @@*/ -#define RADIATIVE_BOUNDARY(lsh, \ - doBC, \ - stencil, \ - coords, \ - offset, \ - rho, \ - var_to, \ - var_from, \ - dim, \ - cctk_type) \ +#define RADIATIVE_BOUNDARY(lsh, stencil, cctk_type) \ { \ /* check the dimensionality */ \ - if (dim != 3) \ + if (gdim != 3) \ { \ CCTK_VWarn (1, __LINE__, __FILE__, CCTK_THORNSTRING, \ "ApplyBndRadiative: variable dimension of %d not supported", \ - dim); \ + gdim); \ return (-5); \ } \ \ /* Lower x-bound */ \ if (doBC[0]) \ { \ - LOWER_RADIATIVE_BOUNDARY_3D (stencil[0], lsh[1], lsh[2], \ - coords[0], coords[MAXDIM], offset[0], rho[0],\ - var_to, var_from, cctk_type); \ + LOWER_RADIATIVE_BOUNDARY_3D (stencil[0], lsh[1], lsh[2], 0, cctk_type); \ } \ \ /* Upper x-bound */ \ if (doBC[1]) \ { \ - UPPER_RADIATIVE_BOUNDARY_3D (lsh[0] - stencil[1], 0, 0, \ - coords[0], coords[MAXDIM], offset[0], rho[0],\ - var_to, var_from, cctk_type); \ + UPPER_RADIATIVE_BOUNDARY_3D (lsh[0] - stencil[1], 0, 0, 0, cctk_type); \ } \ \ /* Lower y-bound */ \ if (doBC[2]) \ { \ - LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], stencil[2], lsh[2], \ - coords[1], coords[MAXDIM], offset[1], rho[1],\ - var_to, var_from, cctk_type); \ + LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], stencil[2], lsh[2], 1, cctk_type); \ } \ \ /* Upper y-bound */ \ if (doBC[3]) \ { \ - UPPER_RADIATIVE_BOUNDARY_3D (0, lsh[1] - stencil[3], 0, \ - coords[1], coords[MAXDIM], offset[1], rho[1],\ - var_to, var_from, cctk_type); \ + UPPER_RADIATIVE_BOUNDARY_3D (0, lsh[1] - stencil[3], 0, 1, cctk_type); \ } \ \ /* Lower z-bound */ \ if (doBC[4]) \ { \ - LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], lsh[1], stencil[4], \ - coords[2], coords[MAXDIM], offset[2], rho[2],\ - var_to, var_from, cctk_type); \ + LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], lsh[1], stencil[4], 2, cctk_type); \ } \ \ /* Upper z-bound */ \ if (doBC[5]) \ { \ - UPPER_RADIATIVE_BOUNDARY_3D (0, 0, lsh[2] - stencil[5], \ - coords[2], coords[MAXDIM], offset[2], rho[2],\ - var_to, var_from, cctk_type); \ + UPPER_RADIATIVE_BOUNDARY_3D (0, 0, lsh[2] - stencil[5], 2, cctk_type); \ } \ } @@ -1734,68 +1612,50 @@ static int ApplyBndRadiative (const cGH *GH, switch (CCTK_VarTypeI (var_to)) { case CCTK_VARIABLE_CHAR: - 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); + RADIATIVE_BOUNDARY (GH->cctk_lsh, widths, CCTK_CHAR); break; case CCTK_VARIABLE_INT: - 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); + RADIATIVE_BOUNDARY (GH->cctk_lsh, widths, CCTK_INT); break; case CCTK_VARIABLE_REAL: - 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); + RADIATIVE_BOUNDARY (GH->cctk_lsh, widths, CCTK_REAL); break; #ifdef CCTK_INT2 case CCTK_VARIABLE_INT2: - 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); + RADIATIVE_BOUNDARY (GH->cctk_lsh, widths, CCTK_INT2); break; #endif #ifdef CCTK_INT4 case CCTK_VARIABLE_INT4: - 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); + RADIATIVE_BOUNDARY (GH->cctk_lsh, widths, CCTK_INT4); break; #endif #ifdef CCTK_INT8 case CCTK_VARIABLE_INT8: - 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); + RADIATIVE_BOUNDARY (GH->cctk_lsh, widths, CCTK_INT8); break; #endif #ifdef CCTK_REAL4 case CCTK_VARIABLE_REAL4: - 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); + RADIATIVE_BOUNDARY (GH->cctk_lsh, widths, CCTK_REAL4); break; #endif #ifdef CCTK_REAL8 case CCTK_VARIABLE_REAL8: - 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); + RADIATIVE_BOUNDARY (GH->cctk_lsh, widths, CCTK_REAL8); break; #endif #ifdef CCTK_REAL16 case CCTK_VARIABLE_REAL16: - 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); + RADIATIVE_BOUNDARY (GH->cctk_lsh, widths, CCTK_REAL16); break; #endif -- cgit v1.2.3