From a9290e911ba86cc363a9e57c13f4a509b758bb79 Mon Sep 17 00:00:00 2001 From: Hamidreza Norouzi Date: Fri, 5 Apr 2024 05:31:09 -0700 Subject: [PATCH] Reflective boundary condition is added and tested. It requires messaging integration when changing velocity. --- .../AdamsBashforth2/AdamsBashforth2.hpp | 2 +- src/phasicFlow/CMakeLists.txt | 1 + .../boundaryField/boundaryField.cpp | 2 +- .../boundaryField/createBoundaryFields.hpp | 18 ++- .../reflectiveBoundaryField.cpp | 31 ++++ .../reflectiveBoundaryField.hpp | 90 +++++++++++ .../pointField/pointField/pointFields.cpp | 67 ++------ .../repository/repository/repository.hpp | 13 +- .../repository/repositoryTemplates.cpp | 71 ++++++--- .../boundaries/boundaryBase/boundaryBase.cpp | 19 ++- .../boundaries/boundaryBase/boundaryBase.hpp | 6 + .../boundaryReflective/boundaryReflective.cpp | 149 ++++++++++++++++++ .../boundaryReflective/boundaryReflective.hpp | 73 +++++++++ 13 files changed, 453 insertions(+), 89 deletions(-) create mode 100644 src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.cpp create mode 100644 src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.hpp create mode 100644 src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.cpp create mode 100644 src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.hpp diff --git a/src/Integration/AdamsBashforth2/AdamsBashforth2.hpp b/src/Integration/AdamsBashforth2/AdamsBashforth2.hpp index 863bd8a1..21e857cb 100644 --- a/src/Integration/AdamsBashforth2/AdamsBashforth2.hpp +++ b/src/Integration/AdamsBashforth2/AdamsBashforth2.hpp @@ -41,11 +41,11 @@ class AdamsBashforth2 { private: - auto& dy1() { return static_cast(*this); } + public: /// Type info diff --git a/src/phasicFlow/CMakeLists.txt b/src/phasicFlow/CMakeLists.txt index 3efe0d7e..49772f57 100644 --- a/src/phasicFlow/CMakeLists.txt +++ b/src/phasicFlow/CMakeLists.txt @@ -78,6 +78,7 @@ structuredData/boundaries/boundaryBase/boundaryBaseKernels.cpp structuredData/boundaries/boundaryExit/boundaryExit.cpp structuredData/boundaries/boundaryNone/boundaryNone.cpp structuredData/boundaries/boundaryPeriodic/boundaryPeriodic.cpp +structuredData/boundaries/boundaryReflective/boundaryReflective.cpp structuredData/boundaries/boundaryList.cpp structuredData/pointStructure/pointStructure/pointStructure.cpp structuredData/pointStructure/selectors/pStructSelector/pStructSelector.cpp diff --git a/src/phasicFlow/containers/pointField/boundaryField/boundaryField/boundaryField.cpp b/src/phasicFlow/containers/pointField/boundaryField/boundaryField/boundaryField.cpp index ecdb524f..2d2cc3bf 100644 --- a/src/phasicFlow/containers/pointField/boundaryField/boundaryField/boundaryField.cpp +++ b/src/phasicFlow/containers/pointField/boundaryField/boundaryField/boundaryField.cpp @@ -68,7 +68,7 @@ pFlow::uniquePtr> printKeys ( fatalError << "Ctor Selector "<< bType << "for type "<< - Yellow_Text(getTypeName()) << " dose not exist.\n" + Yellow_Text(getTypeName()) << " does not exist.\n" <<"Avaiable ones are: \n\n" , boundaryBasevCtorSelector_ diff --git a/src/phasicFlow/containers/pointField/boundaryField/createBoundaryFields.hpp b/src/phasicFlow/containers/pointField/boundaryField/createBoundaryFields.hpp index 8167a6ef..3b00050f 100644 --- a/src/phasicFlow/containers/pointField/boundaryField/createBoundaryFields.hpp +++ b/src/phasicFlow/containers/pointField/boundaryField/createBoundaryFields.hpp @@ -21,15 +21,23 @@ Licence: #ifndef __createBoundaryFields_hpp__ #define __createBoundaryFields_hpp__ -#include "boundaryField.hpp" -#include "exitBoundaryField.hpp" +#include "boundaryField/boundaryField.hpp" +#include "exitBoundaryField/exitBoundaryField.hpp" +#include "periodicBoundaryField/periodicBoundaryField.hpp" +#include "reflectiveBoundaryField/reflectiveBoundaryField.hpp" + +#define createDerivedBoundary(DataType, MemorySpaceType) \ + template class pFlow::exitBoundaryField; \ + template class pFlow::periodicBoundaryField; \ + template class pFlow::reflectiveBoundaryField; #define createBaseBoundary(DataType, MemorySpaceType) \ - template class pFlow::boundaryField; + template class pFlow::boundaryField; -#define createBoundary(DataType, MemorySpaceType, BoundaryType) \ - template class pFlow::BoundaryType##BoundaryField; +#define createBoundaryFields(DataType, MemorySpaceType) \ + createBaseBoundary(DataType, MemorySpaceType); \ + createDerivedBoundary(DataType, MemorySpaceType); #endif //__createBoundaryFields_hpp__ diff --git a/src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.cpp b/src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.cpp new file mode 100644 index 00000000..91bcc7b5 --- /dev/null +++ b/src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.cpp @@ -0,0 +1,31 @@ +/*------------------------------- 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. +-----------------------------------------------------------------------------*/ + +template + pFlow::reflectiveBoundaryField::reflectiveBoundaryField +( + const boundaryBase& boundary, + const pointStructure& pStruct, + InternalFieldType& internal +) +: + BoundaryFieldType(boundary, pStruct, internal) +{ + //this->addEvent(message::BNDR_DELETE); +} \ No newline at end of file diff --git a/src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.hpp b/src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.hpp new file mode 100644 index 00000000..e79f4fe6 --- /dev/null +++ b/src/phasicFlow/containers/pointField/boundaryField/reflectiveBoundaryField/reflectiveBoundaryField.hpp @@ -0,0 +1,90 @@ +/*------------------------------- 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 __reflectiveBoundaryField_hpp__ +#define __reflectiveBoundaryField_hpp__ + +#include "boundaryField.hpp" + +namespace pFlow +{ + +template< class T, class MemorySpace = void> +class reflectiveBoundaryField +: + public boundaryField +{ +public: + + using ReflectiveBoundaryFieldType = reflectiveBoundaryField; + + using BoundaryFieldType = boundaryField; + + using InternalFieldType = typename BoundaryFieldType::InternalFieldType; + + using memory_space = typename BoundaryFieldType::memory_space; + + using execution_space = typename BoundaryFieldType::execution_space; + + + +public: + + TypeInfoTemplate211("boundaryField","reflective", T, memory_space::name()); + + reflectiveBoundaryField( + const boundaryBase& boundary, + const pointStructure& pStruct, + InternalFieldType& internal); + + + add_vCtor + ( + BoundaryFieldType, + ReflectiveBoundaryFieldType, + boundaryBase + ); + + + bool hearChanges + ( + real t, + real dt, + uint32 iter, + const message& msg, + const anyList& varList + ) override + { + BoundaryFieldType::hearChanges(t,dt,iter, msg,varList); + + if(msg.equivalentTo(message::BNDR_DELETE)) + { + // do nothing; + } + + return true; + } + +}; + +} + +#include "reflectiveBoundaryField.cpp" + +#endif diff --git a/src/phasicFlow/containers/pointField/pointField/pointFields.cpp b/src/phasicFlow/containers/pointField/pointField/pointFields.cpp index 6f3cf375..032f6b2d 100644 --- a/src/phasicFlow/containers/pointField/pointField/pointFields.cpp +++ b/src/phasicFlow/containers/pointField/pointField/pointFields.cpp @@ -18,76 +18,37 @@ Licence: -----------------------------------------------------------------------------*/ - #include "pointFields.hpp" #include "createBoundaryFields.hpp" -#include "periodicBoundaryField.hpp" - -#define createAllBoundary(DataType, MemorySpaceType) \ - template class pFlow::exitBoundaryField; \ - template class pFlow::periodicBoundaryField; +#define createPointFields(DataType) \ + template class pFlow::pointField; \ + createBoundaryFields(DataType, pFlow::HostSpace); \ + \ + template class pFlow::pointField; \ + createBoundaryFields(DataType, void); // uint8 -template class pFlow::pointField; -createBaseBoundary(pFlow::uint8, pFlow::HostSpace); -createAllBoundary(pFlow::uint8, pFlow::HostSpace); - -template class pFlow::pointField; -createBaseBoundary(pFlow::uint8, void); -createAllBoundary(pFlow::uint8, void); +createPointFields(pFlow::uint8); /// uint32 -template class pFlow::pointField; -createBaseBoundary(pFlow::uint32, pFlow::HostSpace); -createAllBoundary(pFlow::uint32, pFlow::HostSpace); - -template class pFlow::pointField; -createBaseBoundary(pFlow::uint32, void); -createAllBoundary(pFlow::uint32, void); +createPointFields(pFlow::uint32); /// uint64 -template class pFlow::pointField; -createBaseBoundary(pFlow::uint64, pFlow::HostSpace); -createAllBoundary(pFlow::uint64, pFlow::HostSpace); - -template class pFlow::pointField; -createBaseBoundary(pFlow::uint64, void); -createAllBoundary(pFlow::uint64, void); +createPointFields(pFlow::uint64); /// real -template class pFlow::pointField; -createBaseBoundary(pFlow::real, pFlow::HostSpace); -createAllBoundary(pFlow::real, pFlow::HostSpace); - -template class pFlow::pointField; -createBaseBoundary(pFlow::real, void); -createAllBoundary(pFlow::real, void); +createPointFields(pFlow::real); /// realx3 -template class pFlow::pointField; -createBaseBoundary(pFlow::realx3, pFlow::HostSpace); -createAllBoundary(pFlow::realx3, pFlow::HostSpace); - - -template class pFlow::pointField; -createBaseBoundary(pFlow::realx3, void); -createAllBoundary(pFlow::realx3, void); +createPointFields(pFlow::realx3); /// realx4 -template class pFlow::pointField; -createBaseBoundary(pFlow::realx4, pFlow::HostSpace); -createAllBoundary(pFlow::realx4, pFlow::HostSpace); +createPointFields(pFlow::realx4); -template class pFlow::pointField; -createBaseBoundary(pFlow::realx4, void); -createAllBoundary(pFlow::realx4, void); - - -/// word +/// word, only on host template class pFlow::pointField; -createBaseBoundary(pFlow::word, pFlow::HostSpace); -createAllBoundary(pFlow::word, pFlow::HostSpace); \ No newline at end of file +createBoundaryFields(pFlow::word, pFlow::HostSpace); diff --git a/src/phasicFlow/repository/repository/repository.hpp b/src/phasicFlow/repository/repository/repository.hpp index 93698ab7..5e712c16 100644 --- a/src/phasicFlow/repository/repository/repository.hpp +++ b/src/phasicFlow/repository/repository/repository.hpp @@ -33,7 +33,7 @@ namespace pFlow class repository { -protected: +private: // - repository name word name_; @@ -53,10 +53,11 @@ protected: template - word reportTypeError (IOobject& object); + word reportTypeError (IOobject& object)const; template - bool checkForObjectType(IOobject& object); + static + bool checkForObjectType(IOobject& object); public: @@ -167,9 +168,13 @@ public: } } - // - return a ref to the underlaying data in the object + /// return a ref to the underlaying data in the object template T& lookupObject(const word& name); + + /// return a const ref to the underlaying data in the object + template + const T& lookupObject(const word& name)const; // - search the name and return a ref to repository repository& lookupRepository(const word& name); diff --git a/src/phasicFlow/repository/repository/repositoryTemplates.cpp b/src/phasicFlow/repository/repository/repositoryTemplates.cpp index a4c9ee1e..c1525430 100644 --- a/src/phasicFlow/repository/repository/repositoryTemplates.cpp +++ b/src/phasicFlow/repository/repository/repositoryTemplates.cpp @@ -20,7 +20,7 @@ Licence: template -pFlow::word pFlow::repository::reportTypeError(IOobject& object) +pFlow::word pFlow::repository::reportTypeError(IOobject& object)const { word err; err = "Object " + object.name() + " with type " + Type1::TYPENAME() + @@ -40,28 +40,57 @@ template T& pFlow::repository::lookupObject(const word& name) { if( auto [iter, success] = objects_.findIf(name); success ) - { + { - if( checkType(iter->second) ) - { - return static_cast(*iter->second); - + if( checkType(iter->second) ) + { + return static_cast(*iter->second); + - }else - { - fatalErrorInFunction << - reportTypeError(*iter->second)<(*iter->second); - } - - } - else - { - fatalErrorInFunction << + }else + { + fatalErrorInFunction << + reportTypeError(*iter->second)<(*iter->second); + } + + } + else + { + fatalErrorInFunction << "Object with name " << name << " is not found in repository " << this->name()<(*iter->second); - } + fatalExit; + return static_cast(*iter->second); + } +} + +template +const T& pFlow::repository::lookupObject(const word& name)const +{ + if( auto [iter, success] = objects_.findIf(name); success ) + { + + if( checkType(iter->second) ) + { + return static_cast(*iter->second); + } + else + { + fatalErrorInFunction << + reportTypeError(*iter->second)<(*iter->second); + } + + } + else + { + fatalErrorInFunction << + "Object with name " << name << " is not found in repository " << this->name()<(*iter->second); + } } \ No newline at end of file diff --git a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp index a7f6c96d..7a39276a 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.cpp @@ -24,6 +24,7 @@ Licence: #include "internalPoints.hpp" #include "anyList.hpp" #include "Time.hpp" +#include "pointStructure.hpp" #include "boundaryBaseKernels.hpp" void pFlow::boundaryBase::setSize(uint32 newSize) @@ -114,9 +115,9 @@ bool pFlow::boundaryBase::removeIndices message msgBndry = message::BNDR_RESET; - uint32 iter = internal_.time().currentIter(); - real t = internal_.time().currentTime(); - real dt = internal_.time().dt(); + uint32 iter = time().currentIter(); + real t = time().currentTime(); + real dt = time().dt(); if( !this->notify(iter, t, dt, msgBndry, aList) ) { @@ -197,6 +198,16 @@ pFlow::boundaryBase::boundaryBase unSyncLists(); } +const pFlow::pointStructure &pFlow::boundaryBase::pStruct() const +{ + return boundaries_.pStruct(); +} + +const pFlow::Time &pFlow::boundaryBase::time() const +{ + return boundaries_.pStruct().time(); +} + pFlow::boundaryBase &pFlow::boundaryBase::mirrorBoundary() { return boundaries_[mirrorBoundaryIndex()]; @@ -256,7 +267,7 @@ pFlow::uniquePtr pFlow::boundaryBase::create { printKeys ( - fatalError << "Ctor Selector "<< bType << " dose not exist. \n" + fatalError << "Ctor Selector "<< bType << " does not exist. \n" <<"Avaiable ones are: \n\n" , dictionaryvCtorSelector_ diff --git a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp index 261b7f94..06766b97 100644 --- a/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp +++ b/src/phasicFlow/structuredData/boundaries/boundaryBase/boundaryBase.hpp @@ -37,6 +37,8 @@ namespace pFlow class internalPoints; class dictionary; class boundaryList; +class pointStructure; +class Time; class boundaryBase : @@ -217,6 +219,10 @@ public: return internal_; } + const pointStructure& pStruct()const; + + const Time& time()const; + inline const auto& indexList()const { diff --git a/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.cpp b/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.cpp new file mode 100644 index 00000000..6d818a23 --- /dev/null +++ b/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.cpp @@ -0,0 +1,149 @@ +/*------------------------------- 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 "boundaryReflective.hpp" +#include "pointFields.hpp" +#include "dictionary.hpp" +#include "Time.hpp" + +pFlow::boundaryReflective::boundaryReflective +( + const dictionary &dict, + const plane &bplane, + internalPoints &internal, + boundaryList &bndrs, + uint32 thisIndex +) +: + boundaryBase + ( + dict, + bplane, + internal, + bndrs, + thisIndex + ) +{ + restitution_ = dict.getValOrSet("restitution", restitution_); + velocityName_ = dict.getValOrSet("velocityName",velocityName_); + diameterName_ = dict.getValOrSet("diameterName", diameterName_); +} + +bool pFlow::boundaryReflective::beforeIteration( + uint32 iterNum, + real t, + real dt) +{ + return true; +} + +bool pFlow::boundaryReflective::iterate +( + uint32 iterNum, + real t, + real dt +) +{ + return true; +} + +bool pFlow::boundaryReflective::afterIteration +( + uint32 iterNum, + real t, + real dt +) +{ + if(empty())return true; + + uint32 s = size(); + uint32Vector_D inContactFlags("inContactFlags",s+1, s+1, RESERVE()); + inContactFlags.fill(0u); + auto inContactFlagsD = inContactFlags.deviceViewAll(); + + auto points = thisPoints(); + + auto p = boundaryPlane().infPlane(); + const auto &diam = time().lookupObject(diameterName_); + auto diams = diam.BoundaryField(thisBoundaryIndex()).thisField(); + + uint32 numInContact = 0; + + Kokkos::parallel_reduce + ( + "pFlow::boundaryReflective::afterIteration", + deviceRPolicyStatic(0u,s), + LAMBDA_HD(uint32 i, uint32& nContactToUpdate) + { + if(p.inPositiveDistance(points(i), 0.5*diams(i))) + { + inContactFlagsD(i)=1; + nContactToUpdate++; + } + }, + numInContact + ); + + // no particle in contact + if(numInContact == 0 ) + { + return true; + } + + uint32Vector_D inContactList("inContactList", numInContact); + const auto& inContactListD = inContactList.deviceViewAll(); + + exclusiveScan(inContactFlagsD, 0u, s+1, inContactFlagsD, 0u); + + Kokkos::parallel_for + ( + "pFlow::boundaryReflective::afterIteration", + deviceRPolicyStatic(0, s), + LAMBDA_HD(uint32 i) + { + if(inContactFlagsD(i)!= inContactFlagsD(i+1)) + inContactListD(inContactFlagsD(i)) = points.index(i); + } + ); + Kokkos::fence(); + + const auto& velocity = time().lookupObject(velocityName_); + const auto& velocityD = velocity.deviceViewAll(); + + Kokkos::parallel_for( + "pFlow::boundaryReflective::velocityChange", + deviceRPolicyStatic(0,numInContact), + LAMBDA_HD(uint32 i) + { + auto& vel = velocityD(inContactListD(i)); + real vn = dot(p.normal(), vel); + + if(vn < 0) + { + realx3 vt = vel - vn*p.normal(); + vel = restitution_*(vt - vn*p.normal()); + } + } + ); + Kokkos::fence(); + + // TODO: notify integration for changes in the velocity + + return true; +} diff --git a/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.hpp b/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.hpp new file mode 100644 index 00000000..43e1ce13 --- /dev/null +++ b/src/phasicFlow/structuredData/boundaries/boundaryReflective/boundaryReflective.hpp @@ -0,0 +1,73 @@ +/*------------------------------- 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 __boundaryReflective_hpp__ +#define __boundaryReflective_hpp__ + + +#include "boundaryBase.hpp" + +namespace pFlow +{ + +class boundaryReflective +: + public boundaryBase +{ +private: + + real restitution_ = 0.95; + + word velocityName_{"velocity"}; + + word diameterName_{"diameter"}; + +public: + + TypeInfo("boundary"); + + boundaryReflective( + const dictionary &dict, + const plane &bplane, + internalPoints &internal, + boundaryList &bndrs, + uint32 thisIndex); + + + ~boundaryReflective()override = default; + + add_vCtor + ( + boundaryBase, + boundaryReflective, + dictionary + ); + + bool beforeIteration(uint32 iterNum, real t, real dt) override; + + bool iterate(uint32 iterNum, real t, real dt) override; + + bool afterIteration(uint32 iterNum, real t, real dt) override; + + +}; + +} + +#endif \ No newline at end of file