diff options
Diffstat (limited to 'Carpet/CarpetLib')
-rw-r--r-- | Carpet/CarpetLib/param.ccl | 6 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/data.cc | 8 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/interpolate_3d_2tl.cc | 4 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/interpolate_3d_3tl.cc | 4 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/interpolate_3d_4tl.cc | 4 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/interpolate_3d_5tl.cc | 4 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/interpolate_eno_3d_3tl.cc | 4 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc | 82 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/restrict_3d_cc_o3_rf2.cc | 191 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/restrict_3d_cc_o5_rf2.cc | 65 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/restrict_3d_cc_rf2.cc | 41 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/restrict_3d_rf2.cc | 20 | ||||
-rw-r--r-- | Carpet/CarpetLib/src/restrict_3d_vc_rf2.cc | 59 |
13 files changed, 295 insertions, 197 deletions
diff --git a/Carpet/CarpetLib/param.ccl b/Carpet/CarpetLib/param.ccl index 3af3c401d..541ec2f49 100644 --- a/Carpet/CarpetLib/param.ccl +++ b/Carpet/CarpetLib/param.ccl @@ -39,6 +39,12 @@ BOOLEAN interpolate_from_buffer_zones "Use buffer points for interpolation" STEE +BOOLEAN use_loopcontrol_in_operators "Use LoopControl to parallelize AMR operators" STEERABLE=always +{ +} "no" + + + restricted: BOOLEAN use_higher_order_restriction "Use third order cell centered restriction operators instead of first order" diff --git a/Carpet/CarpetLib/src/data.cc b/Carpet/CarpetLib/src/data.cc index 2aed5d914..fa1b3c92e 100644 --- a/Carpet/CarpetLib/src/data.cc +++ b/Carpet/CarpetLib/src/data.cc @@ -60,6 +60,13 @@ call_operator (void ibbox3 const & restrict dstregbbox, void * const extraargs) { + DECLARE_CCTK_PARAMETERS; + + if (use_loopcontrol_in_operators) { + (* the_operator) + (src, srcpadext, srcext, dst, dstpadext, dstext, srcbbox, dstbbox, + srcregbbox, dstregbbox, extraargs); + } else { #ifndef _OPENMP (* the_operator) (src, srcpadext, srcext, dst, dstpadext, dstext, srcbbox, dstbbox, @@ -111,6 +118,7 @@ call_operator (void assert (alldstregbboxes == ibset (dstregbbox)); # endif #endif + } } template <typename T> diff --git a/Carpet/CarpetLib/src/interpolate_3d_2tl.cc b/Carpet/CarpetLib/src/interpolate_3d_2tl.cc index 33c00d8fd..bace93e5f 100644 --- a/Carpet/CarpetLib/src/interpolate_3d_2tl.cc +++ b/Carpet/CarpetLib/src/interpolate_3d_2tl.cc @@ -1,13 +1,13 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> #include <cstdlib> -#include <loopcontrol.h> - #include "operator_prototypes_3d.hh" #include "typeprops.hh" diff --git a/Carpet/CarpetLib/src/interpolate_3d_3tl.cc b/Carpet/CarpetLib/src/interpolate_3d_3tl.cc index f6dc7790d..6e39691f3 100644 --- a/Carpet/CarpetLib/src/interpolate_3d_3tl.cc +++ b/Carpet/CarpetLib/src/interpolate_3d_3tl.cc @@ -1,13 +1,13 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> #include <cstdlib> -#include <loopcontrol.h> - #include "operator_prototypes_3d.hh" #include "typeprops.hh" diff --git a/Carpet/CarpetLib/src/interpolate_3d_4tl.cc b/Carpet/CarpetLib/src/interpolate_3d_4tl.cc index 4f52f5723..f41d0b12b 100644 --- a/Carpet/CarpetLib/src/interpolate_3d_4tl.cc +++ b/Carpet/CarpetLib/src/interpolate_3d_4tl.cc @@ -1,13 +1,13 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> #include <cstdlib> -#include <loopcontrol.h> - #include "operator_prototypes_3d.hh" #include "typeprops.hh" diff --git a/Carpet/CarpetLib/src/interpolate_3d_5tl.cc b/Carpet/CarpetLib/src/interpolate_3d_5tl.cc index 0cfa144ca..c84f2c793 100644 --- a/Carpet/CarpetLib/src/interpolate_3d_5tl.cc +++ b/Carpet/CarpetLib/src/interpolate_3d_5tl.cc @@ -1,13 +1,13 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> #include <cstdlib> -#include <loopcontrol.h> - #include "operator_prototypes_3d.hh" #include "typeprops.hh" diff --git a/Carpet/CarpetLib/src/interpolate_eno_3d_3tl.cc b/Carpet/CarpetLib/src/interpolate_eno_3d_3tl.cc index ac04121c0..abcee485a 100644 --- a/Carpet/CarpetLib/src/interpolate_eno_3d_3tl.cc +++ b/Carpet/CarpetLib/src/interpolate_eno_3d_3tl.cc @@ -1,13 +1,13 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> #include <cstdlib> -#include <loopcontrol.h> - #include "operator_prototypes_3d.hh" #include "typeprops.hh" diff --git a/Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc b/Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc index 320635fd4..05aaad6ae 100644 --- a/Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc +++ b/Carpet/CarpetLib/src/prolongate_3d_cc_eno_rf2.cc @@ -1,6 +1,8 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> @@ -966,6 +968,7 @@ namespace CarpetLib { ibbox3 const & restrict regbbox, void * extraargs) { + DECLARE_CCTK_PARAMETERS; assert (not extraargs); static_assert (ORDER>=0, "ORDER must be non-negative"); @@ -1091,6 +1094,8 @@ namespace CarpetLib { + if (not use_loopcontrol_in_operators) { + // Loop over fine region // Label scheme: l 8 fk fj fi @@ -1284,9 +1289,86 @@ namespace CarpetLib { // end k loop l9:; + + + + } else { // use_loopcontrol_in_operators + // Loop over fine region +#pragma omp parallel + CCTK_LOOP3(prolongate_3d_cc_eno_rf2, + i,j,k, 0,0,0, regiext,regjext,regkext, + dstipadext,dstjpadext,dstkpadext) + { + const ptrdiff_t is = (srcioff + i) / 2; + const ptrdiff_t js = (srcjoff + j) / 2; + const ptrdiff_t ks = (srckoff + k) / 2; + const ptrdiff_t im = (srcioff + i) % 2; + const ptrdiff_t jm = (srcjoff + j) % 2; + const ptrdiff_t km = (srckoff + k) % 2; + const ptrdiff_t id = dstioff + i; + const ptrdiff_t jd = dstjoff + j; + const ptrdiff_t kd = dstkoff + k; + if (km == 0) { + if (jm == 0) { + if (im == 0) { + check_indices3<T,ORDER,0,0,0> (is,js,ks, srciext, srcjext, srckext); + dst[DSTIND3(id,jd,kd)] = + interp3<T,ORDER,0,0,0> + (& src[SRCIND3(is,js,ks)], srcdi, srcdj, srcdk); + } else { + check_indices3<T,ORDER,1,0,0> (is,js,ks, srciext, srcjext, srckext); + dst[DSTIND3(id,jd,kd)] = + interp3<T,ORDER,1,0,0> + (& src[SRCIND3(is,js,ks)], srcdi, srcdj, srcdk); + } + } else { + if (im == 0) { + check_indices3<T,ORDER,0,1,0> (is,js,ks, srciext, srcjext, srckext); + dst[DSTIND3(id,jd,kd)] = + interp3<T,ORDER,0,1,0> + (& src[SRCIND3(is,js,ks)], srcdi, srcdj, srcdk); + } else { + check_indices3<T,ORDER,1,1,0> (is,js,ks, srciext, srcjext, srckext); + dst[DSTIND3(id,jd,kd)] = + interp3<T,ORDER,1,1,0> + (& src[SRCIND3(is,js,ks)], srcdi, srcdj, srcdk); + } + } + } else { + if (jm == 0) { + if (im == 0) { + check_indices3<T,ORDER,0,0,1> (is,js,ks, srciext, srcjext, srckext); + dst[DSTIND3(id,jd,kd)] = + interp3<T,ORDER,0,0,1> + (& src[SRCIND3(is,js,ks)], srcdi, srcdj, srcdk); + } else { + check_indices3<T,ORDER,1,0,1> (is,js,ks, srciext, srcjext, srckext); + dst[DSTIND3(id,jd,kd)] = + interp3<T,ORDER,1,0,1> + (& src[SRCIND3(is,js,ks)], srcdi, srcdj, srcdk); + } + } else { + if (im == 0) { + check_indices3<T,ORDER,0,1,1> (is,js,ks, srciext, srcjext, srckext); + dst[DSTIND3(id,jd,kd)] = + interp3<T,ORDER,0,1,1> + (& src[SRCIND3(is,js,ks)], srcdi, srcdj, srcdk); + } else { + check_indices3<T,ORDER,1,1,1> (is,js,ks, srciext, srcjext, srckext); + dst[DSTIND3(id,jd,kd)] = + interp3<T,ORDER,1,1,1> + (& src[SRCIND3(is,js,ks)], srcdi, srcdj, srcdk); + } + } + } + } CCTK_ENDLOOP3(prolongate_3d_cc_eno_rf2); + + } // if use_loopcontrol_in_operators + } + // Specialise for complex types diff --git a/Carpet/CarpetLib/src/restrict_3d_cc_o3_rf2.cc b/Carpet/CarpetLib/src/restrict_3d_cc_o3_rf2.cc index 12404ea25..3d1ce9acf 100644 --- a/Carpet/CarpetLib/src/restrict_3d_cc_o3_rf2.cc +++ b/Carpet/CarpetLib/src/restrict_3d_cc_o3_rf2.cc @@ -1,6 +1,8 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> @@ -138,103 +140,102 @@ namespace CarpetLib { // Loop over coarse region -#pragma omp parallel for collapse(3) - for (int k=0; k<regkext; ++k) { - for (int j=0; j<regjext; ++j) { - for (int i=0; i<regiext; ++i) { - -#ifdef CARPET_DEBUG - if(not (2 * k + 2 + srckoff < srckext and - 2 * j + 2 + srcjoff < srcjext and - 2 * i + 2 + srcioff < srciext and - 2 * k - 1 + srckoff >= 0 and - 2 * j - 1 + srcjoff >= 0 and - 2 * i - 1 + srcioff >= 0)) +#pragma omp parallel + CCTK_LOOP3(restrict_3d_cc_o3_rf2, + i,j,k, 0,0,0, regiext,regjext,regkext, + dstipadext,dstjpadext,dstkpadext) { - cout << "restrict_3d_cc_o3_rf2.cc\n"; - cout << "regext " << regext << "\n"; - cout << "srcext " << srcext << "\n"; - cout << "srcbbox=" << srcbbox << "\n"; - cout << "dstbbox=" << dstbbox << "\n"; - cout << "regbbox=" << regbbox << "\n"; - cout << "i,j,k=" << i << " " << j << " " << k << "\n"; - assert(2 * k + 2 + srckoff < srckext); - assert(2 * j + 2 + srcjoff < srcjext); - assert(2 * i + 2 + srcioff < srciext); - assert(2 * k - 1 + srckoff >= 0); - assert(2 * j - 1 + srcjoff >= 0); - assert(2 * i - 1 + srcioff >= 0); - } -#endif - dst [DSTIND3(i, j, k)] = - - f1 * src [SRCIND3(2*i-1, 2*j-1, 2*k-1)] - + f2 * src [SRCIND3(2*i , 2*j-1, 2*k-1)] - + f2 * src [SRCIND3(2*i+1, 2*j-1, 2*k-1)] - - f1 * src [SRCIND3(2*i+2, 2*j-1, 2*k-1)] - + f2 * src [SRCIND3(2*i-1, 2*j , 2*k-1)] - - f3 * src [SRCIND3(2*i , 2*j , 2*k-1)] - - f3 * src [SRCIND3(2*i+1, 2*j , 2*k-1)] - + f2 * src [SRCIND3(2*i+2, 2*j , 2*k-1)] - + f2 * src [SRCIND3(2*i-1, 2*j+1, 2*k-1)] - - f3 * src [SRCIND3(2*i , 2*j+1, 2*k-1)] - - f3 * src [SRCIND3(2*i+1, 2*j+1, 2*k-1)] - + f2 * src [SRCIND3(2*i+2, 2*j+1, 2*k-1)] - - f1 * src [SRCIND3(2*i-1, 2*j+2, 2*k-1)] - + f2 * src [SRCIND3(2*i , 2*j+2, 2*k-1)] - + f2 * src [SRCIND3(2*i+1, 2*j+2, 2*k-1)] - - f1 * src [SRCIND3(2*i+2, 2*j+2, 2*k-1)] - + f2 * src [SRCIND3(2*i-1, 2*j-1, 2*k )] - - f3 * src [SRCIND3(2*i , 2*j-1, 2*k )] - - f3 * src [SRCIND3(2*i+1, 2*j-1, 2*k )] - + f2 * src [SRCIND3(2*i+2, 2*j-1, 2*k )] - - f3 * src [SRCIND3(2*i-1, 2*j , 2*k )] - + f4 * src [SRCIND3(2*i , 2*j , 2*k )] - + f4 * src [SRCIND3(2*i+1, 2*j , 2*k )] - - f3 * src [SRCIND3(2*i+2, 2*j , 2*k )] - - f3 * src [SRCIND3(2*i-1, 2*j+1, 2*k )] - + f4 * src [SRCIND3(2*i , 2*j+1, 2*k )] - + f4 * src [SRCIND3(2*i+1, 2*j+1, 2*k )] - - f3 * src [SRCIND3(2*i+2, 2*j+1, 2*k )] - + f2 * src [SRCIND3(2*i-1, 2*j+2, 2*k )] - - f3 * src [SRCIND3(2*i , 2*j+2, 2*k )] - - f3 * src [SRCIND3(2*i+1, 2*j+2, 2*k )] - + f2 * src [SRCIND3(2*i+2, 2*j+2, 2*k )] - + f2 * src [SRCIND3(2*i-1, 2*j-1, 2*k+1)] - - f3 * src [SRCIND3(2*i , 2*j-1, 2*k+1)] - - f3 * src [SRCIND3(2*i+1, 2*j-1, 2*k+1)] - + f2 * src [SRCIND3(2*i+2, 2*j-1, 2*k+1)] - - f3 * src [SRCIND3(2*i-1, 2*j , 2*k+1)] - + f4 * src [SRCIND3(2*i , 2*j , 2*k+1)] - + f4 * src [SRCIND3(2*i+1, 2*j , 2*k+1)] - - f3 * src [SRCIND3(2*i+2, 2*j , 2*k+1)] - - f3 * src [SRCIND3(2*i-1, 2*j+1, 2*k+1)] - + f4 * src [SRCIND3(2*i , 2*j+1, 2*k+1)] - + f4 * src [SRCIND3(2*i+1, 2*j+1, 2*k+1)] - - f3 * src [SRCIND3(2*i+2, 2*j+1, 2*k+1)] - + f2 * src [SRCIND3(2*i-1, 2*j+2, 2*k+1)] - - f3 * src [SRCIND3(2*i , 2*j+2, 2*k+1)] - - f3 * src [SRCIND3(2*i+1, 2*j+2, 2*k+1)] - + f2 * src [SRCIND3(2*i+2, 2*j+2, 2*k+1)] - - f1 * src [SRCIND3(2*i-1, 2*j-1, 2*k+2)] - + f2 * src [SRCIND3(2*i , 2*j-1, 2*k+2)] - + f2 * src [SRCIND3(2*i+1, 2*j-1, 2*k+2)] - - f1 * src [SRCIND3(2*i+2, 2*j-1, 2*k+2)] - + f2 * src [SRCIND3(2*i-1, 2*j , 2*k+2)] - - f3 * src [SRCIND3(2*i , 2*j , 2*k+2)] - - f3 * src [SRCIND3(2*i+1, 2*j , 2*k+2)] - + f2 * src [SRCIND3(2*i+2, 2*j , 2*k+2)] - + f2 * src [SRCIND3(2*i-1, 2*j+1, 2*k+2)] - - f3 * src [SRCIND3(2*i , 2*j+1, 2*k+2)] - - f3 * src [SRCIND3(2*i+1, 2*j+1, 2*k+2)] - + f2 * src [SRCIND3(2*i+2, 2*j+1, 2*k+2)] - - f1 * src [SRCIND3(2*i-1, 2*j+2, 2*k+2)] - + f2 * src [SRCIND3(2*i , 2*j+2, 2*k+2)] - + f2 * src [SRCIND3(2*i+1, 2*j+2, 2*k+2)] - - f1 * src [SRCIND3(2*i+2, 2*j+2, 2*k+2)]; - - } + +#ifdef CARPET_DEBUG + if(not (2 * k + 2 + srckoff < srckext and + 2 * j + 2 + srcjoff < srcjext and + 2 * i + 2 + srcioff < srciext and + 2 * k - 1 + srckoff >= 0 and + 2 * j - 1 + srcjoff >= 0 and + 2 * i - 1 + srcioff >= 0)) + { + cout << "restrict_3d_cc_o3_rf2.cc\n"; + cout << "regext " << regext << "\n"; + cout << "srcext " << srcext << "\n"; + cout << "srcbbox=" << srcbbox << "\n"; + cout << "dstbbox=" << dstbbox << "\n"; + cout << "regbbox=" << regbbox << "\n"; + cout << "i,j,k=" << i << " " << j << " " << k << "\n"; + assert(2 * k + 2 + srckoff < srckext); + assert(2 * j + 2 + srcjoff < srcjext); + assert(2 * i + 2 + srcioff < srciext); + assert(2 * k - 1 + srckoff >= 0); + assert(2 * j - 1 + srcjoff >= 0); + assert(2 * i - 1 + srcioff >= 0); } - } +#endif + dst [DSTIND3(i, j, k)] = + - f1 * src [SRCIND3(2*i-1, 2*j-1, 2*k-1)] + + f2 * src [SRCIND3(2*i , 2*j-1, 2*k-1)] + + f2 * src [SRCIND3(2*i+1, 2*j-1, 2*k-1)] + - f1 * src [SRCIND3(2*i+2, 2*j-1, 2*k-1)] + + f2 * src [SRCIND3(2*i-1, 2*j , 2*k-1)] + - f3 * src [SRCIND3(2*i , 2*j , 2*k-1)] + - f3 * src [SRCIND3(2*i+1, 2*j , 2*k-1)] + + f2 * src [SRCIND3(2*i+2, 2*j , 2*k-1)] + + f2 * src [SRCIND3(2*i-1, 2*j+1, 2*k-1)] + - f3 * src [SRCIND3(2*i , 2*j+1, 2*k-1)] + - f3 * src [SRCIND3(2*i+1, 2*j+1, 2*k-1)] + + f2 * src [SRCIND3(2*i+2, 2*j+1, 2*k-1)] + - f1 * src [SRCIND3(2*i-1, 2*j+2, 2*k-1)] + + f2 * src [SRCIND3(2*i , 2*j+2, 2*k-1)] + + f2 * src [SRCIND3(2*i+1, 2*j+2, 2*k-1)] + - f1 * src [SRCIND3(2*i+2, 2*j+2, 2*k-1)] + + f2 * src [SRCIND3(2*i-1, 2*j-1, 2*k )] + - f3 * src [SRCIND3(2*i , 2*j-1, 2*k )] + - f3 * src [SRCIND3(2*i+1, 2*j-1, 2*k )] + + f2 * src [SRCIND3(2*i+2, 2*j-1, 2*k )] + - f3 * src [SRCIND3(2*i-1, 2*j , 2*k )] + + f4 * src [SRCIND3(2*i , 2*j , 2*k )] + + f4 * src [SRCIND3(2*i+1, 2*j , 2*k )] + - f3 * src [SRCIND3(2*i+2, 2*j , 2*k )] + - f3 * src [SRCIND3(2*i-1, 2*j+1, 2*k )] + + f4 * src [SRCIND3(2*i , 2*j+1, 2*k )] + + f4 * src [SRCIND3(2*i+1, 2*j+1, 2*k )] + - f3 * src [SRCIND3(2*i+2, 2*j+1, 2*k )] + + f2 * src [SRCIND3(2*i-1, 2*j+2, 2*k )] + - f3 * src [SRCIND3(2*i , 2*j+2, 2*k )] + - f3 * src [SRCIND3(2*i+1, 2*j+2, 2*k )] + + f2 * src [SRCIND3(2*i+2, 2*j+2, 2*k )] + + f2 * src [SRCIND3(2*i-1, 2*j-1, 2*k+1)] + - f3 * src [SRCIND3(2*i , 2*j-1, 2*k+1)] + - f3 * src [SRCIND3(2*i+1, 2*j-1, 2*k+1)] + + f2 * src [SRCIND3(2*i+2, 2*j-1, 2*k+1)] + - f3 * src [SRCIND3(2*i-1, 2*j , 2*k+1)] + + f4 * src [SRCIND3(2*i , 2*j , 2*k+1)] + + f4 * src [SRCIND3(2*i+1, 2*j , 2*k+1)] + - f3 * src [SRCIND3(2*i+2, 2*j , 2*k+1)] + - f3 * src [SRCIND3(2*i-1, 2*j+1, 2*k+1)] + + f4 * src [SRCIND3(2*i , 2*j+1, 2*k+1)] + + f4 * src [SRCIND3(2*i+1, 2*j+1, 2*k+1)] + - f3 * src [SRCIND3(2*i+2, 2*j+1, 2*k+1)] + + f2 * src [SRCIND3(2*i-1, 2*j+2, 2*k+1)] + - f3 * src [SRCIND3(2*i , 2*j+2, 2*k+1)] + - f3 * src [SRCIND3(2*i+1, 2*j+2, 2*k+1)] + + f2 * src [SRCIND3(2*i+2, 2*j+2, 2*k+1)] + - f1 * src [SRCIND3(2*i-1, 2*j-1, 2*k+2)] + + f2 * src [SRCIND3(2*i , 2*j-1, 2*k+2)] + + f2 * src [SRCIND3(2*i+1, 2*j-1, 2*k+2)] + - f1 * src [SRCIND3(2*i+2, 2*j-1, 2*k+2)] + + f2 * src [SRCIND3(2*i-1, 2*j , 2*k+2)] + - f3 * src [SRCIND3(2*i , 2*j , 2*k+2)] + - f3 * src [SRCIND3(2*i+1, 2*j , 2*k+2)] + + f2 * src [SRCIND3(2*i+2, 2*j , 2*k+2)] + + f2 * src [SRCIND3(2*i-1, 2*j+1, 2*k+2)] + - f3 * src [SRCIND3(2*i , 2*j+1, 2*k+2)] + - f3 * src [SRCIND3(2*i+1, 2*j+1, 2*k+2)] + + f2 * src [SRCIND3(2*i+2, 2*j+1, 2*k+2)] + - f1 * src [SRCIND3(2*i-1, 2*j+2, 2*k+2)] + + f2 * src [SRCIND3(2*i , 2*j+2, 2*k+2)] + + f2 * src [SRCIND3(2*i+1, 2*j+2, 2*k+2)] + - f1 * src [SRCIND3(2*i+2, 2*j+2, 2*k+2)]; + + } CCTK_ENDLOOP3(restrict_3d_cc_o3_rf2); } diff --git a/Carpet/CarpetLib/src/restrict_3d_cc_o5_rf2.cc b/Carpet/CarpetLib/src/restrict_3d_cc_o5_rf2.cc index 2cd1b1f1f..5c2f88642 100644 --- a/Carpet/CarpetLib/src/restrict_3d_cc_o5_rf2.cc +++ b/Carpet/CarpetLib/src/restrict_3d_cc_o5_rf2.cc @@ -1,6 +1,8 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> @@ -124,8 +126,6 @@ namespace CarpetLib { DECLARE_CCTK_PARAMETERS; - typedef typename typeprops<T>::real RT; - if (any (srcbbox.stride() >= regbbox.stride() or @@ -205,41 +205,40 @@ namespace CarpetLib { // Loop over coarse region -#pragma omp parallel for collapse(3) - for (int k=0; k<regkext; ++k) { - for (int j=0; j<regjext; ++j) { - for (int i=0; i<regiext; ++i) { +#pragma omp parallel + CCTK_LOOP3(restrict_3d_cc_o5_rf2, + i,j,k, 0,0,0, regiext,regjext,regkext, + dstipadext,dstjpadext,dstkpadext) + { #ifdef CARPET_DEBUG - if(not (2 * k + 3 + srckoff < srckext and - 2 * j + 3 + srcjoff < srcjext and - 2 * i + 3 + srcioff < srciext and - 2 * k - 2 + srckoff >= 0 and - 2 * j - 2 + srcjoff >= 0 and - 2 * i - 2 + srcioff >= 0)) - { - cout << "restrict_3d_cc_o3_rf2.cc\n"; - cout << "regext " << regext << "\n"; - cout << "srcext " << srcext << "\n"; - cout << "srcbbox=" << srcbbox << "\n"; - cout << "dstbbox=" << dstbbox << "\n"; - cout << "regbbox=" << regbbox << "\n"; - cout << "i,j,k=" << i << " " << j << " " << k << "\n"; - assert(2 * k + 2 + srckoff < srckext); - assert(2 * j + 2 + srcjoff < srcjext); - assert(2 * i + 2 + srcioff < srciext); - assert(2 * k - 1 + srckoff >= 0); - assert(2 * j - 1 + srcjoff >= 0); - assert(2 * i - 1 + srcioff >= 0); - } + if(not (2 * k + 3 + srckoff < srckext and + 2 * j + 3 + srcjoff < srcjext and + 2 * i + 3 + srcioff < srciext and + 2 * k - 2 + srckoff >= 0 and + 2 * j - 2 + srcjoff >= 0 and + 2 * i - 2 + srcioff >= 0)) + { + cout << "restrict_3d_cc_o3_rf2.cc\n"; + cout << "regext " << regext << "\n"; + cout << "srcext " << srcext << "\n"; + cout << "srcbbox=" << srcbbox << "\n"; + cout << "dstbbox=" << dstbbox << "\n"; + cout << "regbbox=" << regbbox << "\n"; + cout << "i,j,k=" << i << " " << j << " " << k << "\n"; + assert(2 * k + 2 + srckoff < srckext); + assert(2 * j + 2 + srcjoff < srcjext); + assert(2 * i + 2 + srcioff < srciext); + assert(2 * k - 1 + srckoff >= 0); + assert(2 * j - 1 + srcjoff >= 0); + assert(2 * i - 1 + srcioff >= 0); + } #endif - dst [DSTIND3(i, j, k)] = - CarpetLib::restrict3 - (& src[SRCIND3(2*i, 2*j, 2*k)], srcdi, srcdj, srcdk); + dst [DSTIND3(i, j, k)] = + CarpetLib::restrict3 + (& src[SRCIND3(2*i, 2*j, 2*k)], srcdi, srcdj, srcdk); - } - } - } + } CCTK_ENDLOOP3(restrict_3d_cc_o5_rf2); } diff --git a/Carpet/CarpetLib/src/restrict_3d_cc_rf2.cc b/Carpet/CarpetLib/src/restrict_3d_cc_rf2.cc index ee3fafc18..e186cbced 100644 --- a/Carpet/CarpetLib/src/restrict_3d_cc_rf2.cc +++ b/Carpet/CarpetLib/src/restrict_3d_cc_rf2.cc @@ -1,6 +1,8 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> @@ -127,26 +129,25 @@ namespace CarpetLib { // Loop over coarse region -#pragma omp parallel for collapse(3) - for (int k=0; k<regkext; ++k) { - for (int j=0; j<regjext; ++j) { - for (int i=0; i<regiext; ++i) { - - // TODO: Introduce higher-order restriction operators (but - // don't use these for hydro!) - dst [DSTIND3(i, j, k)] = - + f1*f1*f1 * src [SRCIND3(2*i , 2*j , 2*k )] - + f2*f1*f1 * src [SRCIND3(2*i+1, 2*j , 2*k )] - + f1*f2*f1 * src [SRCIND3(2*i , 2*j+1, 2*k )] - + f2*f2*f1 * src [SRCIND3(2*i+1, 2*j+1, 2*k )] - + f1*f1*f2 * src [SRCIND3(2*i , 2*j , 2*k+1)] - + f2*f1*f2 * src [SRCIND3(2*i+1, 2*j , 2*k+1)] - + f1*f2*f2 * src [SRCIND3(2*i , 2*j+1, 2*k+1)] - + f2*f2*f2 * src [SRCIND3(2*i+1, 2*j+1, 2*k+1)]; - - } - } - } +#pragma omp parallel + CCTK_LOOP3(restrict_3d_cc_rf2, + i,j,k, 0,0,0, regiext,regjext,regkext, + dstipadext,dstjpadext,dstkpadext) + { + + // TODO: Introduce higher-order restriction operators (but + // don't use these for hydro!) + dst [DSTIND3(i, j, k)] = + + f1*f1*f1 * src [SRCIND3(2*i , 2*j , 2*k )] + + f2*f1*f1 * src [SRCIND3(2*i+1, 2*j , 2*k )] + + f1*f2*f1 * src [SRCIND3(2*i , 2*j+1, 2*k )] + + f2*f2*f1 * src [SRCIND3(2*i+1, 2*j+1, 2*k )] + + f1*f1*f2 * src [SRCIND3(2*i , 2*j , 2*k+1)] + + f2*f1*f2 * src [SRCIND3(2*i+1, 2*j , 2*k+1)] + + f1*f2*f2 * src [SRCIND3(2*i , 2*j+1, 2*k+1)] + + f2*f2*f2 * src [SRCIND3(2*i+1, 2*j+1, 2*k+1)]; + + } CCTK_ENDLOOP3(restrict_3d_cc_rf2); } diff --git a/Carpet/CarpetLib/src/restrict_3d_rf2.cc b/Carpet/CarpetLib/src/restrict_3d_rf2.cc index c1dade7fb..9f0936e1f 100644 --- a/Carpet/CarpetLib/src/restrict_3d_rf2.cc +++ b/Carpet/CarpetLib/src/restrict_3d_rf2.cc @@ -1,6 +1,8 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> @@ -113,15 +115,15 @@ namespace CarpetLib { // Loop over coarse region - for (int k=0; k<regkext; ++k) { - for (int j=0; j<regjext; ++j) { - for (int i=0; i<regiext; ++i) { - - dst [DSTIND3(i, j, k)] = src [SRCIND3(2*i, 2*j, 2*k)]; - - } - } - } +#pragma omp parallel + CCTK_LOOP3(restrict_3d_rf2, + i,j,k, 0,0,0, regiext,regjext,regkext, + dstipadext,dstjpadext,dstkpadext) + { + + dst [DSTIND3(i, j, k)] = src [SRCIND3(2*i, 2*j, 2*k)]; + + } CCTK_ENDLOOP3(restrict_3d_rf2); } diff --git a/Carpet/CarpetLib/src/restrict_3d_vc_rf2.cc b/Carpet/CarpetLib/src/restrict_3d_vc_rf2.cc index e1539aa06..62450333c 100644 --- a/Carpet/CarpetLib/src/restrict_3d_vc_rf2.cc +++ b/Carpet/CarpetLib/src/restrict_3d_vc_rf2.cc @@ -1,6 +1,8 @@ #include <cctk.h> #include <cctk_Parameters.h> +#include <loopcontrol.h> + #include <algorithm> #include <cassert> #include <cmath> @@ -152,8 +154,6 @@ namespace CarpetLib { DECLARE_CCTK_PARAMETERS; - typedef typename typeprops<T>::real RT; - // false: vertex centered, true: cell centered ivect const icent (centi, centj, centk); assert (all (icent == 0 or icent == 1)); @@ -250,36 +250,35 @@ namespace CarpetLib { // Loop over coarse region -#pragma omp parallel for collapse(3) - for (int k=0; k<regkext; ++k) { - for (int j=0; j<regjext; ++j) { - for (int i=0; i<regiext; ++i) { -#ifdef CARPET_DEBUG - if(not (2 * k + centk < srckext and - 2 * j + centj < srcjext and - 2 * i + centi < srciext)) +#pragma omp parallel + CCTK_LOOP3(restrict_3d_vc_rf2, + i,j,k, 0,0,0, regiext,regjext,regkext, + dstipadext,dstjpadext,dstkpadext) { - cout << "restrict_3d_vc_rf2.cc\n"; - cout << "regext " << regext << "\n"; - cout << "srcext " << srcext << "\n"; - cout << "srcbbox=" << srcbbox << "\n"; - cout << "dstbbox=" << dstbbox << "\n"; - cout << "regbbox=" << regbbox << "\n"; - cout << "srcregbbox=" << srcregbbox << "\n"; - cout << "icent=" << icent << "\n"; - } - assert(2 * k + centk < srckext and - 2 * j + centj < srcjext and - 2 * i + centi < srciext); -#endif - - dst [DSTIND3(i, j, k)] = - restrict3<T,centi,centj,centk>::call - (& src[SRCIND3(2*i, 2*j, 2*k)], srcdi, srcdj, srcdk); - - } +#ifdef CARPET_DEBUG + if(not (2 * k + centk < srckext and + 2 * j + centj < srcjext and + 2 * i + centi < srciext)) + { + cout << "restrict_3d_vc_rf2.cc\n"; + cout << "regext " << regext << "\n"; + cout << "srcext " << srcext << "\n"; + cout << "srcbbox=" << srcbbox << "\n"; + cout << "dstbbox=" << dstbbox << "\n"; + cout << "regbbox=" << regbbox << "\n"; + cout << "srcregbbox=" << srcregbbox << "\n"; + cout << "icent=" << icent << "\n"; } - } + assert(2 * k + centk < srckext and + 2 * j + centj < srcjext and + 2 * i + centi < srciext); +#endif + + dst [DSTIND3(i, j, k)] = + restrict3<T,centi,centj,centk>::call + (& src[SRCIND3(2*i, 2*j, 2*k)], srcdi, srcdj, srcdk); + + } CCTK_ENDLOOP3(restrict_3d_vc_rf2); } |