easy-iso

commit a7b35095e15c4279f6a3c9f6e172a25f64315ddc
parent c6c54438238981e78bc2203ed40a201c0ed2ed08
Author: Fabian Wermelinger <fabianw@mavt.ethz.ch>
Date:   Thu,  1 Dec 2016 19:15:51 +0100

fixed shared memory issue in wavelet loaders

Diffstat:
M
src/Interpolator.h
 | 
35
+++++++++++++++++++----------------
M
src/InterpolatorMPI.h
 | 
45
++++++++++++++++++++++++---------------------

2 files changed, 43 insertions(+), 37 deletions(-)
diff --git a/src/Interpolator.h b/src/Interpolator.h
@@ -236,25 +236,28 @@ void Interpolator<TKernel>::_load_wavelet(ArgumentParser& p)
     const size_t num_elem = static_cast<size_t>(maxDim[0]) * maxDim[1] * maxDim[2];
     Real* const data = new Real[num_elem];
 
-    Real blockdata[_BLOCKSIZE_][_BLOCKSIZE_][_BLOCKSIZE_];
     const int nBlocks = NBX * NBY * NBZ;
-#pragma omp parallel for
-    for (int i = 0; i < nBlocks; ++i)
+#pragma omp parallel
     {
-        const int ix = i%NBX;
-        const int iy = (i/NBX)%NBY;
-        const int iz = (i/(NBX*NBY))%NBZ;
-        const double zratio = myreader.load_block2(ix, iy, iz, blockdata);
+        Real blockdata[_BLOCKSIZE_][_BLOCKSIZE_][_BLOCKSIZE_];
+#pragma omp for
+        for (int i = 0; i < nBlocks; ++i)
+        {
+            const int ix = i%NBX;
+            const int iy = (i/NBX)%NBY;
+            const int iz = (i/(NBX*NBY))%NBZ;
+            const double zratio = myreader.load_block2(ix, iy, iz, blockdata);
 
-        for (int z=0; z < _BLOCKSIZE_; ++z)
-            for (int y=0; y < _BLOCKSIZE_; ++y)
-            {
-                assert(iy*_BLOCKSIZE_+y < maxDim[1]);
-                assert(iz*_BLOCKSIZE_+z < maxDim[2]);
-                const size_t offset = _BLOCKSIZE_*(static_cast<size_t>(ix) + NBX*(y+iy*_BLOCKSIZE_ + (z+static_cast<size_t>(iz)*_BLOCKSIZE_)*NBY*_BLOCKSIZE_));
-                Real* const dst = data + offset;
-                memcpy(dst, &blockdata[z][y][0], _BLOCKSIZE_*sizeof(Real));
-            }
+            for (int z=0; z < _BLOCKSIZE_; ++z)
+                for (int y=0; y < _BLOCKSIZE_; ++y)
+                {
+                    assert(iy*_BLOCKSIZE_+y < maxDim[1]);
+                    assert(iz*_BLOCKSIZE_+z < maxDim[2]);
+                    const size_t offset = _BLOCKSIZE_*(static_cast<size_t>(ix) + NBX*(y+iy*_BLOCKSIZE_ + (z+static_cast<size_t>(iz)*_BLOCKSIZE_)*NBY*_BLOCKSIZE_));
+                    Real* const dst = data + offset;
+                    memcpy(dst, &blockdata[z][y][0], _BLOCKSIZE_*sizeof(Real));
+                }
+        }
     }
 
     this->m_data = std::move(Matrix_t(maxDim[0], maxDim[1], maxDim[2], data));
diff --git a/src/InterpolatorMPI.h b/src/InterpolatorMPI.h
@@ -384,29 +384,32 @@ void InterpolatorMPI<Tinterp>::_load_wavelet_MPI(ArgumentParser& p)
     const size_t num_elem = static_cast<size_t>(maxDim[0]) * maxDim[1] * maxDim[2];
     Real* const data = new Real[num_elem];
 
-    Real blockdata[_BLOCKSIZE_][_BLOCKSIZE_][_BLOCKSIZE_];
     const int nBlocks = myNBX * myNBY * myNBZ;
-#pragma omp parallel for
-    for (int i = 0; i < nBlocks; ++i)
+#pragma omp parallel
     {
-        const int ix = i%myNBX;
-        const int iy = (i/myNBX)%myNBY;
-        const int iz = (i/(myNBX*myNBY))%myNBZ;
-        const double zratio = myreader.load_block2(
-                ix+myNBX*m_myPEIndex[0],
-                iy+myNBY*m_myPEIndex[1],
-                iz+myNBZ*m_myPEIndex[2],
-                blockdata);
-
-        for (int z=0; z < _BLOCKSIZE_; ++z)
-            for (int y=0; y < _BLOCKSIZE_; ++y)
-            {
-                assert(iy*_BLOCKSIZE_+y < maxDim[1]);
-                assert(iz*_BLOCKSIZE_+z < maxDim[2]);
-                const size_t offset = _BLOCKSIZE_*(static_cast<size_t>(ix) + myNBX*(y+iy*_BLOCKSIZE_ + (z+static_cast<size_t>(iz)*_BLOCKSIZE_)*myNBY*_BLOCKSIZE_));
-                Real* const dst = data + offset;
-                memcpy(dst, &blockdata[z][y][0], _BLOCKSIZE_*sizeof(Real));
-            }
+        Real blockdata[_BLOCKSIZE_][_BLOCKSIZE_][_BLOCKSIZE_];
+#pragma omp for
+        for (int i = 0; i < nBlocks; ++i)
+        {
+            const int ix = i%myNBX;
+            const int iy = (i/myNBX)%myNBY;
+            const int iz = (i/(myNBX*myNBY))%myNBZ;
+            const double zratio = myreader.load_block2(
+                    ix+myNBX*m_myPEIndex[0],
+                    iy+myNBY*m_myPEIndex[1],
+                    iz+myNBZ*m_myPEIndex[2],
+                    blockdata);
+
+            for (int z=0; z < _BLOCKSIZE_; ++z)
+                for (int y=0; y < _BLOCKSIZE_; ++y)
+                {
+                    assert(iy*_BLOCKSIZE_+y < maxDim[1]);
+                    assert(iz*_BLOCKSIZE_+z < maxDim[2]);
+                    const size_t offset = _BLOCKSIZE_*(static_cast<size_t>(ix) + myNBX*(y+iy*_BLOCKSIZE_ + (z+static_cast<size_t>(iz)*_BLOCKSIZE_)*myNBY*_BLOCKSIZE_));
+                    Real* const dst = data + offset;
+                    memcpy(dst, &blockdata[z][y][0], _BLOCKSIZE_*sizeof(Real));
+                }
+        }
     }
 
     this->m_data = std::move(Matrix_t(maxDim[0], maxDim[1], maxDim[2], data));