From c0ee29e39c751c8133f68d6aa67efcd56d6a8420 Mon Sep 17 00:00:00 2001 From: Hamidreza Norouzi Date: Sun, 28 Jan 2024 14:43:10 -0800 Subject: [PATCH] boundaryExit beforeIteration, not tested --- .../sphereParticles/sphereParticles.cpp | 2 +- src/Particles/particles/particles.hpp | 4 +- src/phasicFlow/CMakeLists.txt | 2 +- src/phasicFlow/Kokkos/ViewAlgorithms.hpp | 28 ++--- .../algorithms/kokkosAlgorithms.hpp | 52 ++++---- .../containers/Vector/stdVectorHelper.hpp | 17 +++ .../containers/pointField/pointField.cpp | 43 ++++++- .../containers/pointField/pointField.hpp | 2 +- src/phasicFlow/containers/span/span.hpp | 8 +- .../repository/Time/baseTimeControl.cpp | 7 +- .../repository/Time/baseTimeControl.hpp | 33 ++--- .../boundaries/boundaryBase/boundaryBase.cpp | 2 +- .../boundaries/boundaryBase/boundaryBase.hpp | 43 +++++-- ...eldAccess.hpp => scatteredFieldAccess.hpp} | 53 +++++--- .../boundaries/boundaryExit/boundaryExit.cpp | 114 +++++++++++++++++- .../boundaries/boundaryExit/boundaryExit.hpp | 2 +- .../boundaryList.cpp | 84 ++++++++++--- .../boundaryList.hpp | 29 +++-- .../boundaries/boundaryNone/boundaryNone.cpp | 5 +- .../boundaries/boundaryNone/boundaryNone.hpp | 2 +- src/phasicFlow/structuredData/plane/plane.hpp | 5 + .../pointStructure/internalPoints.cpp | 22 +++- .../pointStructure/internalPoints.hpp | 5 + .../pointStructure/pointFlag.hpp | 4 + .../pointStructure/pointFlagKernels.hpp | 49 ++++++++ .../pointStructure/pointStructure.cpp | 9 +- 26 files changed, 484 insertions(+), 142 deletions(-) rename src/phasicFlow/structuredData/boundaries/boundaryBase/{scatterFieldAccess.hpp => scatteredFieldAccess.hpp} (52%) rename src/phasicFlow/structuredData/{pointStructure => boundaries}/boundaryList.cpp (65%) rename src/phasicFlow/structuredData/{pointStructure => boundaries}/boundaryList.hpp (81%) diff --git a/src/Particles/SphereParticles/sphereParticles/sphereParticles.cpp b/src/Particles/SphereParticles/sphereParticles/sphereParticles.cpp index fcb7f757..517c4e86 100644 --- a/src/Particles/SphereParticles/sphereParticles/sphereParticles.cpp +++ b/src/Particles/SphereParticles/sphereParticles/sphereParticles.cpp @@ -274,7 +274,7 @@ pFlow::sphereParticles::sphereParticles( ( "propertyId", "", - objectFile::READ_ALWAYS, + objectFile::READ_NEVER, objectFile::WRITE_NEVER ), dynPointStruct(), diff --git a/src/Particles/particles/particles.hpp b/src/Particles/particles/particles.hpp index 6b847e63..847b89ac 100644 --- a/src/Particles/particles/particles.hpp +++ b/src/Particles/particles/particles.hpp @@ -87,11 +87,11 @@ protected: return dynPointStruct_.pointPosition(); } - inline + /*inline auto& velocity() { return dynPointStruct_.velocity(); - } + }*/ inline auto& shapeIndex() diff --git a/src/phasicFlow/CMakeLists.txt b/src/phasicFlow/CMakeLists.txt index ae7df3ac..2afe777f 100644 --- a/src/phasicFlow/CMakeLists.txt +++ b/src/phasicFlow/CMakeLists.txt @@ -76,7 +76,7 @@ structuredData/pointStructure/pointStructure.cpp structuredData/boundaries/boundaryBase/boundaryBase.cpp structuredData/boundaries/boundaryExit/boundaryExit.cpp structuredData/boundaries/boundaryNone/boundaryNone.cpp -structuredData/pointStructure/boundaryList.cpp +structuredData/boundaries/boundaryList.cpp structuredData/pointStructure/selectors/pStructSelector/pStructSelector.cpp structuredData/pointStructure/selectors/selectBox/selectBox.cpp structuredData/pointStructure/selectors/selectRange/selectRange.cpp diff --git a/src/phasicFlow/Kokkos/ViewAlgorithms.hpp b/src/phasicFlow/Kokkos/ViewAlgorithms.hpp index 7e8134a9..1ee86959 100644 --- a/src/phasicFlow/Kokkos/ViewAlgorithms.hpp +++ b/src/phasicFlow/Kokkos/ViewAlgorithms.hpp @@ -241,19 +241,19 @@ void copy( } template < - typename dType, - typename sType, + typename Type, typename... sProperties> INLINE_FUNCTION_H void getNth( - dType& dst, - const ViewType1D& src, + Type& dst, + const ViewType1D& src, const uint32 n ) { - range32 span(n,n+1); - auto subV = Kokkos::subview(src, span); - hostViewType1D dstView("getNth",1); + + auto subV = Kokkos::subview(src, Kokkos::make_pair(n,n+1)); + hostViewType1D dstView("getNth",1); + //hostViewTypeScalar Kokkos::deep_copy(dstView,subV); dst = *dstView.data(); } @@ -425,13 +425,12 @@ int32 binarySearch( template< typename Type, typename... properties, - typename dType, typename... dProperties> void exclusiveScan( const ViewType1D& view, uint32 start, uint32 end, - ViewType1D& dView, + ViewType1D& dView, uint32 dStart ) { @@ -439,7 +438,7 @@ void exclusiveScan( ( areAccessible< typename ViewType1D::execution_space, - typename ViewType1D::memory_space>(), + typename ViewType1D::memory_space>(), "In exclusiveScan, view and dView should have the same space" ); @@ -448,7 +447,7 @@ void exclusiveScan( uint32 numElems = end-start; - pFlow::algorithms::KOKKOS::exclusiveScan( + pFlow::algorithms::KOKKOS::exclusiveScan( view.data()+start, dView.data()+dStart, numElems); @@ -458,13 +457,12 @@ void exclusiveScan( template< typename Type, typename... properties, - typename dType, typename... dProperties> void inclusiveScan( const ViewType1D& view, uint32 start, uint32 end, - ViewType1D& dView, + ViewType1D& dView, uint32 dStart) { using ExecutionSpace = typename ViewType1D::execution_space; @@ -473,14 +471,14 @@ void inclusiveScan( ( areAccessible< typename ViewType1D::execution_space, - typename ViewType1D::memory_space>(), + typename ViewType1D::memory_space>(), "In exclusiveScan, view and dView should have the same space" ); uint32 numElems = end-start; - pFlow::algorithms::KOKKOS::inclusiveScan( + pFlow::algorithms::KOKKOS::inclusiveScan( view.data()+start, dView.data()+dStart, numElems); diff --git a/src/phasicFlow/algorithms/kokkosAlgorithms.hpp b/src/phasicFlow/algorithms/kokkosAlgorithms.hpp index 218e6e6d..6102dd7a 100644 --- a/src/phasicFlow/algorithms/kokkosAlgorithms.hpp +++ b/src/phasicFlow/algorithms/kokkosAlgorithms.hpp @@ -30,15 +30,15 @@ namespace pFlow::algorithms::KOKKOS template INLINE_FUNCTION_H -int32 count(const Type* first, int32 numElems, const Type& val) +uint32 count(const Type* first, uint32 numElems, const Type& val) { using policy = Kokkos::RangePolicy< ExecutionSpace, - Kokkos::IndexType >; - int32 num = 0; + Kokkos::IndexType >; + uint32 num = 0; Kokkos::parallel_reduce("count", policy(0, numElems), - LAMBDA_HD(int32 i, int32& updateVal){ + LAMBDA_HD(uint32 i, uint32& updateVal){ if(equal(first[i],val)) updateVal++; }, num); @@ -50,17 +50,17 @@ int32 count(const Type* first, int32 numElems, const Type& val) template INLINE_FUNCTION_H -void fillSequence(Type* first, int32 numElems, const Type& firstVal) +void fillSequence(Type* first, uint32 numElems, const Type& firstVal) { using policy = Kokkos::RangePolicy< ExecutionSpace, - Kokkos::IndexType >; + Kokkos::IndexType >; Kokkos::parallel_for( "fillSequence", policy(0, numElems), - LAMBDA_HD(int32 i){ + LAMBDA_HD(uint32 i){ first[i] = firstVal+i; }); Kokkos::fence(); @@ -68,15 +68,15 @@ void fillSequence(Type* first, int32 numElems, const Type& firstVal) template INLINE_FUNCTION_H -void fillSelected(Type* first, const indexType* indices, const int32 numElems, const Type val) +void fillSelected(Type* first, const indexType* indices, const uint32 numElems, const Type val) { using policy = Kokkos::RangePolicy< ExecutionSpace, - Kokkos::IndexType >; + Kokkos::IndexType >; Kokkos::parallel_for( "fillSelected", policy(0,numElems), - LAMBDA_HD(int32 i){ + LAMBDA_HD(uint32 i){ first[indices[i]]= val; }); Kokkos::fence(); @@ -84,16 +84,16 @@ void fillSelected(Type* first, const indexType* indices, const int32 numElems, c template INLINE_FUNCTION_H -void fillSelected(Type* first, const indexType* indices, const Type* vals, const int32 numElems) +void fillSelected(Type* first, const indexType* indices, const Type* vals, const uint32 numElems) { using policy = Kokkos::RangePolicy< ExecutionSpace, - Kokkos::IndexType >; + Kokkos::IndexType >; Kokkos::parallel_for( "fillSelected", policy(0,numElems), - LAMBDA_HD(int32 i){ + LAMBDA_HD(uint32 i){ first[indices[i]]= vals[i]; }); Kokkos::fence(); @@ -101,17 +101,17 @@ void fillSelected(Type* first, const indexType* indices, const Type* vals, const template INLINE_FUNCTION_H -Type max(const Type* first, int32 numElems) +Type max(const Type* first, uint32 numElems) { using policy = Kokkos::RangePolicy< ExecutionSpace, - Kokkos::IndexType >; + Kokkos::IndexType >; Type maxElement=0; Kokkos::parallel_reduce( "max", policy(0, numElems), - LAMBDA_HD(int32 i, Type& maxUpdate){ + LAMBDA_HD(uint32 i, Type& maxUpdate){ if(maxUpdate(maxElement)); @@ -144,38 +144,38 @@ Type min(const Type* first, int32 numElems) //void sort(Type* first, int32 numElems, CompareFunc compare); //void permuteSort(const Type* first, PermuteType* pFirst, int32 numElems); -template -void exclusiveScan(Type* first, DestType* dFirst, int32 numElems) +template +void exclusiveScan(Type* first, Type* dFirst, uint32 numElems) { using policy = Kokkos::RangePolicy< ExecutionSpace, - Kokkos::IndexType >; + Kokkos::IndexType >; Kokkos::parallel_scan( "exclusiveScan", policy(0, numElems), - LAMBDA_HD(const int32 i, DestType& valToUpdate, const bool final) + LAMBDA_HD(const uint32 i, Type& valToUpdate, const bool final) { - const int32 val = first[i]; + const Type val = first[i]; if(final) dFirst[i] = valToUpdate; valToUpdate += val; }); } -template -void inclusiveScan(Type* first, DestType* dFirst, int32 numElems) +template +void inclusiveScan(Type* first, Type* dFirst, uint32 numElems) { using policy = Kokkos::RangePolicy< ExecutionSpace, - Kokkos::IndexType >; + Kokkos::IndexType >; Kokkos::parallel_scan( "inclusiveScan", policy(0, numElems), - LAMBDA_HD(const int32 i, int32& valToUpdate, const bool final) + LAMBDA_HD(const uint32 i, Type& valToUpdate, const bool final) { - const int32 val = first[i]; + const Type val = first[i]; valToUpdate += val; if(final) dFirst[i] = valToUpdate; diff --git a/src/phasicFlow/containers/Vector/stdVectorHelper.hpp b/src/phasicFlow/containers/Vector/stdVectorHelper.hpp index 37a85a6f..210ea048 100644 --- a/src/phasicFlow/containers/Vector/stdVectorHelper.hpp +++ b/src/phasicFlow/containers/Vector/stdVectorHelper.hpp @@ -48,6 +48,23 @@ public: template using vecAllocator = std::allocator; +template +inline +span makeSpan(std::vector& container) +{ + return span(container.data(), container.size()); +} + +template +inline +span makeSpan(const std::vector& container) +{ + return span( + const_cast(container.data()), + container.size()); +} + + template inline bool writeSpan( diff --git a/src/phasicFlow/containers/pointField/pointField.cpp b/src/phasicFlow/containers/pointField/pointField.cpp index edc42ce1..5b444402 100644 --- a/src/phasicFlow/containers/pointField/pointField.cpp +++ b/src/phasicFlow/containers/pointField/pointField.cpp @@ -1,4 +1,3 @@ -#include "pointField.hpp" /*------------------------------- phasicFlow --------------------------------- O C enter of O O E ngineering and @@ -41,7 +40,25 @@ pFlow::pointField::pointField boundaryFieldList_(pStruct.boundaries(), *this), pStruct_(pStruct), defaultValue_(defVal) -{} +{ + + if(IOobject::implyRead()) + { + REPORT(1)<< "Reading field "<< Green_Text(IOobject::name())<< + " from "< pFlow::pointField::pointField @@ -66,7 +83,23 @@ pFlow::pointField::pointField boundaryFieldList_(pStruct.boundaries(), *this), pStruct_(pStruct), defaultValue_(defVal) -{} +{ + if(IOobject::implyRead()) + { + REPORT(1)<< "Reading field "<< Green_Text(IOobject::name())<< + " from "< pFlow::pointField::pointField @@ -98,7 +131,9 @@ pFlow::pointField::pointField boundaryFieldList_(pStruct.boundaries(), *this), pStruct_(pStruct), defaultValue_(defVal) -{} +{ + +} /* diff --git a/src/phasicFlow/containers/pointField/pointField.hpp b/src/phasicFlow/containers/pointField/pointField.hpp index 9057174c..c6191dc8 100644 --- a/src/phasicFlow/containers/pointField/pointField.hpp +++ b/src/phasicFlow/containers/pointField/pointField.hpp @@ -52,7 +52,7 @@ public: using execution_space = typename InternalFieldType::execution_space; -protected: +private: //// - data members diff --git a/src/phasicFlow/containers/span/span.hpp b/src/phasicFlow/containers/span/span.hpp index 902ad990..91deb71f 100644 --- a/src/phasicFlow/containers/span/span.hpp +++ b/src/phasicFlow/containers/span/span.hpp @@ -157,7 +157,7 @@ public: }; -template class Container> +/*template class Container> size_t makeSpan(Container& container) { return span(container.data(), container.size()); @@ -169,7 +169,7 @@ size_t makeSpan(const Container& container) return span( const_cast(container.data()), container.size()); -} +}*/ template @@ -192,7 +192,7 @@ span charSpan(span s) s.size()*el); } -template class Container> +/*template class Container> span makeSpan(Container& container) { return span(container.data(), container.size()); @@ -204,7 +204,7 @@ span makeSpan(const Container& container) return span( const_cast(container.data()), container.size()); -} +}*/ diff --git a/src/phasicFlow/repository/Time/baseTimeControl.cpp b/src/phasicFlow/repository/Time/baseTimeControl.cpp index 642af864..54d8ea45 100644 --- a/src/phasicFlow/repository/Time/baseTimeControl.cpp +++ b/src/phasicFlow/repository/Time/baseTimeControl.cpp @@ -25,7 +25,8 @@ Licence: pFlow::baseTimeControl::baseTimeControl ( const dictionary &dict, - const word& intervalPrefix + const word& intervalPrefix, + real defStartTime ) { auto tControl = dict.getVal("timeControl"); @@ -47,7 +48,7 @@ pFlow::baseTimeControl::baseTimeControl if(isTimeStep_) { - auto startTime = (dict.getValOrSet("startTime", 0.0)); + auto startTime = (dict.getValOrSet("startTime", defStartTime)); auto endTime = (dict.getValOrSet("endTime", largeValue)); auto interval = dict.getVal(intervalWord); rRange_ = realStridedRange(startTime, endTime, interval); @@ -55,7 +56,7 @@ pFlow::baseTimeControl::baseTimeControl } else { - auto startTime = (dict.getValOrSet("startTime", 0.0)); + auto startTime = (dict.getValOrSet("startTime", 0)); auto endTime = (dict.getValOrSet("endTime", largestPosInt32)); auto interval = dict.getVal(intervalWord); iRange_ = int32StridedRagne(startTime, endTime, interval); diff --git a/src/phasicFlow/repository/Time/baseTimeControl.hpp b/src/phasicFlow/repository/Time/baseTimeControl.hpp index f82fda5f..647cd84e 100644 --- a/src/phasicFlow/repository/Time/baseTimeControl.hpp +++ b/src/phasicFlow/repository/Time/baseTimeControl.hpp @@ -1,20 +1,20 @@ /*------------------------------- phasicFlow --------------------------------- - O C enter of - O O E ngineering and - O O M ultiscale modeling of - OOOOOOO F luid flow + 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 + 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. + 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. + 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 __baseTimeControl_hpp__ @@ -30,8 +30,8 @@ namespace pFlow class baseTimeControl { private: - - bool isTimeStep_; + + bool isTimeStep_; int32StridedRagne iRange_; @@ -40,7 +40,10 @@ private: public: - baseTimeControl(const dictionary& dict, const word& intervalPrefix = ""); + baseTimeControl( + const dictionary& dict, + const word& intervalPrefix = "", + real defStartTime = 0.0); bool timeEvent(uint32 iter, real t, real dt)const; diff --git a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp index 10decdd5..d79b8ee0 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp @@ -50,9 +50,9 @@ pFlow::boundaryBase::boundaryBase : subscriber(dict.name()), boundaryPlane_(bplane), + indexList_(groupNames(dict.name(),"indexList")), neighborLength_(dict.getVal("neighborLength")), internal_(internal), - indexList_(groupNames(dict.name(),"indexList")), mirrorProcessoNo_(dict.getVal("mirrorProcessorNo")), name_(dict.name()), type_(dict.getVal("type")) diff --git a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp index 7f21e18e..45b28cc7 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp @@ -25,7 +25,7 @@ Licence: #include "subscriber.hpp" #include "VectorSingles.hpp" #include "plane.hpp" -#include "scatterFieldAccess.hpp" +#include "scatteredFieldAccess.hpp" #include "streams.hpp" @@ -44,29 +44,27 @@ class boundaryBase public: using pointFieldAccessType = - scatterFieldAccess; + deviceScatteredFieldAccess; - -protected: +private: const plane& boundaryPlane_; + /// list of particles indices on device + uint32Vector_D indexList_; + /// The length defined for creating neighbor list real neighborLength_; /// a reference to internalPoints& internal_; - /// list of particles indices on device - uint32Vector_D indexList_; - uint32 mirrorProcessoNo_; word name_; word type_; - public: TypeInfo("boundaryBase"); @@ -85,7 +83,7 @@ public: boundaryBase& operator=(boundaryBase&&) = default; - virtual ~boundaryBase() = default; + ~boundaryBase() override = default; create_vCtor @@ -117,22 +115,42 @@ public: return name_; } + bool empty()const + { + return indexList_.size()==0; + } + auto size()const { return indexList_.size(); } + const plane& boundaryPlane()const + { + return boundaryPlane_; + } + auto capacity()const { return indexList_.capacity(); } + const internalPoints& internal()const + { + return internal_; + } + + internalPoints& internal() + { + return internal_; + } + /// @brief set the size of indexList /// Always make sure that size+1 <= capacity void setSize(uint32 newSize); virtual - bool beforeIteratoin(uint32 iterNum, real t) = 0 ; + bool beforeIteration(uint32 iterNum, real t, real dt) = 0 ; virtual bool iterate(uint32 iterNum, real t) = 0; @@ -146,6 +164,11 @@ public: return indexList_; } + auto& indexList() + { + return indexList_; + } + pointFieldAccessType thisPoints(); virtual diff --git a/src/phasicFlow/structuredData/boundaries/boundaryBase/scatterFieldAccess.hpp b/src/phasicFlow/structuredData/boundaries/boundaryBase/scatteredFieldAccess.hpp similarity index 52% rename from src/phasicFlow/structuredData/boundaries/boundaryBase/scatterFieldAccess.hpp rename to src/phasicFlow/structuredData/boundaries/boundaryBase/scatteredFieldAccess.hpp index 27c16f2a..19c9caa7 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryBase/scatterFieldAccess.hpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryBase/scatteredFieldAccess.hpp @@ -1,6 +1,6 @@ -#ifndef __scatterFieldAccess_hpp__ -#define __scatterFieldAccess_hpp__ +#ifndef __scatteredFieldAccess_hpp__ +#define __scatteredFieldAccess_hpp__ #include "phasicFlowKokkos.hpp" @@ -9,8 +9,8 @@ namespace pFlow { -template -class scatterFieldAccess +template +class scatteredFieldAccess { public: @@ -22,7 +22,7 @@ public: using execution_space = typename viewType::execution_space; -protected: +private: uint32 size_ = 0; @@ -32,12 +32,12 @@ protected: public: - scatterFieldAccess(): + scatteredFieldAccess(): indices_(), fieldVals_() {} - scatterFieldAccess( + scatteredFieldAccess( uint32 sz, ViewType1D ind, ViewType1D fVals) @@ -47,15 +47,15 @@ public: fieldVals_(fVals) {} - scatterFieldAccess(const scatterFieldAccess&) = default; + scatteredFieldAccess(const scatteredFieldAccess&) = default; - scatterFieldAccess(scatterFieldAccess&&) = default; + scatteredFieldAccess(scatteredFieldAccess&&) = default; - scatterFieldAccess& operator=(const scatterFieldAccess&) = default; + scatteredFieldAccess& operator=(const scatteredFieldAccess&) = default; - scatterFieldAccess& operator=(scatterFieldAccess&&) = default; + scatteredFieldAccess& operator=(scatteredFieldAccess&&) = default; - ~scatterFieldAccess() = default; + ~scatteredFieldAccess() = default; // - Methods @@ -95,10 +95,35 @@ public: return size_ == 0; } + T getFirstCopy()const + { + T val; + uint32 n = indices_(0); + getNth(val, fieldVals_, n); + return val; + } + + T getLastCopy()const + { + T val; + if(empty())return val; + uint32 n = indices_(size()-1); + getNth(val, fieldVals_, n); + return val; + } + }; -} +template + using deviceScatteredFieldAccess = + scatteredFieldAccess; + +template + using hostScatteredFieldAccess = + scatteredFieldAccess; + +} // pFlow -#endif //__scatterFieldAccess_hpp__ +#endif //__scatteredFieldAccess_hpp__ diff --git a/src/phasicFlow/structuredData/boundaries/boundaryExit/boundaryExit.cpp b/src/phasicFlow/structuredData/boundaries/boundaryExit/boundaryExit.cpp index db5d4d5c..e8f005e8 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryExit/boundaryExit.cpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryExit/boundaryExit.cpp @@ -19,8 +19,10 @@ Licence: -----------------------------------------------------------------------------*/ #include "boundaryExit.hpp" +#include "internalPoints.hpp" #include "dictionary.hpp" + pFlow::boundaryExit::boundaryExit ( const dictionary& dict, @@ -36,12 +38,120 @@ pFlow::boundaryExit::boundaryExit dict.getValOrSet("checkForExitInterval", 1), 1); } -bool pFlow::boundaryExit::beforeIteratoin +bool pFlow::boundaryExit::beforeIteration ( uint32 iterNum, - real t + real t, + real dt ) { + // nothing have to be done + if(empty()) + { + return true; + } + + uint32 s = size(); + deviceViewType1D delFlags("delFlags",s+1); + deviceViewType1D keepFlags("keepFlags", s+1); + fill(delFlags, 0, s+1, 0u); + fill(keepFlags, 0, s+1, 0u); + + using policy = Kokkos::RangePolicy< + pFlow::DefaultExecutionSpace, + Kokkos::Schedule, + Kokkos::IndexType>; + + auto points = thisPoints(); + uint32 numDeleted = 0; + auto p = boundaryPlane().infPlane(); + + Kokkos::parallel_reduce + ( + "boundaryExit::beforeIteration", + policy(0,s), + LAMBDA_HD(uint32 i, uint32& delToUpdate) + { + if(p.pointInNegativeSide(points(i))) + { + delFlags(i)=1; + delToUpdate++; + } + else + { + keepFlags(i) = 1; + } + }, + numDeleted + ); + + // no point is deleted + if(numDeleted == 0 ) + { + return true; + } + + exclusiveScan(delFlags, 0u, s+1, delFlags, 0u); + exclusiveScan(keepFlags, 0u, s+1, keepFlags, 0u); + + deviceViewType1D deleteList("deleteList", numDeleted); + + Kokkos::parallel_for + ( + "boundaryExit::parllel_for", + policy(0, size()), + LAMBDA_HD(uint32 i) + { + if(delFlags(i)!= delFlags(i+1)) + deleteList(delFlags(i)) = i; + } + ); + Kokkos::fence(); + + deviceScatteredFieldAccess deleteIndices( + numDeleted, + deleteList, + indexList().deviceVectorAll()); + + // tell internal to remove these points from its list + if(!internal().deletePoints(deleteIndices)) + { + fatalErrorInFunction<< + "error in deleting points from boundary "<< name()< newIndices("newIndices", newSize); + auto oldIndices = indexList().deviceVectorAll(); + + Kokkos::parallel_for( + "fillIndices", + policy(0,s), + LAMBDA_HD(uint32 i) + { + if(keepFlags(i)!= keepFlags(i+1)) + { + newIndices(keepFlags(i)) = oldIndices(i); + } + } + ); + Kokkos::fence(); + copy(oldIndices, newIndices); + indexList().resize(newSize); + } + + // notify your observers + + WARNING<<"notify observers in boundary exit "<(simD.sizeOfBoundaries()), - internal_(internal), - simDomain_(simD), + ListPtr(pStruct.simDomain().sizeOfBoundaries()), + pStruct_(pStruct), timeControl_ ( - simDomain_.subDict("boundaries"), - "update" + pStruct.simDomain().subDict("boundaries"), + "update", + pStruct.currentTime() ) {} @@ -92,16 +91,16 @@ bool pFlow::boundaryList::setLists() { if(listSet_)return true; - for(auto i=0; iset ( i, boundaryBase::create ( - simDomain_.boundaryDict(i), - simDomain_.boundaryPlane(i), - internal_ + pStruct_.simDomain().boundaryDict(i), + pStruct_.simDomain().boundaryPlane(i), + pStruct_ ) ); } @@ -109,3 +108,54 @@ bool pFlow::boundaryList::setLists() return true; } + +bool pFlow::boundaryList::beforeIteration +( + uint32 iter, + real t, + real dt +) +{ + // it is time to update lists + if(timeControl_.timeEvent(iter, t, dt)) + { + if(!updateLists()) + { + fatalErrorInFunction; + return false; + } + + WARNING<<"Maybe notification about the update list "<beforeIteration(iter, t, dt)) + { + fatalErrorInFunction<< + "Error in beforeIteration in boundary "<name()<::operator[](i); + return ListPtr::operator[](i); } + + const baseTimeControl& timeControl()const + { + return timeControl_; + } + + bool beforeIteration(uint32 iter, real t, real dt); + + bool iterate(uint32 iter, real t, real dt); + + bool afterIteration(uint32 iter, real t, real dt); + }; } // pFlow diff --git a/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.cpp b/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.cpp index e3362956..9ac6b2b7 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.cpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.cpp @@ -30,10 +30,11 @@ pFlow::boundaryNone::boundaryNone boundaryBase(dict, bplane, internal) {} -bool pFlow::boundaryNone::beforeIteratoin +bool pFlow::boundaryNone::beforeIteration ( uint32 iterNum, - real t + real t, + real dt ) { return true; diff --git a/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.hpp b/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.hpp index da93977e..332bedd3 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.hpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryNone/boundaryNone.hpp @@ -50,7 +50,7 @@ public: dictionary ); - bool beforeIteratoin(uint32 iterNum, real t) override; + bool beforeIteration(uint32 iterNum, real t, real dt) override; bool iterate(uint32 iterNum, real t) override; diff --git a/src/phasicFlow/structuredData/plane/plane.hpp b/src/phasicFlow/structuredData/plane/plane.hpp index 239b700b..93930a85 100644 --- a/src/phasicFlow/structuredData/plane/plane.hpp +++ b/src/phasicFlow/structuredData/plane/plane.hpp @@ -122,6 +122,11 @@ public: // return the parallel plane to this plane plane parallelPlane(real distance)const; + infinitePlane infPlane()const + { + return infinitePlane(normal(), d()); + } + static bool validPlane4( diff --git a/src/phasicFlow/structuredData/pointStructure/internalPoints.cpp b/src/phasicFlow/structuredData/pointStructure/internalPoints.cpp index 2b092e94..e5766ab2 100644 --- a/src/phasicFlow/structuredData/pointStructure/internalPoints.cpp +++ b/src/phasicFlow/structuredData/pointStructure/internalPoints.cpp @@ -124,7 +124,21 @@ typename pFlow::internalPoints::PointsTypeHost return aPoints; } - +bool pFlow::internalPoints::deletePoints +( + scatteredFieldAccess delPoints +) +{ + if(!pFlagsD_.deletePoints(delPoints)) + { + fatalErrorInFunction<< + "Error in deleting points from internal points"<& dist ) { + pFlagSync_ = false; return pFlagsD_.markPointRegions ( dm, @@ -143,9 +158,7 @@ pFlow::uint32 pFlow::internalPoints::updateFlag dist[3], dist[4], dist[5] - ); - - pFlagSync_ = false; + ); } void pFlow::internalPoints::fillNeighborsLists @@ -166,7 +179,6 @@ void pFlow::internalPoints::fillNeighborsLists rearList, frontList); - pFlagSync_ = false; } FUNCTION_H diff --git a/src/phasicFlow/structuredData/pointStructure/internalPoints.hpp b/src/phasicFlow/structuredData/pointStructure/internalPoints.hpp index 91ed9de3..17413878 100644 --- a/src/phasicFlow/structuredData/pointStructure/internalPoints.hpp +++ b/src/phasicFlow/structuredData/pointStructure/internalPoints.hpp @@ -166,6 +166,11 @@ public: return pFlagsD_.activeRange(); } + ///@brief delete points at indices given in delPoints. + /// The default is that delPoints contains sorted indices + FUNCTION_H + bool deletePoints( + scatteredFieldAccess delPoints); FUNCTION_H uint32 updateFlag( diff --git a/src/phasicFlow/structuredData/pointStructure/pointFlag.hpp b/src/phasicFlow/structuredData/pointStructure/pointFlag.hpp index 83dbcec9..911ea654 100644 --- a/src/phasicFlow/structuredData/pointStructure/pointFlag.hpp +++ b/src/phasicFlow/structuredData/pointStructure/pointFlag.hpp @@ -23,6 +23,7 @@ Licence: #include "phasicFlowKokkos.hpp" #include "domain.hpp" +#include "scatteredFieldAccess.hpp" namespace pFlow { @@ -260,6 +261,9 @@ public: ViewType1D points); + bool deletePoints( + scatteredFieldAccess points); + /// @brief mark points based on their position in the domain. /// This should be the first method to be called when updating /// boundaries (step 1 of 2). diff --git a/src/phasicFlow/structuredData/pointStructure/pointFlagKernels.hpp b/src/phasicFlow/structuredData/pointStructure/pointFlagKernels.hpp index c5eaa78c..46a26828 100644 --- a/src/phasicFlow/structuredData/pointStructure/pointFlagKernels.hpp +++ b/src/phasicFlow/structuredData/pointStructure/pointFlagKernels.hpp @@ -148,6 +148,55 @@ pFlow::uint32 pFlow::pointFlag::scanPointFlag() return numActive; }*/ +template +bool pFlow::pointFlag::deletePoints +( + scatteredFieldAccess points +) +{ + if(points.empty())return true; + + uint32 minIndex = points.getFirstCopy(); + uint32 maxIndex = points.getLastCopy(); + + if(maxIndexactiveRange_.end())return false; + if(minIndex>; + + uint32 numDeleted = 0; + Kokkos::parallel_reduce + ( + "pointFlagKernels::deletePoints", + policy(0u, points.size()), + CLASS_LAMBDA_HD(uint32 i, uint32& valDelUpdate) + { + uint32 n = points(i); + if(isActive(n)) + { + valDelUpdate++; + flags_[n] = Flag::DELETED; + } + }, + numDeleted + ); + + if(numDeleted >= numActive_) + { + activeRange_ = {0, 0}; + numDeleted == numActive_; + } + + numActive_ = numActive_ - numDeleted; + isAllActive_ = + (activeRange_.numElements() == numActive_)&& numActive_>0; + + return true; +} + template pFlow::uint32 pFlow::pointFlag::markPointRegions diff --git a/src/phasicFlow/structuredData/pointStructure/pointStructure.cpp b/src/phasicFlow/structuredData/pointStructure/pointStructure.cpp index 54214a92..dcba001d 100644 --- a/src/phasicFlow/structuredData/pointStructure/pointStructure.cpp +++ b/src/phasicFlow/structuredData/pointStructure/pointStructure.cpp @@ -114,17 +114,17 @@ pFlow::pointStructure::pointStructure ), boundaries_ ( - simulationDomain_(), - *this + *this ) { - REPORT(0)<< "Reading point structure from "<< - IOobject::path()<