diff options
-rw-r--r-- | Carpet/CarpetInterp2/src/fasterp.cc | 5 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc | 31 | ||||
-rw-r--r-- | Carpet/CarpetReduce/src/mask_coords.c | 12 |
3 files changed, 34 insertions, 14 deletions
diff --git a/Carpet/CarpetInterp2/src/fasterp.cc b/Carpet/CarpetInterp2/src/fasterp.cc index 084e3651d..61383fa69 100644 --- a/Carpet/CarpetInterp2/src/fasterp.cc +++ b/Carpet/CarpetInterp2/src/fasterp.cc @@ -598,7 +598,10 @@ namespace CarpetInterp2 { & gsh[0], & lower.AT(m)[0], & upper.AT(m)[0], & delta.AT(m)[0]); assert (not ierr); - delta.AT(m) /= Carpet::maxspacereflevelfact; + //delta.AT(m) /= Carpet::maxspacereflevelfact; + gh const * const hh = Carpet::vhh.AT(m); + ibbox const & baseext = hh->baseextent(Carpet::mglevel, 0); + delta.AT(m) /= baseext.stride(); idelta.AT(m) = 1.0 / delta.AT(m); if (veryverbose) { cout << "GetCoordRange[" << m << "]: lower=" << lower.AT(m) << " upper=" << upper.AT(m) << " delta=" << delta.AT(m) << endl; diff --git a/Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc b/Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc index bbd742f40..82ab24c48 100644 --- a/Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc +++ b/Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc @@ -289,7 +289,7 @@ namespace CarpetLib { } // check that divided differences do not change sign: if so go back to first order! - if (lV*rV < 0) + if (lV*rV <= 0) { // switch back to first order! res = 0; @@ -338,7 +338,7 @@ namespace CarpetLib { } // check that divided differences do not change sign: if so go back to first order! - if (V[0]*V[2] < 0) + if (V[0]*V[2] <= 0) { // switch back to first order! res = 0; @@ -389,6 +389,10 @@ namespace CarpetLib { } } + /*typedef coeffs1d<RT,1,di> coeffs; + for (ptrdiff_t i=coeffs::imin; i<coeffs::imax; ++i) { + res += coeffs::get(i) * interp0<T,1> (p + i*d1); + }*/ return res; } @@ -419,7 +423,6 @@ namespace CarpetLib { { static_assert (dj==0 or dj==1, "dj must be 0 or 1"); typedef typename typeprops<T>::real RT; - typedef coeffs1d<RT,ORDER,dj> coeffs; T res = typeprops<T>::fromreal (0); // get function values needed for the stencil T f[coeffs1d<RT,ORDER,dj>::maxncoeffs]; @@ -446,7 +449,7 @@ namespace CarpetLib { } // check that divided differences do not change sign: if so go back to first order! - if (lV*rV < 0) + if (lV*rV <= 0) { // switch back to first order! res = 0; @@ -493,7 +496,7 @@ namespace CarpetLib { } // check that divided differences do not change sign: if so go back to first order! - if (V[0]*V[2] < 0) + if (V[0]*V[2] <= 0) { // switch back to first order! res = 0; @@ -543,6 +546,11 @@ namespace CarpetLib { res += coeffs1::get(i) * f[i-coeffs1d<RT,ORDER,dj>::minimin]; } } + + /*typedef coeffs1d<RT,1,dj> coeffs; + for (ptrdiff_t i=coeffs::imin; i<coeffs::imax; ++i) { + res += coeffs::get(i) * interp1<T,1,di> (p + i*d2, d1); + }*/ return res; } @@ -575,7 +583,6 @@ namespace CarpetLib { { static_assert (dk==0 or dk==1, "dk must be 0 or 1"); typedef typename typeprops<T>::real RT; - typedef coeffs1d<RT,ORDER,dk> coeffs; T res = typeprops<T>::fromreal (0); // get function values needed for the stencil T f[coeffs1d<RT,ORDER,dk>::maxncoeffs]; @@ -602,7 +609,7 @@ namespace CarpetLib { } // check that divided differences do not change sign: if so go back to first order! - if (lV*rV < 0) + if (lV*rV <= 0) { // switch back to first order! res = 0; @@ -649,7 +656,7 @@ namespace CarpetLib { } // check that divided differences do not change sign: if so go back to first order! - if (V[0]*V[2] < 0) + if (V[0]*V[2] <= 0) { // switch back to first order! res = 0; @@ -690,7 +697,7 @@ namespace CarpetLib { } // check that result is reasonable! - if ((res - f[-coeffs1d<RT,ORDER,dj>::minimin-1+dk]) * (f[-coeffs1d<RT,ORDER,dk>::minimin+dk] - res) < 0) + if ((res - f[-coeffs1d<RT,ORDER,dk>::minimin-1+dk]) * (f[-coeffs1d<RT,ORDER,dk>::minimin+dk] - res) < 0) { res = 0; // switch back to first order @@ -700,6 +707,11 @@ namespace CarpetLib { } } + /*typedef coeffs1d<RT,1,dk> coeffs; + for (ptrdiff_t i=coeffs::imin; i<coeffs::imax; ++i) { + res += coeffs::get(i) * interp2<T,1,di,dj> (p + i*d3, d1, d2); + }*/ + return res; } @@ -1160,6 +1172,7 @@ namespace CarpetLib { #define TYPECASE(N,T) \ \ \ + \ template \ void \ prolongate_3d_cc_eno_rf2<T,2> (T const * restrict const src, \ diff --git a/Carpet/CarpetReduce/src/mask_coords.c b/Carpet/CarpetReduce/src/mask_coords.c index 814f1ecfb..c324671d1 100644 --- a/Carpet/CarpetReduce/src/mask_coords.c +++ b/Carpet/CarpetReduce/src/mask_coords.c @@ -109,10 +109,12 @@ CoordBase_SetupMask (CCTK_ARGUMENTS) } if (f==0) { /* lower face */ - bmax[d] = imin[d] + bnd_points[2*d+f]; + //bmax[d] = imin[d] + bnd_points[2*d+f]; + bmax[d] = imin[d] + cctk_nghostzones[d]; } else { /* upper face */ - bmin[d] = imax[d] - bnd_points[2*d+f]; + //bmin[d] = imax[d] - bnd_points[2*d+f]; + bmin[d] = imax[d] - cctk_nghostzones[d]; } /* Loop over the boundary */ @@ -201,10 +203,12 @@ CoordBase_SetupMask (CCTK_ARGUMENTS) } if (f==0) { /* lower face */ - bmax[d] = imin[d] + bnd_points[2*d+f]; + //bmax[d] = imin[d] + bnd_points[2*d+f]; + bmax[d] = imin[d] + cctk_nghostzones[d]; } else { /* upper face */ - bmin[d] = imax[d] - bnd_points[2*d+f]; + //bmin[d] = imax[d] - bnd_points[2*d+f]; + bmin[d] = imax[d] - cctk_nghostzones[d]; } /* Loop over the boundary */ |