Commit 2cf3f49f authored by Christoph Rettinger's avatar Christoph Rettinger
Browse files

Merge branch 'master' into 'geometry_utils'

# Conflicts:
#   tests/pe_coupling/CMakeLists.txt
parents 9c7f1f48 2d079319
......@@ -111,7 +111,7 @@ void TimingNode<TP>::swap(TimingNode<TP>& tt)
}
}
/// Finds the spezified timer in the timing hierarchy
/// Finds the specified timer in the timing hierarchy
/// \param name timer name which may include more than one hierarchy separated by "."
/// \code findTimer(tn, "firstLevel.secondLevel.thirdLevel.timerName"); \endcode
/// \relates TimingNode
......@@ -130,6 +130,31 @@ const Timer<TP>& findTimer( const TimingNode<TP>& tn, const std::string& name)
}
}
/// Checks if the specified timer exists in the timing hierarchy
/// \param name timer name which may include more than one hierarchy separated by "."
/// \code timerExists(tn, "firstLevel.secondLevel.thirdLevel.timerName"); \endcode
/// \relates TimingNode
template< typename TP > // Timing policy
bool timerExists( const TimingNode<TP>& tn, const std::string& name )
{
auto pos = name.find_first_of(".");
if (pos != std::string::npos)
{
if( tn.tree_.find(name.substr(0, pos)) != tn.tree_.end() )
{
return timerExists( tn.tree_.at(name.substr(0, pos)), name.substr(pos+1, std::string::npos));
}
else
{
return false;
}
}
else
{
return tn.tree_.find(name) != tn.tree_.end();
}
}
/// Resets the timer in the TimingNode structure and all sub timers
/// \relates TimingNode
template< typename TP > // Timing policy
......
......@@ -59,13 +59,13 @@ public:
void swap(TimingTree<TP>& tt);
/// Starts a timer beneath the current hirarchy level
/// Starts a timer beneath the current hierarchy level
void start(const std::string& name);
/// Stops the last started timer and jumps back one hirarchy level
/// Stops the last started timer and jumps back one hierarchy level
void stop(const std::string& name);
/// Checks if specified timer is currently running.
bool isTimerRunning(const std::string& name) const;
/// Resets the the timing hirarchy
/// Resets the the timing hierarchy
void reset();
//** Reduction ******************************************************************************************************
......@@ -81,7 +81,9 @@ public:
/// Returns the raw tree data structure
const TimingNode<TP>& getRawData() const;
const Timer<TP>& operator[](const std::string& name) const;
inline bool timerExists ( const std::string & n ) const;
/// Returns the name of the currently running timer
/// Might be expensive due to value search.
......@@ -89,7 +91,7 @@ public:
private:
/// Tree data structure
TimingNode<TP> root_;
/// Pointer to the current hirarchy level.
/// Pointer to the current hierarchy level.
TimingNode<TP>* current_;
};
......@@ -204,7 +206,7 @@ const TimingNode<TP>& TimingTree<TP>::getRawData() const
return root_;
}
/// Finds the spezified timer in the timing hierarchy
/// Finds the specified timer in the timing hierarchy
/// \param name timer name which may include more than one hierarchy separated by "."
/// \code tt["firstLevel.secondLevel.thirdLevel.timerName"].total(); \endcode
template< typename TP > // Timing policy
......@@ -213,6 +215,15 @@ const Timer<TP>& TimingTree<TP>::operator[](const std::string& name) const
return findTimer(root_, name);
}
/// Checks if the specified timer exists in the timing hierarchy
/// \param name timer name which may include more than one hierarchy separated by "."
/// \code tt.timerExists("firstLevel.secondLevel.thirdLevel.timerName"); \endcode
template< typename TP > // Timing policy
bool TimingTree<TP>::timerExists(const std::string& name) const
{
return walberla::timing::timerExists(root_, name);
}
template< typename TP > // Timing policy
std::string TimingTree<TP>::getCurrentTimerName() const
{
......
......@@ -242,7 +242,7 @@ namespace field {
virtual T * allocateMemory ( uint_t size )
{
void * ptr = aligned_malloc_with_offset(size * sizeof(T), alignment, offset_ % alignment );
void * ptr = aligned_malloc_with_offset(size * sizeof(T) + alignment, alignment, offset_ % alignment );
if(!ptr)
throw std::bad_alloc();
......
......@@ -490,6 +490,65 @@ namespace internal {
}
//===================================================================================================================
//
// Aligned Allocation
//
//===================================================================================================================
template<typename T>
shared_ptr<field::FieldAllocator<T> > getAllocator(uint_t alignment)
{
if( alignment == 0 )
return shared_ptr<field::FieldAllocator<T> >(); // leave to default - auto-detection of alignment
else if ( alignment == 16 )
return make_shared< field::AllocateAligned<T, 16> >();
else if ( alignment == 32 )
return make_shared< field::AllocateAligned<T, 32> >();
else if ( alignment == 64 )
return make_shared< field::AllocateAligned<T, 64> >();
else if ( alignment == 128 )
return make_shared< field::AllocateAligned<T, 128> >();
else {
PyErr_SetString( PyExc_ValueError, "Alignment parameter has to be one of 0, 16, 32, 64, 128." );
throw boost::python::error_already_set();
return shared_ptr<field::FieldAllocator<T> >();
}
}
template< typename GhostLayerField_T >
class GhostLayerFieldDataHandling : public blockforest::AlwaysInitializeBlockDataHandling< GhostLayerField_T >
{
public:
typedef typename GhostLayerField_T::value_type Value_T;
GhostLayerFieldDataHandling( const weak_ptr<StructuredBlockStorage> &blocks, const uint_t nrOfGhostLayers,
const Value_T &initValue, const Layout layout, uint_t alignment = 0 ) :
blocks_( blocks ), nrOfGhostLayers_( nrOfGhostLayers ), initValue_( initValue ), layout_( layout ),
alignment_( alignment ) {}
GhostLayerField_T * initialize( IBlock * const block )
{
auto blocks = blocks_.lock();
WALBERLA_CHECK_NOT_NULLPTR( blocks, "Trying to access 'AlwaysInitializeBlockDataHandling' for a block "
"storage object that doesn't exist anymore" );
GhostLayerField_T * field = new GhostLayerField_T ( blocks->getNumberOfXCells( *block ),
blocks->getNumberOfYCells( *block ),
blocks->getNumberOfZCells( *block ),
nrOfGhostLayers_, initValue_, layout_,
getAllocator<Value_T>(alignment_) );
return field;
}
private:
weak_ptr< StructuredBlockStorage > blocks_;
uint_t nrOfGhostLayers_;
Value_T initValue_;
Layout layout_;
uint_t alignment_;
};
//===================================================================================================================
//
......@@ -843,10 +902,10 @@ namespace internal {
{
public:
CreateFieldExporter( uint_t xs, uint_t ys, uint_t zs, uint_t fs, uint_t gl,
Layout layout, const boost::python::object & type,
const shared_ptr<boost::python::object> & resultPointer )
Layout layout, const boost::python::object & type, uint_t alignment,
const shared_ptr<boost::python::object> & resultPointer )
: xs_( xs ), ys_(ys), zs_(zs), fs_(fs), gl_(gl),
layout_( layout), type_( type ), resultPointer_( resultPointer )
layout_( layout), type_( type ), alignment_(alignment), resultPointer_( resultPointer )
{}
template< typename FieldType>
......@@ -862,7 +921,8 @@ namespace internal {
if( python_coupling::isCppEqualToPythonType<T>( (PyTypeObject *)type_.ptr() ) )
{
T initVal = T(); //extract<T> ( initValue_ );
*resultPointer_ = object( make_shared< GhostLayerField<T,F_SIZE> >( xs_,ys_,zs_, gl_, initVal, layout_ ) );
*resultPointer_ = object( make_shared< GhostLayerField<T,F_SIZE> >( xs_,ys_,zs_, gl_, initVal, layout_,
getAllocator<T>(alignment_)));
}
}
......@@ -874,6 +934,7 @@ namespace internal {
uint_t gl_;
Layout layout_;
boost::python::object type_;
uint_t alignment_;
shared_ptr<boost::python::object> resultPointer_;
};
......@@ -881,7 +942,8 @@ namespace internal {
boost::python::object createPythonField( boost::python::list size,
boost::python::object type,
uint_t ghostLayers,
Layout layout)
Layout layout,
uint_t alignment)
{
using namespace boost::python;
uint_t xSize = extract<uint_t> ( size[0] );
......@@ -898,7 +960,7 @@ namespace internal {
}
auto result = make_shared<boost::python::object>();
CreateFieldExporter exporter( xSize,ySize, zSize, fSize, ghostLayers, layout, type, result );
CreateFieldExporter exporter( xSize,ySize, zSize, fSize, ghostLayers, layout, type, alignment, result );
python_coupling::for_each_noncopyable_type< FieldTypes > ( exporter );
if ( *result == object() )
......@@ -952,12 +1014,13 @@ namespace internal {
class AddToStorageExporter
{
public:
AddToStorageExporter( const shared_ptr<StructuredBlockStorage> & blocks,
AddToStorageExporter(const shared_ptr<StructuredBlockStorage> & blocks,
const std::string & name, uint_t fs, uint_t gl, Layout layout,
const boost::python::object & type,
const boost::python::object & initObj )
const boost::python::object & initObj,
uint_t alignment )
: blocks_( blocks ), name_( name ), fs_( fs ),
gl_(gl),layout_( layout), type_( type ), initObj_( initObj), found_(false)
gl_(gl),layout_( layout), type_( type ), initObj_( initObj), alignment_(alignment), found_(false)
{}
template< typename FieldType>
......@@ -972,11 +1035,15 @@ namespace internal {
if( !found_ && python_coupling::isCppEqualToPythonType<T>( (PyTypeObject *)type_.ptr() ) )
{
if ( initObj_ == object() )
field::addToStorage< GhostLayerField<T,F_SIZE> >( blocks_, name_, T(), layout_, gl_ );
else
field::addToStorage< GhostLayerField<T,F_SIZE> >( blocks_, name_, extract<T>(initObj_), layout_, gl_ );
typedef internal::GhostLayerFieldDataHandling< GhostLayerField<T,F_SIZE > > DataHandling;
if ( initObj_ == object() ) {
auto dataHandling = walberla::make_shared< DataHandling >( blocks_, gl_, T(), layout_, alignment_ );
blocks_->addBlockData( dataHandling, name_ );
}
else {
auto dataHandling = walberla::make_shared< DataHandling >( blocks_, gl_, extract<T>(initObj_), layout_, alignment_ );
blocks_->addBlockData( dataHandling, name_ );
}
found_ = true;
}
}
......@@ -990,12 +1057,14 @@ namespace internal {
Layout layout_;
boost::python::object type_;
boost::python::object initObj_;
uint_t alignment_;
bool found_;
};
template<typename FieldTypes>
void addToStorage( const shared_ptr<StructuredBlockStorage> & blocks, const std::string & name,
boost::python::object type, uint_t fs, uint_t gl, Layout layout, boost::python::object initValue )
boost::python::object type, uint_t fs, uint_t gl, Layout layout, boost::python::object initValue,
uint_t alignment)
{
using namespace boost::python;
......@@ -1005,7 +1074,7 @@ namespace internal {
}
auto result = make_shared<boost::python::object>();
AddToStorageExporter exporter( blocks, name, fs, gl, layout, type, initValue );
AddToStorageExporter exporter( blocks, name, fs, gl, layout, type, initValue, alignment );
python_coupling::for_each_noncopyable_type< FieldTypes > ( std::ref(exporter) );
if ( ! exporter.successful() ) {
......@@ -1252,7 +1321,8 @@ void exportFields()
def( "createField", &internal::createPythonField<FieldTypes>, ( ( arg("size") ),
( arg("type") ),
( arg("ghostLayers") = uint_t(1) ),
( arg("layout") = zyxf ) ) );
( arg("layout") = zyxf ),
( arg("alignment") = 0 )) );
def( "createFlagField", &internal::createPythonFlagField, ( ( arg("size") ),
( arg("nrOfBits") = uint_t(32) ),
......@@ -1264,7 +1334,8 @@ void exportFields()
( arg("fSize") = 1 ),
( arg("ghostLayers") = uint_t(1) ),
( arg("layout") = zyxf ),
( arg("initValue") = object() ) ) );
( arg("initValue") = object() ),
( arg("alignment") = 0 ) ) );
def( "addFlagFieldToStorage",&internal::addFlagFieldToStorage, ( ( arg("blocks") ),
( arg("name") ),
......
......@@ -21,13 +21,13 @@
#pragma once
#include "blockforest/StructuredBlockForest.h"
#include "core/math/Vector3.h"
#include "domain_decomposition/StructuredBlockStorage.h"
#include "pe/rigidbody/BodyIterators.h"
#include "pe/synchronization/SyncForces.h"
#include <map>
#include <vector>
#include <array>
namespace walberla {
namespace pe_coupling {
......@@ -36,9 +36,14 @@ class BodiesForceTorqueContainer
{
public:
BodiesForceTorqueContainer( const shared_ptr<StructuredBlockStorage> & blockStorage, const BlockDataID & bodyStorageID )
: blockStorage_( blockStorage ), bodyStorageID_( bodyStorageID )
{ }
typedef std::map< walberla::id_t, std::array<real_t,6> > ForceTorqueStorage_T;
BodiesForceTorqueContainer( const shared_ptr<StructuredBlockForest> & blockForest, const BlockDataID & bodyStorageID )
: blockForest_( blockForest ), bodyStorageID_( bodyStorageID )
{
// has to be added to the forest (not the storage) to register correctly
bodyForceTorqueStorageID_ = blockForest->addBlockData(make_shared<blockforest::AlwaysCreateBlockDataHandling<ForceTorqueStorage_T> >(), "BodiesForceTorqueContainer");
}
void operator()()
{
......@@ -50,38 +55,18 @@ public:
// clear map
clear();
// sum up all forces/torques from shadow copies on local body (owner)
pe::reduceForces( blockStorage_->getBlockStorage(), bodyStorageID_ );
// send total forces/torques to shadow owners
pe::distributeForces( blockStorage_->getBlockStorage(), bodyStorageID_ );
// (re-)build map
for( auto blockIt = blockStorage_->begin(); blockIt != blockStorage_->end(); ++blockIt )
for( auto blockIt = blockForest_->begin(); blockIt != blockForest_->end(); ++blockIt )
{
auto bodyForceTorqueStorage = blockIt->getData<ForceTorqueStorage_T>(bodyForceTorqueStorageID_);
for( auto bodyIt = pe::BodyIterator::begin(*blockIt, bodyStorageID_); bodyIt != pe::BodyIterator::end(); ++bodyIt )
{
auto & f = bodyForceTorqueMap_[ bodyIt->getSystemID() ];
auto & f = (*bodyForceTorqueStorage)[ bodyIt->getSystemID() ];
// only add if body has not been added already before (from another block)
if( f.empty() )
{
const auto & force = bodyIt->getForce();
f.push_back( force[0] );
f.push_back( force[1] );
f.push_back( force[2] );
const auto & torque = bodyIt->getTorque();
f.push_back( torque[0] );
f.push_back( torque[1] );
f.push_back( torque[2] );
}
// reset of force/torque on remote bodies necessary to erase the multiple occurrences of forces/torques on bodies
// (due to call to distributeForces() before)
if ( bodyIt->isRemote() )
{
bodyIt->resetForceAndTorque();
}
const auto & force = bodyIt->getForce();
const auto & torque = bodyIt->getTorque();
f = {{force[0], force[1], force[2], torque[0], torque[1], torque[2] }};
}
}
......@@ -89,42 +74,52 @@ public:
void setOnBodies()
{
// owning process sets the force/torque on the bodies
for( auto blockIt = blockStorage_->begin(); blockIt != blockStorage_->end(); ++blockIt )
// set the force/torque stored in the block-local map onto all bodies
for( auto blockIt = blockForest_->begin(); blockIt != blockForest_->end(); ++blockIt )
{
for( auto bodyIt = pe::LocalBodyIterator::begin(*blockIt, bodyStorageID_); bodyIt != pe::LocalBodyIterator::end(); ++bodyIt )
auto bodyForceTorqueStorage = blockIt->getData<ForceTorqueStorage_T>(bodyForceTorqueStorageID_);
for( auto bodyIt = pe::BodyIterator::begin(*blockIt, bodyStorageID_); bodyIt != pe::BodyIterator::end(); ++bodyIt )
{
const auto &f = bodyForceTorqueMap_[bodyIt->getSystemID()];
WALBERLA_ASSERT( !f.empty(), "When attempting to set force/torque on local body " << bodyIt->getSystemID() << " at position " << bodyIt->getPosition() << ", body was not found in map!");
bodyIt->addForce ( f[0], f[1], f[2] );
bodyIt->addTorque( f[3], f[4], f[5] );
const auto f = bodyForceTorqueStorage->find( bodyIt->getSystemID() );
if( f != bodyForceTorqueStorage->end() )
{
const auto & ftValues = f->second;
bodyIt->addForce ( ftValues[0], ftValues[1], ftValues[2] );
bodyIt->addTorque( ftValues[3], ftValues[4], ftValues[5] );
}
// else: new body has arrived that was not known before
}
}
}
void clear()
{
bodyForceTorqueMap_.clear();
for( auto blockIt = blockForest_->begin(); blockIt != blockForest_->end(); ++blockIt )
{
auto bodyForceTorqueStorage = blockIt->getData<ForceTorqueStorage_T>(bodyForceTorqueStorageID_);
bodyForceTorqueStorage->clear();
}
}
void swap( BodiesForceTorqueContainer & other )
{
std::swap( this->bodyForceTorqueMap_, other.bodyForceTorqueMap_ );
std::swap( bodyForceTorqueStorageID_, other.bodyForceTorqueStorageID_);
}
private:
shared_ptr<StructuredBlockStorage> blockStorage_;
shared_ptr<StructuredBlockStorage> blockForest_;
const BlockDataID bodyStorageID_;
std::map< walberla::id_t, std::vector<real_t> > bodyForceTorqueMap_;
BlockDataID bodyForceTorqueStorageID_;
};
class BodyContainerSwapper
{
public:
BodyContainerSwapper( const shared_ptr<BodiesForceTorqueContainer> & cont1,
const shared_ptr<BodiesForceTorqueContainer> & cont2 )
BodyContainerSwapper( const shared_ptr<BodiesForceTorqueContainer> & cont1, const shared_ptr<BodiesForceTorqueContainer> & cont2 )
: cont1_( cont1 ), cont2_( cont2 )
{ }
......
//======================================================================================================================
//
// 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 BodySelectorFunctions.cpp
//! \ingroup pe_coupling
//! \author Christoph Rettinger <christoph.rettinger@fau.de>
//
//======================================================================================================================
#include "BodySelectorFunctions.h"
#include "pe/rigidbody/RigidBody.h"
namespace walberla {
namespace pe_coupling {
bool selectAllBodies(pe::BodyID /*bodyID*/)
{
return true;
}
bool selectRegularBodies(pe::BodyID bodyID)
{
return !bodyID->hasInfiniteMass() && !bodyID->isGlobal();
}
bool selectFixedBodies(pe::BodyID bodyID)
{
return bodyID->hasInfiniteMass() && !bodyID->isGlobal();
}
bool selectGlobalBodies(pe::BodyID bodyID)
{
return bodyID->isGlobal();
}
} // namespace pe_coupling
} // namespace walberla
......@@ -21,29 +21,18 @@
#pragma once
#include "pe/Types.h"
#include "pe/rigidbody/RigidBody.h"
namespace walberla {
namespace pe_coupling {
bool selectAllBodies(pe::BodyID /*bodyID*/)
{
return true;
}
bool selectRegularBodies(pe::BodyID bodyID)
{
return !bodyID->hasInfiniteMass() && !bodyID->isGlobal();
}
bool selectFixedBodies(pe::BodyID bodyID)
{
return bodyID->hasInfiniteMass() && !bodyID->isGlobal();
}
bool selectGlobalBodies(pe::BodyID bodyID)
{
return bodyID->isGlobal();
}
bool selectAllBodies(pe::BodyID /*bodyID*/);
bool selectRegularBodies(pe::BodyID bodyID);
bool selectFixedBodies(pe::BodyID bodyID);
bool selectGlobalBodies(pe::BodyID bodyID);
} // namespace pe_coupling
} // namespace walberla
......@@ -35,7 +35,7 @@
#include <functional>
#include <map>
#include <array>
namespace walberla {
namespace pe_coupling {
......@@ -88,43 +88,25 @@ public:
// Since they are reset after the call to collision response, they have to be stored explicitly before.
// Then they are set again after each intermediate step.
// sum up all forces/torques from shadow copies on local body (owner)
pe::reduceForces( blockStorage_->getBlockStorage(), bodyStorageID_ );
// send total forces/torques to shadow owners
pe::distributeForces( blockStorage_->getBlockStorage(), bodyStorageID_ );
// generate map from all known bodies (process local) to total forces/torques
std::map< walberla::id_t, std::vector< real_t > > forceMap;
// this has to be done on a block-local basis, since the same body could reside on several blocks from this process
typedef domain_decomposition::IBlockID::IDType BlockID_T;
std::map< BlockID_T, std::map< walberla::id_t, std::array< real_t, 6 > > > forceTorqueMap;
for( auto blockIt = blockStorage_->begin(); blockIt != blockStorage_->end(); ++blockIt )
{
BlockID_T blockID = blockIt->getId().getID();
auto& blockLocalForceTorqueMap = forceTorqueMap[blockID];
// iterate over local and remote bodies and store force/torque in map
// Remote bodies are required since during the then following collision response time steps,
// bodies can change ownership (i.e. a now remote body could become a local body ).
// Since only the owning process re-sets the force/torque later, remote body values have to be stored as well.
for( auto bodyIt = pe::BodyIterator::begin(*blockIt, bodyStorageID_); bodyIt != pe::BodyIterator::end(); ++bodyIt )
{
auto & f = forceMap[ bodyIt->getSystemID() ];
auto & f = blockLocalForceTorqueMap[ bodyIt->getSystemID() ];
// only add if body has not been added already before (from another block)
if( f.empty() )
{
const auto & force = bodyIt->getForce();
f.push_back( force[0] );
f.push_back( force[1] );
f.push_back( force[2] );
const auto & torque = bodyIt->getTorque();
f.push_back( torque[0] );
f.push_back( torque[1] );
f.push_back( torque[2] );
}
const auto & force = bodyIt->getForce();
const auto & torque = bodyIt->getTorque();