From dcdc72dd6f10779d1d4294a4b2afb2dab218369d Mon Sep 17 00:00:00 2001 From: rideout Date: Mon, 5 May 2003 22:31:49 +0000 Subject: Support for new boundary API, which places an integer boundary_width argument on calls to Select*ForBC*. If a negative value is given then the boundary condition will look into the table for a 2d element integer array specifying the width of each boundary face. A number bug fixes, including the restriction to only execute a boundary condition on a single group at a time. (consider that the next group may have different staggering) Improvements of debugging code. Wrapper functions OldApplyBnd added, which support the old d-element 'stencil width' array API. git-svn-id: http://svn.cactuscode.org/arrangements/CactusBase/Boundary/trunk@221 6a38eb6e-646e-4a02-a296-d141613ad6c4 --- src/Boundary.c | 110 +++++++++++++++---- src/Boundary.h | 17 +-- src/CopyBoundary.c | 240 ++++++++++++++++++++++++++--------------- src/FlatBoundary.c | 239 ++++++++++++++++++++++++++++++----------- src/NoneBoundary.c | 8 +- src/RadiationBoundary.c | 279 ++++++++++++++++++++++++++++++++++++------------ src/RobinBoundary.c | 212 ++++++++++++++++++++++++++++-------- src/ScalarBoundary.c | 238 ++++++++++++++++++++++++++++++++--------- src/StaticBoundary.c | 228 +++++++++++++++++++++++++++++---------- 9 files changed, 1177 insertions(+), 394 deletions(-) diff --git a/src/Boundary.c b/src/Boundary.c index b79aea2..84c2f54 100644 --- a/src/Boundary.c +++ b/src/Boundary.c @@ -39,6 +39,7 @@ struct BCVAR { struct BCVAR *next; /* pointer to next entry in list */ int faces; /* set of faces for this application of bc */ + int width; /* width of the boundary, if it is equal for all faces */ int table; /* table handle holding extra arguments */ int var; /* index of grid variable to which to apply the bc */ }; @@ -80,27 +81,32 @@ CCTK_INT Bdry_Boundary_RegisterPhysicalBC(CCTK_POINTER GH, CCTK_STRING bc_name); CCTK_INT Bdry_Boundary_SelectVarForBC(CCTK_POINTER GH, CCTK_INT faces, + CCTK_INT boundary_width, CCTK_INT table_handle, CCTK_STRING var_name, CCTK_STRING bc_name); CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, CCTK_INT faces, + CCTK_INT boundary_width, CCTK_INT table_handle, CCTK_INT var_index, CCTK_STRING bc_name); CCTK_INT Bdry_Boundary_SelectGroupForBC(CCTK_POINTER GH, CCTK_INT faces, + CCTK_INT boundary_width, CCTK_INT table_handle, CCTK_STRING group_name, CCTK_STRING bc_name); CCTK_INT Bdry_Boundary_SelectGroupForBCI(CCTK_POINTER GH, CCTK_INT faces, + CCTK_INT boundary_width, CCTK_INT table_handle, CCTK_INT group_index, CCTK_STRING bc_name); CCTK_INT Bdry_Boundary_SelectedGVs(CCTK_POINTER GH, CCTK_INT array_size, CCTK_POINTER var_indices, CCTK_POINTER faces, + CCTK_POINTER boundary_widths, CCTK_POINTER table_handles, CCTK_STRING bc_name); @@ -245,6 +251,11 @@ CCTK_INT Bdry_Boundary_RegisterPhysicalBC(CCTK_POINTER GH, @vtype CCTK_INT @vio in @endvar + @var width + @vdesc if >=0, width of boundary in all directions + @vtype CCTK_INT + @vio in + @endvar @var table_handle @vdesc handle of table which holds arguments to be passed to bc @vtype CCTK_INT @@ -269,12 +280,19 @@ CCTK_INT Bdry_Boundary_RegisterPhysicalBC(CCTK_POINTER GH, @@*/ CCTK_INT Bdry_Boundary_SelectVarForBC(CCTK_POINTER GH, CCTK_INT faces, + CCTK_INT width, CCTK_INT table_handle, CCTK_STRING var_name, CCTK_STRING bc_name) { int retval, var_index; +#ifdef DEBUG + printf("Boundary_SelectVarForBC:\n"); + printf(" called with faces=%d, width=%d, table_handle=%d, var_name=%s, bc_name=%s\n", + faces, width, table_handle, var_name, bc_name); +#endif + retval = 0; var_index = CCTK_VarIndex(var_name); @@ -285,7 +303,7 @@ CCTK_INT Bdry_Boundary_SelectVarForBC(CCTK_POINTER GH, retval = -3; } else { - retval = Bdry_Boundary_SelectVarForBCI(GH, faces, table_handle, var_index, + retval = Bdry_Boundary_SelectVarForBCI(GH, faces, width, table_handle, var_index, bc_name); } @@ -314,6 +332,11 @@ CCTK_INT Bdry_Boundary_SelectVarForBC(CCTK_POINTER GH, @vtype CCTK_INT @vio in @endvar + @var width + @vdesc if >=0, width of boundary in all directions + @vtype CCTK_INT + @vio in + @endvar @var table_handle @vdesc handle of table which holds arguments to be passed to bc @vtype CCTK_INT @@ -339,6 +362,7 @@ CCTK_INT Bdry_Boundary_SelectVarForBC(CCTK_POINTER GH, @@*/ CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, CCTK_INT faces, + CCTK_INT width, CCTK_INT table_handle, CCTK_INT var_index, CCTK_STRING bc_name) @@ -358,8 +382,8 @@ CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, #ifdef DEBUG printf("Boundary_SelectVarForBCI:\n"); - printf(" called with table_handle=%d, var_index=%d, bc_name=%s\n", - table_handle, var_index, bc_name); + printf(" called with faces=%d, width=%d, table_handle=%d, var_index=%d, bc_name=%s\n", + faces, width, table_handle, var_index, bc_name); printf(" vi %d corresponds to %s\n", var_index, CCTK_VarName(var_index)); #endif @@ -387,6 +411,7 @@ CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, /* populate new entry with data */ new_entry -> faces = faces; + new_entry -> width = width; new_entry -> table = table_handle; new_entry -> var = var_index; /* new_entry -> next will be filled in later */ @@ -414,10 +439,10 @@ CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, } } - /* If current has not been initialized yet, this is a new bc_name + /* If current_bcdata is NULL, we got to the end of the above loop, this is a new bc_name * that does not appear in the bcdata list. */ - if (!current) /* bc_name was not found in bcdata_list */ + if (!current_bcdata) /* bc_name was not found in bcdata_list */ { #ifdef DEBUG printf("Boundary_SelectVarForBCI: adding new entry to bcdata list\n"); @@ -429,12 +454,6 @@ CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING, "Unable to allocate memory for internal 'bcdata' list"); retval = -4; - } else - { - new_bcdata->var_list = new_entry; /* new_entry will be first element */ - new_entry->next = NULL; /* of var list for this bcdata entry */ - new_bcdata->bc_name = Util_Strdup(bc_name); - new_bcdata->num = 1; } /* Place new entry into bcdata list, maintaining case independent sort. */ @@ -468,6 +487,9 @@ CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, bcdata_list = new_bcdata; new_bcdata->next = NULL; } + + /* Set current_bcdata to new_bcdata, so the new bcdata entry will be filled in below */ + current_bcdata = new_bcdata; } #ifdef DEBUG @@ -475,6 +497,18 @@ CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, " previous is %p\n", current, previous); #endif + if (!current) /* This is the first element in the var_list */ + { +#ifdef DEBUG + printf(" New element of var list.\n"); +#endif + current_bcdata->var_list = new_entry; /* new_entry will be first element */ + new_entry->next = NULL; /* of var list for this bcdata entry */ + current_bcdata->bc_name = Util_Strdup(bc_name); + current_bcdata->num = 1; + } + + /* Enter new_entry into correct location in linked list. * Note that this loop is skipped if new_entry was already inserted as first * element of a new var list above (since in that case current will be @@ -546,6 +580,11 @@ CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, @vtype CCTK_INT @vio in @endvar + @var width + @vdesc if >=0, width of boundary in all directions + @vtype CCTK_INT + @vio in + @endvar @var table_handle @vdesc handle of table which holds arguments to be passed to bc @vtype CCTK_INT @@ -570,6 +609,7 @@ CCTK_INT Bdry_Boundary_SelectVarForBCI(CCTK_POINTER GH, @@*/ CCTK_INT Bdry_Boundary_SelectGroupForBC(CCTK_POINTER GH, CCTK_INT faces, + CCTK_INT width, CCTK_INT table_handle, CCTK_STRING group_name, CCTK_STRING bc_name) @@ -590,7 +630,7 @@ CCTK_INT Bdry_Boundary_SelectGroupForBC(CCTK_POINTER GH, } else { /* call Bdry_Boundary_SelectGroupForBCI() */ - retval = Bdry_Boundary_SelectGroupForBCI(GH, faces, table_handle, gi, + retval = Bdry_Boundary_SelectGroupForBCI(GH, faces, width, table_handle, gi, bc_name); } @@ -618,6 +658,11 @@ CCTK_INT Bdry_Boundary_SelectGroupForBC(CCTK_POINTER GH, @vtype CCTK_INT @vio in @endvar + @var width + @vdesc if >=0, width of boundary in all directions + @vtype CCTK_INT + @vio in + @endvar @var table_handle @vdesc handle of table which holds arguments to be passed to bc @vtype CCTK_INT @@ -643,6 +688,7 @@ CCTK_INT Bdry_Boundary_SelectGroupForBC(CCTK_POINTER GH, @@*/ CCTK_INT Bdry_Boundary_SelectGroupForBCI(CCTK_POINTER GH, CCTK_INT faces, + CCTK_INT width, CCTK_INT table_handle, CCTK_INT group_index, CCTK_STRING bc_name) @@ -661,7 +707,8 @@ CCTK_INT Bdry_Boundary_SelectGroupForBCI(CCTK_POINTER GH, } #ifdef DEBUG - printf("Boundary_SelectGroupForBCI: group %s has %d vars\n", CCTK_GroupName(group_index), num_vars); + printf("Boundary_SelectGroupForBCI: group %s has %d vars\n", + CCTK_GroupName(group_index), num_vars); #endif /* loop over variables in group */ @@ -669,7 +716,7 @@ CCTK_INT Bdry_Boundary_SelectGroupForBCI(CCTK_POINTER GH, max_vi = vi + num_vars; for (; vinext) { #ifdef DEBUG - printf(" looping through bcdata list, at bcdata entry for %s bc\n", current_bcdata->bc_name); + printf(" looping through bcdata list, at bcdata entry for %s bc\n", + current_bcdata->bc_name); #endif if (!bc_name || CCTK_Equals(current_bcdata->bc_name,bc_name)) @@ -776,7 +831,10 @@ CCTK_INT Bdry_Boundary_SelectedGVs(CCTK_POINTER GH, printf(" current->next is %p\n",current->next); #endif + /* One should be able to pass NULL to ignore the values in any of these + arrays... */ ((int *) faces)[i] = current->faces; + ((int *) widths)[i] = current->width; ((int *) table_handles)[i] = current->table; ((int *) var_indices)[i] = current->var; } @@ -825,12 +883,13 @@ void Boundary_ApplyPhysicalBCs(CCTK_ARGUMENTS) { DECLARE_CCTK_ARGUMENTS; phys_bc_fn_ptr bc_fn; - int num_vars, *vars, *faces, *tables, err, max_num_vars; + int num_vars, *vars, *faces, *widths, *tables, err, max_num_vars; struct BCDATA *current_bcdata; max_num_vars = 0; vars = NULL; /* so that it won't be freed if it was never malloced */ faces = NULL; /* avoids a compiler warning */ + widths = NULL; /* avoids a compiler warning */ tables = NULL; /* avoids a compiler warning */ /* Step through each requested physical boundary condition */ @@ -856,17 +915,19 @@ void Boundary_ApplyPhysicalBCs(CCTK_ARGUMENTS) /* reallocate arrays if necessary */ vars = (int *) realloc(vars, num_vars*sizeof(int)); faces = (int *) realloc(faces, num_vars*sizeof(int)); + widths = (int *) realloc(widths, num_vars*sizeof(int)); tables = (int *) realloc(tables, num_vars*sizeof(int)); } else { vars = (int *) malloc(num_vars*sizeof(int)); faces = (int *) malloc(num_vars*sizeof(int)); + widths = (int *) malloc(num_vars*sizeof(int)); tables = (int *) malloc(num_vars*sizeof(int)); } } /* get selected vars for this bc_name */ - err = Bdry_Boundary_SelectedGVs(cctkGH, num_vars, vars, faces, tables, + err = Bdry_Boundary_SelectedGVs(cctkGH, num_vars, vars, faces, widths, tables, current_bcdata->bc_name); if (err<0) /* This is a redundant test for now, Bdry_Boundary_SelectedGVs never returns <0 */ @@ -889,7 +950,7 @@ void Boundary_ApplyPhysicalBCs(CCTK_ARGUMENTS) #endif err = Util_TableGetFnPointer(physbc_table_handle, (CCTK_FPOINTER *)&bc_fn, current_bcdata->bc_name); - if (err<0) + if (err<0) { CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING, "Boundary_ApplyPhysicalBCs: Util_TableGetFnPointer " @@ -903,7 +964,7 @@ void Boundary_ApplyPhysicalBCs(CCTK_ARGUMENTS) " cctkGH %p, num_vars %d, vars, tables\n", (void *) bc_fn, (const void *) cctkGH, num_vars); #endif - err = (*bc_fn)(cctkGH, num_vars, vars, faces, tables); + err = (*bc_fn)(cctkGH, num_vars, vars, faces, widths, tables); if (err<0) { CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING, @@ -917,6 +978,7 @@ void Boundary_ApplyPhysicalBCs(CCTK_ARGUMENTS) { free(vars); free(faces); + free(widths); free(tables); } } @@ -946,7 +1008,7 @@ void Boundary_ClearSelection(void) current_bcdata = current_bcdata->next) { #ifdef DEBUG - printf("Boundary_ClearSelection: freeing %p\n",current_bcdata); + printf("Boundary_ClearSelection: freeing var list rooted at %p\n",current_bcdata); #endif /* Free selections list */ @@ -958,6 +1020,7 @@ void Boundary_ClearSelection(void) } current_bcdata->var_list = NULL; + current_bcdata->num = 0; } } @@ -1042,15 +1105,16 @@ static void print_selections_database(void) struct BCDATA *current_bcdata; struct BCVAR *current; - printf(" Current list of selected vars:\n"); - printf(" bc name vi var name table handle\n"); + printf("Current list of selected vars:\n"); for (current_bcdata = bcdata_list; current_bcdata; current_bcdata = current_bcdata->next) { + printf("%d entries for %s:\n", current_bcdata->num, current_bcdata->bc_name); + printf(" vi var name table handle\n"); for (current = current_bcdata->var_list; current; current = current->next) { - printf(" %9s %3d %12s %2d\n", current_bcdata->bc_name, current->var, + printf("%3d %12s %2d\n", current->var, CCTK_VarName(current->var), current->table); } } diff --git a/src/Boundary.h b/src/Boundary.h index 99de680..1fbe1ff 100644 --- a/src/Boundary.h +++ b/src/Boundary.h @@ -21,18 +21,18 @@ extern "C" condition: */ typedef CCTK_INT (*phys_bc_fn_ptr)(const CCTK_POINTER, const CCTK_INT, const CCTK_INT *, const CCTK_INT *, - const CCTK_INT *); + const CCTK_INT *, const CCTK_INT *); /* prototype for routine registed as providing 'None' boundary condition */ CCTK_INT BndNone(const cGH *GH, CCTK_INT num_vars, CCTK_INT *var_indicies, - CCTK_INT *faces, CCTK_INT *table_handle); + CCTK_INT *faces, CCTK_INT *boundary_widths, CCTK_INT *table_handle); /* Scalar boundaries */ /* prototype for routine registed as providing 'Scalar' boundary condition */ CCTK_INT BndScalar(const cGH *GH, CCTK_INT num_vars, CCTK_INT *var_indicies, - CCTK_INT *faces, CCTK_INT *table_handle); + CCTK_INT *faces, CCTK_INT *boundary_widths, CCTK_INT *table_handle); int BndScalarDirGI (const cGH *GH, int stencil_size, @@ -77,7 +77,7 @@ int BndScalarVN (const cGH *GH, /* prototype for routine registed as providing 'Copy' boundary condition */ CCTK_INT BndCopy(const cGH *GH, CCTK_INT num_vars, CCTK_INT *var_indicies, - CCTK_INT *faces, CCTK_INT *table_handle); + CCTK_INT *faces, CCTK_INT *boundary_widths, CCTK_INT *table_handle); int BndCopyDirGI (const cGH *GH, int stencil_size, @@ -121,7 +121,7 @@ int BndCopyVN (const cGH *GH, /* prototype for routine registed as providing 'Static' boundary condition */ CCTK_INT BndStatic(const cGH *GH, CCTK_INT num_vars, CCTK_INT *var_indicies, - CCTK_INT *faces, CCTK_INT *table_handle); + CCTK_INT *faces, CCTK_INT *boundary_widths, CCTK_INT *table_handle); int BndStaticDirGI (const cGH *GH, int stencil_size, @@ -158,7 +158,8 @@ int BndStaticVN (const cGH *GH, /* prototype for routine registed as providing 'Radiative' boundary conditions */ CCTK_INT BndRadiative(const cGH *GH, CCTK_INT num_vars, CCTK_INT *var_indicies, - CCTK_INT *faces, CCTK_INT *table_handle); + CCTK_INT *faces, CCTK_INT *boundary_widths, + CCTK_INT *table_handle); int BndRadiativeDirGI (const cGH *GH, int stencil_size, @@ -219,7 +220,7 @@ int BndRadiativeVN (const cGH *GH, /* prototype for routine registed as providing 'Robin' boundary condition */ CCTK_INT BndRobin(const cGH *GH, CCTK_INT num_vars, CCTK_INT *var_indicies, - CCTK_INT *faces, CCTK_INT *table_handle); + CCTK_INT *faces, CCTK_INT *boundary_widths, CCTK_INT *table_handle); int BndRobinGI (const cGH *GH, const int *stencil, @@ -247,7 +248,7 @@ int BndRobinVN (const cGH *GH, /* prototype for routine registed as providing 'Flat' boundary condition */ CCTK_INT BndFlat(const cGH *GH, CCTK_INT num_vars, CCTK_INT *var_indicies, - CCTK_INT *faces, int *table_handle); + CCTK_INT *faces, CCTK_INT *boundary_widths, CCTK_INT *table_handle); int BndFlatDirGI (const cGH *GH, int stencil_size, diff --git a/src/CopyBoundary.c b/src/CopyBoundary.c index 8b3f85c..3b91eb5 100644 --- a/src/CopyBoundary.c +++ b/src/CopyBoundary.c @@ -41,7 +41,6 @@ static int ApplyBndCopy (const cGH *GH, /******************************************************************** ******************** External Routines ************************ ********************************************************************/ - /*@@ @routine BndCopy @date 13 Feb 2003 @@ -51,6 +50,14 @@ static int ApplyBndCopy (const cGH *GH, the Copy boundary condition @enddesc @calls ApplyBndCopy + CCTK_GroupDimFromVarI + Util_TableGetIntArray + Util_TableQueryValueInfo + CCTK_VWarn + Util_TableGetString + CCTK_VarIndex + Util_TableGetInt + @var GH @vdesc Pointer to CCTK grid hierarchy @vtype const cGH * @@ -72,6 +79,11 @@ static int ApplyBndCopy (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 @@ -79,32 +91,40 @@ static int ApplyBndCopy (const cGH *GH, @endvar @returntype int @returndesc - return code of @seeroutine ApplyBndScalar + return code of @seeroutine ApplyBndCopy -11 invalid table handle -12 no "COPY FROM" key in table + -21 error reading boundary width array from table + -22 wrong size boundary width array in table @endreturndesc @@*/ -int BndCopy(const cGH *GH, int num_vars, int *vars, int *faces, int *tables) +int BndCopy(const cGH *GH, int num_vars, int *vars, int *faces, int *widths, + int *tables) { - int i, j, k, gdim, max_gdim, err, value_type, value_size, retval; + int i, j, k, gi, gdim, max_gdim, err, value_type, value_size, retval; char *copy_from_name; /* variables to pass to ApplyBndCopy */ - /*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 copy_from; /* variable (index) from which to copy the boundary data */ + int *width_alldirs; /* width of boundary on each face */ + int dir; /* direction in which to apply bc */ + int copy_from; /* variable (index) from which to copy the boundary + data */ - 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 max_gdim) - { - realloc(stencil_alldirs, gdim*sizeof(int)); - max_gdim = gdim; - } - for (k=0; k 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 + { + for (k=0; k<2*gdim; ++k) + { + width_alldirs[k] = widths[i]; } } - if (!retval && (retval = ApplyBndCopy(GH, 0, stencil_alldirs, dir, + /* Apply the boundary condition */ + if (!retval && (retval = ApplyBndCopy(GH, 0, width_alldirs, dir, vars[i], copy_from, j)) < 0) { CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING, "ApplyBndCopy() returned %d", retval); } } - free(stencil_alldirs); + free(width_alldirs); return retval; } @@ -371,13 +409,21 @@ int BndCopyVI (const cGH *GH, int vi_to, int vi_from) { - int retval, num_vars; + int retval, num_vars, dim, i; + int *boundary_widths; /* boundary widths as expected by ApplyBndCopy */ + /* Set up boundary_widths array */ + dim = CCTK_GroupDimFromVarI (vi_to); + boundary_widths = (int *) malloc(2*dim*sizeof(int)); + for (i=0; i<2*dim; i+=2) + { + boundary_widths[i] = stencil[i/2]; + } num_vars = CCTK_NumVars (); if (vi_to >= 0 && vi_to < num_vars && vi_from >= 0 && vi_from < num_vars) { - retval = ApplyBndCopy (GH, -1, stencil, 0, vi_to, vi_from, 1); + retval = ApplyBndCopy (GH, -1, boundary_widths, 0, vi_to, vi_from, 1); } else { @@ -525,14 +571,22 @@ int BndCopyGI (const cGH *GH, int gi_from) { int first_vi_to, first_vi_from, retval; + int i, dim, *boundary_widths; + /* Set up boundary_widths array */ + dim = CCTK_GroupDimI (gi_to); + boundary_widths = (int *) malloc(2*dim*sizeof(int)); + for (i=0; i<2*dim; i+=2) + { + boundary_widths[i] = stencil[i/2]; + } first_vi_to = CCTK_FirstVarIndexI (gi_to); first_vi_from = CCTK_FirstVarIndexI (gi_from); if (first_vi_to >= 0 && first_vi_from >= 0) { - retval = ApplyBndCopy (GH, -1, stencil, 0, first_vi_to, first_vi_from, - CCTK_NumVarsInGroupI (gi_to)); + retval = ApplyBndCopy (GH, -1, boundary_widths, 0, first_vi_to, + first_vi_from, CCTK_NumVarsInGroupI (gi_to)); } else { @@ -683,15 +737,23 @@ int BndCopyGN (const cGH *GH, const char *gname_from) { int gi_to, gi_from, num_groups, retval; - + int i, dim, *boundary_widths; gi_to = CCTK_GroupIndex (gname_to); gi_from = CCTK_GroupIndex (gname_from); num_groups = CCTK_NumGroups (); + /* Set up boundary_widths array */ + dim = CCTK_GroupDimI (gi_to); + boundary_widths = (int *) malloc(2*dim*sizeof(int)); + for (i=0; i<2*dim; i+=2) + { + boundary_widths[i] = stencil[i/2]; + } + if (gi_to >= 0 && gi_to < num_groups && gi_from >= 0 && gi_from < num_groups) { - retval = BndCopyGI (GH, stencil, gi_to, gi_from); + retval = BndCopyGI (GH, boundary_widths, gi_to, gi_from); } else { @@ -844,15 +906,23 @@ int BndCopyVN (const cGH *GH, const char *vname_from) { int vi_to, vi_from, num_vars, retval; - + int i, dim, *boundary_widths; vi_to = CCTK_VarIndex (vname_to); vi_from = CCTK_VarIndex (vname_from); num_vars = CCTK_NumVars (); + /* Set up boundary_widths array */ + dim = CCTK_GroupDimFromVarI(vi_to); + boundary_widths = (int *) malloc(2*dim*sizeof(int)); + for (i=0; i<2*dim; i+=2) + { + boundary_widths[i] = stencil[i/2]; + } + if (vi_to >= 0 && vi_to < num_vars && vi_from >= 0 && vi_from < num_vars) { - retval = BndCopyVI (GH, stencil, vi_to, vi_from); + retval = BndCopyVI (GH, boundary_widths, vi_to, vi_from); } else { @@ -1006,12 +1076,12 @@ void CCTK_FCALL CCTK_FNAME (BndCopyVN) 0 for success -1 if dimension is not supported -2 if direction parameter is invalid - -3 if stencil width array parameter is NULL + -3 if boundary width array parameter is NULL @endreturndesc @@*/ static int ApplyBndCopy (const cGH *GH, - int stencil_dir, - const int *stencil_alldirs, + int width_dir, + const int *in_widths, int dir, int first_var_to, int first_var_from, @@ -1021,7 +1091,7 @@ static int ApplyBndCopy (const cGH *GH, int timelvl_to, timelvl_from; int gindex, gdim; int var_to, var_from, vtypesize; - 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; @@ -1052,15 +1122,17 @@ static int ApplyBndCopy (const cGH *GH, /* set up stencil 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, "ApplyBndCopy: NULL pointer passed for stencil width array"); + CCTK_WARN (1, "ApplyBndCopy: NULL pointer passed for boundary width " + "array"); return (-3); } @@ -1118,23 +1190,23 @@ static int ApplyBndCopy (const cGH *GH, if (gdim > 0) { /* lower x */ - COPY_BOUNDARY (doBC[0], stencil[0], lssh[1], lssh[2], i, j, k); + COPY_BOUNDARY (doBC[0], widths[0], lssh[1], lssh[2], i, j, k); /* upper x */ - COPY_BOUNDARY (doBC[1], stencil[0], lssh[1], lssh[2], lssh[0]-i-1, j, k); + COPY_BOUNDARY (doBC[1], widths[1], lssh[1], lssh[2], lssh[0]-i-1, j, k); } if (gdim > 1) { /* lower y */ - COPY_BOUNDARY (doBC[2], lssh[0], stencil[1], lssh[2], i, j, k); + COPY_BOUNDARY (doBC[2], lssh[0], widths[2], lssh[2], i, j, k); /* upper y */ - COPY_BOUNDARY (doBC[3], lssh[0], stencil[1], lssh[2], i, lssh[1]-j-1, k); + COPY_BOUNDARY (doBC[3], lssh[0], widths[3], lssh[2], i, lssh[1]-j-1, k); } if (gdim > 2) { /* lower z */ - COPY_BOUNDARY (doBC[4], lssh[0], lssh[1], stencil[2], i, j, k); + COPY_BOUNDARY (doBC[4], lssh[0], lssh[1], widths[4], i, j, k); /* upper z */ - COPY_BOUNDARY (doBC[5], lssh[0], lssh[1], stencil[2], i, j, lssh[2]-k-1); + COPY_BOUNDARY (doBC[5], lssh[0], lssh[1], widths[5], i, j, lssh[2]-k-1); } } diff --git a/src/FlatBoundary.c b/src/FlatBoundary.c index 0ac4727..53c2322 100644 --- a/src/FlatBoundary.c +++ b/src/FlatBoundary.c @@ -38,7 +38,12 @@ static int ApplyBndFlat (const cGH *GH, int dir, int first_var, int num_vars); - +static int OldApplyBndFlat (const cGH *GH, + int stencil_dir, + const int *stencil_alldirs, + int dir, + int first_var, + int num_vars); /******************************************************************** ******************** External Routines ************************ @@ -52,6 +57,7 @@ static int ApplyBndFlat (const cGH *GH, the Flat boundary condition @enddesc @calls ApplyBndFlat + @var GH @vdesc Pointer to CCTK grid hierarchy @vtype const cGH * @@ -73,6 +79,11 @@ static int ApplyBndFlat (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 +91,35 @@ static int ApplyBndFlat (const cGH *GH, @endvar @returntype int @returndesc - return code of @seeroutine ApplyBndScalar + return code of @seeroutine ApplyBndFlat + -21 error reading boundary width array from table + -22 wrong size boundary width array in table @endreturndesc @@*/ -int BndFlat(const cGH *GH, int num_vars, int *vars, int *faces, int *tables) +int BndFlat(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 ApplyBndFlat */ - /*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 boundary in all directions */ + int dir; /* direction in which to apply bc */ - 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 max_gdim) { - realloc(stencil_alldirs, gdim*sizeof(int)); + width_alldirs = realloc(width_alldirs, gdim*sizeof(int)); max_gdim = gdim; } - for (k=0; k= 0) { - retval = ApplyBndFlat (GH, -1, stencil, 0, first_vi, + retval = OldApplyBndFlat (GH, -1, stencil, 0, first_vi, CCTK_NumVarsInGroupI (gi)); } else @@ -598,7 +625,7 @@ int BndFlatVI (const cGH *GH, if (vi >= 0 && vi < CCTK_NumVars ()) { - retval = ApplyBndFlat (GH, -1, stencil, 0, vi, 1); + retval = OldApplyBndFlat (GH, -1, stencil, 0, vi, 1); } else { @@ -846,13 +873,13 @@ void CCTK_FCALL CCTK_FNAME (BndFlatVN) @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 @@ -889,19 +916,20 @@ void CCTK_FCALL CCTK_FNAME (BndFlatVN) 0 for success -1 if dimension is not supported -2 if direction parameter is invalid - -3 if stencil width array parameter is NULL + -3 if boundary width array parameter is NULL @endreturndesc @@*/ static int ApplyBndFlat (const cGH *GH, - int stencil_dir, - const int *stencil_alldirs, + int width_dir, + const int *in_widths, int dir, int first_var, int num_vars) { 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], + widths[2*MAXDIM]; SymmetryGHex *sGHex; @@ -929,18 +957,20 @@ static int ApplyBndFlat (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, "ApplyBndFlat: NULL pointer passed for stencil width array"); + CCTK_WARN (1, "ApplyBndFlat: NULL pointer passed for boundary width " + "array"); return (-3); } @@ -1005,11 +1035,11 @@ static int ApplyBndFlat (const cGH *GH, } #endif /* DEBUG_BOUNDARY */ /* lower x */ - FLAT_BOUNDARY (doBC[0], stencil[0], lssh[1], lssh[2], - i, j, k, stencil[0], j, k); + FLAT_BOUNDARY (doBC[0], widths[0], lssh[1], lssh[2], + i, j, k, widths[0], j, k); /* upper x */ - FLAT_BOUNDARY (doBC[1], stencil[0], lssh[1], lssh[2], - lssh[0]-i-1, j, k, lssh[0]-stencil[0]-1, j, k); + FLAT_BOUNDARY (doBC[1], widths[1], lssh[1], lssh[2], + lssh[0]-i-1, j, k, lssh[0]-widths[1]-1, j, k); } if (gdim > 1) @@ -1025,11 +1055,11 @@ static int ApplyBndFlat (const cGH *GH, } #endif /* DEBUG_BOUNDARY */ /* lower y */ - FLAT_BOUNDARY (doBC[2], lssh[0], stencil[1], lssh[2], - i, j, k, i, stencil[1], k); + FLAT_BOUNDARY (doBC[2], lssh[0], widths[2], lssh[2], + i, j, k, i, widths[2], k); /* upper y */ - FLAT_BOUNDARY (doBC[3], lssh[0], stencil[1], lssh[2], - i, lssh[1]-j-1, k, i, lssh[1]-stencil[1]-1, k); + FLAT_BOUNDARY (doBC[3], lssh[0], widths[3], lssh[2], + i, lssh[1]-j-1, k, i, lssh[1]-widths[3]-1, k); } if (gdim > 2) { @@ -1044,13 +1074,96 @@ static int ApplyBndFlat (const cGH *GH, } #endif /* DEBUG_BOUNDARY */ /* lower z */ - FLAT_BOUNDARY (doBC[4], lssh[0], lssh[1], stencil[2], - i, j, k, i, j, stencil[2]); + FLAT_BOUNDARY (doBC[4], lssh[0], lssh[1], widths[4], + i, j, k, i, j, widths[4]); /* upper z */ - FLAT_BOUNDARY (doBC[5], lssh[0], lssh[1], stencil[2], - i, j, lssh[2]-k-1, i, j, lssh[2]-stencil[2]-1); + FLAT_BOUNDARY (doBC[5], lssh[0], lssh[1], widths[5], + i, j, lssh[2]-k-1, i, j, lssh[2]-widths[5]-1); } } return(0); } + +/*@@ + @routine OldApplyBndFlat + @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 stencil width in direction dir + @vtype int + @vio in + @endvar + @var stencil_alldirs + @vdesc stencil 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 first_var + @vdesc index of first variable to apply boundaries to + @vtype int + @vio in + @endvar + @var num_vars + @vdesc number of variables + @vtype int + @vio in + @endvar + + @calls CCTK_GroupIndexFromVarI + CCTK_GroupDimI + ApplyBndFlat + @history + @returntype int + @returndesc + returncode from @seeroutine ApplyBndFlat + @endreturndesc +@@*/ + +int OldApplyBndFlat(const cGH *GH, int stencil_dir, const int *stencil_alldirs, + int dir, 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 = ApplyBndFlat(GH, stencil_dir, boundary_widths, dir, first_var, + num_vars); + + return retval; +} diff --git a/src/NoneBoundary.c b/src/NoneBoundary.c index 59d1267..ab2f294 100644 --- a/src/NoneBoundary.c +++ b/src/NoneBoundary.c @@ -64,6 +64,11 @@ CCTK_FILEVERSION(CactusBase_Boundary_NoneBoundary_c); @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 @@ -75,7 +80,7 @@ CCTK_FILEVERSION(CactusBase_Boundary_NoneBoundary_c); @endreturndesc @@*/ -int BndNone(const cGH *GH, int num_vars, int *var_indices, int *faces, +int BndNone(const cGH *GH, int num_vars, int *var_indices, int *faces, int *widths, int *table_handles) { #ifdef DEBUG @@ -87,6 +92,7 @@ int BndNone(const cGH *GH, int num_vars, int *var_indices, int *faces, num_vars = num_vars; var_indices = var_indices; faces = faces; + widths = widths; table_handles = table_handles; return 0; diff --git a/src/RadiationBoundary.c b/src/RadiationBoundary.c index 4df443b..88d6b1e 100644 --- a/src/RadiationBoundary.c +++ b/src/RadiationBoundary.c @@ -95,6 +95,15 @@ static int ApplyBndRadiative (const cGH *GH, int first_var_to, int first_var_from, int num_vars); +static int OldApplyBndRadiative (const cGH *GH, + int stencil_dir, + const int *stencil_alldirs, + int dir, + CCTK_REAL var0, + CCTK_REAL speed, + int first_var_to, + int first_var_from, + int num_vars); /******************************************************************** ******************** External Routines ************************ @@ -106,9 +115,10 @@ static int ApplyBndRadiative (const cGH *GH, @author David Rideout @desc Top level function which is registered as handling - this boundary condition + the Radiative boundary condition @enddesc @calls ApplyBndRadiative + @var GH @vdesc Pointer to CCTK grid hierarchy @vtype const cGH * @@ -130,6 +140,11 @@ static int ApplyBndRadiative (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 @@ -137,23 +152,23 @@ static int ApplyBndRadiative (const cGH *GH, @endvar @returntype int @returndesc - return code of @seeroutine ApplyBndRadiative
+ return code of @seeroutine ApplyBndRadiative + -21 error reading boundary width array from table + -22 wrong size boundary width array in table @endreturndesc @@*/ -int BndRadiative(const cGH *GH, int num_vars, int *vars, int *faces, - int *tables) +int BndRadiative(const cGH *GH, int num_vars, int *vars, int *faces, + int *widths, int *tables) { - int i, j, k, value_type, value_size, err, - gdim, max_gdim, retval; + int i, j, k, gi, gdim, max_gdim, value_type, value_size, err, retval; char *prev_time_level_name; /* variables to pass to ApplyBndRadiative */ - /*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 boundary in all directions */ + int dir; /* direction in which to apply bc */ CCTK_REAL limit, speed; - int prev_time_level; /* variable (index) which holds the previous time level */ + int prev_time_level; /* variable index which holds the previous time level */ #ifdef DEBUG printf("BndRadiative() got passed: GH=%p, num_vars=%d:\n", (const void *) GH, @@ -168,15 +183,20 @@ int BndRadiative(const cGH *GH, int num_vars, int *vars, int *faces, /* CCTK_WARN(0, "stopping code"); */ #endif - 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 max_gdim) - { - realloc(stencil_alldirs, gdim*sizeof(int)); - max_gdim = gdim; - } - for (k=0; k 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 + { + for (k=0; k<2*gdim; ++k) + { + width_alldirs[k] = widths[i]; + } + } + + /* Apply the boundary condition */ + if ((retval = ApplyBndRadiative(GH, 0, width_alldirs, dir, limit, speed, vars[i], prev_time_level, j)) < 0) { CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING, @@ -270,7 +313,7 @@ int BndRadiative(const cGH *GH, int num_vars, int *vars, int *faces, #ifdef DEBUG printf("BndRadiative(): returning %d\n",retval); #endif - free(stencil_alldirs); + free(width_alldirs); return retval; } @@ -502,7 +545,7 @@ int BndRadiativeGI (const cGH *GH, first_vi_from = CCTK_FirstVarIndexI (gi_from); if (first_vi_to >= 0 && first_vi_from >= 0) { - retval = ApplyBndRadiative (GH, -1, stencil, 0, var0, speed, + retval = OldApplyBndRadiative (GH, -1, stencil, 0, var0, speed, first_vi_to, first_vi_from, CCTK_NumVarsInGroupI (gi_to)); } @@ -872,7 +915,7 @@ int BndRadiativeVI (const cGH *GH, num_vars = CCTK_NumVars (); if (vi_to >= 0 && vi_to < num_vars && vi_from >= 0 && vi_from < num_vars) { - retval = ApplyBndRadiative (GH, -1, stencil, 0, var0, speed, + retval = OldApplyBndRadiative (GH, -1, stencil, 0, var0, speed, vi_to, vi_from, 1); } else @@ -1408,7 +1451,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) /* Upper x-bound */ \ if (doBC[1]) \ { \ - UPPER_RADIATIVE_BOUNDARY_3D (lsh[0] - stencil[0], 0, 0, \ + UPPER_RADIATIVE_BOUNDARY_3D (lsh[0] - stencil[1], 0, 0, \ coords[0], coords[MAXDIM], offset[0], \ var_to, var_from, cctk_type); \ } \ @@ -1416,7 +1459,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) /* Lower y-bound */ \ if (doBC[2]) \ { \ - LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], stencil[1], lsh[2], \ + LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], stencil[2], lsh[2], \ coords[1], coords[MAXDIM], offset[1], \ var_to, var_from, cctk_type); \ } \ @@ -1424,7 +1467,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) /* Upper y-bound */ \ if (doBC[3]) \ { \ - UPPER_RADIATIVE_BOUNDARY_3D (0, lsh[1] - stencil[1], 0, \ + UPPER_RADIATIVE_BOUNDARY_3D (0, lsh[1] - stencil[3], 0, \ coords[1], coords[MAXDIM], offset[1], \ var_to, var_from, cctk_type); \ } \ @@ -1432,7 +1475,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) /* Lower z-bound */ \ if (doBC[4]) \ { \ - LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], lsh[1], stencil[2], \ + LOWER_RADIATIVE_BOUNDARY_3D (lsh[0], lsh[1], stencil[4], \ coords[2], coords[MAXDIM], offset[2], \ var_to, var_from, cctk_type); \ } \ @@ -1440,7 +1483,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) /* Upper z-bound */ \ if (doBC[5]) \ { \ - UPPER_RADIATIVE_BOUNDARY_3D (0, 0, lsh[2] - stencil[2], \ + UPPER_RADIATIVE_BOUNDARY_3D (0, 0, lsh[2] - stencil[5], \ coords[2], coords[MAXDIM], offset[2], \ var_to, var_from, cctk_type); \ } \ @@ -1466,13 +1509,13 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) @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 @@ -1506,7 +1549,7 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) @vtype int @vio in @endvar - CCTK_VarTypeI + @calls CCTK_VarTypeI CCTK_GroupDimFromVarI RADIATIVE_BOUNDARY @history @@ -1521,15 +1564,15 @@ void CCTK_FCALL CCTK_FNAME (BndRadiativeVN) 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 -5 if variable dimension is other than 3D -6 if a coordinate is not found @endreturndesc @@*/ static int ApplyBndRadiative (const cGH *GH, - int stencil_dir, - const int *stencil_alldirs, + int width_dir, + const int *in_widths, int dir, CCTK_REAL var0, CCTK_REAL speed, @@ -1544,10 +1587,9 @@ static int ApplyBndRadiative (const cGH *GH, SymmetryGHex *sGHex; char coord_system_name[10]; CCTK_REAL dxyz[MAXDIM], rho[MAXDIM], *coords[MAXDIM+1]; - int doBC[2*MAXDIM], stencil[MAXDIM], offset[MAXDIM]; + int doBC[2*MAXDIM], widths[2*MAXDIM], offset[MAXDIM]; CCTK_REAL dtv, dtvh, dtvvar0, dtvvar0H; - /* check the direction parameter */ if (abs (dir) > MAXDIM) { @@ -1569,19 +1611,20 @@ static int ApplyBndRadiative (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, "ApplyBndRadiative: NULL pointer passed " - "for stencil width array"); + "for boundary width array"); return (-3); } @@ -1674,26 +1717,26 @@ static int ApplyBndRadiative (const cGH *GH, switch (CCTK_VarTypeI (var_to)) { case CCTK_VARIABLE_CHAR: - RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, stencil, coords, offset, + RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, 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, stencil, coords, offset, + RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, 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, stencil, coords, offset, + RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, 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, stencil, coords, offset, + RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, GH->data[var_to][timelvl_to], GH->data[var_from][timelvl_from], gdim, CCTK_INT2); break; @@ -1701,7 +1744,7 @@ static int ApplyBndRadiative (const cGH *GH, #ifdef CCTK_INT4 case CCTK_VARIABLE_INT4: - RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, stencil, coords, offset, + RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, GH->data[var_to][timelvl_to], GH->data[var_from][timelvl_from], gdim, CCTK_INT4); break; @@ -1709,7 +1752,7 @@ static int ApplyBndRadiative (const cGH *GH, #ifdef CCTK_INT8 case CCTK_VARIABLE_INT8: - RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, stencil, coords, offset, + RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, GH->data[var_to][timelvl_to], GH->data[var_from][timelvl_from], gdim, CCTK_INT8); break; @@ -1717,7 +1760,7 @@ static int ApplyBndRadiative (const cGH *GH, #ifdef CCTK_REAL4 case CCTK_VARIABLE_REAL4: - RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, stencil, coords, offset, + RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, GH->data[var_to][timelvl_to], GH->data[var_from][timelvl_from], gdim, CCTK_REAL4); break; @@ -1725,7 +1768,7 @@ static int ApplyBndRadiative (const cGH *GH, #ifdef CCTK_REAL8 case CCTK_VARIABLE_REAL8: - RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, stencil, coords, offset, + RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, GH->data[var_to][timelvl_to], GH->data[var_from][timelvl_from], gdim, CCTK_REAL8); break; @@ -1733,7 +1776,7 @@ static int ApplyBndRadiative (const cGH *GH, #ifdef CCTK_REAL16 case CCTK_VARIABLE_REAL16: - RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, stencil, coords, offset, + RADIATIVE_BOUNDARY (GH->cctk_lsh, doBC, widths, coords, offset, GH->data[var_to][timelvl_to], GH->data[var_from][timelvl_from], gdim,CCTK_REAL16); break; @@ -1749,3 +1792,107 @@ static int ApplyBndRadiative (const cGH *GH, return (0); } + + +/*@@ + @routine OldApplyBndRadiative + @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 width_dir + @vdesc boundary width in direction dir + @vtype int + @vio in + @endvar + @var in_widths + @vdesc boundary widths for all directions + @vtype int [ dimension of variable(s) ] + @vio in + @endvar + @var dir + @vdesc direction to copy boundaries (0 for copying all directions) + @vtype int + @vio in + @endvar + @var var0 + @vdesc asymptotic value of function at infinity + @vtype CCTK_REAL + @vio in + @endvar + @var speed + @vdesc wave speed + @vtype CCTK_REAL + @vio in + @endvar + @var first_var_to + @vdesc index of first variable to copy boundaries to + @vtype int + @vio in + @endvar + @var first_var_from + @vdesc index of first variable to copy boundaries from + @vtype int + @vio in + @endvar + @var num_vars + @vdesc number of variables + @vtype int + @vio in + @endvar + + @calls CCTK_GroupIndexFromVarI + ApplyBndRadiative + @returntype int + @returndesc + returncode from @seeroutine ApplyBndRadiative + @endreturndesc +@@*/ + +int OldApplyBndRadiative(const cGH *GH, + int width_dir, + const int *stencil_alldirs, + int dir, + CCTK_REAL var0, + CCTK_REAL speed, + int first_var_to, + int first_var_from, + int num_vars) +{ + int retval, *boundary_widths, dim, i; + static int warned; + + /* Convert stencil_alldirs to new format */ + dim = CCTK_GroupDimFromVarI(first_var_to); + 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 = ApplyBndRadiative(GH, width_dir, boundary_widths, dir, var0, speed, + first_var_to, first_var_from, num_vars); + + return retval; +} diff --git a/src/RobinBoundary.c b/src/RobinBoundary.c index 3d2d67e..7c8020c 100644 --- a/src/RobinBoundary.c +++ b/src/RobinBoundary.c @@ -37,6 +37,12 @@ static int ApplyBndRobin (const cGH *GH, int npow, int first_var, int num_vars); +static int OldApplyBndRobin (const cGH *GH, + const int *stencil, + CCTK_REAL finf, + int npow, + int first_var, + int num_vars); /******************************************************************** @@ -72,6 +78,11 @@ static int ApplyBndRobin (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 @@ -79,32 +90,40 @@ static int ApplyBndRobin (const cGH *GH, @endvar @returntype int @returndesc - return code of @seeroutine ApplyBndRobin
+ return code of @seeroutine ApplyBndRobin + -21 error reading boundary width array from table + -22 wrong size boundary width array in table @endreturndesc @@*/ -int BndRobin(const cGH *GH, int num_vars, int *vars, int *faces, int *tables) +int BndRobin(const cGH *GH, int num_vars, int *vars, int *faces, int *widths, + int *tables) { - int i, j, k, err, gdim, max_gdim, retval; + int i, j, k, gi, err, gdim, max_gdim, retval; /* variables to pass to ApplyBndRobin */ - int *stencil; /* width of stencil in all directions */ - CCTK_REAL finf; /* value of function at infinity */ - int npow; /* decay rate */ + int *width_alldirs; /* width of boundary in all directions */ + CCTK_REAL finf; /* value of function at infinity */ + int npow; /* decay rate */ #ifdef DEBUG printf("BndRobin(): got passed GH=%p, num_vars=%d, var_indices[0]=%d, table_handles[0]=%d\n", (const void *) GH, num_vars, var_indices[0], table_handles[0]); #endif - retval = 0; stencil = NULL; max_gdim = 0; + retval = 0; width_alldirs = NULL; max_gdim = 0; /* loop through variables, j at a time */ for (i=0; i max_gdim) - { - realloc(stencil, gdim*sizeof(int)); - max_gdim = gdim; - } - for (k=0; k 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 + { + for (k=0; k<2*gdim; ++k) + { + width_alldirs[k] = widths[i]; + } + } + + /* Apply the boundary condition */ + if ((retval = ApplyBndRobin(GH, width_alldirs, finf, npow, vars[i], j)) < 0) { CCTK_VWarn(1, __LINE__, __FILE__, CCTK_THORNSTRING, "ApplyBndRobin() returned %d", retval); @@ -169,7 +210,7 @@ int BndRobin(const cGH *GH, int num_vars, int *vars, int *faces, int *tables) #ifdef DEBUG printf("BndRobin(): returning %d\n",retval); #endif - free(stencil); + free(width_alldirs); return retval; } @@ -263,7 +304,7 @@ int BndRobinGI (const cGH *GH, first_vi = CCTK_FirstVarIndexI (gi); if (first_vi >= 0) { - retval = ApplyBndRobin (GH, stencil, finf, npow, first_vi, + retval = OldApplyBndRobin (GH, stencil, finf, npow, first_vi, CCTK_NumVarsInGroupI (gi)); } else @@ -419,7 +460,7 @@ int BndRobinVI (const cGH *GH, if (vi >= 0 && vi < CCTK_NumVars ()) { - retval = ApplyBndRobin (GH, stencil, finf, npow, vi, 1); + retval = OldApplyBndRobin (GH, stencil, finf, npow, vi, 1); } else { @@ -712,8 +753,8 @@ void CCTK_FCALL CCTK_FNAME (BndRobinVN) @vtype const cGH * @vio in @endvar - @var stencil - @vdesc stencil width array + @var in_widths + @vdesc boundary width array @vtype int [ dimension of variable ] @vio in @endvar @@ -752,7 +793,7 @@ void CCTK_FCALL CCTK_FNAME (BndRobinVN) @returndesc 0 for success -1 if variable dimension is not supported - -2 if NULL pointer passed as stencil width array + -2 if NULL pointer passed as boundary width array -3 if stencil width is other than 1 -4 if variable type is not supported -5 if variable dimension is other than 3D @@ -760,7 +801,7 @@ void CCTK_FCALL CCTK_FNAME (BndRobinVN) @endreturndesc @@*/ static int ApplyBndRobin (const cGH *GH, - const int *stencil, + const int *in_widths, CCTK_REAL finf, int npow, int first_var, @@ -787,16 +828,17 @@ static int ApplyBndRobin (const cGH *GH, return (-1); } - /* check the stencil width */ - if (! stencil) + /* check the boundary width */ + if (! in_widths) { - CCTK_WARN (1, "ApplyBndRobin: NULL pointer passed for stencil width array"); + CCTK_WARN (1, "ApplyBndRobin: NULL pointer passed for boundary width " + "array"); return (-2); } - for (dim = 0; dim < gdim; dim++) + for (dim = 0; dim < 2*gdim; dim++) { - if (stencil[dim] != 1) + if (in_widths[dim] != 1) { CCTK_WARN (1, "ApplyBndRobin: Stencil width must be 1 " "for Robin boundary conditions"); @@ -917,3 +959,87 @@ static int ApplyBndRobin (const cGH *GH, return(0); } + + +/*@@ + @routine OldApplyBndRobin + @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 in_widths + @vdesc boundary width array + @vtype int [ dimension of variable ] + @vio in + @endvar + @var finf + @vdesc value of f at infinity + @vtype CCTK_REAL + @vio in + @endvar + @var npow + @vdesc power of decay rate + @vtype int + @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 + ApplyBndRobin + @returntype int + @returndesc + returncode from @seeroutine ApplyBndRobin + @endreturndesc +@@*/ +static int OldApplyBndRobin (const cGH *GH, + const int *in_widths, + CCTK_REAL finf, + int npow, + 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] = in_widths[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 = ApplyBndRobin(GH, boundary_widths, finf, npow, first_var, num_vars); + + return retval; +} 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 max_gdim) - { - realloc(stencil_alldirs, gdim*sizeof(int)); - max_gdim = gdim; - } - for (k=0; k 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; +} diff --git a/src/StaticBoundary.c b/src/StaticBoundary.c index 79027c5..e88415e 100644 --- a/src/StaticBoundary.c +++ b/src/StaticBoundary.c @@ -30,7 +30,12 @@ static int ApplyBndStatic (const cGH *GH, int dir, int first_var, int num_vars); - +static int OldApplyBndStatic (const cGH *GH, + int stencil_dir, + const int *stencil_alldirs, + int dir, + int first_var, + int num_vars); /******************************************************************** ******************** External Routines ************************ @@ -65,6 +70,11 @@ static int ApplyBndStatic (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 @@ -72,32 +82,38 @@ static int ApplyBndStatic (const cGH *GH, @endvar @returntype int @returndesc - return code of @seeroutine ApplyBndStatic
+ return code of @seeroutine ApplyBndStatic + -21 error reading boundary width array from table + -22 wrong size boundary width array in table @endreturndesc @@*/ -int BndStatic(const cGH *GH, int num_vars, int *vars, int *faces, int *tables) +int BndStatic(const cGH *GH, int num_vars, int *vars, int *faces, int *widths, + int *tables) { - int i, j, k, err, gdim, max_gdim, retval; + int i, j, k, gi, err, gdim, max_gdim, retval; /* variables to pass to ApplyBndStatic */ - /*int stencil_dir;*/ /* width of stencil in direction dir, not needed */ - int *stencil; /* width of stencil in all directions */ - int dir; /* direction in which to apply bc */ + int *width_alldirs; /* width of boundary in all directions */ + int dir; /* direction in which to apply bc */ #ifdef DEBUG printf("BndStatic(): got passed GH=%p, num_vars=%d, var_indices[0]=%d, table_handles[0]=%d\n", (const void *) GH, num_vars, var_indices[0], table_handles[0]); #endif - retval = 0; stencil = NULL; max_gdim = 0; + retval = 0; width_alldirs = NULL; max_gdim = 0; /* loop through variables, j at a time */ for (i=0; i max_gdim) { - realloc(stencil, gdim*sizeof(int)); + width_alldirs = realloc(width_alldirs, gdim*sizeof(int)); max_gdim = gdim; } - for (k=0; k= 0 && vi < num_vars) { - retval = ApplyBndStatic (GH, -1, stencil, 0, vi, 1); + retval = OldApplyBndStatic (GH, -1, stencil, 0, vi, 1); } else { @@ -461,7 +488,7 @@ int BndStaticGI (const cGH *GH, first_vi = CCTK_FirstVarIndexI (gi); if (first_vi >= 0) { - retval = ApplyBndStatic (GH, -1, stencil, 0, first_vi, + retval = OldApplyBndStatic (GH, -1, stencil, 0, first_vi, CCTK_NumVarsInGroupI (gi)); } else @@ -854,13 +881,13 @@ void CCTK_FCALL CCTK_FNAME (BndStaticVN) @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 @@ -902,8 +929,8 @@ void CCTK_FCALL CCTK_FNAME (BndStaticVN) @endreturndesc @@*/ static int ApplyBndStatic (const cGH *GH, - int stencil_dir, - const int *stencil_alldirs, + int width_dir, + const int *in_widths, int dir, int first_var, int num_vars) @@ -912,7 +939,8 @@ static int ApplyBndStatic (const cGH *GH, int timelvl_to, timelvl_from; int gindex, gdim; int var, vtypesize; - 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; /* Only apply boundary condition if more than one timelevel */ @@ -945,18 +973,20 @@ static int ApplyBndStatic (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, "ApplyBndStatic: NULL pointer passed for stencil width array"); + CCTK_WARN (1, "ApplyBndStatic: NULL pointer passed for boundary width " + "array"); return (-3); } @@ -1014,25 +1044,113 @@ static int ApplyBndStatic (const cGH *GH, if (gdim > 0) { /* lower x */ - STATIC_BOUNDARY (doBC[0], stencil[0], lssh[1], lssh[2], i, j, k); + STATIC_BOUNDARY (doBC[0], widths[0], lssh[1], lssh[2], i, j, k); /* upper x */ - STATIC_BOUNDARY (doBC[1], stencil[0], lssh[1], lssh[2], lssh[0]-i-1, j, k); + STATIC_BOUNDARY (doBC[1], widths[1], lssh[1], lssh[2], lssh[0]-i-1, j, + k); } if (gdim > 1) { /* lower y */ - STATIC_BOUNDARY (doBC[2], lssh[0], stencil[1], lssh[2], i, j, k); + STATIC_BOUNDARY (doBC[2], lssh[0], widths[2], lssh[2], i, j, k); /* upper y */ - STATIC_BOUNDARY (doBC[3], lssh[0], stencil[1], lssh[2], i, lssh[1]-j-1, k); + STATIC_BOUNDARY (doBC[3], lssh[0], widths[3], lssh[2], i, lssh[1]-j-1, + k); } if (gdim > 2) { /* lower z */ - STATIC_BOUNDARY (doBC[4], lssh[0], lssh[1], stencil[2], i, j, k); + STATIC_BOUNDARY (doBC[4], lssh[0], lssh[1], widths[4], i, j, k); /* upper z */ - STATIC_BOUNDARY (doBC[5], lssh[0], lssh[1], stencil[2], i, j, lssh[2]-k-1); + STATIC_BOUNDARY (doBC[5], lssh[0], lssh[1], widths[5], i, j, + lssh[2]-k-1); } } return(0); } + + +/*@@ + @routine OldApplyBndStatic + @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 for static boundaries (0 for copying all directions) + @vtype int + @vio in + @endvar + @var first_var + @vdesc index of first variable to apply static boundaries 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 OldApplyBndStatic (const cGH *GH, + int stencil_dir, + const int *stencil_alldirs, + int dir, + 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 = ApplyBndStatic(GH, stencil_dir, boundary_widths, dir, first_var, num_vars); + + return retval; +} -- cgit v1.2.3