diff --git a/DEMSystems/CMakeLists.txt b/DEMSystems/CMakeLists.txt index fba75be5..6cc134e2 100644 --- a/DEMSystems/CMakeLists.txt +++ b/DEMSystems/CMakeLists.txt @@ -4,8 +4,8 @@ domainDistribute/domainDistribute.cpp DEMSystem/DEMSystem.cpp sphereDEMSystem/sphereFluidParticles.cpp sphereDEMSystem/sphereDEMSystem.cpp -# - +grainDEMSystem/grainFluidParticles.cpp +grainDEMSystem/grainDEMSystem.cpp ) set(link_libs Kokkos::kokkos phasicFlow Particles Geometry Property Interaction Interaction Utilities) diff --git a/DEMSystems/DEMSystem/DEMSystem.hpp b/DEMSystems/DEMSystem/DEMSystem.hpp index 2ad2e23d..23814d3e 100644 --- a/DEMSystems/DEMSystem/DEMSystem.hpp +++ b/DEMSystems/DEMSystem/DEMSystem.hpp @@ -121,6 +121,9 @@ public: virtual span diameter() = 0; + virtual + span courseGrainFactor() = 0; + virtual span acceleration()=0; diff --git a/DEMSystems/grainDEMSystem/grainDEMSystem.cpp b/DEMSystems/grainDEMSystem/grainDEMSystem.cpp new file mode 100644 index 00000000..0673012f --- /dev/null +++ b/DEMSystems/grainDEMSystem/grainDEMSystem.cpp @@ -0,0 +1,268 @@ +/*------------------------------- 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 "grainDEMSystem.hpp" +#include "vocabs.hpp" + +bool pFlow::grainDEMSystem::loop() +{ + + do + { + // + if(! insertion_().insertParticles( + Control().time().currentIter(), + Control().time().currentTime(), + Control().time().dt() ) ) + { + fatalError<< + "particle insertion failed in grainDEMSystem.\n"; + return false; + } + + geometry_->beforeIteration(); + + interaction_->beforeIteration(); + + particles_->beforeIteration(); + + + interaction_->iterate(); + + particles_->iterate(); + + geometry_->iterate(); + + particles_->afterIteration(); + + geometry_->afterIteration(); + + }while(Control()++); + + return true; +} + +pFlow::grainDEMSystem::grainDEMSystem( + word demSystemName, + const std::vector& domains, + int argc, + char* argv[], + bool requireRVel) +: + DEMSystem(demSystemName, domains, argc, argv), + requireRVel_(requireRVel) +{ + + REPORT(0)<<"\nReading proprties . . . "<( + propertyFile__, + Control().caseSetup().path()); + + REPORT(0)<< "\nCreating surface geometry for grainDEMSystem . . . "<(Control(), Property()); + + + insertion_ = makeUnique( + particles_(), + particles_().grains()); + + REPORT(0)<<"\nCreating interaction model for grain-grain contact and grain-wall contact . . ."<(domains, maxD); + +} + + +pFlow::grainDEMSystem::~grainDEMSystem() +{ + +} + +bool pFlow::grainDEMSystem::updateParticleDistribution( + real extentFraction, + const std::vector domains) +{ + + if(!particleDistribution_->changeDomainsSize( + extentFraction, + maxBounndingSphereSize(), + domains)) + { + fatalErrorInFunction<< + "cannot change the domain size"<locateParticles( + positionHost_, + particles_->pStruct().activePointsMaskHost())) + { + fatalErrorInFunction<< + "error in locating particles among sub-domains"< + pFlow::grainDEMSystem::numParInDomains()const +{ + const auto& distribute = particleDistribution_(); + int32 numDomains = distribute.numDomains(); + std::vector nums(numDomains); + for(int32 i=0; i +pFlow::grainDEMSystem::parIndexInDomain(int32 di)const +{ + return particleDistribution_->particlesInDomain(di); +} + +pFlow::span pFlow::grainDEMSystem::diameter() +{ + return span(diameterHost_.data(), diameterHost_.size()); +} + +pFlow::span pFlow::grainDEMSystem::courseGrainFactor() +{ + return span(particles_->courseGrainFactorHost().data(), particles_->courseGrainFactorHost().size()); +} + +pFlow::span pFlow::grainDEMSystem::acceleration() +{ + return span(nullptr, 0); +} + +pFlow::span pFlow::grainDEMSystem::velocity() +{ + return span(velocityHost_.data(), velocityHost_.size()); +} + +pFlow::span pFlow::grainDEMSystem::position() +{ + return span(positionHost_.data(), positionHost_.size()); +} + +pFlow::span pFlow::grainDEMSystem::rAcceleration() +{ + return span(nullptr, 0); +} + +pFlow::span pFlow::grainDEMSystem::rVelocity() +{ + return span(rVelocityHost_.data(), rVelocityHost_.size()); +} + +pFlow::span pFlow::grainDEMSystem::rPosition() +{ + return span(nullptr, 0); +} + +pFlow::span pFlow::grainDEMSystem::parFluidForce() +{ + auto& hVec = particles_->fluidForceHost(); + + return span(hVec.data(), hVec.size()); +} + +pFlow::span pFlow::grainDEMSystem::parFluidTorque() +{ + auto& hVec = particles_->fluidTorqueHost(); + return span(hVec.data(), hVec.size()); +} + +bool pFlow::grainDEMSystem::sendFluidForceToDEM() +{ + particles_->fluidForceHostUpdatedSync(); + return true; +} + +bool pFlow::grainDEMSystem::sendFluidTorqueToDEM() +{ + particles_->fluidTorqueHostUpdatedSync(); + return true; +} + +bool pFlow::grainDEMSystem::beforeIteration() +{ + velocityHost_ = std::as_const(particles_()).velocity().hostView(); + positionHost_ = std::as_const(particles_()).pointPosition().hostView(); + diameterHost_ = particles_->diameter().hostView(); + + if(requireRVel_) + rVelocityHost_ = std::as_const(particles_()).rVelocity().hostView(); + + + return true; +} + + +bool pFlow::grainDEMSystem::iterate( + real upToTime, + real timeToWrite, + word timeName) +{ + + Control().time().setStopAt(upToTime); + Control().time().setOutputToFile(timeToWrite, timeName); + + return loop(); + + return true; +} + +bool pFlow::grainDEMSystem::iterate(real upToTime) +{ + Control().time().setStopAt(upToTime); + return loop(); + return true; +} + +pFlow::real + pFlow::grainDEMSystem::maxBounndingSphereSize()const +{ + real minD, maxD; + particles_->boundingSphereMinMax(minD, maxD); + + return maxD; +} \ No newline at end of file diff --git a/DEMSystems/grainDEMSystem/grainDEMSystem.hpp b/DEMSystems/grainDEMSystem/grainDEMSystem.hpp new file mode 100644 index 00000000..07451942 --- /dev/null +++ b/DEMSystems/grainDEMSystem/grainDEMSystem.hpp @@ -0,0 +1,163 @@ +/*------------------------------- 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 __grainDEMSystem_hpp__ +#define __grainDEMSystem_hpp__ + +#include "DEMSystem.hpp" +#include "property.hpp" +#include "uniquePtr.hpp" +#include "geometry.hpp" +#include "grainFluidParticles.hpp" +#include "interaction.hpp" +#include "Insertions.hpp" +#include "domainDistribute.hpp" + + +namespace pFlow +{ + +class grainDEMSystem +: + public DEMSystem +{ +protected: + + // protected members + + uniquePtr property_ = nullptr; + + uniquePtr geometry_ = nullptr; + + uniquePtr particles_ = nullptr; + + uniquePtr insertion_ = nullptr; + + uniquePtr interaction_ = nullptr; + + uniquePtr particleDistribution_=nullptr; + + // to be used for CPU communications + ViewType1D velocityHost_; + + ViewType1D positionHost_; + + ViewType1D diameterHost_; + + bool requireRVel_ = false; + + ViewType1D rVelocityHost_; + + +// protected member functions + auto& Property() + { + return property_(); + } + + auto& Geometry() + { + return geometry_(); + } + + auto& Particles() + { + return particles_(); + } + + auto& Interaction() + { + return interaction_(); + } + + bool loop(); + +public: + + TypeInfo("grainDEMSystem"); + + grainDEMSystem( + word demSystemName, + const std::vector& domains, + int argc, + char* argv[], + bool requireRVel=false); + + virtual ~grainDEMSystem(); + + grainDEMSystem(const grainDEMSystem&)=delete; + + grainDEMSystem& operator = (const grainDEMSystem&)=delete; + + add_vCtor( + DEMSystem, + grainDEMSystem, + word); + + + bool updateParticleDistribution(real extentFraction, const std::vector domains) override; + + int32 numParInDomain(int32 di)const override; + + std::vector numParInDomains()const override; + + span parIndexInDomain(int32 di)const override; + + span diameter() override; + + span courseGrainFactor() override; + + span acceleration() override; + + span velocity() override; + + span position() override; + + span rAcceleration() override; + + span rVelocity() override; + + span rPosition() override; + + span parFluidForce() override; + + span parFluidTorque() override; + + bool sendFluidForceToDEM() override; + + bool sendFluidTorqueToDEM() override; + + bool beforeIteration() override; + + bool iterate( + real upToTime, + real timeToWrite, + word timeName) override; + + bool iterate(real upToTime) override; + + real maxBounndingSphereSize()const override; + +}; + + +} // pFlow + +#endif \ No newline at end of file diff --git a/DEMSystems/grainDEMSystem/grainFluidParticles.cpp b/DEMSystems/grainDEMSystem/grainFluidParticles.cpp new file mode 100644 index 00000000..f6710d71 --- /dev/null +++ b/DEMSystems/grainDEMSystem/grainFluidParticles.cpp @@ -0,0 +1,103 @@ +/*------------------------------- 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 "grainFluidParticles.hpp" +#include "grainFluidParticlesKernels.hpp" + + +pFlow::grainFluidParticles::grainFluidParticles( + systemControl &control, + const property& prop +) +: + grainParticles(control, prop), + fluidForce_( + objectFile( + "fluidForce", + "", + objectFile::READ_IF_PRESENT, + objectFile::WRITE_ALWAYS + ), + dynPointStruct(), + zero3 + ), + fluidTorque_( + objectFile( + "fluidTorque", + "", + objectFile::READ_IF_PRESENT, + objectFile::WRITE_NEVER + ), + dynPointStruct(), + zero3 + ) +{} + +bool pFlow::grainFluidParticles::beforeIteration() +{ + grainParticles::beforeIteration(); + + // copy the data (if required) from device to host + courseGrainFactorHost_ = coarseGrainFactor().hostView(); + + return true; +} + + +bool pFlow::grainFluidParticles::iterate() +{ + + accelerationTimer().start(); + pFlow::grainFluidParticlesKernels::acceleration( + control().g(), + mass().deviceViewAll(), + contactForce().deviceViewAll(), + fluidForce_.deviceViewAll(), + I().deviceViewAll(), + contactTorque().deviceViewAll(), + fluidTorque_.deviceViewAll(), + pStruct().activePointsMaskDevice(), + accelertion().deviceViewAll(), + rAcceleration().deviceViewAll() + ); + accelerationTimer().end(); + + intCorrectTimer().start(); + + dynPointStruct().correct(this->dt(), accelertion()); + + rVelIntegration().correct(this->dt(), rVelocity(), rAcceleration()); + + intCorrectTimer().end(); + + return true; +} + +void pFlow::grainFluidParticles::fluidForceHostUpdatedSync() +{ + copy(fluidForce_.deviceView(), fluidForceHost_); + return; +} + +void pFlow::grainFluidParticles::fluidTorqueHostUpdatedSync() +{ + copy(fluidTorque_.deviceView(), fluidTorqueHost_); + return; +} \ No newline at end of file diff --git a/DEMSystems/grainDEMSystem/grainFluidParticles.hpp b/DEMSystems/grainDEMSystem/grainFluidParticles.hpp new file mode 100644 index 00000000..f8adcee1 --- /dev/null +++ b/DEMSystems/grainDEMSystem/grainFluidParticles.hpp @@ -0,0 +1,102 @@ +/*------------------------------- 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. + +-----------------------------------------------------------------------------*/ + +/*! + @class pFlow::grainFluidParticles + + @brief Class for managing grain particles with fluid interactions + + This is a top-level class that contains the essential components for + defining grain prticles with fluid interactions in a DEM simulation. + */ + +#ifndef __grainFluidParticles_hpp__ +#define __grainFluidParticles_hpp__ + +#include "grainParticles.hpp" + + +namespace pFlow +{ + +class grainFluidParticles +: + public grainParticles +{ +protected: + + /// pointField of rotational acceleration of particles on device + realx3PointField_D fluidForce_; + + hostViewType1D fluidForceHost_; + + realx3PointField_D fluidTorque_; + + hostViewType1D fluidTorqueHost_; + + hostViewType1D courseGrainFactorHost_; + +public: + + /// construct from systemControl and property + grainFluidParticles(systemControl &control, const property& prop); + + ~grainFluidParticles() override = default; + + /// before iteration step + bool beforeIteration() override; + + /// iterate particles + bool iterate() override; + + auto& fluidForce() + { + return fluidForce_; + } + + auto& fluidTorque() + { + return fluidTorque_; + } + + auto& fluidForceHost() + { + return fluidForceHost_; + } + + auto& fluidTorqueHost() + { + return fluidTorqueHost_; + } + + auto& courseGrainFactorHost() + { + return courseGrainFactorHost_; + } + + void fluidForceHostUpdatedSync(); + + void fluidTorqueHostUpdatedSync(); + +}; //grainFluidParticles + +} // pFlow + +#endif //__sphereFluidParticles_hpp__ diff --git a/DEMSystems/grainDEMSystem/grainFluidParticlesKernels.hpp b/DEMSystems/grainDEMSystem/grainFluidParticlesKernels.hpp new file mode 100644 index 00000000..da82e5ed --- /dev/null +++ b/DEMSystems/grainDEMSystem/grainFluidParticlesKernels.hpp @@ -0,0 +1,79 @@ +/*------------------------------- 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 __grainFluidParticlesKernels_hpp__ +#define __grainFluidParticlesKernels_hpp__ + +namespace pFlow::grainFluidParticlesKernels +{ + +using rpAcceleration = Kokkos::RangePolicy< + DefaultExecutionSpace, + Kokkos::Schedule, + Kokkos::IndexType + >; + +template +void acceleration( + realx3 g, + deviceViewType1D mass, + deviceViewType1D force, + deviceViewType1D fluidForce, + deviceViewType1D I, + deviceViewType1D torque, + deviceViewType1D fluidTorque, + IncludeFunctionType incld, + deviceViewType1D lAcc, + deviceViewType1D rAcc + ) +{ + + auto activeRange = incld.activeRange(); + if(incld.isAllActive()) + { + Kokkos::parallel_for( + "pFlow::grainParticlesKernels::acceleration", + rpAcceleration(activeRange.first, activeRange.second), + LAMBDA_HD(int32 i){ + lAcc[i] = (force[i]+fluidForce[i])/mass[i] + g; + rAcc[i] = (torque[i]+fluidTorque[i])/I[i]; + }); + } + else + { + Kokkos::parallel_for( + "pFlow::grainParticlesKernels::acceleration", + rpAcceleration(activeRange.first, activeRange.second), + LAMBDA_HD(int32 i){ + if(incld(i)) + { + lAcc[i] = (force[i]+fluidForce[i])/mass[i] + g; + rAcc[i] = (torque[i]+fluidTorque[i])/I[i]; + } + }); + + } + + Kokkos::fence(); +} + +} + +#endif diff --git a/DEMSystems/sphereDEMSystem/sphereDEMSystem.cpp b/DEMSystems/sphereDEMSystem/sphereDEMSystem.cpp index ddf21f9e..29e2fff4 100644 --- a/DEMSystems/sphereDEMSystem/sphereDEMSystem.cpp +++ b/DEMSystems/sphereDEMSystem/sphereDEMSystem.cpp @@ -164,6 +164,11 @@ pFlow::span pFlow::sphereDEMSystem::diameter() return span(diameterHost_.data(), diameterHost_.size()); } +pFlow::span pFlow::sphereDEMSystem::courseGrainFactor() +{ + return span(nullptr, 0); +} + pFlow::span pFlow::sphereDEMSystem::acceleration() { return span(nullptr, 0); diff --git a/DEMSystems/sphereDEMSystem/sphereDEMSystem.hpp b/DEMSystems/sphereDEMSystem/sphereDEMSystem.hpp index 029ec46d..b2092490 100644 --- a/DEMSystems/sphereDEMSystem/sphereDEMSystem.hpp +++ b/DEMSystems/sphereDEMSystem/sphereDEMSystem.hpp @@ -122,6 +122,8 @@ public: span diameter() override; + span courseGrainFactor() override; + span acceleration() override; span velocity() override; diff --git a/src/Particles/GrainParticles/grainParticles/grainParticles.hpp b/src/Particles/GrainParticles/grainParticles/grainParticles.hpp index 9a0a35ac..b6c7a802 100644 --- a/src/Particles/GrainParticles/grainParticles/grainParticles.hpp +++ b/src/Particles/GrainParticles/grainParticles/grainParticles.hpp @@ -100,15 +100,23 @@ private: realVector& Is, uint32Vector& shIndex ); - /*bool initializeParticles(); + +protected: - bool insertSphereParticles( - const wordVector& names, - const int32IndexContainer& indices, - bool setId = true); + Timer& accelerationTimer() + { + return accelerationTimer_; + } - virtual uniquePtr> getFieldObjectList()const override; - */ + Timer& intCorrectTimer() + { + return intCorrectTimer_; + } + + integration& rVelIntegration() + { + return rVelIntegration_(); + } public: @@ -185,7 +193,12 @@ public: return grainDiameter_; } - const realPointField_D& coarseGrainFactor() const + const realPointField_D& grainDiameter()const + { + return grainDiameter_; + } + + const realPointField_D& coarseGrainFactor() const { return coarseGrainFactor_; }