diff --git a/solvers/iterateGeometry/iterateGeometry.cpp b/solvers/iterateGeometry/iterateGeometry.cpp index d541a5bf..1279d1c6 100755 --- a/solvers/iterateGeometry/iterateGeometry.cpp +++ b/solvers/iterateGeometry/iterateGeometry.cpp @@ -35,12 +35,10 @@ Licence: #include "commandLine.hpp" //#include "readControlDict.hpp" -using namespace pFlow; - int main( int argc, char* argv[] ) { -commandLine cmds( +pFlow::commandLine cmds( "iterateGeometry", "Performs simulation without particles, only geometry is solved"); @@ -55,8 +53,8 @@ commandLine cmds( // this should be palced in each main -processors::initProcessors(argc, argv); -initialize_pFlowProcessors(); +pFlow::processors::initProcessors(argc, argv); +pFlow::initialize_pFlowProcessors(); #include "initialize_Control.hpp" #include "setProperty.hpp" @@ -71,7 +69,7 @@ initialize_pFlowProcessors(); // this should be palced in each main #include "finalize.hpp" -processors::finalizeProcessors(); +pFlow::processors::finalizeProcessors(); } diff --git a/solvers/iterateSphereParticles/createDEMComponents.hpp b/solvers/iterateSphereParticles/createDEMComponents.hpp index 85cb798a..e9975f44 100755 --- a/solvers/iterateSphereParticles/createDEMComponents.hpp +++ b/solvers/iterateSphereParticles/createDEMComponents.hpp @@ -20,6 +20,6 @@ Licence: // REPORT(0)<<"\nReading sphere particles . . ."<numSurfaces() != motionComponentName_.size() ) { fatalErrorInFunction<< - "Number of surfaces is not equal to number of motion component names"<numSurfaces() << + ") is not equal to number of motion component names("<< + motionComponentName_.size()<<")"<::findMotionIndex() return true; } +namespace pFlow::GMotion +{ + template + void moveGeometry( + real dt, + uint32 nPoints, + const ModelInterface& mModel, + const deviceViewType1D& pointMIndexD, + const deviceViewType1D& pointsD + ) + { + Kokkos::parallel_for( + "geometryMotion::movePoints", + deviceRPolicyStatic(0, nPoints), + LAMBDA_HD(uint32 i){ + auto newPos = mModel.transferPoint(pointMIndexD[i], pointsD[i], dt); + pointsD[i] = newPos; + }); + + Kokkos::fence(); + } +} template bool pFlow::geometryMotion::moveGeometry() @@ -84,8 +105,15 @@ template auto& pointMIndexD= pointMotionIndex_.deviceViewAll(); auto& pointsD = points().deviceViewAll(); + pFlow::GMotion::moveGeometry( + dt, + numPoints(), + motionModel_.getModelInterface(iter, t, dt), + pointMotionIndex_.deviceViewAll(), + points().deviceViewAll() + ); - Kokkos::parallel_for( + /*Kokkos::parallel_for( "geometryMotion::movePoints", deviceRPolicyStatic(0, numPoints()), LAMBDA_HD(uint32 i){ @@ -93,7 +121,7 @@ template pointsD[i] = newPos; }); - Kokkos::fence(); + Kokkos::fence();*/ // move the motion components motionModel_.move(iter, t,dt); diff --git a/src/Interaction/CMakeLists.txt b/src/Interaction/CMakeLists.txt index c0265808..96197f6a 100644 --- a/src/Interaction/CMakeLists.txt +++ b/src/Interaction/CMakeLists.txt @@ -7,8 +7,8 @@ contactSearch/methods/cellBased/NBS/NBS.cpp contactSearch/methods/cellBased/NBS/cellsWallLevel0.cpp contactSearch/boundaries/boundaryContactSearch/boundaryContactSearch.cpp -contactSearch/boundaries/twoPartContactSearch/twoPartContactSearchKernels.cpp -contactSearch/boundaries/twoPartContactSearch/twoPartContactSearch.cpp +#contactSearch/boundaries/twoPartContactSearch/twoPartContactSearchKernels.cpp +#contactSearch/boundaries/twoPartContactSearch/twoPartContactSearch.cpp contactSearch/boundaries/periodicBoundaryContactSearch/ppwBndryContactSearchKernels.cpp contactSearch/boundaries/periodicBoundaryContactSearch/ppwBndryContactSearch.cpp contactSearch/boundaries/periodicBoundaryContactSearch/wallBoundaryContactSearch.cpp diff --git a/src/Interaction/contactLists/sortedContactList.hpp b/src/Interaction/contactLists/sortedContactList.hpp index 73f0ce4b..4f41caac 100644 --- a/src/Interaction/contactLists/sortedContactList.hpp +++ b/src/Interaction/contactLists/sortedContactList.hpp @@ -138,7 +138,7 @@ public: start, end, newPair); - idx0!=-1) + idx0!=static_cast(-1)) { values_[idx] = values0_[idx0]; } @@ -147,7 +147,7 @@ public: start, end, newPair); - idx0!=-1) + idx0!= static_cast(-1) ) { values_[idx] = values0_[idx0]; diff --git a/src/Interaction/contactLists/unsortedContactList.hpp b/src/Interaction/contactLists/unsortedContactList.hpp index 5922c4dc..7aa98237 100644 --- a/src/Interaction/contactLists/unsortedContactList.hpp +++ b/src/Interaction/contactLists/unsortedContactList.hpp @@ -95,7 +95,8 @@ public: // swap conainer and values swapViews(values0_, values_); swapViews(container0_, this->container_); - return UnsortedPairs::beforeBroadSearch(); + UnsortedPairs::beforeBroadSearch(); + return true; } bool afterBroadSearch() @@ -110,7 +111,7 @@ public: rpFillPairs(0,this->capacity()), *this); Kokkos::fence(); - + return true; } @@ -123,7 +124,7 @@ public: INLINE_FUNCTION_HD bool getValue(const PairType& p, ValueType& val)const { - if(auto idx = this->find(p); idx!=-1) + if(auto idx = this->find(p); idx!=static_cast(-1)) { val = getValue(idx); return true; @@ -140,7 +141,7 @@ public: INLINE_FUNCTION_HD bool setValue(const PairType& p, const ValueType& val)const { - if(uint32 idx = this->find(p); idx!=-1) + if(uint32 idx = this->find(p); idx!=static_cast(-1)) { setValue(idx, val); return true;; @@ -155,7 +156,7 @@ public: { if( uint32 idx0 = container0_.find(this->getPair(idx)); - idx0!=-1 ) + idx0!= static_cast(-1) ) { values_[idx] = values0_[idx0]; } diff --git a/src/Interaction/contactLists/unsortedPairs.hpp b/src/Interaction/contactLists/unsortedPairs.hpp index fd4f117d..fb622bdb 100644 --- a/src/Interaction/contactLists/unsortedPairs.hpp +++ b/src/Interaction/contactLists/unsortedPairs.hpp @@ -105,7 +105,7 @@ public: uint32 insert(idType i, idType j)const { if(auto insertResult = container_.insert(PairType(i,j)); insertResult.failed()) - return -1; + return static_cast(-1); else return insertResult.index(); @@ -115,7 +115,7 @@ public: uint32 insert(const PairType& p)const { if(auto insertResult = container_.insert(p); insertResult.failed()) - return -1; + return static_cast(-1); else return insertResult.index(); @@ -154,7 +154,7 @@ public: idx != Kokkos::UnorderedMapInvalidIndex ) return idx; else - return -1; + return static_cast(-1); } INLINE_FUNCTION_HD diff --git a/src/Interaction/contactSearch/ContactSearch/ContactSearch.hpp b/src/Interaction/contactSearch/ContactSearch/ContactSearch.hpp index d80813ca..f2dda414 100644 --- a/src/Interaction/contactSearch/ContactSearch/ContactSearch.hpp +++ b/src/Interaction/contactSearch/ContactSearch/ContactSearch.hpp @@ -160,7 +160,8 @@ public: csPairContainerType& pwPairs, bool force = false)override { - Particles().boundingSphere().updateBoundaries(DataDirection::SlaveToMaster); + if(i==0u) + Particles().boundingSphere().updateBoundaries(DataDirection::SlaveToMaster); return csBoundaries_[i].broadSearch( iter, t, @@ -183,6 +184,23 @@ public: { return ppwContactSearch_().performedSearch(); } + + + uint32 updateInterval()const override + { + return ppwContactSearch_().updateInterval(); + } + + real sizeRatio()const override + { + return ppwContactSearch_().sizeRatio(); + } + + + real cellExtent()const override + { + return ppwContactSearch_().cellExtent(); + } }; diff --git a/src/Interaction/contactSearch/boundaries/boundaryContactSearch/boundaryContactSearch.hpp b/src/Interaction/contactSearch/boundaries/boundaryContactSearch/boundaryContactSearch.hpp index 78c61bb8..55ecb47e 100644 --- a/src/Interaction/contactSearch/boundaries/boundaryContactSearch/boundaryContactSearch.hpp +++ b/src/Interaction/contactSearch/boundaries/boundaryContactSearch/boundaryContactSearch.hpp @@ -72,11 +72,6 @@ public: return contactSearch_; } - void fill(const std::any &val) override - { - return; - } - bool hearChanges( real t, real dt, diff --git a/src/Interaction/contactSearch/boundaries/periodicBoundaryContactSearch/ppwBndryContactSearchKernels.cpp b/src/Interaction/contactSearch/boundaries/periodicBoundaryContactSearch/ppwBndryContactSearchKernels.cpp index ea266e78..a0d47e78 100644 --- a/src/Interaction/contactSearch/boundaries/periodicBoundaryContactSearch/ppwBndryContactSearchKernels.cpp +++ b/src/Interaction/contactSearch/boundaries/periodicBoundaryContactSearch/ppwBndryContactSearchKernels.cpp @@ -80,14 +80,14 @@ pFlow::uint32 pFlow::pweBndryContactSearchKernels::broadSearchPP if(!searchCells.inCellRange(ind))continue; uint32 thisI = head(ind.x(),ind.y(),ind.z()); - while (thisI!=-1) + while (thisI!=static_cast(-1)) { auto d_n = sizeRatio*diams[thisI]; // first item is for this boundary and second itme, for mirror if(sphereSphereCheckB(p_m, points[thisI], d_m, d_n)&& - ppPairs.insert(thisI,mrrI) == -1) + ppPairs.insert(thisI,mrrI) == static_cast(-1)) { getFullUpdate++; } diff --git a/src/Interaction/contactSearch/boundaries/periodicBoundaryContactSearch/wallBoundaryContactSearch.cpp b/src/Interaction/contactSearch/boundaries/periodicBoundaryContactSearch/wallBoundaryContactSearch.cpp index 99bbb8ac..acfad226 100644 --- a/src/Interaction/contactSearch/boundaries/periodicBoundaryContactSearch/wallBoundaryContactSearch.cpp +++ b/src/Interaction/contactSearch/boundaries/periodicBoundaryContactSearch/wallBoundaryContactSearch.cpp @@ -113,7 +113,10 @@ pFlow::uint32 pFlow::wallBoundaryContactSearch::findPairsElementRangeCount uint32 nNotInserted = 0; uint32 nThis = pPoints.size(); - + const auto& numElements = numElements_; + const auto& elementBox = elementBox_; + const auto& validBox = validBox_; + Kokkos::parallel_reduce( "pFlow::wallBoundaryContactSearch::findPairsElementRangeCount", deviceRPolicyDynamic(0,nThis), @@ -123,11 +126,11 @@ pFlow::uint32 pFlow::wallBoundaryContactSearch::findPairsElementRangeCount int32x3 ind; if( searchCells.pointIndexInDomain(p, ind) ) { - for(uint32 nTri=0; nTri(-1)) { notInsertedUpdate++; } diff --git a/src/Interaction/contactSearch/contactSearch/contactSearch.hpp b/src/Interaction/contactSearch/contactSearch/contactSearch.hpp index 235403c4..865ec291 100644 --- a/src/Interaction/contactSearch/contactSearch/contactSearch.hpp +++ b/src/Interaction/contactSearch/contactSearch/contactSearch.hpp @@ -140,6 +140,15 @@ public: virtual bool performedBroadSearch()const = 0; + virtual + uint32 updateInterval()const = 0; + + virtual + real sizeRatio()const = 0; + + virtual + real cellExtent()const = 0; + static uniquePtr create( const dictionary& dict, diff --git a/src/Interaction/contactSearch/methods/cellBased/NBS/NBS.hpp b/src/Interaction/contactSearch/methods/cellBased/NBS/NBS.hpp index 11d2e108..9cd5ab6f 100644 --- a/src/Interaction/contactSearch/methods/cellBased/NBS/NBS.hpp +++ b/src/Interaction/contactSearch/methods/cellBased/NBS/NBS.hpp @@ -134,6 +134,16 @@ public: { return 1; } + + real sizeRatio()const + { + return sizeRatio_; + } + + real cellExtent()const + { + return cellExtent_; + } auto getCellIterator([[maybe_unused]] uint32 lvl)const { diff --git a/src/Interaction/contactSearch/methods/cellBased/NBS/NBSLoop.hpp b/src/Interaction/contactSearch/methods/cellBased/NBS/NBSLoop.hpp index ab3eeea7..ecbc63e5 100644 --- a/src/Interaction/contactSearch/methods/cellBased/NBS/NBSLoop.hpp +++ b/src/Interaction/contactSearch/methods/cellBased/NBS/NBSLoop.hpp @@ -43,7 +43,7 @@ while( m != mapperNBS::NoPos) auto lm = m; if(lm>ln) Swap(lm,ln); - if( pairs.insert(lm,ln) == -1) + if( pairs.insert(lm,ln) == static_cast(-1)) { getFullUpdate++; } @@ -86,7 +86,7 @@ while( m != mapperNBS::NoPos) auto ln = n; auto lm = m; if(lm>ln) Swap(lm,ln); - if( pairs.insert(lm,ln) == -1) + if( pairs.insert(lm,ln) == static_cast(-1)) { getFullUpdate++; } diff --git a/src/Interaction/contactSearch/methods/cellBased/NBS/cellsWallLevel0.cpp b/src/Interaction/contactSearch/methods/cellBased/NBS/cellsWallLevel0.cpp index 00932281..e1b6e3a7 100644 --- a/src/Interaction/contactSearch/methods/cellBased/NBS/cellsWallLevel0.cpp +++ b/src/Interaction/contactSearch/methods/cellBased/NBS/cellsWallLevel0.cpp @@ -85,21 +85,26 @@ bool pFlow::cellsWallLevel0::broadSearch bool pFlow::cellsWallLevel0::build(const cells & searchBox) { + const auto& points = points_; + const auto& vertices = vertices_; + const auto& elementBox = elementBox_; + const auto cellExtent = cellExtent_; + Kokkos::parallel_for( "pFlow::cellsWallLevel0::build", deviceRPolicyStatic(0,numElements_), - CLASS_LAMBDA_HD(uint32 i) + LAMBDA_HD(uint32 i) { - auto v = vertices_[i]; - auto p1 = points_[v.x()]; - auto p2 = points_[v.y()]; - auto p3 = points_[v.z()]; + auto v = vertices[i]; + auto p1 = points[v.x()]; + auto p2 = points[v.y()]; + auto p3 = points[v.z()]; realx3 minP; realx3 maxP; - searchBox.extendBox(p1, p2, p3, cellExtent_, minP, maxP); - elementBox_[i] = iBoxType(searchBox.pointIndex(minP), searchBox.pointIndex(maxP)); + searchBox.extendBox(p1, p2, p3, cellExtent, minP, maxP); + elementBox[i] = iBoxType(searchBox.pointIndex(minP), searchBox.pointIndex(maxP)); }); Kokkos::fence(); @@ -153,7 +158,12 @@ pFlow::int32 pFlow::cellsWallLevel0::findPairsElementRangeCount { uint32 getFull =0; - + const auto& elementBox = elementBox_; + const auto& normals = normals_; + const auto& points = points_; + const auto& vertices = vertices_; + const auto cellExtent = cellExtent_; + Kokkos::parallel_reduce( "pFlow::cellsWallLevel0::findPairsElementRangeCount", tpPWContactSearch(numElements_, Kokkos::AUTO), @@ -163,10 +173,10 @@ pFlow::int32 pFlow::cellsWallLevel0::findPairsElementRangeCount const uint32 iTri = teamMember.league_rank(); - const auto triBox = elementBox_[iTri]; + const auto triBox = elementBox[iTri]; const auto triPlane = infinitePlane( - normals_[iTri], - points_[vertices_[iTri].x()]); + normals[iTri], + points[vertices[iTri].x()]); uint32 getFull2 = 0; @@ -186,11 +196,12 @@ pFlow::int32 pFlow::cellsWallLevel0::findPairsElementRangeCount while( n != particleMap.NoPos) { // id is wall id the pair is (particle id, wall id) - if( abs(triPlane.pointFromPlane(pPoints[n]))< pDiams[n]*sizeRatio*cellExtent_) + if( abs(triPlane.pointFromPlane(pPoints[n]))< pDiams[n]*sizeRatio*cellExtent) { if( pairs.insert( static_cast(n), - static_cast(iTri) ) == -1 ) + static_cast(iTri) ) == static_cast(-1) + ) innerUpdate++; } n = particleMap.next(n); diff --git a/src/Interaction/contactSearch/methods/particleWallContactSearchs.hpp b/src/Interaction/contactSearch/methods/particleWallContactSearchs.hpp index d72f9404..1bf5dab3 100644 --- a/src/Interaction/contactSearch/methods/particleWallContactSearchs.hpp +++ b/src/Interaction/contactSearch/methods/particleWallContactSearchs.hpp @@ -113,6 +113,21 @@ public: return false; } + uint32 updateInterval()const + { + return updateInterval_; + } + + real sizeRatio()const + { + return getMethod().sizeRatio(); + } + + real cellExtent()const + { + return getMethod().cellExtent(); + } + }; diff --git a/src/Interaction/interaction/interaction.cpp b/src/Interaction/interaction/interaction.cpp index c6099405..03abbb4b 100644 --- a/src/Interaction/interaction/interaction.cpp +++ b/src/Interaction/interaction/interaction.cpp @@ -72,7 +72,7 @@ pFlow::uniquePtr pFlow::interaction::create gSettings::sleepMiliSeconds( - 1000*(pFlowProcessors().localSize()-pFlowProcessors().localRank()-1)); + 100*(pFlowProcessors().localSize()-pFlowProcessors().localRank()-1)); pOutput.space(2)<<"Creating interaction "< +void pFlow::boundarySphereInteraction::allocatePPPairs() +{ + ppPairs_.reset(nullptr); + ppPairs_ = makeUnique(1); +} + +template +void pFlow::boundarySphereInteraction::allocatePWPairs() +{ + pwPairs_.reset(nullptr); + pwPairs_ = makeUnique(1); +} + template pFlow::boundarySphereInteraction::boundarySphereInteraction( @@ -28,8 +43,6 @@ pFlow::boundarySphereInteraction::boundarySphereInteraction( geometryMotion_(geomMotion), sphParticles_(sphPrtcls) { - ppPairs_ = makeUnique(1); - pwPairs_ = makeUnique(1); } template diff --git a/src/Interaction/sphereInteraction/boundaries/boundarySphereInteraction/boundarySphereInteraction.hpp b/src/Interaction/sphereInteraction/boundaries/boundarySphereInteraction/boundarySphereInteraction.hpp index 56058b6d..708fb61e 100644 --- a/src/Interaction/sphereInteraction/boundaries/boundarySphereInteraction/boundarySphereInteraction.hpp +++ b/src/Interaction/sphereInteraction/boundaries/boundarySphereInteraction/boundarySphereInteraction.hpp @@ -22,7 +22,7 @@ Licence: #include "virtualConstructor.hpp" #include "generalBoundary.hpp" -#include "unsortedContactList.hpp" +#include "sortedContactList.hpp" #include "sphereParticles.hpp" namespace pFlow @@ -51,7 +51,7 @@ public: using IndexType = uint32; using ContactListType = - unsortedContactList; + sortedContactList; private: @@ -60,9 +60,15 @@ private: /// const reference to sphere particles const sphereParticles& sphParticles_; - uniquePtr ppPairs_; + uniquePtr ppPairs_ = nullptr; - uniquePtr pwPairs_; + uniquePtr pwPairs_ = nullptr; + +protected: + + void allocatePPPairs(); + + void allocatePWPairs(); public: @@ -124,15 +130,30 @@ public: return pwPairs_(); } + bool ppPairsAllocated()const + { + if( ppPairs_)return true; + return false; + } + + bool pwPairsAllocated()const + { + if( pwPairs_)return true; + return false; + } + virtual bool sphereSphereInteraction( real dt, - const ContactForceModel& cfModel) + const ContactForceModel& cfModel, + uint32 step) { // for default boundary, no thing to be done - return true; + return false; } + + bool hearChanges ( real t, @@ -149,11 +170,6 @@ public: return true; } - void fill(const std::any& val)override - { - notImplementedFunction; - } - static uniquePtr create( const boundaryBase& boundary, diff --git a/src/Interaction/sphereInteraction/boundaries/periodicBoundarySphereInteraction/periodicBoundarySphereInteraction.cpp b/src/Interaction/sphereInteraction/boundaries/periodicBoundarySphereInteraction/periodicBoundarySphereInteraction.cpp index 38a4336b..3200c4e7 100644 --- a/src/Interaction/sphereInteraction/boundaries/periodicBoundarySphereInteraction/periodicBoundarySphereInteraction.cpp +++ b/src/Interaction/sphereInteraction/boundaries/periodicBoundarySphereInteraction/periodicBoundarySphereInteraction.cpp @@ -30,7 +30,10 @@ pFlow::periodicBoundarySphereInteraction::periodicBoundarySphereIntera { if(boundary.thisBoundaryIndex()%2==1) { - masterInteraction_ = true; + masterInteraction_ = true; + this->allocatePPPairs(); + this->allocatePWPairs(); + } else { @@ -42,10 +45,11 @@ template bool pFlow::periodicBoundarySphereInteraction::sphereSphereInteraction ( real dt, - const ContactForceModel &cfModel + const ContactForceModel &cfModel, + uint32 step ) { - if(!masterInteraction_) return true; + if(!masterInteraction_) return false; pFlow::periodicBoundarySIKernels::sphereSphereInteraction( dt, @@ -61,5 +65,5 @@ bool pFlow::periodicBoundarySphereInteraction::sphereSphereInteraction this->sphParticles().contactForce().deviceViewAll(), this->sphParticles().contactTorque().deviceViewAll()); - return true; + return false; } \ No newline at end of file diff --git a/src/Interaction/sphereInteraction/boundaries/periodicBoundarySphereInteraction/periodicBoundarySphereInteraction.hpp b/src/Interaction/sphereInteraction/boundaries/periodicBoundarySphereInteraction/periodicBoundarySphereInteraction.hpp index e2887787..5f1e7a45 100644 --- a/src/Interaction/sphereInteraction/boundaries/periodicBoundarySphereInteraction/periodicBoundarySphereInteraction.hpp +++ b/src/Interaction/sphereInteraction/boundaries/periodicBoundarySphereInteraction/periodicBoundarySphereInteraction.hpp @@ -83,7 +83,8 @@ public: bool sphereSphereInteraction( real dt, - const ContactForceModel& cfModel)override; + const ContactForceModel& cfModel, + uint32 step)override; }; diff --git a/src/Interaction/sphereInteraction/sphereInteraction/sphereInteraction.cpp b/src/Interaction/sphereInteraction/sphereInteraction/sphereInteraction.cpp index b4546521..2aa10ef8 100644 --- a/src/Interaction/sphereInteraction/sphereInteraction/sphereInteraction.cpp +++ b/src/Interaction/sphereInteraction/sphereInteraction/sphereInteraction.cpp @@ -137,7 +137,9 @@ pFlow::sphereInteraction::sphereInteraction ppInteractionTimer_("sphere-sphere interaction", &this->timers()), pwInteractionTimer_("sphere-wall interaction", &this->timers()), contactListMangementTimer_("contact-list management", &this->timers()), - boundaryInteractionTimer_("interaction for boundary", &this->timers()) + boundaryContactSearchTimer_("contact search for boundary", &this->timers()), + boundaryInteractionTimer_("interaction for boundary", &this->timers()), + contactListBoundaryTimer_("contact-list management for boundary", &this->timers()) { if(!createSphereInteraction()) @@ -179,11 +181,14 @@ bool pFlow::sphereInteraction::iterate() contactListMangementTimer_.pause(); } + contactListBoundaryTimer_.start(); for(uint32 i=0; i<6u; i++) { - boundaryInteraction_[i].ppPairs().beforeBroadSearch(); - boundaryInteraction_[i].pwPairs().beforeBroadSearch(); + auto& BI = boundaryInteraction_[i]; + if(BI.ppPairsAllocated()) BI.ppPairs().beforeBroadSearch(); + if(BI.pwPairsAllocated()) BI.pwPairs().beforeBroadSearch(); } + contactListBoundaryTimer_.pause(); if( sphParticles_.numActive()<=0)return true; @@ -200,21 +205,27 @@ bool pFlow::sphereInteraction::iterate() fatalExit; } + boundaryContactSearchTimer_.start(); for(uint32 i=0; i<6u; i++) { - if( !contactSearch_().boundaryBroadSearch( - i, - iter, - t, - dt, - boundaryInteraction_[i].ppPairs(), - boundaryInteraction_[i].pwPairs())) + auto& BI =boundaryInteraction_[i]; + if(BI.ppPairsAllocated()) { - fatalErrorInFunction<< - "failed to perform broadSearch for boundary index "<::iterate() contactListMangementTimer_.end(); } + contactListBoundaryTimer_.resume(); for(uint32 i=0; i<6u; i++) { - boundaryInteraction_[i].ppPairs().afterBroadSearch(); - boundaryInteraction_[i].pwPairs().afterBroadSearch(); + auto& BI = boundaryInteraction_[i]; + if(BI.ppPairsAllocated()) BI.ppPairs().afterBroadSearch(); + if(BI.pwPairsAllocated()) BI.pwPairs().afterBroadSearch(); } + contactListBoundaryTimer_.end(); + + { + boundaryInteractionTimer_.start(); + std::array requireStep{ + boundaryInteraction_[0].isBoundaryMaster(), + boundaryInteraction_[1].isBoundaryMaster(), + boundaryInteraction_[2].isBoundaryMaster(), + boundaryInteraction_[3].isBoundaryMaster(), + boundaryInteraction_[4].isBoundaryMaster(), + boundaryInteraction_[5].isBoundaryMaster()}; + int step = 1; + const auto& cfModel = this->forceModel_(); + while( std::any_of(requireStep.begin(), requireStep.end(), [](bool v) { return v==true; })) + { + for(uint32 i=0; i<6u; i++) + { + if(step==1u || requireStep[i] ) + { + requireStep[i] = boundaryInteraction_[i].sphereSphereInteraction( + dt, + this->forceModel_(), + step + ); + } + } + step++; + } + boundaryInteractionTimer_.pause(); + } ppInteractionTimer_.start(); sphereSphereInteraction(); ppInteractionTimer_.end(); @@ -239,14 +282,36 @@ bool pFlow::sphereInteraction::iterate() sphereWallInteraction(); pwInteractionTimer_.end(); - boundaryInteractionTimer_.start(); - for(uint32 i=0; i<6u; i++) { - boundaryInteraction_[i].sphereSphereInteraction( - dt, - this->forceModel_()); + boundaryInteractionTimer_.resume(); + std::array requireStep{ + !boundaryInteraction_[0].isBoundaryMaster(), + !boundaryInteraction_[1].isBoundaryMaster(), + !boundaryInteraction_[2].isBoundaryMaster(), + !boundaryInteraction_[3].isBoundaryMaster(), + !boundaryInteraction_[4].isBoundaryMaster(), + !boundaryInteraction_[5].isBoundaryMaster()}; + + int step = 2; + const auto& cfModel = this->forceModel_(); + while(std::any_of(requireStep.begin(), requireStep.end(), [](bool v) { return v==true; })) + { + for(uint32 i=0; i<6u; i++) + { + if(requireStep[i]) + { + requireStep[i] = boundaryInteraction_[i].sphereSphereInteraction( + dt, + this->forceModel_(), + step + ); + } + } + step++; } boundaryInteractionTimer_.end(); + } + return true; } diff --git a/src/Interaction/sphereInteraction/sphereInteraction/sphereInteraction.hpp b/src/Interaction/sphereInteraction/sphereInteraction/sphereInteraction.hpp index 4ccb9836..74cb2442 100644 --- a/src/Interaction/sphereInteraction/sphereInteraction/sphereInteraction.hpp +++ b/src/Interaction/sphereInteraction/sphereInteraction/sphereInteraction.hpp @@ -96,9 +96,12 @@ private: /// timer for managing contact lists (only inernal points) Timer contactListMangementTimer_; + Timer boundaryContactSearchTimer_; /// timer for boundary interaction time Timer boundaryInteractionTimer_; + Timer contactListBoundaryTimer_; + bool createSphereInteraction(); diff --git a/src/Interaction/sphereInteraction/sphereInteraction/sphereInteractionKernels.hpp b/src/Interaction/sphereInteraction/sphereInteraction/sphereInteractionKernels.hpp index c581d598..751c0fec 100644 --- a/src/Interaction/sphereInteraction/sphereInteraction/sphereInteractionKernels.hpp +++ b/src/Interaction/sphereInteraction/sphereInteraction/sphereInteractionKernels.hpp @@ -261,7 +261,7 @@ struct pwInteractionFunctor int32 propId_i = propId_[i]; int32 wPropId_j = wPropId_[tj]; - realx3 FCn, FCt, Mri, Mrj, Mij, Mji; + realx3 FCn, FCt, Mri, Mrj, Mij; //output<< "before "<(-1)) { if( ((position_[n]-p).length() - 0.5*(diameters_[n]+d )) <= 0.0 ) { @@ -85,7 +85,7 @@ pFlow::collisionCheck::mapLastAddedParticle() "size mismatch of next and position"<(-1)); const auto& p = position_[n]; if(!searchBox_.isInside(p)) diff --git a/src/Particles/Insertion/insertion/insertion.cpp b/src/Particles/Insertion/insertion/insertion.cpp index d8d0881d..bb551a47 100644 --- a/src/Particles/Insertion/insertion/insertion.cpp +++ b/src/Particles/Insertion/insertion/insertion.cpp @@ -74,9 +74,6 @@ pFlow::insertion::pStruct() const return particles_.pStruct(); } - - - bool pFlow::insertion::readInsertionDict() { diff --git a/src/Particles/Insertion/insertion/insertion.hpp b/src/Particles/Insertion/insertion/insertion.hpp index f8419aa3..f39f5c35 100644 --- a/src/Particles/Insertion/insertion/insertion.hpp +++ b/src/Particles/Insertion/insertion/insertion.hpp @@ -32,7 +32,9 @@ class pointStructure; /** * Base class for particle insertion */ -class insertion : public fileDictionary +class insertion +: + public fileDictionary { private: @@ -118,6 +120,8 @@ public: /*/// read from iIstream virtual bool read(iIstream& is) = 0;*/ + using fileDictionary::write; + /// write to iOstream bool write(iOstream& os, const IOPattern& iop)const override ; }; diff --git a/src/Particles/SphereParticles/sphereParticles/sphereParticles.cpp b/src/Particles/SphereParticles/sphereParticles/sphereParticles.cpp index ea29aee8..df7b98b8 100644 --- a/src/Particles/SphereParticles/sphereParticles/sphereParticles.cpp +++ b/src/Particles/SphereParticles/sphereParticles/sphereParticles.cpp @@ -389,7 +389,9 @@ pFlow::sphereParticles::sphereParticles( intPredictTimer_( "Integration-predict", &this->timers() ), intCorrectTimer_( - "Integration-correct", &this->timers() ) + "Integration-correct", &this->timers() ), + fieldUpdateTimer_( + "fieldUpdate", &this->timers() ) { auto intMethod = control.settingsDict().getVal("integrationMethod"); @@ -514,13 +516,14 @@ bool pFlow::sphereParticles::beforeIteration() rVelIntegration_().predict(dt,rVelocity_, rAcceleration_); intPredictTimer_.end(); + fieldUpdateTimer_.start(); propertyId_.updateBoundariesSlaveToMasterIfRequested(); diameter_.updateBoundariesSlaveToMasterIfRequested(); mass_.updateBoundariesSlaveToMasterIfRequested(); I_.updateBoundariesSlaveToMasterIfRequested(); rVelocity_.updateBoundariesSlaveToMasterIfRequested(); rAcceleration_.updateBoundariesSlaveToMasterIfRequested(); - + fieldUpdateTimer_.end(); return true; } diff --git a/src/Particles/SphereParticles/sphereParticles/sphereParticles.hpp b/src/Particles/SphereParticles/sphereParticles/sphereParticles.hpp index fc0aa82c..d71e751a 100644 --- a/src/Particles/SphereParticles/sphereParticles/sphereParticles.hpp +++ b/src/Particles/SphereParticles/sphereParticles/sphereParticles.hpp @@ -79,10 +79,10 @@ private: /// timer for integration computations (correction step) Timer intCorrectTimer_; - + Timer fieldUpdateTimer_; private: - bool initializeParticles(); + bool getParticlesInfoFromShape( const wordVector& shapeNames, @@ -119,11 +119,14 @@ public: */ /*bool insertParticles ( - const realx3Vector& position, + const realx3Vector& position, const wordVector& shapes, const setFieldList& setField ) override ;*/ + // TODO: make this method private later + bool initializeParticles(); + /// const reference to shapes object const auto& spheres() const { diff --git a/src/Particles/dynamicPointStructure/dynamicPointStructure.cpp b/src/Particles/dynamicPointStructure/dynamicPointStructure.cpp index c57475ec..6ad16556 100644 --- a/src/Particles/dynamicPointStructure/dynamicPointStructure.cpp +++ b/src/Particles/dynamicPointStructure/dynamicPointStructure.cpp @@ -38,6 +38,7 @@ pFlow::dynamicPointStructure::dynamicPointStructure *this, zero3 ), + velocityUpdateTimer_("velocity boundary update", &timers()), integrationMethod_ ( control.settingsDict().getVal("integrationMethod") @@ -81,11 +82,11 @@ pFlow::dynamicPointStructure::dynamicPointStructure bool pFlow::dynamicPointStructure::beforeIteration() { - return pointStructure::beforeIteration(); - /*real dt = this->dt(); - - auto& acc = time().lookupObject("acceleration"); - return predict(dt, acc);*/ + if(!pointStructure::beforeIteration())return false; + velocityUpdateTimer_.start(); + velocity_.updateBoundariesSlaveToMasterIfRequested(); + velocityUpdateTimer_.end(); + return true; } bool pFlow::dynamicPointStructure::iterate() diff --git a/src/Particles/dynamicPointStructure/dynamicPointStructure.hpp b/src/Particles/dynamicPointStructure/dynamicPointStructure.hpp index 0e6a9fcd..f82f1dca 100644 --- a/src/Particles/dynamicPointStructure/dynamicPointStructure.hpp +++ b/src/Particles/dynamicPointStructure/dynamicPointStructure.hpp @@ -44,6 +44,8 @@ private: uniquePtr integrationVel_ = nullptr; + Timer velocityUpdateTimer_; + /// @brief integration method for velocity and position word integrationMethod_; diff --git a/src/Particles/particles/particles.cpp b/src/Particles/particles/particles.cpp index 187f45da..a8fc970a 100644 --- a/src/Particles/particles/particles.cpp +++ b/src/Particles/particles/particles.cpp @@ -74,7 +74,8 @@ pFlow::particles::particles(systemControl& control) dynPointStruct_, zero3 ), - idHandler_(particleIdHandler::create(dynPointStruct_)) + idHandler_(particleIdHandler::create(dynPointStruct_)), + baseFieldBoundaryUpdateTimer_("baseFieldBoundaryUpdate",&timers()) { this->addToSubscriber(dynPointStruct_); @@ -84,18 +85,18 @@ pFlow::particles::particles(systemControl& control) bool pFlow::particles::beforeIteration() { - zeroForce(); - zeroTorque(); - if( !dynPointStruct_.beforeIteration()) { return false; } + zeroForce(); + zeroTorque(); + baseFieldBoundaryUpdateTimer_.start(); shapeIndex_.updateBoundariesSlaveToMasterIfRequested(); accelertion_.updateBoundariesSlaveToMasterIfRequested(); idHandler_().updateBoundariesSlaveToMasterIfRequested(); - + baseFieldBoundaryUpdateTimer_.end(); return true; } diff --git a/src/Particles/particles/particles.hpp b/src/Particles/particles/particles.hpp index 1f3c7ab1..12e3da8c 100644 --- a/src/Particles/particles/particles.hpp +++ b/src/Particles/particles/particles.hpp @@ -54,6 +54,8 @@ private: /// handling new ids for new particles uniquePtr idHandler_ = nullptr; + Timer baseFieldBoundaryUpdateTimer_; + /// messages for this objects static inline const message defaultMessage_{ message::DEFAULT }; diff --git a/src/Particles/particles/regularParticleIdHandler/regularParticleIdHandler.cpp b/src/Particles/particles/regularParticleIdHandler/regularParticleIdHandler.cpp index 1916f430..b8a55aba 100644 --- a/src/Particles/particles/regularParticleIdHandler/regularParticleIdHandler.cpp +++ b/src/Particles/particles/regularParticleIdHandler/regularParticleIdHandler.cpp @@ -15,8 +15,11 @@ pFlow::regularParticleIdHandler::regularParticleIdHandler pFlow::Pair pFlow::regularParticleIdHandler::getIdRange(uint32 nNewParticles) { + + if(nNewParticles==0) return {0,0}; + uint32 startId; - if(maxId_==-1) + if(maxId_== static_cast(-1)) { startId = 0; } @@ -37,7 +40,7 @@ bool pFlow::regularParticleIdHandler::initialIdCheck() uint32 maxId = max( *this ); /// particles should get ids from 0 to size-1 - if(maxId == -1) + if(maxId == static_cast(-1)) { fillSequence(*this,0u); maxId_ = size()-1; diff --git a/src/Particles/particles/regularParticleIdHandler/regularParticleIdHandler.hpp b/src/Particles/particles/regularParticleIdHandler/regularParticleIdHandler.hpp index a237e0da..3e4f8b41 100644 --- a/src/Particles/particles/regularParticleIdHandler/regularParticleIdHandler.hpp +++ b/src/Particles/particles/regularParticleIdHandler/regularParticleIdHandler.hpp @@ -31,7 +31,7 @@ class regularParticleIdHandler { private: - uint32 maxId_ = -1; + uint32 maxId_ = static_cast(-1); bool initialIdCheck()override; public: diff --git a/src/Particles/particles/shape/baseShapeNames.hpp b/src/Particles/particles/shape/baseShapeNames.hpp index bc513b06..7e0a00c3 100644 --- a/src/Particles/particles/shape/baseShapeNames.hpp +++ b/src/Particles/particles/shape/baseShapeNames.hpp @@ -132,7 +132,7 @@ public: return true; } } - idx = -1; + idx = static_cast(-1); return false; } @@ -144,6 +144,8 @@ public: // - IO + using fileDictionary::write; + bool write(iOstream& os)const override; }; diff --git a/src/phasicFlow/CMakeLists.txt b/src/phasicFlow/CMakeLists.txt index de526e20..7de9ceda 100644 --- a/src/phasicFlow/CMakeLists.txt +++ b/src/phasicFlow/CMakeLists.txt @@ -59,6 +59,7 @@ Timer/Timers.cpp containers/Vector/Vectors.cpp +containers/VectorHD/wordVectorHost.cpp containers/VectorHD/VectorSingles.cpp containers/Field/Fields.cpp containers/symArrayHD/symArrays.cpp @@ -94,6 +95,8 @@ structuredData/pointStructure/selectors/selectorStridedRange/selectorStridedRang structuredData/pointStructure/selectors/selectorRandomPoints/selectorRandomPoints.cpp structuredData/pointStructure/selectors/selectBox/selectBox.cpp structuredData/pointStructure/selectors/selectorGeometric/selectorGeometrics.cpp +structuredData/pointStructure/pointStructure/pointSorting/mortonIndexing.cpp +structuredData/pointStructure/pointStructure/pointSorting/pointSorting.cpp triSurface/subSurface.cpp triSurface/triSurface.cpp diff --git a/src/phasicFlow/Kokkos/ViewAlgorithms.hpp b/src/phasicFlow/Kokkos/ViewAlgorithms.hpp index db071f0b..28e32abd 100644 --- a/src/phasicFlow/Kokkos/ViewAlgorithms.hpp +++ b/src/phasicFlow/Kokkos/ViewAlgorithms.hpp @@ -421,10 +421,10 @@ binarySearch( ) { if (end <= start) - return -1; + return static_cast(-1); if (auto res = binarySearch_(view.data() + start, end - start, val); - res != -1) + res != static_cast(-1)) { return res + start; } diff --git a/src/phasicFlow/commandLine/CLI/TypeTools.hpp b/src/phasicFlow/commandLine/CLI/TypeTools.hpp index 2b87ec60..667d3eac 100644 --- a/src/phasicFlow/commandLine/CLI/TypeTools.hpp +++ b/src/phasicFlow/commandLine/CLI/TypeTools.hpp @@ -146,11 +146,11 @@ template class is_direct_constructible { static auto test(int, std::true_type) -> decltype( // NVCC warns about narrowing conversions here #ifdef __CUDACC__ -#pragma diag_suppress 2361 +#pragma nv_diag_suppress 2361 #endif TT { std::declval() } #ifdef __CUDACC__ -#pragma diag_default 2361 +#pragma nv_diag_default 2361 #endif , std::is_move_assignable()); diff --git a/src/phasicFlow/containers/Field/Field.hpp b/src/phasicFlow/containers/Field/Field.hpp index 887b994d..738811bc 100644 --- a/src/phasicFlow/containers/Field/Field.hpp +++ b/src/phasicFlow/containers/Field/Field.hpp @@ -24,6 +24,7 @@ Licence: #include "types.hpp" #include "VectorSingle.hpp" +#include "wordVectorHost.hpp" #include "Vector.hpp" #include "streams.hpp" diff --git a/src/phasicFlow/containers/Field/Fields.cpp b/src/phasicFlow/containers/Field/Fields.cpp index b0cbad75..e76daa11 100644 --- a/src/phasicFlow/containers/Field/Fields.cpp +++ b/src/phasicFlow/containers/Field/Fields.cpp @@ -31,5 +31,4 @@ template class pFlow::Field; template class pFlow::Field; - - +template class pFlow::Field; diff --git a/src/phasicFlow/containers/List/ListPtr/ListPtrI.hpp b/src/phasicFlow/containers/List/ListPtr/ListPtrI.hpp index 77ef5a99..861239de 100644 --- a/src/phasicFlow/containers/List/ListPtr/ListPtrI.hpp +++ b/src/phasicFlow/containers/List/ListPtr/ListPtrI.hpp @@ -39,6 +39,7 @@ inline bool pFlow::ListPtr::copy(const ListPtrType& src) } template +inline T* pFlow::ListPtr::ptr(size_t i) { @@ -51,6 +52,7 @@ T* pFlow::ListPtr::ptr(size_t i) } template +inline const T* pFlow::ListPtr::ptr ( size_t i @@ -66,6 +68,7 @@ const T* pFlow::ListPtr::ptr } template +inline auto pFlow::ListPtr::pos ( size_t i @@ -84,6 +87,7 @@ auto pFlow::ListPtr::pos } template +inline auto pFlow::ListPtr::pos ( size_t i @@ -102,6 +106,7 @@ auto pFlow::ListPtr::pos } template +inline pFlow::ListPtr::ListPtr ( const ListPtrType& src @@ -119,6 +124,7 @@ pFlow::ListPtr::ListPtr template +inline pFlow::ListPtr& pFlow::ListPtr::operator= ( const ListPtrType& rhs @@ -144,6 +150,7 @@ pFlow::ListPtr& pFlow::ListPtr::operator= } template +inline pFlow::uniquePtr pFlow::ListPtr::set ( size_t i, T* ptr @@ -155,6 +162,7 @@ pFlow::uniquePtr pFlow::ListPtr::set template +inline pFlow::uniquePtr pFlow::ListPtr::set ( size_t i, @@ -179,6 +187,7 @@ pFlow::uniquePtr pFlow::ListPtr::set template template +inline pFlow::uniquePtr pFlow::ListPtr::setSafe ( size_t i, @@ -190,6 +199,7 @@ pFlow::uniquePtr pFlow::ListPtr::setSafe } template +inline void pFlow::ListPtr::push_back ( T* ptr @@ -199,6 +209,7 @@ void pFlow::ListPtr::push_back } template +inline void pFlow::ListPtr::push_back(uniquePtr&& ptr) { list_.push_back( ptr.release() ); @@ -206,13 +217,15 @@ void pFlow::ListPtr::push_back(uniquePtr&& ptr) template template +inline void pFlow::ListPtr::push_backSafe(Args&&... args) { - auto ptr=makeUnique(std::forward(args)...) ; + uniquePtr ptr = makeUnique(std::forward(args)...) ; push_back(ptr); } template +inline T& pFlow::ListPtr::operator[] ( size_t i @@ -231,6 +244,7 @@ T& pFlow::ListPtr::operator[] } template +inline const T& pFlow::ListPtr::operator[] ( size_t i @@ -248,18 +262,21 @@ const T& pFlow::ListPtr::operator[] } template +inline size_t pFlow::ListPtr::size()const { return list_.size(); } template +inline auto pFlow::ListPtr::empty() const { return list_.emtpy(); } template +inline pFlow::uniquePtr pFlow::ListPtr::release ( size_t i @@ -273,15 +290,14 @@ pFlow::uniquePtr pFlow::ListPtr::release template +inline void pFlow::ListPtr::clear() { - int i =0; for( auto iter = list_.begin(); iter != list_.end(); ++iter ) { if(*iter != nullptr) - { - + { delete *iter; *iter = nullptr; } @@ -291,6 +307,7 @@ void pFlow::ListPtr::clear() } template +inline void pFlow::ListPtr::clear ( size_t i diff --git a/src/phasicFlow/containers/VectorHD/VectorSingle.cpp b/src/phasicFlow/containers/VectorHD/VectorSingle.cpp index 24248594..0e6f02b2 100644 --- a/src/phasicFlow/containers/VectorHD/VectorSingle.cpp +++ b/src/phasicFlow/containers/VectorHD/VectorSingle.cpp @@ -55,30 +55,14 @@ void pFlow::VectorSingle::changeCapacity bool withInit ) { - if constexpr( isTriviallyCopyable_ ) + + if(withInit) { - if(withInit) - { - resizeInit(view_, actualCap); - } - else - { - resizeNoInit(view_, actualCap); - } - } - else if constexpr( isHostAccessible_ ) - { - viewType newView(view_.label(), actualCap); - - for(auto i=0u; i::reallocateCapacitySize uint32 s ) { - if constexpr (isTriviallyCopyable_) - { - reallocNoInit(view_, cap); - } - else - { - viewType newView(view_.label(), cap); - view_ = newView; - } - - return setSize(s); + reallocNoInit(view_, cap); + return setSize(s); } template @@ -209,21 +184,19 @@ pFlow::VectorSingle::VectorSingle : VectorSingle(name, src.capacity(), src.size(), RESERVE()) { - if constexpr(isTriviallyCopyable_) - { - copy(deviceView(), src.deviceView()); - } - else if constexpr( isHostAccessible_) - { - for(auto i=0u; i +pFlow::VectorSingle::VectorSingle +( + const word& name, + const ViewType1D& src +) +: + VectorSingle(name, src.size(), src.size(), RESERVE()) +{ + copy(deviceView(), src); } template @@ -299,18 +272,7 @@ INLINE_FUNCTION_H auto pFlow::VectorSingle::hostViewAll()const { auto hView = Kokkos::create_mirror_view(view_); - if constexpr(isTriviallyCopyable_) - { - copy(hView, view_); - } - else if constexpr( isHostAccessible_ ) - { - // nothing to be done, since it is already a host memory - } - else - { - static_assert("hostViewAll is not valid for non-trivially copyable data type on device memory"); - } + copy(hView, view_); return hView; } @@ -319,19 +281,7 @@ INLINE_FUNCTION_H auto pFlow::VectorSingle::hostView()const { auto hView = Kokkos::create_mirror_view(deviceView()); - if constexpr(isTriviallyCopyable_) - { - copy(hView, deviceView()); - } - else if constexpr( isHostAccessible_ ) - { - // nothing to be done, since it is already a host memory - } - else - { - static_assert("hostView is not valid for non-trivially copyable data type on device memory"); - } - + copy(hView, deviceView()); return hView; } @@ -367,6 +317,7 @@ template INLINE_FUNCTION_H void pFlow::VectorSingle::reserve(uint32 cap) { + if(cap == capacity() ) return; changeCapacity(cap); } @@ -470,23 +421,9 @@ void pFlow::VectorSingle::assign changeSize(srcSize); } - if constexpr( isTriviallyCopyable_ ) - { - // - unmanaged view in the host - hostViewType1D temp(src.data(), srcSize ); - copy(deviceView(), temp); - } - else if constexpr( isHostAccessible_) - { - for(auto i=0u; i temp(src.data(), srcSize ); + copy(deviceView(), temp); } @@ -516,21 +453,8 @@ void pFlow::VectorSingle::assignFromHost(const VectorTypeHost& sr changeSize(srcSize); } - if constexpr(isTriviallyCopyable_) - { - copy(deviceView(), src.hostView()); - } - else if constexpr( isHostAccessible_) - { - for(auto i=0u; i @@ -553,22 +477,48 @@ void pFlow::VectorSingle::assign changeSize(srcSize); } + copy(deviceView(), src.deviceView()); +} - if constexpr(isTriviallyCopyable_) - { - copy(deviceView(), src.deviceView()); + +template +template +INLINE_FUNCTION_H +void pFlow::VectorSingle::assignFromDevice( + const VectorSingle& src, + bool srcCapacity +) +{ + uint32 srcSize = src.size(); + uint32 srcCap = src.capacity(); + + if(srcCapacity && srcCap != capacity()){ + reallocateCapacitySize(srcCap, srcSize); } - else if constexpr( isHostAccessible_) - { - for(auto i=0u; i +INLINE_FUNCTION_H +void pFlow::VectorSingle::append(const ViewType1D& appVec) +{ + uint32 appSize = appVec.size(); + if(appSize == 0) return; + + uint32 oldS = size(); + uint32 newSize = oldS + appSize; + + changeSize(newSize); + + auto appendView = Kokkos::subview( + view_, + Kokkos::make_pair(oldS, newSize)); + + copy(appendView, appVec); + } template @@ -589,22 +539,7 @@ void pFlow::VectorSingle::append hostViewType1D temp(appVec.data(), srcSize ); auto dest = Kokkos::subview(view_, Kokkos::make_pair(oldSize,newSize)); - if constexpr( isTriviallyCopyable_) - { - copy(dest, temp); - } - else if constexpr( isHostAccessible_) - { - for(auto i=0u; i @@ -626,21 +561,8 @@ void pFlow::VectorSingle::append view_, Kokkos::make_pair(oldS, newSize)); - if constexpr( isTriviallyCopyable_) - { - copy(appendView, appVec.deviceView()); - } - else if constexpr( isHostAccessible_) - { - for(auto i=0u; i @@ -813,7 +735,7 @@ bool pFlow::VectorSingle::insertSetElement template INLINE_FUNCTION_H -bool pFlow::VectorSingle::reorderItems(uint32IndexContainer indices) +bool pFlow::VectorSingle::reorderItems(const uint32IndexContainer& indices) { if(indices.size() == 0) { @@ -831,10 +753,8 @@ bool pFlow::VectorSingle::reorderItems(uint32IndexContainer indic return false; } - uint32 newSize = indices.size(); + uint32 newSize = indices.size(); - setSize(newSize); - viewType sortedView(this->name(), newSize); using policy = Kokkos::RangePolicy< execution_space,Kokkos::IndexType>; @@ -875,7 +795,10 @@ bool pFlow::VectorSingle::reorderItems(uint32IndexContainer indic Kokkos::fence(); } - copy(deviceView(), sortedView); - + setSize(newSize); + + copy(deviceView(), sortedView); + + return true; } \ No newline at end of file diff --git a/src/phasicFlow/containers/VectorHD/VectorSingle.hpp b/src/phasicFlow/containers/VectorHD/VectorSingle.hpp index 256ccb6c..f64f78b1 100644 --- a/src/phasicFlow/containers/VectorHD/VectorSingle.hpp +++ b/src/phasicFlow/containers/VectorHD/VectorSingle.hpp @@ -39,9 +39,6 @@ Licence: namespace pFlow { -//- Forward -template -class VectorSingle; template class VectorSingle @@ -101,6 +98,8 @@ private: static constexpr bool isTriviallyCopyable_ = std::is_trivially_copyable_v; + static_assert(isTriviallyCopyable_, "This type is not trivially copyable"); + /// Evaluate capacity based on the input size static INLINE_FUNCTION_H uint32 evalCapacity(uint32 n) { @@ -158,6 +157,9 @@ public: /// Copy construct with a new name (perform deep copy) VectorSingle(const word& name, const VectorSingle& src); + + /// Copy construct with a new name (perform deep copy) + VectorSingle(const word& name, const ViewType1D& src); /// Copy assignment (perform deep copy from rhs to *this) VectorSingle& operator = (const VectorSingle& rhs) ; @@ -287,25 +289,10 @@ public: template INLINE_FUNCTION_H - void assignFromDevice(const VectorSingle& src, bool srcCapacity = true) - { - uint32 srcSize = src.size(); - uint32 srcCap = src.capacity(); + void assignFromDevice(const VectorSingle& src, bool srcCapacity = true); - if(srcCapacity && srcCap != capacity()){ - reallocateCapacitySize(srcCap, srcSize); - } - else { - changeSize(srcSize); - } - - if constexpr(isTriviallyCopyable_){ - copy(deviceView(), src.deviceView()); - } - else{ - static_assert("Not a valid operation for this data type "); - } - } + INLINE_FUNCTION_H + void append(const ViewType1D& appVec); INLINE_FUNCTION_H void append(const std::vector& appVec); @@ -331,7 +318,7 @@ public: const ViewType1D vals); INLINE_FUNCTION_H - bool reorderItems(uint32IndexContainer indices); + bool reorderItems(const uint32IndexContainer& indices); /// @brief push a new element at the end (host call only) /// resize if necessary and works on host accessible vector. diff --git a/src/phasicFlow/containers/VectorHD/VectorSingles.hpp b/src/phasicFlow/containers/VectorHD/VectorSingles.hpp index d0c3acdf..2479715a 100644 --- a/src/phasicFlow/containers/VectorHD/VectorSingles.hpp +++ b/src/phasicFlow/containers/VectorHD/VectorSingles.hpp @@ -25,6 +25,7 @@ Licence: #include "types.hpp" #include "VectorSingle.hpp" +#include "wordVectorHost.hpp" namespace pFlow { @@ -77,6 +78,8 @@ typedef VectorSingle realx3x3Vector_D; typedef VectorSingle realx3x3Vector_H; +typedef VectorSingle wordVector_H; + } #endif \ No newline at end of file diff --git a/src/phasicFlow/containers/VectorHD/wordVectorHost.cpp b/src/phasicFlow/containers/VectorHD/wordVectorHost.cpp new file mode 100644 index 00000000..1941a239 --- /dev/null +++ b/src/phasicFlow/containers/VectorHD/wordVectorHost.cpp @@ -0,0 +1,23 @@ +/*------------------------------- phasicFlow --------------------------------- + O C enter of + O O E ngineering and + O O M ultiscale modeling of + OOOOOOO F luid flow +------------------------------------------------------------------------------ + Copyright (C): www.cemf.ir + email: hamid.r.norouzi AT gmail.com +------------------------------------------------------------------------------ +Licence: + This file is part of phasicFlow code. It is a free software for simulating + granular and multiphase flows. You can redistribute it and/or modify it under + the terms of GNU General Public License v3 or any other later versions. + + phasicFlow is distributed to help others in their research in the field of + granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the + implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + +-----------------------------------------------------------------------------*/ + + +#include "wordVectorHost.hpp" + diff --git a/src/phasicFlow/containers/VectorHD/wordVectorHost.hpp b/src/phasicFlow/containers/VectorHD/wordVectorHost.hpp new file mode 100644 index 00000000..e638bd27 --- /dev/null +++ b/src/phasicFlow/containers/VectorHD/wordVectorHost.hpp @@ -0,0 +1,558 @@ +/*------------------------------- phasicFlow --------------------------------- + O C enter of + O O E ngineering and + O O M ultiscale modeling of + OOOOOOO F luid flow +------------------------------------------------------------------------------ + Copyright (C): www.cemf.ir + email: hamid.r.norouzi AT gmail.com +------------------------------------------------------------------------------ +Licence: + This file is part of phasicFlow code. It is a free software for simulating + granular and multiphase flows. You can redistribute it and/or modify it under + the terms of GNU General Public License v3 or any other later versions. + + phasicFlow is distributed to help others in their research in the field of + granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the + implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + +-----------------------------------------------------------------------------*/ +#ifndef __wordVectorHost_hpp__ +#define __wordVectorHost_hpp__ + +#include "VectorSingle.hpp" +#include "Vector.hpp" + +namespace pFlow +{ + + +template<> +class VectorSingle +{ +public: + + //- typedefs for accessing data + + using VectorType = VectorSingle; + + using VectorTypeHost = VectorSingle; + + using iterator = word*; + + using const_iterator = const word*; + + using reference = word&; + + using const_reference = const word&; + + using value_type = word; + + using pointer = word*; + + using const_pointer = const word*; + + //- typedefs related to memory management + + using viewType = ViewType1D; + + using device_type = typename viewType::device_type; + + using memory_space = typename viewType::memory_space; + + using execution_space = typename viewType::execution_space; + +private: + + // - Data members + + Vector container_; + + mutable viewType unManagedView_; + + // - protected members and methods + + /// Is the memory of this vector accessible from Host + static constexpr + bool isHostAccessible_ = true; + + /// Is the memory of this vector accessiple from Divce + static constexpr + bool isDeviceAccessible_ = false; + + static constexpr + bool isTriviallyCopyable_ = false; + + +public: + + /// Type info + TypeInfoTemplateNV111("VectorSingle", word , memoerySpaceName()); + + //// - Constructors + + /// Empty vector + VectorSingle() = default; + + /// Empty vector with a name (capacity = 2) + explicit VectorSingle(const word& name) + : + container_(name) + {} + + /// Vector with name and size n + VectorSingle(const word& name, uint32 n) + : + container_(name, n) + {} + + /// Vector with name, size and value + VectorSingle(const word& name, uint32 n, const word& val) + : + container_(name, n, val) + {} + + /// Vector with name, size (n) and reserved capacity + VectorSingle(const word& name, uint32 cap, uint32 n, const RESERVE& r ) + : + container_(name, cap, n, r) + {} + + /// Construct with a name and form std::vector (host memory) + VectorSingle(const word& name, const std::vector & src) + : + container_(name, src) + {} + + /// Construct with a name and form std::vector (host memory) and with a desired capacity. + VectorSingle(const word& name, const std::vector & src, uint32 cap) + : + container_(name, src, cap) + {} + + /// Copy construct (performs deep copy) + VectorSingle(const VectorSingle& src)=default; + + /// Copy construct with a new name (perform deep copy) + VectorSingle(const word& name, const VectorSingle& src) + : + container_(name, src.container_) + {} + + + /// Copy assignment (perform deep copy from rhs to *this) + VectorSingle& operator = (const VectorSingle& rhs)=default; + + /// Move construct + VectorSingle(VectorSingle&&) = default; + + /// Move assignment + VectorSingle& operator= (VectorSingle&&) = default; + + /// @brief Descructor + /// This may not destroy the underlying memory, sice view is + /// shared_ptr and maybe referenced by another object too + ~VectorSingle() = default; + + /// Clone as a uniquePtr (perform deep copy) + INLINE_FUNCTION_H + uniquePtr clone() const + { + return makeUnique(*this); + } + + + //// - Methods + + /// Return *this + INLINE_FUNCTION_H + VectorType& VectorField() + { + return *this; + } + + /// Return *this + INLINE_FUNCTION_H + const VectorType& VectorField()const + { + return *this; + } + + /// Device view range [0,capcity) + INLINE_FUNCTION_H + auto& deviceViewAll() + { + // return un-managed view + unManagedView_ = viewType(container_.data(), container_.capacity()); + return unManagedView_; + } + + /// Device view range [0,capcity) + INLINE_FUNCTION_H + const auto& deviceViewAll() const + { + // return un-managed view + unManagedView_ = viewType(const_cast(container_.data()), container_.capacity()); + return unManagedView_; + } + + /// Device view range [0, size) + INLINE_FUNCTION_H + auto deviceView()const + { + return viewType(const_cast(container_.data()), container_.size()); + } + + /// Return a view accessible on Host in range [0,capacity) + INLINE_FUNCTION_H + auto hostViewAll()const + { + return viewType(const_cast(container_.data()), container_.capacity()); + } + + + /// Return a view accessible on Host in range [0,size) + INLINE_FUNCTION_H + auto hostView()const + { + return viewType(const_cast(container_.data()), container_.size()); + } + + /// Name of the vector + INLINE_FUNCTION_H + word name()const + { + return container_.name(); + } + + + /// Size of the vector + INLINE_FUNCTION_H + uint32 size()const + { + return container_.size(); + } + + + // Capcity of the vector + INLINE_FUNCTION_H + uint32 capacity()const + { + return container_.capacity(); + } + + /// If vector is empty + INLINE_FUNCTION_H + bool empty()const + { + return container_.size()==0uL; + } + + + /// Reserve capacity for vector + /// Preserve the content. + INLINE_FUNCTION_H + void reserve(uint32 cap) + { + container_.reserve(cap); + } + + /// Reallocate memory to new cap and set size to 0. + /*INLINE_FUNCTION_H + void reallocate(uint32 cap); + + /// Reallocate memory to new cap and set size to newSize. + /// Do not preserve the content + INLINE_FUNCTION_H + void reallocate(uint32 cap, uint32 newSize);*/ + + /// Resize the vector and preserve the content + INLINE_FUNCTION_H + void resize(uint32 n) + { + container_.resize(n); + } + + /// Resize the vector and assign the value to it. + INLINE_FUNCTION_H + void resize(uint32 n, const word& val) + { + container_.resize(n, val); + } + + /// Clear the vector, but keep the allocated memory unchanged + INLINE_FUNCTION_H + void clear() + { + container_.clear(); + } + + /// Fill the range [0,size) with val + INLINE_FUNCTION_H + void fill(const word& val) + { + container_.fill(val); + } + + INLINE_FUNCTION_H + void fill(rangeU32 r, const word& val) + { + container_.fill(r.start(), r.end(), val); + } + + /// Change size of the vector and assign val to vector and + INLINE_FUNCTION_H + void assign(size_t n, const word& val) + { + container_.assign(n, val); + } + + /// Assign source vector with specified capacity. + /// The size of *this becomes the size of src. + INLINE_FUNCTION_H + void assign(const std::vector& src, uint32 cap) + { + container_.reserve(cap); + this->assign(src); + } + + + /// Assign source vector. + /// The size of *this becomes the size of src. + /// The capacity of *this becomes the capacity of src. + INLINE_FUNCTION_H + void assign(const std::vector& src) + { + container_.assign(src.begin(), src.end()); + } + + /// Assign source vector from host side. + /// The size of *this becomes the size of src. + /// The capacity of *this becomes the capacity of src. + INLINE_FUNCTION_H + void assignFromHost(const VectorTypeHost& src) + { + notImplementedFunction; + fatalExit; + } + + INLINE_FUNCTION_H + void assign(const VectorType& src, bool srcCapacity = true) + { + notImplementedFunction; + fatalExit; + } + + template + INLINE_FUNCTION_H + void assignFromDevice(const VectorSingle& src, bool srcCapacity = true) + { + notImplementedFunction; + fatalExit; + } + + /*INLINE_FUNCTION_H + void append(const ViewType1D& appVec); + + INLINE_FUNCTION_H + void append(const std::vector& appVec); + + INLINE_FUNCTION_H + void append(const VectorType& appVec);*/ + + INLINE_FUNCTION_H + auto getSpan() + { + return span(container_.data(), container_.size()); + } + + INLINE_FUNCTION_H + auto getSpan()const + { + return span(const_cast(container_.data()), container_.size()); + } + + INLINE_FUNCTION_H + bool insertSetElement(const uint32IndexContainer& indices, const word& val) + { + notImplementedFunction; + return false; + } + + INLINE_FUNCTION_H + bool insertSetElement(const uint32IndexContainer& indices, const std::vector& vals) + { + notImplementedFunction; + return false; + } + + INLINE_FUNCTION_H + bool insertSetElement( + const uint32IndexContainer& indices, + const ViewType1D vals) + { + notImplementedFunction; + return false; + } + + INLINE_FUNCTION_H + bool reorderItems(const uint32IndexContainer& indices) + { + notImplementedFunction; + return false; + } + + /// @brief push a new element at the end (host call only) + /// resize if necessary and works on host accessible vector. + + void push_back(const word& val) + { + container_.push_back(val); + } + + INLINE_FUNCTION_H + pointer data(){ + return container_.data(); + } + + INLINE_FUNCTION_H + const_pointer data()const{ + return container_.data(); + } + + /// Return begin iterator. It works when devices is host accessible. + auto + begin(){ + return container_.begin(); + } + + /// Return begin iterator. it works when host is accessible. + const auto + begin()const { + return container_.begin(); + } + + auto end(){ + return container_.end(); + } + + /// Return end iterator. it works when host is accessible. + const auto end()const{ + return container_.end(); + } + + /// Return reference to element i. it works when host is accessible. + word& operator[](size_t i){ + return container_[i]; + } + + const word& operator[](size_t i)const{ + return container_[i]; + } + + //// - IO operations + + /// Read vector from stream + FUNCTION_H + bool read(iIstream& is) + { + return container_.read(is); + } + + /// Read vector from stream + FUNCTION_H + bool read(iIstream& is, const IOPattern& iop) + { + return container_.read(is, iop); + } + + /// Write the vector to os + FUNCTION_H + bool write(iOstream& os, const IOPattern& iop)const + { + return container_.write(os, iop); + } + + FUNCTION_H + bool write(iOstream& os)const + { + return container_.write(os); + } + + template + FUNCTION_H + bool write(iOstream& os, const IOPattern& iop, const HostMask& mask)const + { + notImplementedFunction; + return false; + } + + + /// Name of the memory space + static + constexpr const char* memoerySpaceName() + { + return memory_space::name(); + } + +}; // class wordVectorHost + +inline iIstream& operator >> (iIstream & is, VectorSingle & ivec ) +{ + if( !ivec.read(is) ) + { + ioErrorInFile (is.name(), is.lineNumber()); + fatalExit; + } + return is; +} + +inline iOstream& operator << (iOstream& os, const VectorSingle& ovec ) +{ + + if( !ovec.write(os) ) + { + ioErrorInFile(os.name(), os.lineNumber()); + fatalExit; + } + + return os; +} + + +} // - pFlow + + + +#endif + + +/* +INLINE_FUNCTION_H + bool append(const deviceViewType1D& dVec, size_t numElems) + { + + if(numElems == 0 )return true; + auto oldSize = size_; + auto newSize = size_ + numElems; + + if(this->empty() || newSize > capacity() ) + { + resize(newSize); + } + else + { + size_ = size_+numElems; + } + + auto dSubView = Kokkos::subview(view_, Kokkos::make_pair(oldSize, newSize)); + copy(dSubView, dVec); + + return true; + } + + INLINE_FUNCTION_H + bool append(const VectorSingle& Vec) + { + return append(Vec.deviceView(), Vec.size()); + }*/ \ No newline at end of file diff --git a/src/phasicFlow/containers/indexContainer/indexContainer.hpp b/src/phasicFlow/containers/indexContainer/indexContainer.hpp index 72118609..4486cd2c 100644 --- a/src/phasicFlow/containers/indexContainer/indexContainer.hpp +++ b/src/phasicFlow/containers/indexContainer/indexContainer.hpp @@ -143,6 +143,16 @@ public: indexContainer(ind.data(), ind.size()) {} + indexContainer(const DeviceViewType& ind): + size_(ind.size()), + views_("indexContainer", size_) + { + copy(views_.h_view, ind); + copy(views_.d_view, ind); + min_ = pFlow::min(views_.d_view, 0, size_); + max_ = pFlow::max(views_.d_view, 0, size_); + } + /// Copy indexContainer(const indexContainer&) = default; @@ -240,13 +250,13 @@ public: void syncViews() { bool findMinMax = false; - if(views_.template need_sync()) + if(views_.template need_sync()) { Kokkos::deep_copy(views_.d_view, views_.h_view); views_.clear_sync_state(); findMinMax = true; } - else if(views_.template need_sync()) + else if(views_.template need_sync()) { Kokkos::deep_copy(views_.h_view, views_.d_view); views_.clear_sync_state(); @@ -260,6 +270,12 @@ public: } } + void syncViews(uint32 newSize) + { + size_ = newSize; + syncViews(); + } + }; diff --git a/src/phasicFlow/containers/pointField/boundaryField/boundaryField/boundaryField.hpp b/src/phasicFlow/containers/pointField/boundaryField/boundaryField/boundaryField.hpp index d30b3f18..20267f33 100644 --- a/src/phasicFlow/containers/pointField/boundaryField/boundaryField/boundaryField.hpp +++ b/src/phasicFlow/containers/pointField/boundaryField/boundaryField/boundaryField.hpp @@ -121,21 +121,45 @@ public: return true; } - + inline + InternalFieldType& internal() + { + return internal_; + } + + inline + const InternalFieldType& internal()const + { + return internal_; + } + FieldAccessType thisField()const { - return FieldAccessType( - this->size(), - this->indexList().deviceViewAll(), - internal_.deviceViewAll()); + if constexpr(isDeviceAccessible()) + return FieldAccessType( + this->size(), + this->indexList().deviceViewAll(), + internal_.deviceViewAll()); + else + return FieldAccessType( + this->size(), + this->boundary().indexListHost().deviceViewAll(), + internal_.deviceViewAll()); } FieldAccessType mirrorField()const { - return FieldAccessType( - this->mirrorBoundary().size(), - this->mirrorBoundary().indexList().deviceViewAll(), - internal_.deviceViewAll()); + if constexpr(isDeviceAccessible()) + return FieldAccessType( + this->mirrorBoundary().size(), + this->mirrorBoundary().indexList().deviceViewAll(), + internal_.deviceViewAll()); + else + return FieldAccessType( + this->mirrorBoundary().size(), + this->mirrorBoundary().indexListHost().deviceViewAll(), + internal_.deviceViewAll()); + } virtual @@ -144,11 +168,6 @@ public: virtual const ProcVectorType& neighborProcField()const; - void fill(const std::any& val)override - { - return; - } - virtual void fill(const T& val) { diff --git a/src/phasicFlow/containers/pointField/boundaryField/boundaryFieldList.hpp b/src/phasicFlow/containers/pointField/boundaryField/boundaryFieldList.hpp index 31a13661..9ef6fcb1 100644 --- a/src/phasicFlow/containers/pointField/boundaryField/boundaryFieldList.hpp +++ b/src/phasicFlow/containers/pointField/boundaryField/boundaryFieldList.hpp @@ -43,7 +43,7 @@ protected: const boundaryList& boundaries_; - uint32 slaveToMasterUpdateIter_ = -1; + uint32 slaveToMasterUpdateIter_ = static_cast(-1); public: @@ -95,7 +95,7 @@ public: bool slaveToMasterUpdateRequested()const { - return slaveToMasterUpdateIter_ != -1; + return slaveToMasterUpdateIter_ != static_cast(-1); } diff --git a/src/phasicFlow/containers/pointField/boundaryField/generalBoundary/generalBoundary.hpp b/src/phasicFlow/containers/pointField/boundaryField/generalBoundary/generalBoundary.hpp index c55a7f71..6a1f6808 100644 --- a/src/phasicFlow/containers/pointField/boundaryField/generalBoundary/generalBoundary.hpp +++ b/src/phasicFlow/containers/pointField/boundaryField/generalBoundary/generalBoundary.hpp @@ -147,10 +147,6 @@ public: } const Time& time()const; - - virtual - void fill(const std::any& val)=0; - }; diff --git a/src/phasicFlow/containers/pointField/internalField/internalField.cpp b/src/phasicFlow/containers/pointField/internalField/internalField.cpp index c59faf3b..4b7c01e9 100644 --- a/src/phasicFlow/containers/pointField/internalField/internalField.cpp +++ b/src/phasicFlow/containers/pointField/internalField/internalField.cpp @@ -52,9 +52,29 @@ bool pFlow::internalField::insert(const anyList& varList) } return true; - } +template +bool pFlow::internalField::rearrange(const anyList& varList) +{ + const word eventName = message::eventName(message::ITEM_REARRANGE); + + const auto& indices = varList.getObject( + eventName); + + field_.reserve( internalPoints_.capacity()); + field_.resize(internalPoints_.size()); + if(!field_.reorderItems(indices)) + { + fatalErrorInFunction<< + "cannot reorder items in field "<< name()< pFlow::internalField::internalField ( @@ -166,15 +186,15 @@ bool pFlow::internalField:: hearChanges { // do nothing } - if(msg.equivalentTo(message::ITEM_REARRANGE)) - { - notImplementedFunction; - return false; - } if(msg.equivalentTo(message::ITEM_INSERT)) { return insert(varList); } + if(msg.equivalentTo(message::ITEM_REARRANGE)) + { + return rearrange(varList); + + } return true; } diff --git a/src/phasicFlow/containers/pointField/internalField/internalField.hpp b/src/phasicFlow/containers/pointField/internalField/internalField.hpp index 423f1ef7..ab26c663 100644 --- a/src/phasicFlow/containers/pointField/internalField/internalField.hpp +++ b/src/phasicFlow/containers/pointField/internalField/internalField.hpp @@ -73,6 +73,8 @@ protected: bool insert(const anyList& varList); + bool rearrange(const anyList& varList); + public: internalField( diff --git a/src/phasicFlow/dictionary/dictionary.hpp b/src/phasicFlow/dictionary/dictionary.hpp index 9baba109..1684987c 100644 --- a/src/phasicFlow/dictionary/dictionary.hpp +++ b/src/phasicFlow/dictionary/dictionary.hpp @@ -302,10 +302,10 @@ public: //// - IO operations /// read from stream - virtual bool read(iIstream& is); + bool read(iIstream& is) override; /// write to stream - virtual bool write(iOstream& os) const; + bool write(iOstream& os) const override; }; diff --git a/src/phasicFlow/dictionary/fileDictionary.hpp b/src/phasicFlow/dictionary/fileDictionary.hpp index b8685dfd..a74c989e 100644 --- a/src/phasicFlow/dictionary/fileDictionary.hpp +++ b/src/phasicFlow/dictionary/fileDictionary.hpp @@ -46,6 +46,9 @@ public: const dictionary& dict, repository* owner=nullptr); + using dictionary::read; + using dictionary::write; + /// read from stream bool read(iIstream& is, const IOPattern& iop) override; diff --git a/src/phasicFlow/eventManagement/message.hpp b/src/phasicFlow/eventManagement/message.hpp index 3ece7411..0593221a 100644 --- a/src/phasicFlow/eventManagement/message.hpp +++ b/src/phasicFlow/eventManagement/message.hpp @@ -48,14 +48,16 @@ public: BNDR_RESET = 10, // boundary indices reset entirely BNDR_DELETE = 11, // boundary indices deleted BNDR_APPEND = 12, // - BNDR_PROCTRANS1 = 13, // transfer of data between processors step 1 - BNDR_PROCTRANS2 = 14 // transfer of data between processors step 2 + BNDR_PROCTRANSFER_SEND = 13, // transfer of data between processors step 1 + BNDR_PROCTRANSFER_RECIEVE = 14, // transfer of data between processors step 2 + BNDR_PROCTRANSFER_WAITFILL = 15, // wait for data and fill the field + BNDR_PROC_SIZE_CHANGED = 16 }; private: - static constexpr size_t numberOfEvents_ = 15; + static constexpr size_t numberOfEvents_ = 17; std::bitset events_{0x0000}; @@ -76,7 +78,9 @@ private: "deletedIndices", "appendedIndices", "transferredIndices", - "numberRecieved" + "numToRecieve", + "insertedIndices", + "size" }; public: diff --git a/src/phasicFlow/globals/pFlowMacros.hpp b/src/phasicFlow/globals/pFlowMacros.hpp index 544fc150..8973b6ba 100755 --- a/src/phasicFlow/globals/pFlowMacros.hpp +++ b/src/phasicFlow/globals/pFlowMacros.hpp @@ -37,6 +37,10 @@ Licence: #define CONSUME_PARAM(x) (void)(x); +#if defined(pFlow_Build_Cuda) && !defined(__CUDACC__) +#define __CUDACC__ +#endif + #ifdef __CUDACC__ #define INLINE_FUNCTION_HD inline __host__ __device__ #define INLINE_FUNCTION_D inline __device__ diff --git a/src/phasicFlow/processors/localProcessors.hpp b/src/phasicFlow/processors/localProcessors.hpp index eab18c08..1e1d6d3e 100644 --- a/src/phasicFlow/processors/localProcessors.hpp +++ b/src/phasicFlow/processors/localProcessors.hpp @@ -145,6 +145,15 @@ public: } #endif + static + bool builtForMPI() + { + #ifdef pFlow_Build_MPI + return true; + #else + return false; + #endif + } }; } diff --git a/src/phasicFlow/processors/processors.cpp b/src/phasicFlow/processors/processors.cpp index f5eb066a..078c9a5d 100644 --- a/src/phasicFlow/processors/processors.cpp +++ b/src/phasicFlow/processors/processors.cpp @@ -29,7 +29,6 @@ Licence: #include "processors.hpp" #include "streams.hpp" -static int numVarsInitialized__ = 0; #ifdef pFlow_Build_MPI diff --git a/src/phasicFlow/repository/IOobject/IOfileHeader.cpp b/src/phasicFlow/repository/IOobject/IOfileHeader.cpp index 96a49f67..f67eec28 100644 --- a/src/phasicFlow/repository/IOobject/IOfileHeader.cpp +++ b/src/phasicFlow/repository/IOobject/IOfileHeader.cpp @@ -42,12 +42,22 @@ pFlow::uniquePtr pFlow::IOfileHeader::outStream()const return osPtr; } -pFlow::IOfileHeader::IOfileHeader -( - const objectFile& objf -) -: - objectFile(objf) +pFlow::uniquePtr pFlow::IOfileHeader::dummyOutStream() const +{ + auto osPtr = makeUnique( CWD()+word("dummyFile") , outFileBinary()); + + if(osPtr && owner()) + { + auto outPrecision = owner()->outFilePrecision(); + osPtr->precision(static_cast(outPrecision)); + } + + return osPtr; +} + +pFlow::IOfileHeader::IOfileHeader( + const objectFile &objf) + : objectFile(objf) {} pFlow::fileSystem pFlow::IOfileHeader::path() const diff --git a/src/phasicFlow/repository/IOobject/IOfileHeader.hpp b/src/phasicFlow/repository/IOobject/IOfileHeader.hpp index 1d2a1e4c..d4549c01 100644 --- a/src/phasicFlow/repository/IOobject/IOfileHeader.hpp +++ b/src/phasicFlow/repository/IOobject/IOfileHeader.hpp @@ -58,6 +58,8 @@ protected: // - ouput file stream uniquePtr outStream()const; + uniquePtr dummyOutStream()const; + public: // with owner diff --git a/src/phasicFlow/repository/IOobject/IOobject.cpp b/src/phasicFlow/repository/IOobject/IOobject.cpp index b5857ea1..6fae2f0a 100644 --- a/src/phasicFlow/repository/IOobject/IOobject.cpp +++ b/src/phasicFlow/repository/IOobject/IOobject.cpp @@ -122,18 +122,35 @@ bool pFlow::IOobject::writeObject() const if(ioPattern().thisCallWrite()) { - - if(auto ptrOS = outStream(); ptrOS ) + if( ioPattern().thisProcWriteData()) { - return writeObject(ptrOS()); + if(auto ptrOS = outStream(); ptrOS ) + { + return writeObject(ptrOS()); + } + else + { + warningInFunction<< + "error in opening file "<< path() <= endTime_ ) return true; - if( abs(currentTime_-endTime_) < 0.5*dt_ )return true; + if( std::abs(currentTime_-endTime_) < 0.5*dt_ )return true; return false; } bool pFlow::timeControl::reachedStopAt()const { if( currentTime_ >= stopAt_ ) return true; - if( abs(currentTime_-stopAt_) < 0.5*dt_ )return true; + if( std::abs(currentTime_-stopAt_) < 0.5*dt_ )return true; return false; } @@ -154,7 +154,7 @@ void pFlow::timeControl::checkForOutputToFile() bool save = false; if(managedExternaly_) { - if( abs(currentTime_-writeTime_) < 0.5*dt_) + if( std::abs(currentTime_-writeTime_) < 0.5*dt_) { save = true; lastSaved_ = currentTime_; @@ -162,12 +162,12 @@ void pFlow::timeControl::checkForOutputToFile() } else { - if ( abs(currentTime_ - lastSaved_ - saveInterval_) < 0.5 * dt_ ) + if ( std::abs(currentTime_ - lastSaved_ - saveInterval_) < 0.5 * dt_ ) { lastSaved_ = currentTime_; save = true; } - else if( abs(currentTime_ - lastSaved_) < min( pow(10.0,-1.0*timePrecision_), 0.5 *dt_) ) + else if( std::abs(currentTime_ - lastSaved_) < std::min( pow(10.0,-1.0*timePrecision_), 0.5 *dt_) ) { lastSaved_ = currentTime_; save = true; diff --git a/src/phasicFlow/smartPointers/uniquePtr.hpp b/src/phasicFlow/smartPointers/uniquePtr.hpp index c1363cde..7bb6a6b9 100644 --- a/src/phasicFlow/smartPointers/uniquePtr.hpp +++ b/src/phasicFlow/smartPointers/uniquePtr.hpp @@ -37,16 +37,15 @@ namespace pFlow template< - typename T, - typename Deleter = std::default_delete + typename T > class uniquePtr : - public std::unique_ptr + public std::unique_ptr { public: - using uniquePtrType = std::unique_ptr; + using uniquePtrType = std::unique_ptr; // using base constructors using uniquePtrType::unique_ptr; diff --git a/src/phasicFlow/streams/TStream/iTstream.cpp b/src/phasicFlow/streams/TStream/iTstream.cpp index 49b80cee..ac8e503e 100755 --- a/src/phasicFlow/streams/TStream/iTstream.cpp +++ b/src/phasicFlow/streams/TStream/iTstream.cpp @@ -325,7 +325,7 @@ void pFlow::iTstream::reset() size_t pFlow::iTstream::tell() { notImplementedFunction; - return -1; + return static_cast(-1); } const pFlow::tokenList& pFlow::iTstream::tokens()const diff --git a/src/phasicFlow/streams/dataIO/dataIO.cpp b/src/phasicFlow/streams/dataIO/dataIO.cpp index 39fbd9e0..851ec4a5 100644 --- a/src/phasicFlow/streams/dataIO/dataIO.cpp +++ b/src/phasicFlow/streams/dataIO/dataIO.cpp @@ -22,7 +22,7 @@ template bool pFlow::writeDataAsciiBinary(iOstream& os, span data) { - if( os.isBinary() ) + if( std::is_trivially_copyable_v && os.isBinary() ) { // first write the number of data uint64 numData = data.size(); @@ -83,7 +83,7 @@ bool pFlow::readDataAsciiBinary ) { - if(is.isBinary()) + if(std::is_trivially_copyable_v && is.isBinary()) { data.clear(); // read length of data diff --git a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp index f8307478..906d5fd9 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp @@ -42,6 +42,7 @@ void pFlow::boundaryBase::setNewIndices unSyncLists(); } + bool pFlow::boundaryBase::appendNewIndices ( const uint32Vector_D& newIndices @@ -144,7 +145,7 @@ bool pFlow::boundaryBase::setRemoveKeepIndices return true; } -bool pFlow::boundaryBase::transferPoints +bool pFlow::boundaryBase::transferPointsToMirror ( uint32 numTransfer, const uint32Vector_D& transferMask, @@ -203,6 +204,8 @@ pFlow::boundaryBase::boundaryBase indexList_(groupNames("indexList", dict.name())), indexListHost_(groupNames("hostIndexList",dict.name())), neighborLength_(dict.getVal("neighborLength")), + updateInetrval_(dict.getVal("updateInterval")), + boundaryExtntionLengthRatio_(dict.getVal("boundaryExtntionLengthRatio")), internal_(internal), boundaries_(bndrs), thisBoundaryIndex_(thisIndex), diff --git a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp index 9b63c9df..978f674b 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp @@ -68,6 +68,12 @@ private: /// The length defined for creating neighbor list real neighborLength_; + /// time steps between successive update of boundary lists + uint32 updateInetrval_; + + /// the extra boundary extension beyound actual limits of boundary + real boundaryExtntionLengthRatio_; + /// a reference to internal points internalPoints& internal_; @@ -106,7 +112,7 @@ protected: const uint32Vector_D& removeIndices, const uint32Vector_D& keepIndices); - bool transferPoints( + bool transferPointsToMirror( uint32 numTransfer, const uint32Vector_D& transferMask, uint32 transferBoundaryIndex, @@ -126,6 +132,12 @@ protected: } } + /// Is this iter the right time for updating bounday list + bool boundaryListUpdate(uint32 iter)const + { + return iter%updateInetrval_ == 0u || iter == 0u; + } + /// Update this boundary data in two steps (1 and 2). /// This is called after calling beforeIteration for /// all boundaries, so any particle addition, deletion, @@ -133,7 +145,7 @@ protected: /// This two-step update help to have a flexible mechanism /// for data transfer, mostly for MPI related jobs. virtual - bool updataBoundary(int step) + bool updataBoundaryData(int step) { return true; } @@ -145,11 +157,13 @@ protected: /// to complete. false: if the operation is complete and no need for /// additional step in operation. virtual - bool transferData(int step) + bool transferData(uint32 iter, int step) { return false; } + + public: TypeInfo("boundaryBase"); @@ -189,20 +203,21 @@ public: /// The length from boundary plane into the domain /// where beyond that distance internal points exist. /// By conventions is it always equal to neighborLength_ + inline real neighborLengthIntoInternal()const { return neighborLength_; } /// The distance length from boundary plane - /// where neighbor particles exist in that distance. + /// where neighbor particles still exist in that distance. /// This length may be modified in each boundary type /// as required. In this case the boundaryExtensionLength /// method should also be modified accordingly. virtual real neighborLength()const { - return neighborLength_; + return (1+boundaryExtntionLengthRatio_)*neighborLength_; } /// The extention length (in vector form) for the boundary @@ -213,7 +228,7 @@ public: virtual realx3 boundaryExtensionLength()const { - return {0,0,0}; + return -boundaryExtntionLengthRatio_*neighborLength_ * boundaryPlane_.normal(); } inline @@ -351,6 +366,18 @@ public: virtual const realx3Vector_D& neighborProcPoints()const; + virtual + uint32 numToTransfer()const + { + return 0u; + } + + virtual + uint32 numToRecieve()const + { + return 0u; + } + /// - static create static uniquePtr create diff --git a/src/phasicFlow/structuredData/boundaries/boundaryBase/scatteredFieldAccess.hpp b/src/phasicFlow/structuredData/boundaries/boundaryBase/scatteredFieldAccess.hpp index f26438db..c1967bc8 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryBase/scatteredFieldAccess.hpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryBase/scatteredFieldAccess.hpp @@ -114,6 +114,17 @@ public: return fieldVals_; } + auto& indices() + { + return indices_; + } + + const auto& indices()const + { + return indices_; + } + + INLINE_FUNCTION_HD uint32 index(uint32 i)const { return indices_[i]; diff --git a/src/phasicFlow/structuredData/boundaries/boundaryExit/boundaryExit.cpp b/src/phasicFlow/structuredData/boundaries/boundaryExit/boundaryExit.cpp index 3b753877..9937433c 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryExit/boundaryExit.cpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryExit/boundaryExit.cpp @@ -47,6 +47,9 @@ bool pFlow::boundaryExit::beforeIteration real dt ) { + + if( !boundaryListUpdate(iterNum))return true; + // nothing have to be done if(empty()) { diff --git a/src/phasicFlow/structuredData/boundaries/boundaryList.cpp b/src/phasicFlow/structuredData/boundaries/boundaryList.cpp index 6cccb194..473e523d 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryList.cpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryList.cpp @@ -91,18 +91,18 @@ pFlow::boundaryList::updateNeighborLists() pFlow::boundaryList::boundaryList(pointStructure& pStruct) : ListPtr(pStruct.simDomain().sizeOfBoundaries()), pStruct_(pStruct), - timeControl_( - pStruct.simDomain().subDict("boundaries"), - "update", - pStruct.currentTime() - ) + neighborListUpdateInterval_( + pStruct.simDomain().subDict("boundaries").getVal( + "neighborListUpdateInterval" + ) + ) { } bool -pFlow::boundaryList::updateNeighborLists(uint32 iter, real t, real dt) +pFlow::boundaryList::updateNeighborLists(uint32 iter, bool force) { - if (timeControl_.timeEvent(iter, t, dt)) + if(iter%neighborListUpdateInterval_==0u || iter == 0u || force) { return updateNeighborLists(); } @@ -155,15 +155,18 @@ pFlow::boundaryList::internalDomainBox() const } bool -pFlow::boundaryList::beforeIteration(uint32 iter, real t, real dt) +pFlow::boundaryList::beforeIteration(uint32 iter, real t, real dt, bool force) { // it is time to update lists - if (timeControl_.timeEvent(iter, t, dt) && !updateNeighborLists()) + if(!updateNeighborLists(iter , force) ) { fatalErrorInFunction; return false; } + // this forces performing updating the list on each boundary + if(force) iter = 0; + for (auto bdry : *this) { if (!bdry->beforeIteration(iter, t, dt)) @@ -176,12 +179,12 @@ pFlow::boundaryList::beforeIteration(uint32 iter, real t, real dt) for (auto bdry : *this) { - bdry->updataBoundary(1); + bdry->updataBoundaryData(1); } for (auto bdry : *this) { - bdry->updataBoundary(2); + bdry->updataBoundaryData(2); } return true; @@ -215,20 +218,11 @@ pFlow::boundaryList::afterIteration(uint32 iter, real t, real dt) { if(requireStep[i]) { - requireStep[i] = this->operator[](i).transferData(step); + requireStep[i] = this->operator[](i).transferData(iter, step); } } step++; } - /*for (auto& bdry : *this) - { - if (!bdry->afterIteration(iter, t, dt)) - { - fatalErrorInFunction << "Error in afterIteration in boundary " - << bdry->name() << endl; - return false; - } - }*/ return true; } diff --git a/src/phasicFlow/structuredData/boundaries/boundaryList.hpp b/src/phasicFlow/structuredData/boundaries/boundaryList.hpp index 5eba6399..8d4c8710 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryList.hpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryList.hpp @@ -17,15 +17,12 @@ Licence: implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. -----------------------------------------------------------------------------*/ - - #ifndef __boundaryList_hpp__ #define __boundaryList_hpp__ #include "domain.hpp" #include "boundaryBase.hpp" #include "ListPtr.hpp" -#include "baseTimeControl.hpp" namespace pFlow @@ -42,7 +39,7 @@ private: //// - data members pointStructure& pStruct_; - baseTimeControl timeControl_; + uint32 neighborListUpdateInterval_; domain extendedDomain_; @@ -72,7 +69,7 @@ public: /// @brief update neighbor list of boundaries based on /// the time intervals - bool updateNeighborLists(uint32 iter, real t, real dt); + bool updateNeighborLists(uint32 iter, bool force = false); bool createBoundaries(); @@ -94,11 +91,6 @@ public: return ListPtr::operator[](i); } - inline - const baseTimeControl& timeControl()const - { - return timeControl_; - } inline const auto& extendedDomain()const @@ -114,7 +106,7 @@ public: box internalDomainBox()const; - bool beforeIteration(uint32 iter, real t, real dt); + bool beforeIteration(uint32 iter, real t, real dt, bool force = false); bool iterate(uint32 iter, real t, real dt); diff --git a/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.hpp b/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.hpp index 372a1c28..8e103d39 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.hpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.hpp @@ -43,7 +43,7 @@ public: boundaryList &bndrs, uint32 thisIndex); - ~boundaryNone() override= default; + ~boundaryNone() final= default; add_vCtor ( @@ -52,11 +52,11 @@ public: dictionary ); - bool beforeIteration(uint32 iterNum, real t, real dt) override; + bool beforeIteration(uint32 iterNum, real t, real dt) final; - bool iterate(uint32 iterNum, real t, real dt) override; + bool iterate(uint32 iterNum, real t, real dt) final; - bool afterIteration(uint32 iterNum, real t, real dt) override; + bool afterIteration(uint32 iterNum, real t, real dt) final; }; diff --git a/src/phasicFlow/structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.cpp b/src/phasicFlow/structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.cpp index 30fdc039..cbdecf71 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.cpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.cpp @@ -33,17 +33,20 @@ pFlow::boundaryPeriodic::boundaryPeriodic ) : boundaryBase(dict, bplane, internal, bndrs, thisIndex), - mirrorBoundaryIndex_(dict.getVal("mirrorBoundaryIndex")) -{} + mirrorBoundaryIndex_(dict.getVal("mirrorBoundaryIndex")), + extensionLength_(dict.getVal("boundaryExtntionLengthRatio")) +{ + extensionLength_ = max(extensionLength_, static_cast(0.1)); +} pFlow::real pFlow::boundaryPeriodic::neighborLength() const { - return (1+extensionLength_)*boundaryBase::neighborLength(); + return (1+extensionLength_)*neighborLengthIntoInternal(); } pFlow::realx3 pFlow::boundaryPeriodic::boundaryExtensionLength() const { - return -extensionLength_*neighborLength()*boundaryBase::boundaryPlane().normal(); + return -extensionLength_*neighborLengthIntoInternal()*boundaryBase::boundaryPlane().normal(); } @@ -57,6 +60,8 @@ bool pFlow::boundaryPeriodic::beforeIteration( { return true; } + + if( !boundaryListUpdate(iterNum))return true; uint32 s = size(); uint32Vector_D transferFlags("transferFlags",s+1, s+1, RESERVE()); @@ -92,7 +97,7 @@ bool pFlow::boundaryPeriodic::beforeIteration( // to obtain the transfer vector realx3 transferVec = displacementVectroToMirror(); - return transferPoints + return transferPointsToMirror ( numTransfered, transferFlags, diff --git a/src/phasicFlow/structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.hpp b/src/phasicFlow/structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.hpp index b4746a5e..88983fa3 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.hpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.hpp @@ -35,8 +35,7 @@ private: uint32 mirrorBoundaryIndex_; - static - inline const real extensionLength_ = 0.1; + real extensionLength_ = 0.1; public: diff --git a/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.cpp b/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.cpp index 6d818a23..74ee4d75 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.cpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.cpp @@ -125,6 +125,8 @@ bool pFlow::boundaryReflective::afterIteration const auto& velocity = time().lookupObject(velocityName_); const auto& velocityD = velocity.deviceViewAll(); + + const auto restitution = restitution_; Kokkos::parallel_for( "pFlow::boundaryReflective::velocityChange", @@ -137,7 +139,7 @@ bool pFlow::boundaryReflective::afterIteration if(vn < 0) { realx3 vt = vel - vn*p.normal(); - vel = restitution_*(vt - vn*p.normal()); + vel = restitution*(vt - vn*p.normal()); } } ); diff --git a/src/phasicFlow/structuredData/cells/cells.hpp b/src/phasicFlow/structuredData/cells/cells.hpp index c77b1a69..c4b50813 100644 --- a/src/phasicFlow/structuredData/cells/cells.hpp +++ b/src/phasicFlow/structuredData/cells/cells.hpp @@ -28,58 +28,40 @@ Licence: namespace pFlow { -template class cells { -public: - - using CellType = triple; - -protected: +private: // - domain - box domain_{realx3(0.0), realx3(1.0)}; + box domainBox_{realx3(0.0), realx3(1.0)}; // - cell size - realx3 cellSize_{1,1,1}; - - CellType numCells_{1,1,1}; + real celldx_{1}; + int32x3 numCells_{1,1,1}; // - protected methods INLINE_FUNCTION_H void calculate() { - numCells_ = (domain_.maxPoint()-domain_.minPoint())/cellSize_ + realx3(1.0); - numCells_ = max( numCells_ , CellType(static_cast(1)) ); + numCells_ = (domainBox_.maxPoint()-domainBox_.minPoint())/celldx_ + realx3(1.0); + numCells_ = max( numCells_ , int32x3(1) ); } public: INLINE_FUNCTION_HD - cells() - {} + cells() = default; INLINE_FUNCTION_H cells(const box& domain, real cellSize) : - domain_(domain), - cellSize_(cellSize) + domainBox_(domain), + celldx_(cellSize) { calculate(); } - - INLINE_FUNCTION_H - cells(const box& domain, int32 nx, int32 ny, int32 nz) - : - domain_(domain), - cellSize_( - (domain_.maxPoint() - domain_.minPoint())/realx3(nx, ny, nz) - ), - numCells_(nx, ny, nz) - {} - INLINE_FUNCTION_HD cells(const cells&) = default; @@ -100,43 +82,36 @@ public: INLINE_FUNCTION_H void setCellSize(real cellSize) { - cellSize_ = cellSize; - calculate(); - } - - INLINE_FUNCTION_H - void setCellSize(realx3 cellSize) - { - cellSize_ = cellSize; + celldx_ = cellSize; calculate(); } INLINE_FUNCTION_HD - realx3 cellSize()const + real cellSize()const { - return cellSize_; + return celldx_; } INLINE_FUNCTION_HD - const CellType& numCells()const + const int32x3& numCells()const { return numCells_; } INLINE_FUNCTION_HD - indexType nx()const + int32 nx()const { return numCells_.x(); } INLINE_FUNCTION_HD - indexType ny()const + int32 ny()const { return numCells_.y(); } INLINE_FUNCTION_HD - indexType nz()const + int32 nz()const { return numCells_.z(); } @@ -149,22 +124,21 @@ public: static_cast(numCells_.z()); } - const auto& domain()const + const auto& domainBox()const { - return domain_; + return domainBox_; } INLINE_FUNCTION_HD - CellType pointIndex(const realx3& p)const + int32x3 pointIndex(const realx3& p)const { - return CellType( (p - domain_.minPoint())/cellSize_ ); + return int32x3( (p - domainBox_.minPoint())/celldx_ ); } INLINE_FUNCTION_HD - bool pointIndexInDomain(const realx3 p, CellType& index)const + bool pointIndexInDomain(const realx3 p, int32x3& index)const { - if( !domain_.isInside(p) ) return false; - + if(!inDomain(p))return false; index = this->pointIndex(p); return true; } @@ -172,11 +146,11 @@ public: INLINE_FUNCTION_HD bool inDomain(const realx3& p)const { - return domain_.isInside(p); + return domainBox_.isInside(p); } INLINE_FUNCTION_HD - bool isInRange(const CellType& cell)const + bool inCellRange(const int32x3& cell)const { if(cell.x()<0)return false; if(cell.y()<0)return false; @@ -188,7 +162,7 @@ public: } INLINE_FUNCTION_HD - bool isInRange(indexType i, indexType j, indexType k)const + bool inCellRange(int32 i, int32 j, int32 k)const { if(i<0)return false; if(j<0)return false; @@ -199,22 +173,7 @@ public: return true; } - INLINE_FUNCTION_HD - void extendBox( - const CellType& p1, - const CellType& p2, - const CellType& p3, - indexType extent, - CellType& minP, - CellType& maxP)const - { - minP = min( min( p1, p2), p3)-extent; - maxP = max( max( p1, p2), p3)+extent; - - minP = bound(minP); - maxP = bound(maxP); - } - + INLINE_FUNCTION_HD void extendBox( const realx3& p1, @@ -224,17 +183,17 @@ public: realx3& minP, realx3& maxP)const { - minP = min(min(p1,p2),p3) - extent*cellSize_ ; - maxP = max(max(p1,p2),p3) + extent*cellSize_ ; + minP = min(min(p1,p2),p3) - extent*celldx_ ; + maxP = max(max(p1,p2),p3) + extent*celldx_ ; minP = bound(minP); maxP = bound(maxP); } INLINE_FUNCTION_HD - CellType bound(CellType p)const + int32x3 bound(int32x3 p)const { - return CellType( + return int32x3( min( numCells_.x()-1, max(0,p.x())), min( numCells_.y()-1, max(0,p.y())), min( numCells_.z()-1, max(0,p.z())) @@ -245,9 +204,9 @@ public: realx3 bound(realx3 p)const { return realx3( - min( domain_.maxPoint().x(), max(domain_.minPoint().x(),p.x())), - min( domain_.maxPoint().y(), max(domain_.minPoint().y(),p.y())), - min( domain_.maxPoint().z(), max(domain_.minPoint().z(),p.z())) + min( domainBox_.maxPoint().x(), max(domainBox_.minPoint().x(),p.x())), + min( domainBox_.maxPoint().y(), max(domainBox_.minPoint().y(),p.y())), + min( domainBox_.maxPoint().z(), max(domainBox_.minPoint().z(),p.z())) ); } }; diff --git a/src/phasicFlow/structuredData/cylinder/cylinder.cpp b/src/phasicFlow/structuredData/cylinder/cylinder.cpp index 8c706d3e..014462c2 100644 --- a/src/phasicFlow/structuredData/cylinder/cylinder.cpp +++ b/src/phasicFlow/structuredData/cylinder/cylinder.cpp @@ -18,7 +18,7 @@ Licence: -----------------------------------------------------------------------------*/ - +#include "math.hpp" #include "cylinder.hpp" #include "zAxis.hpp" #include "streams.hpp" diff --git a/src/phasicFlow/structuredData/domain/regularSimulationDomain.cpp b/src/phasicFlow/structuredData/domain/regularSimulationDomain.cpp index e94d111e..a84325a9 100644 --- a/src/phasicFlow/structuredData/domain/regularSimulationDomain.cpp +++ b/src/phasicFlow/structuredData/domain/regularSimulationDomain.cpp @@ -31,7 +31,10 @@ bool pFlow::regularSimulationDomain::createBoundaryDicts() this->addDict("regularBoundaries", boundaries); auto& rbBoundaries = this->subDict("regularBoundaries"); - real neighborLength = boundaries.getVal("neighborLength"); + auto neighborLength = boundaries.getVal("neighborLength"); + auto boundaryExtntionLengthRatio = + boundaries.getValOrSet("boundaryExtntionLengthRatio", 0.1); + auto updateIntercal = boundaries.getValOrSet("updateInterval", 1u); for(uint32 i=0; i emptySpots ) @@ -337,8 +351,20 @@ pFlow::internalPoints::insertPoints( // we should fill the scattered empty spots else { - notImplementedFunction; - return false; + auto newIndices = pFlagsD_.getEmptyPoints(numNew); + if(numNew != newIndices.size()) + { + fatalErrorInFunction<<"not enough empty points in pointFlag"<< + numNew<< " "<( + msg.addAndName(message::ITEM_INSERT), + newIndices + ); } const auto& indices = varList.getObject( @@ -398,6 +424,108 @@ pFlow::internalPoints::insertPoints( return true; } +bool pFlow::internalPoints::insertPointsOnly( + const realx3Vector_D &points, + message& msg, + anyList &varList +) +{ + + uint32 numNew = static_cast(points.size()); + + auto aRange = pFlagsD_.activeRange(); + uint32 emptySpots = pFlagsD_.capacity() - pFlagsD_.numActive(); + + if(emptySpots!= 0) emptySpots--; + + if( numNew > emptySpots ) + { + // increase the capacity to hold new points + aRange = pFlagsD_.activeRange(); + uint32 newCap = pFlagsD_.changeCapacity(numNew); + unSyncFlag(); + varList.emplaceBack( + msg.addAndName(message::CAP_CHANGED), + newCap); + } + + + // first check if it is possible to add to the beggining of list + if(numNew <= aRange.start()) + { + varList.emplaceBack( + msg.addAndName(message::ITEM_INSERT), + 0u, numNew); + }// check if it is possible to add to the end of the list + else if( numNew <= pFlagsD_.capacity() - aRange.end() ) + { + varList.emplaceBack( + msg.addAndName(message::ITEM_INSERT), + aRange.end(), aRange.end()+numNew); + }// we should fill the scattered empty spots + else + { + + auto newIndices = pFlagsD_.getEmptyPoints(numNew); + if(numNew != newIndices.size()) + { + fatalErrorInFunction<<"not enough empty points in pointFlag"<< + numNew<< " "<( + msg.addAndName(message::ITEM_INSERT), + newIndices + ); + } + + const auto& indices = varList.getObject( + message::eventName(message::ITEM_INSERT) + ); + + auto nAdded = pFlagsD_.addInternalPoints(indices.deviceView()); + unSyncFlag(); + + if(nAdded != numNew ) + { + fatalErrorInFunction; + return false; + } + + pointPosition_.reserve( pFlagsD_.capacity() ); + + if(!pointPosition_.insertSetElement(indices, points.deviceView())) + { + fatalErrorInFunction<< + "Error in inserting new positions into pointPosition"<< + " internalPoints field"<(Flag::DELETED)); + } + else + { + fill(flags_, end, cap, static_cast(Flag::DELETED)); + } + + activeRange_ = {start, end}; + numActive_ = activeRange_.numElements(); + isAllActive_ = true; + fill(flags_, activeRange_, static_cast(Flag::INTERNAL)); + + nLeft_ = nRight_ = 0; + nBottom_ = nTop_ = 0; + nRear_ = nFront_ = 0; + } + public: - friend class internalPoints; pointFlag() = default; @@ -283,6 +305,9 @@ public: } ViewType1D getActivePoints(); + + + ViewType1D getEmptyPoints(uint32 numToGet)const; /// @brief Loop over the active points and mark those out of the box diff --git a/src/phasicFlow/structuredData/pointStructure/internalPoints/pointFlagKernels.hpp b/src/phasicFlow/structuredData/pointStructure/internalPoints/pointFlagKernels.hpp index 91764cf6..e76533b9 100644 --- a/src/phasicFlow/structuredData/pointStructure/internalPoints/pointFlagKernels.hpp +++ b/src/phasicFlow/structuredData/pointStructure/internalPoints/pointFlagKernels.hpp @@ -50,6 +50,52 @@ pFlow::ViewType1D::memo return aPoints; } +template +pFlow::ViewType1D::memory_space> + pFlow::pointFlag::getEmptyPoints(uint32 numToGet)const +{ + + uint32 cap = capacity(); + using rpAPoints = Kokkos::RangePolicy>; + + ViewType1D indices("indices", cap+1); + + ViewType1D emptyPoints("emptyPoints", numToGet); + + auto flags = flags_; + + Kokkos::parallel_for( + "getEmptyPoints", + rPolicy(0, cap), + LAMBDA_HD(uint32 i){ + indices(i) = flags[i] == DELETED; + }); + Kokkos::fence(); + + exclusiveScan(indices, 0u, cap, indices, 0u); + uint32 numFound = 0; + Kokkos::parallel_reduce( + "fillEmptyPoints", + rpAPoints(0, cap-1), + LAMBDA_HD(uint32 i, uint32& numFoundUpdate){ + if(indices(i)!= indices(i+1) && indices(i)< numToGet) + { + emptyPoints(indices(i)) = i; + numFoundUpdate++; + } + }, + numFound + ); + + if(numFound < numToGet) + { + return Kokkos::subview(emptyPoints, Kokkos::make_pair(0u, numFound)); + } + return emptyPoints; + +} + template pFlow::uint32 pFlow::pointFlag::markOutOfBoxDelete @@ -106,7 +152,7 @@ pFlow::uint32 pFlow::pointFlag::markOutOfBoxDelete { minRange = 0; maxRange = 0; - numDeleted == numActive_; + numDeleted = numActive_; } else { @@ -135,7 +181,8 @@ pFlow::pointFlag::addInternalPoints uint32 maxRange; uint32 numAdded = 0; - + const auto& flag = flags_; + Kokkos::parallel_reduce( "pointFlagKernels::addInternalPoints", deviceRPolicyStatic(0, points.extent(0)), @@ -146,10 +193,10 @@ pFlow::pointFlag::addInternalPoints uint32& addToUpdate) { uint32 idx = points(i); - if( flags_[idx] <= DELETED) addToUpdate ++; + if( flag[idx] <= DELETED) addToUpdate ++; minUpdate = min(minUpdate,idx); maxUpdate = max(maxUpdate,idx); - flags_[idx] = INTERNAL; + flag[idx] = INTERNAL; }, Kokkos::Min(minRange), Kokkos::Max(maxRange), @@ -202,7 +249,7 @@ bool pFlow::pointFlag::deletePoints if(numDeleted >= numActive_) { activeRange_ = {0, 0}; - numDeleted == numActive_; + numDeleted = numActive_; } numActive_ = numActive_ - numDeleted; @@ -245,7 +292,7 @@ bool pFlow::pointFlag::deletePoints if(numDeleted >= numActive_) { activeRange_ = {0, 0}; - numDeleted == numActive_; + numDeleted = numActive_; } numActive_ = numActive_ - numDeleted; @@ -264,13 +311,14 @@ bool pFlow::pointFlag::changeFlags ) { auto flg = getBoundaryFlag(boundaryIndex); + const auto& flags = flags_; Kokkos::parallel_for ( "pointFlag::changeFlags", rPolicy(0, changePoints.size()), LAMBDA_HD(uint32 i) { - flags_[changePoints(i)] = flg; + flags[changePoints(i)] = flg; } ); Kokkos::fence(); diff --git a/src/phasicFlow/structuredData/pointStructure/pointStructure/mortonIndexing.cpp b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/mortonIndexing.cpp similarity index 59% rename from src/phasicFlow/structuredData/pointStructure/pointStructure/mortonIndexing.cpp rename to src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/mortonIndexing.cpp index 63969747..2405bca4 100644 --- a/src/phasicFlow/structuredData/pointStructure/pointStructure/mortonIndexing.cpp +++ b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/mortonIndexing.cpp @@ -17,66 +17,60 @@ Licence: implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. -----------------------------------------------------------------------------*/ - #include "mortonIndexing.hpp" #include "cells.hpp" -#include "streams.hpp" -bool pFlow::getSortedIndex( +pFlow::uint32IndexContainer pFlow::getSortedIndices( box boundingBox, - real dx, - range activeRange, - ViewType1D pos, - ViewType1D flag, - int32IndexContainer& sortedIndex) + real dx, + const ViewType1D& pos, + const pFlagTypeDevice& flag +) { + if(flag.numActive() == 0u)return uint32IndexContainer(); // obtain the morton code of the particles - cells allCells( boundingBox, dx); - int32IndexContainer index(activeRange.first, activeRange.second); + cells allCells( boundingBox, dx); + auto aRange = flag.activeRange(); - ViewType1D mortonCode("mortonCode", activeRange.second); + uint32IndexContainer sortedIndex(aRange.start(), aRange.end()); - output<<"before first kernel"< mortonCode("mortonCode", aRange.end()); - using rpMorton = - Kokkos::RangePolicy>; - int32 numActive = 0; - Kokkos::parallel_reduce + + Kokkos::parallel_for ( "mortonIndexing::getIndex::morton", - rpMorton(activeRange.first, activeRange.second), - LAMBDA_HD(int32 i, int32& sumToUpdate){ - if( flag[i] == 1 ) // active point + deviceRPolicyStatic(aRange.start(), aRange.end()), + LAMBDA_HD(uint32 i){ + if( flag.isActive(i)) // active point { auto cellInd = allCells.pointIndex(pos[i]); mortonCode[i] = xyzToMortonCode64(cellInd.x(), cellInd.y(), cellInd.z()); - sumToUpdate++; }else { mortonCode[i] = xyzToMortonCode64 ( - static_cast(-1), - static_cast(-1), - static_cast(-1) + largestPosInt32, + largestPosInt32, + largestPosInt32 ); } - }, - numActive + } ); + Kokkos::fence(); + permuteSort( mortonCode, - activeRange.first, - activeRange.second, - index.deviceView(), + aRange.start(), + aRange.end(), + sortedIndex.deviceView(), 0 ); - index.modifyOnDevice(); - index.setSize(numActive); - index.syncViews(); - sortedIndex = index; + sortedIndex.modifyOnDevice(); + sortedIndex.syncViews(flag.numActive()); - return true; + return sortedIndex; } \ No newline at end of file diff --git a/src/phasicFlow/structuredData/pointStructure/pointStructure/mortonIndexing.hpp b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/mortonIndexing.hpp similarity index 94% rename from src/phasicFlow/structuredData/pointStructure/pointStructure/mortonIndexing.hpp rename to src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/mortonIndexing.hpp index 9ccdc83e..0d8c58d7 100644 --- a/src/phasicFlow/structuredData/pointStructure/pointStructure/mortonIndexing.hpp +++ b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/mortonIndexing.hpp @@ -17,24 +17,22 @@ Licence: implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. -----------------------------------------------------------------------------*/ - #ifndef __mortonIndexing_hpp__ #define __mortonIndexing_hpp__ #include "types.hpp" #include "box.hpp" #include "indexContainer.hpp" +#include "pointFlag.hpp" namespace pFlow { -bool getSortedIndex( +uint32IndexContainer getSortedIndices( box boundingBox, - real dx, - range activeRange, - ViewType1D pos, - ViewType1D flag, - int32IndexContainer& sortedIndex); + real dx, + const ViewType1D& pos, + const pFlagTypeDevice& flag); INLINE_FUNCTION_HD diff --git a/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/pointSorting.cpp b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/pointSorting.cpp new file mode 100644 index 00000000..c8dda774 --- /dev/null +++ b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/pointSorting.cpp @@ -0,0 +1,52 @@ +/*------------------------------- phasicFlow --------------------------------- + O C enter of + O O E ngineering and + O O M ultiscale modeling of + OOOOOOO F luid flow +------------------------------------------------------------------------------ + Copyright (C): www.cemf.ir + email: hamid.r.norouzi AT gmail.com +------------------------------------------------------------------------------ +Licence: + This file is part of phasicFlow code. It is a free software for simulating + granular and multiphase flows. You can redistribute it and/or modify it under + the terms of GNU General Public License v3 or any other later versions. + + phasicFlow is distributed to help others in their research in the field of + granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the + implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + +-----------------------------------------------------------------------------*/ + +#include "pointSorting.hpp" +#include "mortonIndexing.hpp" + +pFlow::pointSorting::pointSorting(const dictionary & dict) +: + performSorting_(dict.getValOrSet("active", Logical(false))), + timeControl_( + performSorting_()? + baseTimeControl(dict, "sorting"): + baseTimeControl(0,1,1, "sorting") + ), + dx_( + performSorting_()? + dict.getVal("dx"): + 1.0 + ) +{ + if( performSorting_() ) + REPORT(1)<<"Point sorting is "< &pos, + const pFlagTypeDevice &flag +) const +{ + return pFlow::getSortedIndices(boundingBox, dx_, pos, flag); +} diff --git a/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/pointSorting.hpp b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/pointSorting.hpp new file mode 100644 index 00000000..d8d9f195 --- /dev/null +++ b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointSorting/pointSorting.hpp @@ -0,0 +1,68 @@ +/*------------------------------- phasicFlow --------------------------------- + O C enter of + O O E ngineering and + O O M ultiscale modeling of + OOOOOOO F luid flow +------------------------------------------------------------------------------ + Copyright (C): www.cemf.ir + email: hamid.r.norouzi AT gmail.com +------------------------------------------------------------------------------ +Licence: + This file is part of phasicFlow code. It is a free software for simulating + granular and multiphase flows. You can redistribute it and/or modify it under + the terms of GNU General Public License v3 or any other later versions. + + phasicFlow is distributed to help others in their research in the field of + granular and multiphase flows, but WITHOUT ANY WARRANTY; without even the + implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. + +-----------------------------------------------------------------------------*/ +#ifndef __pointSorting_hpp__ +#define __pointSorting_hpp__ + +#include "baseTimeControl.hpp" +#include "indexContainer.hpp" +#include "pointFlag.hpp" + + +namespace pFlow +{ + +class box; + + +class pointSorting +{ +private: + + Logical performSorting_; + + baseTimeControl timeControl_; + + real dx_; + +public: + + explicit pointSorting(const dictionary& dict); + + bool performSorting()const + { + return performSorting_(); + } + + bool sortTime(uint32 iter, real t, real dt)const + { + return performSorting_() && timeControl_.timeEvent(iter, t, dt); + } + + uint32IndexContainer getSortedIndices( + const box& boundingBox, + const ViewType1D& pos, + const pFlagTypeDevice& flag + )const; + +}; + +} + +#endif \ No newline at end of file diff --git a/src/phasicFlow/structuredData/pointStructure/pointStructure/pointStructure.cpp b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointStructure.cpp index 9a20395c..bc0e7976 100644 --- a/src/phasicFlow/structuredData/pointStructure/pointStructure/pointStructure.cpp +++ b/src/phasicFlow/structuredData/pointStructure/pointStructure/pointStructure.cpp @@ -22,6 +22,7 @@ Licence: #include "pointStructure.hpp" #include "systemControl.hpp" #include "vocabs.hpp" +#include "anyList.hpp" bool pFlow::pointStructure::setupPointStructure(const realx3Vector& points) { @@ -73,6 +74,7 @@ bool pFlow::pointStructure::setupPointStructure(const PointsTypeHost &points) } boundaries_.createBoundaries(); + boundaries_.updateNeighborLists(0u, true); return true; } @@ -111,10 +113,13 @@ pFlow::pointStructure::pointStructure ( simulationDomain::create(control) ), + pointSorting_(simulationDomain_->subDictOrCreate("pointSorting")), boundaries_ ( *this - ) + ), + boundaryUpdateTimer_("boundaryUpdate", &timers()), + boundaryDataTransferTimer_("boundaryDataTransferTimer", &timers()) { REPORT(1)<< "Reading "<< Green_Text("pointStructure")<< " from "<subDictOrCreate("pointSorting")), boundaries_ ( *this @@ -165,12 +171,58 @@ pFlow::pointStructure::pointStructure( bool pFlow::pointStructure::beforeIteration() { - if( !boundaries_.beforeIteration(currentIter(), currentTime(), dt()) ) + uint32 iter = currentIter(); + real t = currentTime(); + real deltat = dt(); + + if(pointSorting_.sortTime(iter, t, deltat)) { - fatalErrorInFunction<< - "Unable to perform beforeIteration for boundaries"< simulationDomain_ = nullptr; + pointSorting pointSorting_; + boundaryList boundaries_; - + + Timer boundaryUpdateTimer_; + + Timer boundaryDataTransferTimer_; + + + bool setupPointStructure(const realx3Vector& points); bool setupPointStructure(const PointsTypeHost& points); diff --git a/src/phasicFlow/structuredData/zAxis/array2D.hpp b/src/phasicFlow/structuredData/zAxis/array2D.hpp index f76e6080..3a1bb059 100644 --- a/src/phasicFlow/structuredData/zAxis/array2D.hpp +++ b/src/phasicFlow/structuredData/zAxis/array2D.hpp @@ -33,12 +33,12 @@ struct array2D constexpr size_t nCols()const noexcept { - return nCols; + return nCol; } constexpr size_t nRows()const noexcept { - return nRows; + return nRow; } }; diff --git a/src/phasicFlow/triSurface/triangleFunctions.hpp b/src/phasicFlow/triSurface/triangleFunctions.hpp index 7345c07f..c6812d2f 100644 --- a/src/phasicFlow/triSurface/triangleFunctions.hpp +++ b/src/phasicFlow/triSurface/triangleFunctions.hpp @@ -39,7 +39,7 @@ realx3 normal(const realx3& p1, const realx3& p2, const realx3& p3) { auto n = cross(p2-p1, p3-p1); if( equal(n.length(), 0.0) ) - return zero3; + return realx3(0); else return normalize(n); } diff --git a/src/phasicFlow/types/basicTypes/math.hpp b/src/phasicFlow/types/basicTypes/math.hpp index ee140058..fbe3a6f1 100755 --- a/src/phasicFlow/types/basicTypes/math.hpp +++ b/src/phasicFlow/types/basicTypes/math.hpp @@ -21,15 +21,16 @@ Licence: #ifndef __math_hpp__ #define __math_hpp__ + +#include "builtinTypes.hpp" +#include "pFlowMacros.hpp" + #ifdef __CUDACC__ #include "math.h" #else #include #endif -#include "builtinTypes.hpp" -#include "pFlowMacros.hpp" - //* * * * * * * * * * * List of functinos * * * * * * * * // // abs, mod, exp, log, log10, pow, sqrt, cbrt // sin, cos, tan, asin, acos, atan, atan2 @@ -51,22 +52,26 @@ abs(real x) #endif } -#ifndef __CUDACC__ INLINE_FUNCTION_HD int64 abs(int64 x) { +#ifdef __CUDACC__ + return ::abs(x); +#else return std::abs(x); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD int32 abs(int32 x) { +#ifdef __CUDACC__ + return ::abs(x); +#else return std::abs(x); -} #endif +} INLINE_FUNCTION_HD real mod(real x, real y) @@ -78,6 +83,7 @@ mod(real x, real y) #endif } + INLINE_FUNCTION_HD int64 mod(int64 x, int64 y) { @@ -102,147 +108,188 @@ mod(uint32 x, uint32 y) return x % y; } -#ifndef __CUDACC__ INLINE_FUNCTION_HD real remainder(real x, real y) { +#ifdef __CUDACC__ + return ::remainder(x,y); +#else return std::remainder(x, y); -} #endif +} -#ifndef __CUDACC__ -INLINE_FUNCTION_H + +INLINE_FUNCTION_HD real exp(real x) { +#ifdef __CUDACC__ + return ::exp(x); +#else return std::exp(x); -} #endif +} + -#ifndef __CUDACC__ INLINE_FUNCTION_HD real log(real x) { +#ifdef __CUDACC__ + return ::log(x); +#else return std::log(x); -} #endif +} + -#ifndef __CUDACC__ INLINE_FUNCTION_HD real log10(real x) { +#ifdef __CUDACC__ + return ::log10(x); +#else return std::log10(x); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD real pow(real x, real y) { +#ifdef __CUDACC__ + return ::pow(x,y); +#else return std::pow(x, y); -} #endif +} + -#ifndef __CUDACC__ INLINE_FUNCTION_HD real sqrt(real x) { +#ifdef __CUDACC__ + return ::sqrt(x); +#else return std::sqrt(x); -} #endif +} + + -#ifndef __CUDACC__ INLINE_FUNCTION_HD real cbrt(real x) { +#ifdef __CUDACC__ + return ::cbrt(x); +#else return std::cbrt(x); +#endif } -#endif -#ifndef __CUDACC__ INLINE_FUNCTION_HD real sin(real x) { +#ifdef __CUDACC__ + return ::sin(x); +#else return std::sin(x); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD real cos(real x) { +#ifdef __CUDACC__ + return ::cos(x); +#else return std::cos(x); -} #endif +} + + -#ifndef __CUDACC__ INLINE_FUNCTION_HD real tan(real x) { +#ifdef __CUDACC__ + return ::tan(x); +#else return std::tan(x); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD real asin(real x) { +#ifdef __CUDACC__ + return ::asin(x); +#else return std::asin(x); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD real acos(real x) { +#ifdef __CUDACC__ + return ::acos(x); +#else return std::acos(x); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD real atan(real x) { +#ifdef __CUDACC__ + return ::atan(x); +#else return std::atan(x); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD real atan2(real y, real x) { +#ifdef __CUDACC__ + return ::atan2(y,x); +#else return std::atan2(y, x); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD real sinh(real x) { +#ifdef __CUDACC__ + return ::sinh(x); +#else return std::sinh(x); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD real cosh(real x) { +#ifdef __CUDACC__ + return ::cosh(x); +#else return std::cosh(x); -} #endif +} + INLINE_FUNCTION_HD real tanh(real x) @@ -274,14 +321,17 @@ acosh(real x) #endif } -#ifndef __CUDACC__ INLINE_FUNCTION_HD real atanh(real x) { +#ifdef __CUDACC__ + return ::atanh(x); +#else return std::atanh(x); -} #endif +} + INLINE_FUNCTION_HD real min(real x, real y) @@ -293,40 +343,54 @@ min(real x, real y) #endif } -#ifndef __CUDACC__ INLINE_FUNCTION_HD int64 min(int32 x, int32 y) { +#ifdef __CUDACC__ + return ::min(x,y); +#else return std::min(x, y); -} #endif +} + -#ifndef __CUDACC__ INLINE_FUNCTION_HD int64 min(int64 x, int64 y) { +#ifdef __CUDACC__ + return ::min(x,y); +#else return std::min(x, y); -} #endif +} + + -#ifndef __CUDACC__ INLINE_FUNCTION_HD uint64 min(uint64 x, uint64 y) { +#ifdef __CUDACC__ + return ::min(x,y); +#else return std::min(x, y); -} #endif +} + + -#ifndef __CUDACC__ INLINE_FUNCTION_HD uint32 min(uint32 x, uint32 y) { +#ifdef __CUDACC__ + return ::min(x,y); +#else return std::min(x, y); -} #endif +} + INLINE_FUNCTION_HD real max(real x, real y) @@ -338,37 +402,48 @@ max(real x, real y) #endif } -#ifndef __CUDACC__ INLINE_FUNCTION_HD int64 max(int64 x, int64 y) { +#ifdef __CUDACC__ + return ::max(x, y); +#else return std::max(x, y); -} #endif +} + + -#ifndef __CUDACC__ INLINE_FUNCTION_HD int32 max(int32 x, int32 y) { +#ifdef __CUDACC__ + return ::max(x, y); +#else return std::max(x, y); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD uint64 max(uint64 x, uint64 y) { +#ifdef __CUDACC__ + return ::max(x, y); +#else return std::max(x, y); -} #endif +} -#ifndef __CUDACC__ INLINE_FUNCTION_HD uint32 max(uint32 x, uint32 y) { +#ifdef __CUDACC__ + return ::max(x, y); +#else return std::max(x, y); -} #endif +} + } // pFlow diff --git a/src/setHelpers/setProperty.hpp b/src/setHelpers/setProperty.hpp index 43ec5ef2..978662d4 100644 --- a/src/setHelpers/setProperty.hpp +++ b/src/setHelpers/setProperty.hpp @@ -23,6 +23,6 @@ Licence: REPORT(0)<<"\nReading proprties . . . "<(vertecies)); solidNames_.push_back(name); } @@ -283,7 +283,7 @@ void pFlow::stlFile::addSolid realx3x3Vector&& vertecies ) { - solids_.push_backSafe(std::move(vertecies)); + solids_.push_back(makeUnique(vertecies)); solidNames_.push_back(name); } @@ -325,7 +325,7 @@ bool pFlow::stlFile::write()const { oFstream os(file_); os.precision(8); - for(label i=0; i= size() ) @@ -369,7 +369,7 @@ const pFlow::realx3x3Vector& pFlow::stlFile::solid const pFlow::word& pFlow::stlFile::name ( - label i + size_t i )const { if(i >= size() ) diff --git a/src/phasicFlow/triSurface/stlFile.hpp b/utilities/Utilities/geometryPhasicFlow/stlWall/stlFile.hpp similarity index 97% rename from src/phasicFlow/triSurface/stlFile.hpp rename to utilities/Utilities/geometryPhasicFlow/stlWall/stlFile.hpp index 6f47fbd5..eea2f928 100644 --- a/src/phasicFlow/triSurface/stlFile.hpp +++ b/utilities/Utilities/geometryPhasicFlow/stlWall/stlFile.hpp @@ -106,10 +106,10 @@ public: size_t size()const; // - vertecies of ith solid - const realx3x3Vector& solid(uint64 i)const; + const realx3x3Vector& solid(size_t i)const; // - name of ith solid - const word& name(uint64 i)const; + const word& name(size_t i)const; }; diff --git a/utilities/checkPhasicFlow/checkPhasicFlow.cpp b/utilities/checkPhasicFlow/checkPhasicFlow.cpp index 3ffe4168..cbd36d03 100755 --- a/utilities/checkPhasicFlow/checkPhasicFlow.cpp +++ b/utilities/checkPhasicFlow/checkPhasicFlow.cpp @@ -20,6 +20,7 @@ Licence: #include "KokkosTypes.hpp" #include "systemControl.hpp" +#include "localProcessors.hpp" #include "commandLine.hpp" @@ -40,9 +41,11 @@ if(!cmds.parse(argc, argv)) return 0; #include "initialize.hpp" output< pStructPtr = nullptr; + pFlow::uniquePtr pStructPtr = nullptr; if(!setOnly) @@ -84,13 +82,13 @@ int main( int argc, char* argv[] ) // position particles based on the dict content REPORT(0)<< "Positioning points . . . \n"<(Control, finalPos); + pStructPtr = pFlow::makeUnique(Control, finalPos); REPORT(1)<< "Created pStruct with "<< pStructPtr().size() << " points and capacity "<< @@ -100,11 +98,11 @@ int main( int argc, char* argv[] ) else { // read the content of pStruct from 0/pStructure - pStructPtr = makeUnique(Control); + pStructPtr = pFlow::makeUnique(Control); } - List> allObjects; + pFlow::List> allObjects; if(!positionOnly) { @@ -113,7 +111,7 @@ int main( int argc, char* argv[] ) auto& sfDict = cpDict.subDict("setFields"); - setFieldList defValueList(sfDict.subDict("defaultValue")); + pFlow::setFieldList defValueList(sfDict.subDict("defaultValue")); for(auto& sfEntry: defValueList) { @@ -128,7 +126,7 @@ int main( int argc, char* argv[] ) } } - output<("shapeName"); + auto& shapeName = Control.time().template lookupObject("shapeName"); REPORT(0)<< "Converting shapeName field to shapeIndex field"<