diff --git a/src/Integration/AdamsBashforth3/AdamsBashforth3.H b/src/Integration/AdamsBashforth3/AdamsBashforth3.H index c432926e..6b9d9e94 100644 --- a/src/Integration/AdamsBashforth3/AdamsBashforth3.H +++ b/src/Integration/AdamsBashforth3/AdamsBashforth3.H @@ -78,7 +78,7 @@ protected: //realx3PointField_D& dy2_; // this is a device - HistoryFieldType history_; + HistoryFieldType& history_; using rpIntegration = Kokkos::RangePolicy< DefaultExecutionSpace, diff --git a/src/Integration/AdamsBashforth4/AdamsBashforth4.C b/src/Integration/AdamsBashforth4/AdamsBashforth4.C new file mode 100644 index 00000000..13d25435 --- /dev/null +++ b/src/Integration/AdamsBashforth4/AdamsBashforth4.C @@ -0,0 +1,108 @@ +/*------------------------------- 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 "AdamsBashforth4.H" + + + +pFlow::AdamsBashforth4::AdamsBashforth4 +( + const word& baseName, + repository& owner, + const pointStructure& pStruct, + const word& method +) +: + integration(baseName, owner, pStruct, method), + history_( + owner.emplaceObject( + objectFile( + groupNames(baseName,"AB4History"), + "", + objectFile::READ_IF_PRESENT, + objectFile::WRITE_ALWAYS), + pStruct, + AB4History({zero3,zero3, zero3}))) + +{ + +} + +bool pFlow::AdamsBashforth4::predict +( + real UNUSED(dt), + realx3Vector_D& UNUSED(y), + realx3Vector_D& UNUSED(dy) +) +{ + + return true; +} + +bool pFlow::AdamsBashforth4::correct +( + real dt, + realx3Vector_D& y, + realx3Vector_D& dy +) +{ + + if(this->pStruct().allActive()) + { + return intAll(dt, y, dy, this->pStruct().activeRange()); + } + else + { + return intRange(dt, y, dy, this->pStruct().activePointsMaskD()); + } + + return true; +} + +bool pFlow::AdamsBashforth4::intAll( + real dt, + realx3Vector_D& y, + realx3Vector_D& dy, + range activeRng) +{ + auto d_dy = dy.deviceVectorAll(); + auto d_y = y.deviceVectorAll(); + auto d_history = history_.deviceVectorAll(); + + Kokkos::parallel_for( + "AdamsBashforth4::correct", + rpIntegration (activeRng.first, activeRng.second), + LAMBDA_HD(int32 i){ + d_y[i] += dt*( + static_cast(55.0 / 24.0) * d_dy[i] + - static_cast(59.0 / 24.0) * d_history[i].dy1_ + + static_cast(37.0 / 24.0) * d_history[i].dy2_ + - static_cast( 9.0 / 24.0) * d_history[i].dy3_ + ); + d_history[i].dy3_ = d_history[i].dy2_; + d_history[i].dy2_ = d_history[i].dy1_; + d_history[i].dy1_ = d_dy[i]; + + + }); + Kokkos::fence(); + + return true; +} \ No newline at end of file diff --git a/src/Integration/AdamsBashforth4/AdamsBashforth4.H b/src/Integration/AdamsBashforth4/AdamsBashforth4.H new file mode 100644 index 00000000..00b7c752 --- /dev/null +++ b/src/Integration/AdamsBashforth4/AdamsBashforth4.H @@ -0,0 +1,176 @@ +/*------------------------------- 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 __AdamsBashforth4_H__ +#define __AdamsBashforth4_H__ + + +#include "integration.H" +#include "pointFields.H" + +namespace pFlow +{ + +struct AB4History +{ + realx3 dy1_={0,0,0}; + realx3 dy2_={0,0,0}; + realx3 dy3_={0,0,0}; + + TypeNameNV("AB4History"); +}; + + +INLINE_FUNCTION +iIstream& operator>>(iIstream& str, AB4History& ab4) +{ + str.readBegin("AB4History"); + + str >> ab4.dy1_; + str >> ab4.dy2_; + str >> ab4.dy3_; + + str.readEnd("AB4History"); + + str.check(FUNCTION_NAME); + + return str; + +} + +INLINE_FUNCTION +iOstream& operator<<(iOstream& str, const AB4History& ab4) +{ + str << token::BEGIN_LIST << ab4.dy1_ + << token::SPACE << ab4.dy2_ + << token::SPACE << ab4.dy3_ + << token::END_LIST; + + str.check(FUNCTION_NAME); + + return str; +} + + +class AdamsBashforth4 +: + public integration +{ +protected: + + using HistoryFieldType = pointField; + //realx3PointField_D& dy1_; + + //realx3PointField_D& dy2_; + + + + // this is a device + HistoryFieldType& history_; + + using rpIntegration = Kokkos::RangePolicy< + DefaultExecutionSpace, + Kokkos::Schedule, + Kokkos::IndexType + >; + +public: + + // type info + TypeName("AdamsBashforth4"); + + //// - Constructors + AdamsBashforth4( + const word& baseName, + repository& owner, + const pointStructure& pStruct, + const word& method); + + virtual ~AdamsBashforth4()=default; + + // - add a virtual constructor + add_vCtor( + integration, + AdamsBashforth4, + word); + + + //// - Methods + bool predict( + real UNUSED(dt), + realx3Vector_D & UNUSED(y), + realx3Vector_D& UNUSED(dy)) override; + + bool correct(real dt, + realx3Vector_D & y, + realx3Vector_D& dy) override; + + bool intAll(real dt, + realx3Vector_D& y, + realx3Vector_D& dy, + range activeRng); + + template + bool intRange(real dt, + realx3Vector_D& y, + realx3Vector_D& dy, + activeFunctor activeP ); + +}; + + +template +bool pFlow::AdamsBashforth4::intRange( + real dt, + realx3Vector_D& y, + realx3Vector_D& dy, + activeFunctor activeP ) +{ + auto d_dy = dy.deviceVectorAll(); + auto d_y = y.deviceVectorAll(); + auto d_history = history_.deviceVectorAll(); + auto activeRng = activeP.activeRange(); + + Kokkos::parallel_for( + "AdamsBashforth4::correct", + rpIntegration (activeRng.first, activeRng.second), + LAMBDA_HD(int32 i){ + if( activeP(i)) + { + + d_y[i] += dt*( + static_cast(55.0 / 24.0) * d_dy[i] + - static_cast(59.0 / 24.0) * d_history[i].dy1_ + + static_cast(37.0 / 24.0) * d_history[i].dy2_ + - static_cast( 9.0 / 24.0) * d_history[i].dy3_ + ); + d_history[i].dy3_ = d_history[i].dy2_; + d_history[i].dy2_ = d_history[i].dy1_; + d_history[i].dy1_ = d_dy[i]; + } + }); + Kokkos::fence(); + + return true; +} + +} // pFlow + +#endif //__integration_H__ diff --git a/src/Integration/AdamsBashforth5/AdamsBashforth5.C b/src/Integration/AdamsBashforth5/AdamsBashforth5.C new file mode 100644 index 00000000..19971468 --- /dev/null +++ b/src/Integration/AdamsBashforth5/AdamsBashforth5.C @@ -0,0 +1,105 @@ +/*------------------------------- 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 "AdamsBashforth5.H" + + + +pFlow::AdamsBashforth5::AdamsBashforth5 +( + const word& baseName, + repository& owner, + const pointStructure& pStruct, + const word& method +) +: + integration(baseName, owner, pStruct, method), + history_( + owner.emplaceObject( + objectFile( + groupNames(baseName,"AB5History"), + "", + objectFile::READ_IF_PRESENT, + objectFile::WRITE_ALWAYS), + pStruct, + AB5History({zero3,zero3, zero3}))) + +{ + +} + +bool pFlow::AdamsBashforth5::predict +( + real UNUSED(dt), + realx3Vector_D& UNUSED(y), + realx3Vector_D& UNUSED(dy) +) +{ + + return true; +} + +bool pFlow::AdamsBashforth5::correct +( + real dt, + realx3Vector_D& y, + realx3Vector_D& dy +) +{ + + if(this->pStruct().allActive()) + { + return intAll(dt, y, dy, this->pStruct().activeRange()); + } + else + { + return intRange(dt, y, dy, this->pStruct().activePointsMaskD()); + } + + return true; +} + +bool pFlow::AdamsBashforth5::intAll( + real dt, + realx3Vector_D& y, + realx3Vector_D& dy, + range activeRng) +{ + auto d_dy = dy.deviceVectorAll(); + auto d_y = y.deviceVectorAll(); + auto d_history = history_.deviceVectorAll(); + + Kokkos::parallel_for( + "AdamsBashforth5::correct", + rpIntegration (activeRng.first, activeRng.second), + LAMBDA_HD(int32 i){ + d_y[i] += dt*( + static_cast(1901.0 / 720.0) * d_dy[i] + - static_cast(2774.0 / 720.0) * d_history[i].dy1_ + + static_cast(2616.0 / 720.0) * d_history[i].dy2_ + - static_cast(1274.0 / 720.0) * d_history[i].dy3_ + + static_cast( 251.0 / 720.0) * d_history[i].dy4_ + ); + d_history[i] = {d_dy[i] ,d_history[i].dy1_, d_history[i].dy2_, d_history[i].dy3_}; + }); + Kokkos::fence(); + + return true; +} \ No newline at end of file diff --git a/src/Integration/AdamsBashforth5/AdamsBashforth5.H b/src/Integration/AdamsBashforth5/AdamsBashforth5.H new file mode 100644 index 00000000..c389e8c4 --- /dev/null +++ b/src/Integration/AdamsBashforth5/AdamsBashforth5.H @@ -0,0 +1,173 @@ +/*------------------------------- 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 __AdamsBashforth5_H__ +#define __AdamsBashforth5_H__ + + +#include "integration.H" +#include "pointFields.H" + +namespace pFlow +{ + +struct AB5History +{ + realx3 dy1_={0,0,0}; + realx3 dy2_={0,0,0}; + realx3 dy3_={0,0,0}; + realx3 dy4_={0,0,0}; + + TypeNameNV("AB5History"); +}; + + +INLINE_FUNCTION +iIstream& operator>>(iIstream& str, AB5History& ab5) +{ + str.readBegin("AB5History"); + + str >> ab5.dy1_; + str >> ab5.dy2_; + str >> ab5.dy3_; + str >> ab5.dy4_; + + str.readEnd("AB5History"); + + str.check(FUNCTION_NAME); + + return str; + +} + +INLINE_FUNCTION +iOstream& operator<<(iOstream& str, const AB5History& ab5) +{ + str << token::BEGIN_LIST << ab5.dy1_ + << token::SPACE << ab5.dy2_ + << token::SPACE << ab5.dy3_ + << token::SPACE << ab5.dy4_ + << token::END_LIST; + + str.check(FUNCTION_NAME); + + return str; +} + + +class AdamsBashforth5 +: + public integration +{ +protected: + + using HistoryFieldType = pointField; + + // this is a device + HistoryFieldType& history_; + + using rpIntegration = Kokkos::RangePolicy< + DefaultExecutionSpace, + Kokkos::Schedule, + Kokkos::IndexType + >; + +public: + + // type info + TypeName("AdamsBashforth5"); + + //// - Constructors + AdamsBashforth5( + const word& baseName, + repository& owner, + const pointStructure& pStruct, + const word& method); + + virtual ~AdamsBashforth5()=default; + + // - add a virtual constructor + add_vCtor( + integration, + AdamsBashforth5, + word); + + + //// - Methods + bool predict( + real UNUSED(dt), + realx3Vector_D & UNUSED(y), + realx3Vector_D& UNUSED(dy)) override; + + bool correct(real dt, + realx3Vector_D & y, + realx3Vector_D& dy) override; + + bool intAll(real dt, + realx3Vector_D& y, + realx3Vector_D& dy, + range activeRng); + + template + bool intRange(real dt, + realx3Vector_D& y, + realx3Vector_D& dy, + activeFunctor activeP ); + +}; + + +template +bool pFlow::AdamsBashforth5::intRange( + real dt, + realx3Vector_D& y, + realx3Vector_D& dy, + activeFunctor activeP ) +{ + auto d_dy = dy.deviceVectorAll(); + auto d_y = y.deviceVectorAll(); + auto d_history = history_.deviceVectorAll(); + auto activeRng = activeP.activeRange(); + + Kokkos::parallel_for( + "AdamsBashforth5::correct", + rpIntegration (activeRng.first, activeRng.second), + LAMBDA_HD(int32 i){ + if( activeP(i)) + { + + d_y[i] += dt*( + static_cast(1901.0 / 720.0) * d_dy[i] + - static_cast(2774.0 / 720.0) * d_history[i].dy1_ + + static_cast(2616.0 / 720.0) * d_history[i].dy2_ + - static_cast(1274.0 / 720.0) * d_history[i].dy3_ + + static_cast( 251.0 / 720.0) * d_history[i].dy4_ + ); + d_history[i] = {d_dy[i] ,d_history[i].dy1_, d_history[i].dy2_, d_history[i].dy3_}; + } + }); + Kokkos::fence(); + + return true; +} + +} // pFlow + +#endif //__integration_H__ diff --git a/src/Integration/CMakeLists.txt b/src/Integration/CMakeLists.txt index 6c74df58..c96de103 100644 --- a/src/Integration/CMakeLists.txt +++ b/src/Integration/CMakeLists.txt @@ -1,6 +1,8 @@ list(APPEND SourceFiles integration/integration.C +AdamsBashforth5/AdamsBashforth5.C +AdamsBashforth4/AdamsBashforth4.C AdamsBashforth3/AdamsBashforth3.C AdamsBashforth2/AdamsBashforth2.C )