From e68518078f4e0dec380c1bc1f8f7a836d52da17f Mon Sep 17 00:00:00 2001 From: eschnett Date: Thu, 23 Dec 2010 00:26:45 +0000 Subject: Allow direct calls to the low-level Periodic API. This can be used to apply periodic conditions only on selected grid variables. This is implemented by moving the "logic" (handling of errors and parameters) into the function Periodic_ApplyBC to make BndPeriodicVI "thorn agnostic". [Patch by David Radice] git-svn-id: http://svn.cactuscode.org/arrangements/CactusNumerical/Periodic/trunk@23 1bf05452-ddb3-4880-bfa1-00436340132b --- interface.ccl | 39 +++++++++++++++++++-- src/periodic.c | 105 ++++++++++++++++++++++++++++++++++++++------------------- src/periodic.h | 37 -------------------- 3 files changed, 108 insertions(+), 73 deletions(-) diff --git a/interface.ccl b/interface.ccl index cc228db..1b5497b 100644 --- a/interface.ccl +++ b/interface.ccl @@ -3,8 +3,6 @@ IMPLEMENTS: Periodic -INCLUDES HEADER: periodic.h IN Periodic.h - USES INCLUDE HEADER: Slab.h @@ -28,6 +26,43 @@ CCTK_INT FUNCTION Boundary_SelectedGVs(CCTK_POINTER_TO_CONST IN GH, \ CCTK_INT ARRAY OUT table_handles, CCTK_STRING IN bc_name) REQUIRES FUNCTION Boundary_SelectedGVs + +CCTK_INT FUNCTION Periodic_ApplyVI( \ + CCTK_POINTER_TO_CONST IN cctkGH, \ + CCTK_INT IN size, \ + CCTK_INT IN ARRAY stencil, \ + CCTK_INT IN ARRAY do_periodic, \ + CCTK_INT IN var_index \ + ) +PROVIDES FUNCTION Periodic_ApplyVI WITH BndPeriodicVI LANGUAGE C + +CCTK_INT FUNCTION Periodic_ApplyVN( \ + CCTK_POINTER_TO_CONST IN cctkGH, \ + CCTK_INT IN size, \ + CCTK_INT IN ARRAY stencil, \ + CCTK_INT IN ARRAY do_periodic, \ + CCTK_STRING IN var_name \ + ) +PROVIDES FUNCTION Periodic_ApplyVN WITH BndPeriodicVN LANGUAGE C + +CCTK_INT FUNCTION Periodic_ApplyGI( \ + CCTK_POINTER_TO_CONST IN cctkGH, \ + CCTK_INT IN size, \ + CCTK_INT IN ARRAY stencil, \ + CCTK_INT IN ARRAY do_periodic, \ + CCTK_INT IN group_index \ + ) +PROVIDES FUNCTION Periodic_ApplyGI WITH BndPeriodicGI LANGUAGE C + +CCTK_INT FUNCTION Periodic_ApplyGN( \ + CCTK_POINTER_TO_CONST IN cctkGH, \ + CCTK_INT IN size, \ + CCTK_INT IN ARRAY stencil, \ + CCTK_INT IN ARRAY do_periodic, \ + CCTK_STRING IN group_name \ + ) +PROVIDES FUNCTION Periodic_ApplyGN WITH BndPeriodicGN LANGUAGE C + CCTK_INT FUNCTION GetBoundarySpecification \ (CCTK_INT IN size, \ CCTK_INT OUT ARRAY nboundaryzones, \ diff --git a/src/periodic.c b/src/periodic.c index 1cfb274..b04f45a 100644 --- a/src/periodic.c +++ b/src/periodic.c @@ -6,7 +6,6 @@ #include "cctk.h" #include "cctk_Parameters.h" #include "Slab.h" -#include "periodic.h" static const char * restrict const rcsid = "$Header$"; CCTK_FILEVERSION(TAT_Periodic_periodic_c); @@ -14,18 +13,18 @@ CCTK_FILEVERSION(TAT_Periodic_periodic_c); int -BndPeriodicVI (cGH const * restrict const cctkGH, +BndPeriodicVI (CCTK_POINTER_TO_CONST _GH, + int size, int const * restrict const stencil, + int const do_periodic[3], int const vi) { - DECLARE_CCTK_PARAMETERS; - + cGH const * restrict const cctkGH = _GH; + cGroup group; cGroupDynamicData data; - char * restrict fullname; void * restrict varptr; struct xferinfo * restrict xferinfo; - int do_periodic[3]; int global_bbox[6]; int global_lbnd[3], global_ubnd[3]; int fake_bbox[6]; @@ -39,15 +38,6 @@ BndPeriodicVI (cGH const * restrict const cctkGH, assert (stencil); assert (vi>=0 && vi=0 && gi= 0); if (data.gsh[dir] < 2*stencil[dir]+1) { - CCTK_VWarn (0, __LINE__, __FILE__, CCTK_THORNSTRING, - "The group \"%s\" has in the %c-direction only %d grid points. " - "This is not large enough for a periodic boundary that is %d grid points wide. " - "The group needs to have at least %d grid points in that direction.", - CCTK_GroupNameFromVarI(vi), "xyz"[dir], data.gsh[dir], - stencil[dir], - 2*stencil[dir]+1); + return dir+1; } /* Loop over faces */ @@ -214,22 +199,30 @@ BndPeriodicVI (cGH const * restrict const cctkGH, int -BndPeriodicVN (cGH const * restrict const cctkGH, +BndPeriodicVN (CCTK_POINTER_TO_CONST _GH, + int size, int const * restrict const stencil, + int const do_periodic[3], char const * restrict const vn) { + cGH const * restrict const cctkGH = _GH; + int const vi = CCTK_VarIndex (vn); assert (vi>=0 && vi=0 && v1=0 && gi 0 && ierr < 4) { + int dir = ierr-1; + + int gi; + cGroup group; + cGroupDynamicData data; + + gi = CCTK_GroupIndexFromVarI (vi); + ierr = CCTK_GroupData (gi, &group); + assert (!ierr); + + ierr = CCTK_GroupDynamicData (cctkGH, gi, &data); + assert (!ierr); + + CCTK_VWarn (0, __LINE__, __FILE__, CCTK_THORNSTRING, + "The group \"%s\" has in the %c-direction only %d grid points. " + "This is not large enough for a periodic boundary that is %d grid points wide. " + "The group needs to have at least %d grid points in that direction.", + CCTK_GroupNameFromVarI(vi), "xyz"[dir], data.gsh[dir], + stencil[dir], + 2*stencil[dir]+1); + } + assert (!ierr); free (stencil); diff --git a/src/periodic.h b/src/periodic.h index 32a77e5..e69de29 100644 --- a/src/periodic.h +++ b/src/periodic.h @@ -1,37 +0,0 @@ -/* $Header$ */ - -#ifndef PERIODIC_H -#define PERIODIC_H - -#include "cctk.h" - -int -BndPeriodicVI (cGH const * restrict const cctkGH, - int const * restrict const stencil, - int const vi); - -int -BndPeriodicVN (cGH const * restrict const cctkGH, - int const * restrict const stencil, - char const * restrict const vn); - -int -BndPeriodicGI (cGH const * restrict const cctkGH, - int const * restrict const stencil, - int const gi); - -int -BndPeriodicGN (cGH const * restrict const cctkGH, - int const * restrict const stencil, - char const * restrict const gn); - - - -void -Periodic_RegisterBC (cGH * restrict const cctkGH); - -void -Periodic_ApplyBC (cGH const * restrict const cctkGH); - -#endif /* defined PERIODIC_H */ - -- cgit v1.2.3