diff --git a/src/Geometry/geometry/geometry.cpp b/src/Geometry/geometry/geometry.cpp index 0875dbd6..00855156 100644 --- a/src/Geometry/geometry/geometry.cpp +++ b/src/Geometry/geometry/geometry.cpp @@ -34,8 +34,8 @@ bool pFlow::geometry::createPropertyId() uint32Vector propId( "propId", - surface().capacity(), - surface().size(), + capacity(), + size(), RESERVE()); ForAll(i, materialName_) @@ -61,6 +61,11 @@ bool pFlow::geometry::createPropertyId() } +void pFlow::geometry::zeroForce() +{ + contactForceWall_.fill(zero3); +} + pFlow::geometry::geometry ( systemControl& control, @@ -84,16 +89,6 @@ pFlow::geometry::geometry control ), wallProperty_(prop), - motionComponentName_ - ( - "motionComponentName", - "motionComponentName" - ), - materialName_ - ( - "materialName", - "materialName" - ), propertyId_ ( objectFile @@ -197,16 +192,6 @@ pFlow::geometry::geometry ( prop ), - motionComponentName_ - ( - "motionComponentName", - "motionComponentName" - ), - materialName_ - ( - "materialName", - "materialName" - ), propertyId_ ( objectFile @@ -274,9 +259,47 @@ pFlow::geometry::geometry } } +bool pFlow::geometry::beforeIteration() +{ + zeroForce(); + return true; +} + +bool pFlow::geometry::iterate() +{ + return true; +} + +bool pFlow::geometry::afterIteration() +{ + auto numTris = size(); + auto& normalsD = normals().deviceViewAll(); + auto& areaD = area().deviceViewAll(); + auto& cForceD = contactForceWall_.deviceViewAll(); + auto& sStressD = shearStressWall_.deviceViewAll(); + auto& nStressD = normalStressWall_.deviceViewAll(); + + Kokkos::parallel_for( + "pFlow::geometry::afterIteration", + deviceRPolicyStatic(0, numTris), + LAMBDA_HD(uint32 i) + { + realx3 n = normalsD[i]; + real A = max(areaD[i],smallValue); + realx3 nF = dot(cForceD[i],n)*n; + realx3 sF = cForceD[i] - nF; + nStressD[i] = nF/A; + sStressD[i] = sF/A; + }); + Kokkos::fence(); + + return true; +} + bool pFlow::geometry::read(iIstream &is, const IOPattern &iop) { motionComponentName_.read(is, iop); + materialName_.read(is, iop); if( readWholeObject_ ) @@ -305,6 +328,93 @@ bool pFlow::geometry::write(iOstream &os, const IOPattern &iop) const return multiTriSurface::write(os,iop); } + +pFlow::uniquePtr + pFlow::geometry::create +( + systemControl& control, + const property& prop +) +{ + + // + fileDictionary dict( motionModelFile__, control.time().geometry().path()); + + word model = dict.getVal("motionModel"); + + auto geomModel = angleBracketsNames("geometry", model); + + REPORT(1)<< "Selecting geometry model . . ."< + pFlow::geometry::create + ( + systemControl& control, + const property& prop, + multiTriSurface& surf, + const wordVector& motionCompName, + const wordVector& materialName, + const dictionary& motionDic + ) +{ + + word model = motionDic.getVal("motionModel"); + + auto geomModel = angleBracketsNames("geometry", model); + + REPORT(1)<< "Selecting geometry model . . ."< - pFlow::geometry::create -( - systemControl& control, - const property& prop -) -{ - //motionModelFile__ - auto motionDictPtr = IOobject::make - ( - objectFile - ( - motionModelFile__, - control.geometry().path(), - objectFile::READ_ALWAYS, - objectFile::WRITE_NEVER - ), - motionModelFile__, - true - ); - word model = motionDictPtr().getObject().getVal("motionModel"); - - auto geomModel = angleBracketsNames("geometry", model); - - REPORT(1)<< "Selecting geometry model . . ."< - pFlow::geometry::create - ( - systemControl& control, - const property& prop, - multiTriSurface& surf, - const wordVector& motionCompName, - const wordVector& materialName, - const dictionary& motionDic - ) -{ - - word model = motionDic.getVal("motionModel"); - - auto geomModel = angleBracketsNames("geometry", model); - - REPORT(1)<< "Selecting geometry model . . ."<(*this); - } - - /// Surface - inline - const auto& surface()const - { - return static_cast(*this); - } - inline const auto& motionComponentName()const { @@ -201,91 +159,54 @@ public: } /// Access to contact force - /*inline - realx3TriSurfaceField_D& contactForceWall() + inline + auto& contactForceWall() { return contactForceWall_; } /// Access to contact force inline - const realx3TriSurfaceField_D& contactForceWall() const + const auto& contactForceWall() const { return contactForceWall_; } + /// Property ide of triangles + inline + const auto& propertyId()const + { + return propertyId_; + } + /// Access to property inline const auto& wallProperty()const { return wallProperty_; - }*/ - - /// Owner repository - /*inline - const repository& owner()const - { - return geometryRepository_; - }*/ - - /// Owner repository - /*inline - repository& owner() - { - return geometryRepository_; - }*/ - - /// Path to the repository folder - /*inline auto path() - { - return owner().path(); - }*/ + } /// The name of motion model - /*virtual - word motionModelTypeName()const = 0;*/ + virtual + word motionModelTypeName()const = 0; /// Motion model index of triangles - /*virtual - const int8Vector_HD& triMotionIndex() const =0; + virtual + const uint32Field_D& triMotionIndex() const =0; /// Motion model index of points virtual - const int8Vector_HD& pointMotionIndex()const = 0; - - /// Property ide of triangles - const int8TriSurfaceField_D& propertyId() const - { - return propertyId_; - }*/ - - /// Operations before each iteration - //bool beforeIteration() override; - - /// Operations after each iteration - //bool afterIteration() override; + const uint32Field_D& pointMotionIndex()const = 0; - bool beforeIteration() override - { - notImplementedFunction; - return true; - } - + bool beforeIteration() override; + /// This is called in time loop. Perform the main calculations /// when the component should evolve along time. - bool iterate() override - { - notImplementedFunction; - return true; - } - + bool iterate() override; + /// This is called in time loop, after iterate. - bool afterIteration() override - { - notImplementedFunction; - return true; - } - + bool afterIteration() override; + //- IO bool read(iIstream& is, const IOPattern& iop) override; @@ -297,8 +218,10 @@ public: //- Static members - /*static - uniquePtr create(systemControl& control, const property& prop);*/ + static + uniquePtr create( + systemControl& control, + const property& prop); static uniquePtr create( diff --git a/src/Geometry/geometryMotion/geometryMotion.cpp b/src/Geometry/geometryMotion/geometryMotion.cpp index 53f3509d..6420dadb 100644 --- a/src/Geometry/geometryMotion/geometryMotion.cpp +++ b/src/Geometry/geometryMotion/geometryMotion.cpp @@ -45,17 +45,18 @@ bool pFlow::geometryMotion::findMotionIndex() fatalErrorInFunction<< mName<< " does not exist in the list of motion names -> "<< motionModel_.componentNames(); + return false; } surfMotionIndex.push_back(mInd); - auto surfRange = this->surface().subSurfaceRange(surfI); + auto surfRange = subSurfaceRange(surfI); for(uint32 i=0; isurface().subSurfacePointRange(surfI); + auto pointRange = subSurfacePointRange(surfI); for(uint32 n=0; n::findMotionIndex() } +template + bool pFlow::geometryMotion::moveGeometry() + { + + uint32 iter = this->currentIter(); + real t = this->currentTime(); + real dt = this->dt(); + + auto mModel = motionModel_.getModelInterface(iter, t, dt); + + auto& pointMIndexD= pointMotionIndex_.deviceViewAll(); + auto& pointsD = points().deviceViewAll(); + + + Kokkos::parallel_for( + "geometryMotion::movePoints", + deviceRPolicyStatic(0, numPoints()), + LAMBDA_HD(uint32 i){ + auto newPos = mModel.transferPoint(pointMIndexD[i], pointsD[i], dt); + pointsD[i] = newPos; + }); + + Kokkos::fence(); + + // move the motion components + motionModel_.move(iter, t,dt); + + // end of calculations + + + return true; + } + template pFlow::geometryMotion::geometryMotion ( @@ -137,129 +171,16 @@ pFlow::geometryMotion::geometryMotion } } -/*template - pFlow::geometryMotion::geometryMotion - ( - systemControl& control, - const property& prop, - const multiTriSurface& triSurface, - const wordVector& motionCompName, - const wordVector& propName, - const MotionModel& motionModel - ) - : - geometry( - control, - prop, - triSurface, - motionCompName, - propName - ), - motionModel_( - this->owner().template emplaceObject( - objectFile( - motionModelFile__, - "", - objectFile::READ_NEVER, - objectFile::WRITE_ALWAYS - ), - motionModel - ) - ), - moveGeomTimer_("move geometry", &this->timers()) - { - findMotionIndex(); - } - - template - pFlow::geometryMotion::geometryMotion - ( - systemControl& control, - const property& prop, - const dictionary& dict, - const multiTriSurface& triSurface, - const wordVector& motionCompName, - const wordVector& propName - ) - : - geometry( - control, - prop, - dict, - triSurface, - motionCompName, - propName - ), - motionModel_( - this->owner().template emplaceObject( - objectFile( - motionModelFile__, - "", - objectFile::READ_NEVER, - objectFile::WRITE_ALWAYS - ), - dict - ) - ), - moveGeomTimer_("move geometry", &this->timers()) - { - findMotionIndex(); - } - - template - bool pFlow::geometryMotion::beforeIteration() - { - geometry::beforeIteration(); - return true; - } template bool pFlow::geometryMotion::iterate() { if( motionModel_.isMoving() ) { - moveGeomTimer_.start(); - moveGeometry(); - moveGeomTimer_.end(); + moveGeomTimer_.start(); + moveGeometry(); + this->calculateNormals(); + moveGeomTimer_.end(); } return true; } - - template - bool pFlow::geometryMotion::afterIteration() - { - geometry::afterIteration(); - return true; - } - - template - bool pFlow::geometryMotion::moveGeometry() - { - - real dt = this->dt(); - real t = this->currentTime(); - - auto pointMIndex= pointMotionIndex_.deviceVector(); - auto mModel = motionModel_.getModel(t); - realx3* points = triSurface_.pointsData_D(); - auto numPoints = triSurface_.numPoints(); - - - Kokkos::parallel_for( - "geometryMotion::movePoints", - numPoints, - LAMBDA_HD(int32 i){ - auto newPos = mModel.transferPoint(pointMIndex[i], points[i], dt); - points[i] = newPos; - }); - - Kokkos::fence(); - - // move the motion components - motionModel_.move(t,dt); - - // end of calculations - moveGeomTimer_.end(); - - return true; - }*/ \ No newline at end of file diff --git a/src/Geometry/geometryMotion/geometryMotion.hpp b/src/Geometry/geometryMotion/geometryMotion.hpp index a7b0c3bc..157c05e0 100644 --- a/src/Geometry/geometryMotion/geometryMotion.hpp +++ b/src/Geometry/geometryMotion/geometryMotion.hpp @@ -37,6 +37,8 @@ class geometryMotion public geometry { public: + + using MotionModel = MotionModelType; using ModelComponent = typename MotionModelType::ModelComponent; @@ -60,6 +62,9 @@ private: /// determine the motion index of each triangle bool findMotionIndex(); + /// Move geometry + bool moveGeometry(); + public: /// Type info @@ -77,31 +82,13 @@ public: const wordVector& materialName, const dictionary& motionDict); - - /*geometryMotion( - systemControl& control, - const property& prop, - const multiTriSurface& triSurface, - const wordVector& motionCompName, - const wordVector& propName, - const MotionModel& motionModel);*/ - - /// Construct from components and dictionary that contains - /// motionModel - /*geometryMotion(systemControl& control, - const property& prop, - const dictionary& dict, - const multiTriSurface& triSurface, - const wordVector& motionCompName, - const wordVector& propName);*/ - /// Add virtual constructor - /*add_vCtor + add_vCtor ( geometry, geometryMotion, systemControl - );*/ + ); /// Add virtual constructor add_vCtor @@ -114,40 +101,33 @@ public: // - Methods /// Obtain motion model at time t - /*auto getModel(real t)const + auto getModel(uint32 iter, real t, real dt)const { - return motionModel_.getModel(t); - }*/ + return motionModel_.getModelInterface(iter, t, dt); + } /// TypeName / TypeInfo of motion model - /*word motionModelTypeName()const override + word motionModelTypeName()const override { return motionModel_.typeName(); - }*/ + } /// Access to motion model index of triangles - /*const auto& triMotionIndex()const override + const uint32Field_D& triMotionIndex()const override { return triMotionIndex_; } /// Access to motion model index of points - const int8Vector_HD& pointMotionIndex()const override + const uint32Field_D& pointMotionIndex()const override { return pointMotionIndex_; - }*/ - - /// Operations before each iteration - /*bool beforeIteration() override; + } /// Iterate geometry one time step bool iterate() override ; - /// Operations after each iteration - bool afterIteration() override;*/ - - /// Move geometry - //bool moveGeometry(); + }; @@ -155,9 +135,6 @@ public: #include "geometryMotion.cpp" -#ifndef BUILD_SHARED_LIBS - #include "geometryMotionsInstantiate.cpp" -#endif #endif //__geometryMotion_hpp__ diff --git a/src/Geometry/geometryMotion/geometryMotions.cpp b/src/Geometry/geometryMotion/geometryMotions.cpp index 4eaa179a..8b1d1d55 100644 --- a/src/Geometry/geometryMotion/geometryMotions.cpp +++ b/src/Geometry/geometryMotion/geometryMotions.cpp @@ -17,9 +17,13 @@ Licence: implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. -----------------------------------------------------------------------------*/ - #include "geometryMotions.hpp" -#ifdef BUILD_SHARED_LIBS -#include "geometryMotionsInstantiate.cpp" -#endif + +template class pFlow::geometryMotion; + +template class pFlow::geometryMotion; + +template class pFlow::geometryMotion; + +//template class pFlow::geometryMotion; diff --git a/src/Geometry/geometryMotion/geometryMotions.hpp b/src/Geometry/geometryMotion/geometryMotions.hpp index e0c109a1..1bb078d8 100644 --- a/src/Geometry/geometryMotion/geometryMotions.hpp +++ b/src/Geometry/geometryMotion/geometryMotions.hpp @@ -35,9 +35,10 @@ using vibratingMotionGeometry = geometryMotion; using rotationAxisMotionGeometry = geometryMotion; +using stationaryGeometry = geometryMotion; + //typedef geometryMotion multiRotationAxisMotionGeometry; -typedef geometryMotion stationaryGeometry; diff --git a/src/Geometry/geometryMotion/geometryMotionsInstantiate.cpp b/src/Geometry/geometryMotion/geometryMotionsInstantiate.cpp deleted file mode 100644 index 334e2aad..00000000 --- a/src/Geometry/geometryMotion/geometryMotionsInstantiate.cpp +++ /dev/null @@ -1,32 +0,0 @@ -/*------------------------------- 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 "stationaryWall.hpp" -#include "rotatingAxisMotion.hpp" -//#include "multiRotatingAxisMotion.hpp" -#include "vibratingMotion.hpp" - -template class pFlow::geometryMotion; - -template class pFlow::geometryMotion; - -//template class pFlow::geometryMotion; - -template class pFlow::geometryMotion;