easy-iso

Iso-surface extraction from volume data
git clone https://git.0xfab.ch/easy-iso.git
Log | Files | Refs | Submodules | README | LICENSE

commit a4b27de4a47a267aec48d62c6ce8d414d289b8bd
parent bbcc969a7b81ae81862f540da06339df5c4c6f6e
Author: Fabian Wermelinger <fabianw@mavt.ethz.ch>
Date:   Thu,  1 Dec 2016 09:29:55 +0100

MPI function signatures modified for opem-mpi 1.6.5

Diffstat:
Msrc/InterpolatorMPI.h | 14+++++++-------
Msrc/IsoExtractorMPI.h | 2+-
2 files changed, 8 insertions(+), 8 deletions(-)

diff --git a/src/InterpolatorMPI.h b/src/InterpolatorMPI.h @@ -68,7 +68,7 @@ public: m_PESize{p("xpesize").asInt(1),p("ypesize").asInt(1),p("zpesize").asInt(1)} { this->m_extent = p("extent").asDouble(1.0); - const int periodic[3] = {false}; + int periodic[3] = {false}; int world_size; MPI_Comm_size(comm, &world_size); @@ -153,7 +153,7 @@ public: private: const MPI_Comm m_comm; - const int m_PESize[3]; + int m_PESize[3]; int m_myPEIndex[3]; bool m_isroot; @@ -450,7 +450,7 @@ template <typename Tinterp> void InterpolatorMPI<Tinterp>::_send_faces(const MPI_Comm comm, const MPInbr& nbrs, const Faces& f, std::vector<MPI_Request>& f_req) const { const int dst[6] = {nbrs(-1,0,0), nbrs(1,0,0), nbrs(0,-1,0), nbrs(0,1,0), nbrs(0,0,-1), nbrs(0,0,1)}; - const Real* const data[6] = {f.x0.data, f.x1.data, f.y0.data, f.y1.data, f.z0.data, f.z1.data}; + Real* const data[6] = {f.x0.data, f.x1.data, f.y0.data, f.y1.data, f.z0.data, f.z1.data}; const int Nelem[6] = {f.x0._Nelem, f.x1._Nelem, f.y0._Nelem, f.y1._Nelem, f.z0._Nelem, f.z1._Nelem}; const int tagbase = 0; for (int i = 0; i < 6; ++i) @@ -494,9 +494,9 @@ void InterpolatorMPI<Tinterp>::_send_edges(const MPI_Comm comm, const MPInbr& nb const int dstX[4] = {nbrs(0,-1,-1), nbrs(0,1,-1), nbrs(0,1,1), nbrs(0,-1,1)}; const int dstY[4] = {nbrs(-1,0,-1), nbrs(-1,0,1), nbrs(1,0,1), nbrs(1,0,-1)}; const int dstZ[4] = {nbrs(-1,-1,0), nbrs(1,-1,0), nbrs(1,1,0), nbrs(-1,1,0)}; - const Real* const dataX[4] = {e.x0.data, e.x1.data, e.x2.data, e.x3.data}; - const Real* const dataY[4] = {e.y0.data, e.y1.data, e.y2.data, e.y3.data}; - const Real* const dataZ[4] = {e.z0.data, e.z1.data, e.z2.data, e.z3.data}; + Real* const dataX[4] = {e.x0.data, e.x1.data, e.x2.data, e.x3.data}; + Real* const dataY[4] = {e.y0.data, e.y1.data, e.y2.data, e.y3.data}; + Real* const dataZ[4] = {e.z0.data, e.z1.data, e.z2.data, e.z3.data}; const int NelemX[4] = {e.x0._Nelem, e.x1._Nelem, e.x2._Nelem, e.x3._Nelem}; const int NelemY[4] = {e.y0._Nelem, e.y1._Nelem, e.y2._Nelem, e.y3._Nelem}; const int NelemZ[4] = {e.z0._Nelem, e.z1._Nelem, e.z2._Nelem, e.z3._Nelem}; @@ -568,7 +568,7 @@ template <typename Tinterp> void InterpolatorMPI<Tinterp>::_send_corners(const MPI_Comm comm, const MPInbr& nbrs, const Corners& c, std::vector<MPI_Request>& c_req) const { const int dst[8] = {nbrs(-1,-1,-1), nbrs(-1,1,-1), nbrs(-1,1,1), nbrs(-1,-1,1), nbrs(1,-1,-1), nbrs(1,1,-1), nbrs(1,1,1), nbrs(1,-1,1)}; - const Real* const data[8] = {c.x00.data, c.x01.data, c.x02.data, c.x03.data, c.x10.data, c.x11.data, c.x12.data, c.x13.data}; + Real* const data[8] = {c.x00.data, c.x01.data, c.x02.data, c.x03.data, c.x10.data, c.x11.data, c.x12.data, c.x13.data}; const int Nelem[8] = {c.x00._Nelem, c.x01._Nelem, c.x02._Nelem, c.x03._Nelem, c.x10._Nelem, c.x11._Nelem, c.x12._Nelem, c.x13._Nelem}; const int tagbase = 18; for (int i = 0; i < 8; ++i) diff --git a/src/IsoExtractorMPI.h b/src/IsoExtractorMPI.h @@ -64,7 +64,7 @@ private: void _writePLY_MPI(const std::string basename, const std::vector<Triangle>& tList) { const MPI_Comm comm = this->m_interp.getComm(); - const unsigned int mytriangles = tList.size(); + unsigned long long mytriangles = tList.size(); unsigned long long alltriangles; MPI_Reduce(&mytriangles, &alltriangles, 1, MPI_UNSIGNED_LONG_LONG, MPI_SUM, 0, comm);