From f094848d3eb13e236144a7fbed7fe255fb4fcf93 Mon Sep 17 00:00:00 2001 From: korobkin Date: Thu, 1 May 2008 08:20:48 +0000 Subject: (*) CCTK_FNAME(get_grid_offsets) -> get_grid_offsets (+) Added get_grid_offsets to thorn interface as GetFDGridOffsets (+) Added another thorn interface function GetLSHIndexRanges to inquire actual local index ranges (in Fortran notation), which account for both the ghostzones and shiftout values git-svn-id: https://svn.cct.lsu.edu/repos/numrel/LSUThorns/SummationByParts/trunk@109 f69c4107-0314-4c4f-9ad4-17e986b73f4a --- interface.ccl | 10 ++++++++++ src/call_derivs.c | 4 ++-- src/dissipation.c | 4 ++-- src/get_offset.c | 51 +++++++++++++++++++++++++++++++++++++++++++++++++-- 4 files changed, 63 insertions(+), 6 deletions(-) diff --git a/interface.ccl b/interface.ccl index d619180..381ac63 100644 --- a/interface.ccl +++ b/interface.ccl @@ -71,6 +71,16 @@ PROVIDES FUNCTION GetScalProdDiag WITH SBP_GetScalProdDiag LANGUAGE Fortran CCTK_REAL FUNCTION GetScalProdCoeff () PROVIDES FUNCTION GetScalProdCoeff WITH GetCoeff LANGUAGE Fortran +SUBROUTINE GetFDGridOffsets ( CCTK_INT OUT ARRAY offset ) +PROVIDES FUNCTION GetFDGridOffsets WITH get_grid_offsets LANGUAGE C + +SUBROUTINE GetLSHIndexRanges ( \ + CCTK_POINTER_TO_CONST IN cctkGH, \ + CCTK_INT OUT imn, CCTK_INT OUT imx, \ + CCTK_INT OUT jmn, CCTK_INT OUT jmx, \ + CCTK_INT OUT kmn, CCTK_INT OUT kmx ) +PROVIDES FUNCTION GetLSHIndexRanges WITH get_lsh_iranges LANGUAGE C + CCTK_INT FUNCTION GetDomainSpecification \ (CCTK_INT IN size, \ CCTK_REAL OUT ARRAY physical_min, \ diff --git a/src/call_derivs.c b/src/call_derivs.c index c92dd34..8b764fb 100644 --- a/src/call_derivs.c +++ b/src/call_derivs.c @@ -126,7 +126,7 @@ void DiffGv ( const CCTK_POINTER_TO_CONST cctkGH_, const CCTK_INT dir, const CCTK_REAL *delta, CCTK_REAL *dvar); - void CCTK_FCALL CCTK_FNAME(get_grid_offsets) (CCTK_INT *offset); + void get_grid_offsets (CCTK_INT *offset); ni = cctk_lsh[0]; nj = cctk_lsh[1]; nk = cctk_lsh[2]; @@ -162,7 +162,7 @@ void DiffGv ( const CCTK_POINTER_TO_CONST cctkGH_, const CCTK_INT dir, } if ( use_shiftout ) { /* get values of boundary_shiftout_* from CoordBase */ - CCTK_FNAME(get_grid_offsets) (offset); + get_grid_offsets (offset); } diff --git a/src/dissipation.c b/src/dissipation.c index 87f303d..a3b2bb6 100644 --- a/src/dissipation.c +++ b/src/dissipation.c @@ -145,7 +145,7 @@ void CCTK_FCALL CCTK_FNAME(dissipation_8_4_alt) (const CCTK_REAL *var, const CCTK_REAL *dx, const CCTK_REAL *epsdis, CCTK_REAL *rhs); -void CCTK_FCALL CCTK_FNAME(get_grid_offsets) (CCTK_INT *offset); +void get_grid_offsets (CCTK_INT *offset); static void apply (int const varindex, char const * const optstring, void * const arg); @@ -193,7 +193,7 @@ apply (int const varindex, char const * const optstring, void * const arg) } if ( use_shiftout ) { /* get values of boundary_shiftout_* from CoordBase */ - CCTK_FNAME(get_grid_offsets) (offset); + get_grid_offsets (offset); } SBP_determine_onesided_stencil (cctkGH, onesided); diff --git a/src/get_offset.c b/src/get_offset.c index d1a35d9..a3b8b90 100644 --- a/src/get_offset.c +++ b/src/get_offset.c @@ -3,7 +3,13 @@ #include "cctk_Functions.h" #include "cctk_Parameters.h" -void CCTK_FCALL CCTK_FNAME(get_grid_offsets) (CCTK_INT *offset) { +/*** + * This is a utility function which reads values of shiftout parameters from + * thorn CoordBase into an integer array offset[]. It is declared at the thorn + * interface as GetFDGridOffsets and can be called from Fortran routines. + * This is necessary since CCTK_ParameterGet is not available in Fortran. + */ +void get_grid_offsets (CCTK_INT *offset) { DECLARE_CCTK_PARAMETERS; @@ -34,5 +40,46 @@ void CCTK_FCALL CCTK_FNAME(get_grid_offsets) (CCTK_INT *offset) { offset[5] = *shiftout_z_upper; } - } + +/*** + * This function returns the effective values of local index ranges, + * taking into account ghost zones and boundary_shiftout_* values. + * If the use_shiftout=no, boundary offsets are set to zero. + * Index ranges are set in Fortran convention (starting from 1): + * - imin <= i <= imax + * - jmin <= j <= jmax + * - kmin <= k <= kmax + * This function is declared in interface.ccl as GetLSHIndexRanges + * to be used by other thorns, which need to know the actual index + * ranges for the region where SBP thorn will apply derivative / + * dissipation operators. + */ +void get_lsh_iranges ( const CCTK_POINTER_TO_CONST cctkGH_, + CCTK_INT *imin, CCTK_INT *imax, + CCTK_INT *jmin, CCTK_INT *jmax, + CCTK_INT *kmin, CCTK_INT *kmax) +{ + cGH const * restrict const cctkGH = cctkGH_; + DECLARE_CCTK_PARAMETERS + DECLARE_CCTK_ARGUMENTS + + int i, offset[6]; + + for (i=0;i<6;i++) offset[i] = 0; + if (use_shiftout) get_grid_offsets (offset); + + *imin = 1 + ((cctk_bbox[0]) ? offset[0] : cctk_nghostzones[0]); + *jmin = 1 + ((cctk_bbox[2]) ? offset[2] : cctk_nghostzones[1]); + *kmin = 1 + ((cctk_bbox[4]) ? offset[4] : cctk_nghostzones[2]); + + *imax = cctk_lsh[0] - ((cctk_bbox[1]) ? offset[1] : cctk_nghostzones[0] ); + *jmax = cctk_lsh[1] - ((cctk_bbox[3]) ? offset[3] : cctk_nghostzones[1] ); + *kmax = cctk_lsh[2] - ((cctk_bbox[5]) ? offset[5] : cctk_nghostzones[2] ); + + /*** debug ***/ + CCTK_VInfo (CCTK_THORNSTRING, "imn-kmx = {%d, %d, %d, %d, %d, %d}", + *imin, *imax, *jmin, *jmax, *kmin, *kmax); +} + + -- cgit v1.2.3