reflevel = GetRefinementLevel(cctkGH) map = MultiPatch_GetMap(cctkGH) if (reflevel /= 0 .or. map == 0) return do j = 1, cctk_lsh(2) do i = 1, cctk_lsh(1) do k = cctk_lsh(3)-cctk_nghostzones(3)+1, cctk_lsh(3) if (Sn(i,j,k) == -2) then [my BCs go here] end if end do end do end doUhmm, probably not quite. The code should use the cctk_bbox(6) array which has entries 0 or 1 depending on whether a face is an outer (or symmetry) boundary or not.
The indexing is i = 1 + 2*idir + iface where idir = 0,1,2 and iface = 0,1 (lower, upper face).
The loops are kind of inefficient as you should have the "k" loop as the outer one (as usual). Using cctk_nghostzones(3) as the boundary width is not quite correct as technically you should use the boundary width returned by thorn CoordBase's GetBoundarySpecification aliased function (see its interface.ccl).
oh, but won't this be handled by the
if (Sn(i,j,k) == -2)
line above? from the description of the Sn function in the Interpolate2 thorn:
CCTK_INT source_patch TYPE=gf TAGS='Checkpoint="no" Prolongation="none"' { Sn } "source patch number, -1 for interior points, -2 for outer boundary points, -3 for inter- processor ghost points"
and also from the way this is done in LlamaWaveToy, i was assuming that Sn(i,j,k) == -2 would guarantee that the points belong to the outer boundary (ie, points where i want to specify a physical boundary condition). is this assumption not correct?
thanks, Miguel