diff --git a/src/pe_coupling/all.h b/src/pe_coupling/all.h index 2dfc25e1272ad3e361b9d8bd094aa0f09fb4ce8d..53ba7b855d0848c8010f66bd34afa4ee29ba9d78 100644 --- a/src/pe_coupling/all.h +++ b/src/pe_coupling/all.h @@ -24,6 +24,7 @@ #pragma once +#include "amr/all.h" #include "discrete_particle_methods/all.h" #include "mapping/all.h" #include "geometry/all.h" diff --git a/src/pe_coupling/amr/BlockInfo.h b/src/pe_coupling/amr/BlockInfo.h new file mode 100644 index 0000000000000000000000000000000000000000..6d4873852ce3678479da0b1bcdf66fe51a12b124 --- /dev/null +++ b/src/pe_coupling/amr/BlockInfo.h @@ -0,0 +1,88 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file BlockInfo.h +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#pragma once + +#include <core/mpi/RecvBuffer.h> +#include <core/mpi/SendBuffer.h> + +#include <ostream> + +namespace walberla { +namespace pe_coupling { + +struct BlockInfo { + // lbm quantities + uint_t numberOfCells; + uint_t numberOfFluidCells; + uint_t numberOfNearBoundaryCells; + // pe quantities + uint_t numberOfLocalBodies; + uint_t numberOfShadowBodies; + uint_t numberOfContacts; + // coupling quantities + uint_t numberOfPeSubCycles; + + BlockInfo() + : numberOfCells(0), numberOfFluidCells(0), numberOfNearBoundaryCells(0), + numberOfLocalBodies(0), numberOfShadowBodies(0), numberOfContacts(0), + numberOfPeSubCycles(0) {} + + BlockInfo(const uint_t numCells, const uint_t numFluidCells, const uint_t numNearBoundaryCells, + const uint_t numLocalBodies, const uint_t numShadowBodies, const uint_t numContacts, + const uint_t numPeSubCycles) + : numberOfCells(numCells), numberOfFluidCells(numFluidCells), numberOfNearBoundaryCells(numNearBoundaryCells), + numberOfLocalBodies(numLocalBodies), numberOfShadowBodies(numShadowBodies), numberOfContacts(numContacts), + numberOfPeSubCycles(numPeSubCycles) {} +}; + + +inline +std::ostream& operator<<( std::ostream& os, const BlockInfo& bi ) +{ + os << bi.numberOfCells << " / " << bi.numberOfFluidCells << " / " << bi.numberOfNearBoundaryCells << " / " + << bi.numberOfLocalBodies << " / " << bi.numberOfShadowBodies << " / " << bi.numberOfContacts << " / " + << bi.numberOfPeSubCycles; + return os; +} + +template< typename T, // Element type of SendBuffer + typename G> // Growth policy of SendBuffer +mpi::GenericSendBuffer<T,G>& operator<<( mpi::GenericSendBuffer<T,G> & buf, const BlockInfo& info ) +{ + buf.addDebugMarker( "pca" ); + buf << info.numberOfCells << info.numberOfFluidCells << info.numberOfNearBoundaryCells + << info.numberOfLocalBodies << info.numberOfShadowBodies << info.numberOfContacts + << info.numberOfPeSubCycles; + return buf; +} + +template< typename T> // Element type of SendBuffer +mpi::GenericRecvBuffer<T>& operator>>( mpi::GenericRecvBuffer<T> & buf, BlockInfo& info ) +{ + buf.readDebugMarker( "pca" ); + buf >> info.numberOfCells >> info.numberOfFluidCells >> info.numberOfNearBoundaryCells + >> info.numberOfLocalBodies >> info.numberOfShadowBodies >> info.numberOfContacts + >> info.numberOfPeSubCycles; + return buf; +} + +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/InfoCollection.cpp b/src/pe_coupling/amr/InfoCollection.cpp new file mode 100644 index 0000000000000000000000000000000000000000..19238621785912cc6ee8bcfdaf3466c93a5b058b --- /dev/null +++ b/src/pe_coupling/amr/InfoCollection.cpp @@ -0,0 +1,100 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file InfoCollection.h +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + + +#include "InfoCollection.h" + +namespace walberla { +namespace pe_coupling { + + +void getBlockInfoFromInfoCollection( const PhantomBlock * block, const shared_ptr<InfoCollection>& ic, + BlockInfo & blockInfo ) +{ + WALBERLA_ASSERT_NOT_NULLPTR(block); + + if (block->sourceBlockIsLarger()) + { + // block is a result of refinement -> BlockInfo object only available for the father block + // there should be no particles on the block (otherwise it would not have been refined) + // and refinement in LBM does not change the number of cells + // we assume that the number of fluid and near boundary cells also stays the same + // (ATTENTION: not true for blocks intersecting with a boundary!) + // -> we can use the information of the father block for weight assignment + + auto infoIt = ic->find( block->getId().getFatherId() ); + WALBERLA_CHECK_UNEQUAL( infoIt, ic->end(), "Father block with ID " << block->getId().getFatherId() << " not found in info collection!" ); + + // check the above mentioned assumptions + WALBERLA_ASSERT_EQUAL(infoIt->second.numberOfLocalBodies, uint_t(0)); + WALBERLA_ASSERT_EQUAL(infoIt->second.numberOfShadowBodies, uint_t(0)); + WALBERLA_ASSERT_EQUAL(infoIt->second.numberOfContacts, uint_t(0)); + + blockInfo = infoIt->second; + } + else if (block->sourceBlockHasTheSameSize()) + { + auto infoIt = ic->find( block->getId() ); + WALBERLA_CHECK_UNEQUAL( infoIt, ic->end(), "Block with ID " << block->getId() << " not found in info collection!" ); + blockInfo = infoIt->second; + } + else + { + // source block of block is smaller + + // block is a result of coarsening -> BlockInfo object is available on all 8 child blocks + // there should be no particles on the block (otherwise it would not have been coarsened) + // and refinement in LBM does not change the number of cells + // we assume that the number of fluid and near boundary cells will be the average of all 8 child blocks + // -> we can use the information of the child blocks for weight assignment + + blockforest::BlockID childIdForInit(block->getId(), 0); + auto childForInitIt = ic->find( childIdForInit ); + WALBERLA_CHECK_UNEQUAL( childForInitIt, ic->end(), "Child block with ID " << childIdForInit << " not found in info collection!" ); + BlockInfo combinedInfo = childForInitIt->second; + uint_t numFluidCells(0), numNearBoundaryCells(0); + for (uint_t child = 0; child < 8; ++child) + { + blockforest::BlockID childId(block->getId(), child); + auto childIt = ic->find( childId ); + WALBERLA_CHECK_UNEQUAL( childIt, ic->end(), "Child block with ID " << childId << " not found in info collection!" ); + numFluidCells += childIt->second.numberOfFluidCells; + numNearBoundaryCells += childIt->second.numberOfNearBoundaryCells; + + // check above mentioned assumptions + WALBERLA_ASSERT_EQUAL(childIt->second.numberOfLocalBodies, uint_t(0)); + WALBERLA_ASSERT_EQUAL(childIt->second.numberOfShadowBodies, uint_t(0)); + WALBERLA_ASSERT_EQUAL(childIt->second.numberOfContacts, uint_t(0)); + } + // total number of cells remains unchanged + combinedInfo.numberOfFluidCells = uint_c(numFluidCells / uint_t(8)); //average + combinedInfo.numberOfNearBoundaryCells = uint_c( numNearBoundaryCells / uint_t(8) ); //average + combinedInfo.numberOfLocalBodies = uint_t(0); + combinedInfo.numberOfShadowBodies = uint_t(0); + combinedInfo.numberOfContacts = uint_t(0); //sum + // number of pe sub cycles stays the same + + blockInfo = combinedInfo; + } + +} + +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/InfoCollection.h b/src/pe_coupling/amr/InfoCollection.h new file mode 100644 index 0000000000000000000000000000000000000000..b1f1be21df72d5656d1e6a49252991a9dbf47d6d --- /dev/null +++ b/src/pe_coupling/amr/InfoCollection.h @@ -0,0 +1,121 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file InfoCollection.h +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#pragma once + +#include "BlockInfo.h" + +#include "blockforest/BlockForest.h" +#include "blockforest/BlockID.h" +#include "core/mpi/BufferSystem.h" + +#include "pe/ccd/ICCD.h" +#include "pe/fcd/IFCD.h" +#include "pe/rigidbody/BodyStorage.h" + +#include <map> + +namespace walberla { +namespace pe_coupling { + +typedef std::map<blockforest::BlockID, BlockInfo> InfoCollection; +typedef std::pair<blockforest::BlockID, BlockInfo> InfoCollectionPair; + +template <typename BoundaryHandling_T> +void createWithNeighborhood(BlockForest& bf, const BlockDataID boundaryHandlingID, + const BlockDataID bodyStorageID, const BlockDataID ccdID, const BlockDataID fcdID, + const uint_t numberOfPeSubCycles, + InfoCollection& ic ) +{ + ic.clear(); + + mpi::BufferSystem bs( MPIManager::instance()->comm(), 856 ); + + for (auto blockIt = bf.begin(); blockIt != bf.end(); ++blockIt) + { + blockforest::Block* block = static_cast<blockforest::Block*> (&(*blockIt)); + + // evaluate LBM quantities + BoundaryHandling_T * boundaryHandling = blockIt->getData< BoundaryHandling_T >( boundaryHandlingID ); + auto xyzSize = boundaryHandling->getFlagField()->xyzSize(); + const uint_t numCells = xyzSize.numCells(); + uint_t numFluidCells(0), numNearBoundaryCells(0); + for( auto cellIt = xyzSize.begin(); cellIt != xyzSize.end(); ++cellIt) + { + if( boundaryHandling->isDomain(*cellIt) ) + { + ++numFluidCells; + } + if( boundaryHandling->isNearBoundary(*cellIt)) + { + ++numNearBoundaryCells; + } + } + + // evaluate PE quantities + pe::Storage * bodyStorage = block->getData<pe::Storage>(bodyStorageID); + pe::BodyStorage const & localStorage = (*bodyStorage)[pe::StorageType::LOCAL]; + pe::BodyStorage const & shadowStorage = (*bodyStorage)[pe::StorageType::SHADOW]; + const uint_t numLocalParticles = localStorage.size(); + const uint_t numShadowParticles = shadowStorage.size(); + + pe::ccd::ICCD * ccd = block->getData<pe::ccd::ICCD>(ccdID); + pe::fcd::IFCD * fcd = block->getData<pe::fcd::IFCD>(fcdID); + ccd->generatePossibleContacts(); + pe::Contacts& contacts = fcd->generateContacts( ccd->getPossibleContacts() ); + const uint_t numContacts = contacts.size(); + + BlockInfo blockInfo(numCells, numFluidCells, numNearBoundaryCells, numLocalParticles, numShadowParticles, numContacts, numberOfPeSubCycles); + InfoCollectionPair infoCollectionEntry(block->getId(), blockInfo); + + ic.insert( infoCollectionEntry ); + + for( uint_t nb = uint_t(0); nb < block->getNeighborhoodSize(); ++nb ) + { + bs.sendBuffer( block->getNeighborProcess(nb) ) << infoCollectionEntry; + } + + //note: is it necessary to add child blocks already into the info collection? + // here, we still have full geometrical information and can probably determine number of fluid and near boundary cells more easily + // however, the interesting (and most costly) blocks are never refined and thus their child infos is never needed + // see pe/amr/InfoCollection.cpp for an example + + } + + // size of buffer is unknown and changes with each send + bs.setReceiverInfoFromSendBufferState(false, true); + bs.sendAll(); + + for( auto recvIt = bs.begin(); recvIt != bs.end(); ++recvIt ) + { + while( !recvIt.buffer().isEmpty() ) + { + InfoCollectionPair val; + recvIt.buffer() >> val; + ic.insert(val); + } + } +} + +void getBlockInfoFromInfoCollection( const PhantomBlock * block, const shared_ptr<InfoCollection>& ic, BlockInfo & blockInfo ); + + +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/all.h b/src/pe_coupling/amr/all.h new file mode 100644 index 0000000000000000000000000000000000000000..da4f4edfa80ea7004ca8903c15585a416bc807f2 --- /dev/null +++ b/src/pe_coupling/amr/all.h @@ -0,0 +1,29 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file all.h +//! \ingroup pe_coupling +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +//! \brief Collective header file for module pe_coupling +// +//====================================================================================================================== + +#pragma once + +#include "BlockInfo.h" +#include "InfoCollection.h" + +#include "level_determination/all.h" +#include "weight_assignment/all.h" \ No newline at end of file diff --git a/src/pe_coupling/amr/level_determination/BodyPresenceLevelDetermination.cpp b/src/pe_coupling/amr/level_determination/BodyPresenceLevelDetermination.cpp new file mode 100644 index 0000000000000000000000000000000000000000..4227c7cfacf27ffdf24c61ab5509a971a405b1b9 --- /dev/null +++ b/src/pe_coupling/amr/level_determination/BodyPresenceLevelDetermination.cpp @@ -0,0 +1,82 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file BodyPresenceLevelDetermination.cpp +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#include "BodyPresenceLevelDetermination.h" + +namespace walberla { +namespace pe_coupling { +namespace amr { + + +void BodyPresenceLevelDetermination::operator()( std::vector< std::pair< const Block *, uint_t > > & minTargetLevels, + std::vector< const Block * > &, const BlockForest & /*forest*/ ) +{ + for( auto it = minTargetLevels.begin(); it != minTargetLevels.end(); ++it ) + { + uint_t currentLevelOfBlock = it->first->getLevel(); + + const uint_t numberOfParticlesInDirectNeighborhood = getNumberOfLocalAndShadowBodiesInNeighborhood(it->first); + + uint_t targetLevelOfBlock = currentLevelOfBlock; //keep everything as it is + if ( numberOfParticlesInDirectNeighborhood > uint_t(0) ) + { + // set block to finest level if there are bodies nearby + targetLevelOfBlock = finestLevel_; + } + else + { + // block could coarsen since there are no bodies nearby + if( currentLevelOfBlock > uint_t(0) ) + targetLevelOfBlock = currentLevelOfBlock - uint_t(1); + } + + WALBERLA_CHECK_LESS_EQUAL(std::abs(int_c(targetLevelOfBlock) - int_c(currentLevelOfBlock)), uint_t(1), "Only level difference of maximum 1 allowed!"); + it->second = targetLevelOfBlock; + } +} + +uint_t BodyPresenceLevelDetermination::getNumberOfLocalAndShadowBodiesInNeighborhood(const Block * block) +{ + uint_t numBodies = uint_t(0); + + // add bodies of current block + const auto infoIt = infoCollection_->find(block->getId()); + WALBERLA_CHECK_UNEQUAL(infoIt, infoCollection_->end(), "Block with ID " << block->getId() << " not found in info collection!"); + + numBodies += infoIt->second.numberOfLocalBodies; + numBodies += infoIt->second.numberOfShadowBodies; + + // add bodies of all neighboring blocks + for(uint_t i = 0; i < block->getNeighborhoodSize(); ++i) + { + BlockID neighborBlockID = block->getNeighborId(i); + const auto infoItNeighbor = infoCollection_->find(neighborBlockID); + WALBERLA_CHECK_UNEQUAL(infoItNeighbor, infoCollection_->end(), "Neighbor block with ID " << neighborBlockID << " not found in info collection!"); + + numBodies += infoItNeighbor->second.numberOfLocalBodies; + numBodies += infoItNeighbor->second.numberOfShadowBodies; + } + return numBodies; +} + + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/level_determination/BodyPresenceLevelDetermination.h b/src/pe_coupling/amr/level_determination/BodyPresenceLevelDetermination.h new file mode 100644 index 0000000000000000000000000000000000000000..e57d318f01045172a6f0a43b9479dd19d5a1c0b4 --- /dev/null +++ b/src/pe_coupling/amr/level_determination/BodyPresenceLevelDetermination.h @@ -0,0 +1,62 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file BodyPresenceLevelDetermination.h +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#pragma once + +#include "blockforest/BlockForest.h" +#include "pe_coupling/amr/InfoCollection.h" + +namespace walberla { +namespace pe_coupling { +namespace amr { + +/* + * Class to determine the minimum level a block can be. + * For coupled LBM-PE simulations the following rules apply: + * - a moving body will always remain on the finest block + * - a moving body is not allowed to extend into an area with a coarser block + * - if no moving body is present, the level can be as coarse as possible (restricted by the 2:1 rule) + * Therefore, if a body, local or remote (due to bodies that are larger than a block), is present on any of the + * neighboring blocks of a certain block, this block's target level is the finest level. + * This, together with a refinement checking frequency that depends on the maximum translational body velocity, + * ensures the above given requirements. + */ +class BodyPresenceLevelDetermination +{ +public: + + BodyPresenceLevelDetermination( const shared_ptr<pe_coupling::InfoCollection> & infoCollection, uint_t finestLevel) : + infoCollection_( infoCollection ), finestLevel_( finestLevel) + {} + + void operator()( std::vector< std::pair< const Block *, uint_t > > & minTargetLevels, + std::vector< const Block * > &, const BlockForest & /*forest*/ ); + +private: + + uint_t getNumberOfLocalAndShadowBodiesInNeighborhood(const Block * block); + + shared_ptr<pe_coupling::InfoCollection> infoCollection_; + uint_t finestLevel_; +}; + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/level_determination/GlobalBodyPresenceLevelDetermination.cpp b/src/pe_coupling/amr/level_determination/GlobalBodyPresenceLevelDetermination.cpp new file mode 100644 index 0000000000000000000000000000000000000000..7556523fa6c170b56e1f0bfdf99220697b8be3cb --- /dev/null +++ b/src/pe_coupling/amr/level_determination/GlobalBodyPresenceLevelDetermination.cpp @@ -0,0 +1,90 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file GlobalBodyPresenceLevelDetermination.cpp +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#include "GlobalBodyPresenceLevelDetermination.h" + +#include "pe_coupling/geometry/PeOverlapFraction.h" + +namespace walberla { +namespace pe_coupling { +namespace amr { + +void GlobalBodyPresenceLevelDetermination::operator()( std::vector< std::pair< const Block *, uint_t > > & minTargetLevels, + std::vector< const Block * > &, const BlockForest & /*forest*/ ) +{ + for( auto it = minTargetLevels.begin(); it != minTargetLevels.end(); ++it ) + { + uint_t currentLevelOfBlock = it->first->getLevel(); + + auto blockExtendedAABB = it->first->getAABB().getExtended(blockExtensionLength_); + bool blockPartiallyOverlapsWithGlobalBodies = checkForPartialOverlapWithGlobalBodies(blockExtendedAABB); + + uint_t targetLevelOfBlock = currentLevelOfBlock; //keep everything as it is + if ( blockPartiallyOverlapsWithGlobalBodies ) + { + // set block to finest level since overlap with at least one global body is present + targetLevelOfBlock = finestLevel_; + } + else + { + // block could coarsen since there are no overlaps with global bodies + if( currentLevelOfBlock > uint_t(0) ) + targetLevelOfBlock = currentLevelOfBlock - uint_t(1); + } + + WALBERLA_ASSERT_LESS_EQUAL(std::abs(int_c(targetLevelOfBlock) - int_c(currentLevelOfBlock)), uint_t(1), "Only level difference of maximum 1 allowed!"); + it->second = targetLevelOfBlock; + } +} + +bool GlobalBodyPresenceLevelDetermination::checkForPartialOverlapWithGlobalBodies(const AABB& box) +{ + const Vector3<real_t> boxMidPoint( box.min() + real_t(0.5) * box.sizes()); + const Vector3<real_t> dxVec( box.sizes() ); + const uint_t maxDepthSuperSampling = uint_t(2); + + bool partialOverlapWithAllBodies = false; + + for( auto bodyIt = globalBodyStorage_->begin(); bodyIt != globalBodyStorage_->end(); ++bodyIt ) + { + auto bodyID = bodyIt.getBodyID(); + if( globalBodySelectorFct_(bodyID)) + { + real_t overlapFraction = pe_coupling::overlapFractionPe(*bodyID, boxMidPoint, dxVec, maxDepthSuperSampling); + + // check for partial overlap + if( overlapFraction > real_t(0) && overlapFraction < real_t(1)) + { + partialOverlapWithAllBodies = true; + } + // if fully contained in at least one body, no partial overlap possible + if( floatIsEqual(overlapFraction, real_t(1) )) + { + partialOverlapWithAllBodies = false; + break; + } + } + } + return partialOverlapWithAllBodies; +} + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/level_determination/GlobalBodyPresenceLevelDetermination.h b/src/pe_coupling/amr/level_determination/GlobalBodyPresenceLevelDetermination.h new file mode 100644 index 0000000000000000000000000000000000000000..1acdd271febd925de73daefd21657ac642d9c0aa --- /dev/null +++ b/src/pe_coupling/amr/level_determination/GlobalBodyPresenceLevelDetermination.h @@ -0,0 +1,66 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file GlobalBodyPresenceLevelDetermination.h +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#pragma once + +#include "blockforest/BlockForest.h" +#include "pe/rigidbody/BodyStorage.h" +#include "pe_coupling/utility/BodySelectorFunctions.h" + +namespace walberla { +namespace pe_coupling { +namespace amr { + +/* + * Class to determine the minimum level a block can be based on presence of global bodies. + * Only the global bodies that are given by the body selection function are considered. + * If a global body partially overlaps with a block this is block is determined to be on the finest grid. + * The block can be extended by the given extension length (e.g. the number of ghost layers * dx). + * This ensures correctness of the body mapping across block borders. + * Note: Blocks that are fully contained inside a global body can have any level. + */ +class GlobalBodyPresenceLevelDetermination +{ +public: + + GlobalBodyPresenceLevelDetermination( const shared_ptr<pe::BodyStorage> & globalBodyStorage, + uint_t finestLevel, real_t blockExtensionLength = real_t(0), + const std::function<bool(pe::BodyID)> & globalBodySelectorFct = selectAllBodies) : + globalBodyStorage_( globalBodyStorage ), finestLevel_( finestLevel ), + blockExtensionLength_( blockExtensionLength ), + globalBodySelectorFct_( globalBodySelectorFct ) + {} + + void operator()( std::vector< std::pair< const Block *, uint_t > > & minTargetLevels, + std::vector< const Block * > &, const BlockForest & /*forest*/ ); + +private: + + bool checkForPartialOverlapWithGlobalBodies(const AABB& box); + + shared_ptr<pe::BodyStorage> globalBodyStorage_; + uint_t finestLevel_; + real_t blockExtensionLength_; + std::function<bool(pe::BodyID)> globalBodySelectorFct_; +}; + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/level_determination/all.h b/src/pe_coupling/amr/level_determination/all.h new file mode 100644 index 0000000000000000000000000000000000000000..de0bbcaf5468a48b0f74d59535372f7d64f34d0b --- /dev/null +++ b/src/pe_coupling/amr/level_determination/all.h @@ -0,0 +1,26 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file all.h +//! \ingroup pe_coupling +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +//! \brief Collective header file for module pe_coupling +// +//====================================================================================================================== + +#pragma once + +#include "BodyPresenceLevelDetermination.h" +#include "GlobalBodyPresenceLevelDetermination.h" diff --git a/src/pe_coupling/amr/weight_assignment/MetisAssignmentFunctor.cpp b/src/pe_coupling/amr/weight_assignment/MetisAssignmentFunctor.cpp new file mode 100644 index 0000000000000000000000000000000000000000..374648d3b7971a24e31a10328d929b7027362a9f --- /dev/null +++ b/src/pe_coupling/amr/weight_assignment/MetisAssignmentFunctor.cpp @@ -0,0 +1,101 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file MetisAssignmentFunctor.cpp +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#include "MetisAssignmentFunctor.h" +#include "stencil/D3Q27.h" + +namespace walberla { +namespace pe_coupling { +namespace amr { + +void MetisAssignmentFunctor::operator()( std::vector< std::pair< const PhantomBlock *, walberla::any > > & blockData, const PhantomBlockForest & /*phantomBlockForest*/) +{ + for( auto it = blockData.begin(); it != blockData.end(); ++it ) + { + const PhantomBlock * block = it->first; + //only change of one level is supported! + WALBERLA_ASSERT_LESS( std::abs(int_c(block->getLevel()) - int_c(block->getSourceLevel())), 2 ); + + BlockInfo blockInfo; + pe_coupling::getBlockInfoFromInfoCollection(block, ic_, blockInfo); + + + std::vector<int64_t> metisVertexWeights(ncon_); + + for( uint_t con = uint_t(0); con < ncon_; ++con ) + { + real_t vertexWeight = std::max(weightEvaluationFct_[con](blockInfo), blockBaseWeight_); + + int64_t metisVertexWeight = int64_c( vertexWeight ); + + WALBERLA_ASSERT_GREATER(metisVertexWeight, int64_t(0)); + metisVertexWeights[con] = metisVertexWeight; + } + + blockforest::DynamicParMetisBlockInfo info( metisVertexWeights ); + + info.setVertexCoords( it->first->getAABB().center() ); + + real_t blockVolume = it->first->getAABB().volume(); + real_t approximateEdgeLength = std::cbrt( blockVolume ); + + int64_t faceNeighborWeight = int64_c(approximateEdgeLength * approximateEdgeLength ); //common face + int64_t edgeNeighborWeight = int64_c(approximateEdgeLength); //common edge + int64_t cornerNeighborWeight = int64_c( 1 ); //common corner + + int64_t vertexSize = int64_c(blockVolume); + info.setVertexSize( vertexSize ); + + for( const uint_t idx : blockforest::getFaceNeighborhoodSectionIndices() ) + { + for( uint_t nb = uint_t(0); nb < it->first->getNeighborhoodSectionSize(idx); ++nb ) + { + auto neighborBlockID = it->first->getNeighborId(idx,nb); + info.setEdgeWeight(neighborBlockID, faceNeighborWeight ); + } + } + + for( const uint_t idx : blockforest::getEdgeNeighborhoodSectionIndices() ) + { + for( uint_t nb = uint_t(0); nb < it->first->getNeighborhoodSectionSize(idx); ++nb ) + { + auto neighborBlockID = it->first->getNeighborId(idx,nb); + info.setEdgeWeight(neighborBlockID, edgeNeighborWeight ); + } + } + + for( const uint_t idx : blockforest::getCornerNeighborhoodSectionIndices() ) + { + for( uint_t nb = uint_t(0); nb < it->first->getNeighborhoodSectionSize(idx); ++nb ) + { + auto neighborBlockID = it->first->getNeighborId(idx,nb); + info.setEdgeWeight(neighborBlockID, cornerNeighborWeight ); + } + } + + it->second = info; + + } +} + + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/weight_assignment/MetisAssignmentFunctor.h b/src/pe_coupling/amr/weight_assignment/MetisAssignmentFunctor.h new file mode 100644 index 0000000000000000000000000000000000000000..384746aa6c3777add20c1283fde6c2ca06523a09 --- /dev/null +++ b/src/pe_coupling/amr/weight_assignment/MetisAssignmentFunctor.h @@ -0,0 +1,68 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file MetisAssignmentFunctor.h +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#pragma once + +#include "pe_coupling/amr/InfoCollection.h" + +#include "blockforest/loadbalancing/DynamicParMetis.h" + +#include <functional> +#include <vector> + +namespace walberla { +namespace pe_coupling { +namespace amr { + +class MetisAssignmentFunctor +{ +public: + + MetisAssignmentFunctor( shared_ptr<InfoCollection>& ic, + const std::function<real_t(const BlockInfo&)> & weightEvaluationFct, + const uint_t ncon = 1) + : ic_(ic), ncon_(ncon), weightEvaluationFct_(ncon, weightEvaluationFct) {} + + MetisAssignmentFunctor( shared_ptr<InfoCollection>& ic, + const std::vector< std::function<real_t(const BlockInfo&)> > & weightEvaluationFct ) + : ic_(ic), ncon_(weightEvaluationFct.size()), weightEvaluationFct_(weightEvaluationFct) {} + + void operator()( std::vector< std::pair< const PhantomBlock *, walberla::any > > & blockData, const PhantomBlockForest & phantomBlockForest); + + inline void setWeightEvaluationFct( const std::function<real_t(const BlockInfo &)> & weightEvaluationFct, const uint_t con = 0 ) { weightEvaluationFct_[con] = weightEvaluationFct; } + inline void setWeightEvaluationFcts( const std::vector< std::function<real_t(const BlockInfo &)> > & weightEvaluationFct) { weightEvaluationFct_ = weightEvaluationFct; } + + uint_t getNcon() const { return ncon_; } + + inline void setBlockBaseWeight( const real_t blockBaseWeight ){ blockBaseWeight_ = blockBaseWeight; } + inline real_t getBlockBaseWeight() const { return blockBaseWeight_; } + +private: + + shared_ptr<InfoCollection> ic_; + uint_t ncon_; + std::vector< std::function<real_t(const BlockInfo&)> > weightEvaluationFct_; + real_t blockBaseWeight_ = real_t(1); +}; + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla + diff --git a/src/pe_coupling/amr/weight_assignment/WeightAssignmentFunctor.cpp b/src/pe_coupling/amr/weight_assignment/WeightAssignmentFunctor.cpp new file mode 100644 index 0000000000000000000000000000000000000000..ce64d8fbe78b507a9b15bfe8811f50dcdb2172ed --- /dev/null +++ b/src/pe_coupling/amr/weight_assignment/WeightAssignmentFunctor.cpp @@ -0,0 +1,47 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file WeightAssignmentFunctor.cpp +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#include "WeightAssignmentFunctor.h" + +namespace walberla { +namespace pe_coupling { +namespace amr { + +void WeightAssignmentFunctor::operator()( std::vector< std::pair< const PhantomBlock *, walberla::any > > & blockData, const PhantomBlockForest & ) +{ + for( auto it = blockData.begin(); it != blockData.end(); ++it ) + { + const PhantomBlock * block = it->first; + //only change of one level is supported! + WALBERLA_ASSERT_LESS( std::abs(int_c(block->getLevel()) - int_c(block->getSourceLevel())), 2 ); + + BlockInfo blockInfo; + pe_coupling::getBlockInfoFromInfoCollection(block, ic_, blockInfo); + + real_t blockWeight = std::max(weightEvaluationFct_(blockInfo), blockBaseWeight_); + + it->second = PhantomBlockWeight( double_c( blockWeight ) ); + + } +} + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/weight_assignment/WeightAssignmentFunctor.h b/src/pe_coupling/amr/weight_assignment/WeightAssignmentFunctor.h new file mode 100644 index 0000000000000000000000000000000000000000..5fa776362e591d01a5286eb7922b179875bfbe06 --- /dev/null +++ b/src/pe_coupling/amr/weight_assignment/WeightAssignmentFunctor.h @@ -0,0 +1,58 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file WeightAssignmentFunctor.h +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#pragma once + +#include "pe_coupling/amr/BlockInfo.h" +#include "pe_coupling/amr/InfoCollection.h" + +#include "blockforest/loadbalancing/PODPhantomData.h" + +#include <functional> + +namespace walberla { +namespace pe_coupling { +namespace amr { + +class WeightAssignmentFunctor +{ +public: + typedef walberla::blockforest::PODPhantomWeight<double> PhantomBlockWeight; + + WeightAssignmentFunctor( const shared_ptr<InfoCollection>& ic, + const std::function<real_t(const BlockInfo&)> & weightEvaluationFct ) : + ic_(ic), weightEvaluationFct_(weightEvaluationFct) {} + + void operator()( std::vector< std::pair< const PhantomBlock *, walberla::any > > & blockData, const PhantomBlockForest & ); + + inline void setWeightEvaluationFct( const std::function<real_t(const BlockInfo &)> & weightEvaluationFct ) { weightEvaluationFct_ = weightEvaluationFct;} + + inline void setBlockBaseWeight( const real_t blockBaseWeight ){blockBaseWeight_ = blockBaseWeight;} + inline real_t getBlockBaseWeight() const { return blockBaseWeight_;} + +private: + shared_ptr<InfoCollection> ic_; + std::function<real_t(const BlockInfo&)> weightEvaluationFct_; + real_t blockBaseWeight_ = real_t(1); +}; + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/weight_assignment/WeightEvaluationFunctions.cpp b/src/pe_coupling/amr/weight_assignment/WeightEvaluationFunctions.cpp new file mode 100644 index 0000000000000000000000000000000000000000..81e1fa0c46f4be2dc3698230fc9cf19f084657cc --- /dev/null +++ b/src/pe_coupling/amr/weight_assignment/WeightEvaluationFunctions.cpp @@ -0,0 +1,35 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file WeightEvaluationFunctions.cpp +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#include "pe_coupling/amr/weight_assignment/WeightEvaluationFunctions.h" + +namespace walberla { +namespace pe_coupling { +namespace amr { + + +real_t defaultWeightEvaluationFunction(const BlockInfo& blockInfo) +{ + return real_c(blockInfo.numberOfCells); +} + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/weight_assignment/WeightEvaluationFunctions.h b/src/pe_coupling/amr/weight_assignment/WeightEvaluationFunctions.h new file mode 100644 index 0000000000000000000000000000000000000000..958151ea39697c6c8f26e6bcc4e5034fa5fc9d5b --- /dev/null +++ b/src/pe_coupling/amr/weight_assignment/WeightEvaluationFunctions.h @@ -0,0 +1,39 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file WeightEvaluationFunctions.h +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +// +//====================================================================================================================== + +#pragma once + +#include "pe_coupling/amr/BlockInfo.h" + +namespace walberla { +namespace pe_coupling { +namespace amr { + + +/* + * Examples of weight evaluation functions, useful for coupled LBM-PE simulations: + * - defaultWeightEvaluationFunction: weight is just the number of cells (thus constant on each block) + */ + +real_t defaultWeightEvaluationFunction(const BlockInfo& blockInfo); + +} // namespace amr +} // namespace pe_coupling +} // namespace walberla diff --git a/src/pe_coupling/amr/weight_assignment/all.h b/src/pe_coupling/amr/weight_assignment/all.h new file mode 100644 index 0000000000000000000000000000000000000000..8f9cf169df8cecc0567c3c89706222b28f3285ab --- /dev/null +++ b/src/pe_coupling/amr/weight_assignment/all.h @@ -0,0 +1,27 @@ +//====================================================================================================================== +// +// This file is part of waLBerla. waLBerla is free software: you can +// redistribute it and/or modify it under the terms of the GNU General Public +// License as published by the Free Software Foundation, either version 3 of +// the License, or (at your option) any later version. +// +// waLBerla is distributed in the hope that it will be useful, but WITHOUT +// ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or +// FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License +// for more details. +// +// You should have received a copy of the GNU General Public License along +// with waLBerla (see COPYING.txt). If not, see <http://www.gnu.org/licenses/>. +// +//! \file all.h +//! \ingroup pe_coupling +//! \author Christoph Rettinger <christoph.rettinger@fau.de> +//! \brief Collective header file for module pe_coupling +// +//====================================================================================================================== + +#pragma once + +#include "MetisAssignmentFunctor.h" +#include "WeightAssignmentFunctor.h" +#include "WeightEvaluationFunctions.h"