Merge branch 'Union' into 'master'

New PE Body: Union

See merge request !4

src/core/math/GenericAABB.h

@@ -179,6 +179,8 @@ public:
    template< typename Engine >
    inline vector_type randomPoint( Engine & engine ) const;
 
+   bool checkInvariant() const;
+
    /**
     * \brief Deserializes an GenericAABB from a mpi::GenericRecvBuffer
     *
@@ -206,8 +208,6 @@ protected:
    inline GenericAABB( const GenericAABB & lhs, const GenericAABB & rhs );
 
 private:
-   bool checkInvariant() const;
-
    vector_type minCorner_; /// minimal values
    vector_type maxCorner_; /// maximal values
 };

src/core/math/Quaternion.h

@@ -277,7 +277,8 @@ inline Quaternion<Type>::Quaternion( Vector3<Axis> axis, Type angle )
 {
    static_assert(boost::is_floating_point<Axis>::value, "Axis has to be floating point!" );
 
-   if( math::equal(std::fabs(angle), real_c(0))  ) {
+   auto axisLength = axis.length();
+   if( (floatIsEqual(axisLength, 0)) || (math::equal(std::fabs(angle), real_c(0)))  ) {
       reset();
       return;
    }
@@ -287,7 +288,8 @@ inline Quaternion<Type>::Quaternion( Vector3<Axis> axis, Type angle )
    const Type sina( std::sin( angle*Type(0.5) ) );
    const Type cosa( std::cos( angle*Type(0.5) ) );
 
-   axis = axis.getNormalized();
+   auto invAxisLength = real_c(1) / axisLength;
+   axis *= invAxisLength;
 
    v_[0] = cosa;
    v_[1] = sina * axis[0];

src/pe/Types.h

@@ -84,6 +84,7 @@ class SliderJoint;
 class Sphere;
 class Spring;
 class TriangleMesh;
+template <typename BodyTypeTuple>
 class Union;
 
 typedef RigidBody             BodyType;            //!< Type of the rigid bodies.
@@ -117,10 +118,6 @@ typedef TriangleMesh          MeshType;             //!< Type of the triangle me
 typedef TriangleMesh*         MeshID;               //!< Handle for a triangle mesh primitive.
 typedef const TriangleMesh*   ConstMeshID;          //!< Handle for a constant triangle mesh primitive.
 
-typedef Union                 UnionType;           //!< Type of the union compound geometry.
-typedef Union*                UnionID;             //!< Handle for a union.
-typedef const Union*          ConstUnionID;        //!< Handle for a constant union.
-
 typedef Attachable            AttachableType;      //!< Type of the attachables.
 typedef Attachable*           AttachableID;        //!< Handle for an attachable.
 typedef const Attachable*     ConstAttachableID;   //!< Handle for a constant attachable.

src/pe/basic.h

@@ -43,6 +43,7 @@
 #include "pe/rigidbody/CapsuleFactory.h"
 #include "pe/rigidbody/PlaneFactory.h"
 #include "pe/rigidbody/SphereFactory.h"
+#include "pe/rigidbody/UnionFactory.h"
 
 #include "pe/synchronization/SyncNextNeighbors.h"
 #include "pe/synchronization/SyncShadowOwners.h"

src/pe/collision/Collide.h

@@ -20,11 +20,14 @@
 //
 //======================================================================================================================
 
+#pragma once
+
 #include "pe/Types.h"
 #include "pe/rigidbody/Box.h"
 #include "pe/rigidbody/Capsule.h"
 #include "pe/rigidbody/Plane.h"
 #include "pe/rigidbody/Sphere.h"
+#include "pe/rigidbody/Union.h"
 #include "pe/contact/Contact.h"
 
 #include "GJKEPAHelper.h"
@@ -34,9 +37,17 @@
 #include "core/math/Shims.h"
 #include "geometry/GeometricalFunctions.h"
 
+#include <boost/tuple/tuple.hpp>
+
 namespace walberla {
 namespace pe {
 
+namespace fcd {
+//Forward Declaration
+template < typename BodyTypeTuple, typename BaseT >
+struct DoubleDispatch;
+}
+
 //*************************************************************************************************
 /*!\brief Contact generation between two colliding rigid bodies.
  * \ingroup contact_generation
@@ -74,9 +85,6 @@ bool collide( SphereID s1, SphereID s2, Container& container )
 {
    WALBERLA_ASSERT_UNEQUAL(s1->getSystemID(), s2->getSystemID(), "colliding with itself!\n" << s1 << "\n" << s2);
 
-   WALBERLA_CHECK_EQUAL(s1->getSuperBody(), s1);
-   WALBERLA_CHECK_EQUAL(s2->getSuperBody(), s2);
-
    Vec3 contactNormal = ( s1->getPosition() - s2->getPosition() );
    real_t penetrationDepth = ( contactNormal.length() - s1->getRadius() - s2->getRadius() );
 
@@ -1994,5 +2002,39 @@ bool collide( CapsuleID c, BoxID b, Container& container )
    return collide(b, c, container);
 }
 
+template <typename BodyTypeTuple, typename BodyB, typename Container>
+inline
+bool collide( Union<BodyTypeTuple>* bd1, BodyB* bd2, Container& container )
+{
+   bool collision = false;
+   for( auto it=bd1->begin(); it!=bd1->end(); ++it )
+   {
+      collision |= fcd::DoubleDispatch< BodyTypeTuple, boost::tuple<BodyB> >::execute(*it, bd2, container);
+   }
+   return collision;
+}
+
+template <typename BodyA, typename BodyTypeTuple, typename Container>
+inline
+bool collide( BodyA* bd1, Union<BodyTypeTuple>* bd2, Container& container )
+{
+   return collide (bd2, bd1, container);
+}
+
+template <typename BodyTypeTupleA, typename BodyTypeTupleB, typename Container>
+inline
+bool collide( Union<BodyTypeTupleA>* bd1, Union<BodyTypeTupleB>* bd2, Container& container )
+{
+   bool collision = false;
+   for( auto it1=bd1->begin(); it1!=bd1->end(); ++it1 )
+   {
+      for( auto it2=bd2->begin(); it2!=bd2->end(); ++it2 )
+      {
+         collision |= fcd::DoubleDispatch< BodyTypeTupleA, BodyTypeTupleB >::execute(*it1, *it2, container);
+      }
+   }
+   return collision;
+}
+
 }
 }

src/pe/communication/DynamicMarshalling.h

@@ -33,6 +33,7 @@
 #include "pe/communication/rigidbody/Box.h"
 #include "pe/communication/rigidbody/Capsule.h"
 #include "pe/communication/rigidbody/Sphere.h"
+#include "pe/communication/rigidbody/Union.h"
 
 #include "core/Abort.h"
 
@@ -95,7 +96,8 @@ struct UnmarshalDynamically{
       if (BodyTypeTuple::head_type::getStaticTypeID() == typeID)
       {
          typedef typename BodyTypeTuple::head_type BodyT;
-         return instantiate<BodyT>( buffer, domain, block );
+         BodyT* newBody;
+         return instantiate( buffer, domain, block, newBody );
       } else
       {
          return UnmarshalDynamically<typename BodyTypeTuple::tail_type>::execute(buffer, typeID, domain, block);

src/pe/communication/Instantiate.h

@@ -28,6 +28,7 @@
 // Includes
 //*************************************************************************************************
 
+#include "core/Abort.h"
 #include "core/math/AABB.h"
 #include "core/math/Vector3.h"
 
@@ -55,9 +56,9 @@ void correctBodyPosition(const math::AABB& domain, const Vec3& center, Vec3& pos
 }
 
 template < class BodyT >
-BodyID instantiate( mpi::RecvBuffer& /*buffer*/, const math::AABB& /*domain*/, const math::AABB& /*block*/ )
+BodyT* instantiate( mpi::RecvBuffer& /*buffer*/, const math::AABB& /*domain*/, const math::AABB& /*block*/, BodyT*& /*newBody*/ )
 {
-   throw std::runtime_error( "Body instantiation not implemented!" );
+   WALBERLA_ABORT( "Body instantiation not implemented! (" << demangle(typeid(BodyT).name()) << ")" );
 }
 
 }  // namespace communication

src/pe/communication/Marshalling.cpp

@@ -47,7 +47,7 @@ void marshal( mpi::SendBuffer& buffer, const MPIRigidBodyTrait& obj ) {
  * \param hasSuperBody False if body is not part of a union. Subordinate bodies in unions do not encode velocities but encode relative positions.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, MPIRigidBodyTraitParameter& objparam, bool /*hasSuperBody*/ ) {
+void unmarshal( mpi::RecvBuffer& buffer, MPIRigidBodyTraitParameter& objparam ) {
    buffer >> objparam.owner_;
 }
 
@@ -66,6 +66,7 @@ void marshal( mpi::SendBuffer& buffer, const RigidBody& obj ) {
    buffer << obj.isCommunicating();
    buffer << obj.hasInfiniteMass();
    buffer << obj.getPosition();
+   buffer << obj.hasSuperBody();
    if( obj.hasSuperBody() )
    {
       buffer <<  obj.getRelPosition();
@@ -87,23 +88,24 @@ void marshal( mpi::SendBuffer& buffer, const RigidBody& obj ) {
  * \param hasSuperBody False if body is not part of a union. Subordinate bodies in unions do not encode velocities but encode relative positions.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, RigidBodyParameters& objparam, bool hasSuperBody ) {
-   unmarshal( buffer, objparam.mpiTrait_, hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, RigidBodyParameters& objparam ) {
+   unmarshal( buffer, objparam.mpiTrait_ );
 
    buffer >> objparam.sid_;
    buffer >> objparam.uid_;
    buffer >> objparam.communicating_;
    buffer >> objparam.infiniteMass_;
    buffer >> objparam.gpos_;
+   buffer >> objparam.hasSuperBody_;
 
-   if( hasSuperBody )
+   if( objparam.hasSuperBody_ )
    {
       buffer >> objparam.rpos_;
    }
 
    buffer >> objparam.q_;
 
-   if( !hasSuperBody )
+   if( !objparam.hasSuperBody_ )
    {
       buffer >> objparam.v_;
       buffer >> objparam.w_;
@@ -139,8 +141,8 @@ void marshal( mpi::SendBuffer& buffer, const GeomPrimitive& obj ) {
  * \param hasSuperBody False if body is not part of a union. Passed on to rigid body unmarshalling.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, GeomPrimitiveParameters& objparam, bool hasSuperBody ) {
-   unmarshal( buffer, static_cast<RigidBodyParameters&>( objparam ), hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, GeomPrimitiveParameters& objparam ) {
+   unmarshal( buffer, static_cast<RigidBodyParameters&>( objparam ) );
    buffer >> objparam.material_;
 }
 

src/pe/communication/Marshalling.h

@@ -63,7 +63,7 @@ void marshal( mpi::SendBuffer& buffer, const MPIRigidBodyTrait& obj );
  * \param hasSuperBody False if body is not part of a union. Subordinate bodies in unions do not encode velocities but encode relative positions.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, MPIRigidBodyTraitParameter& objparam, bool hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, MPIRigidBodyTraitParameter& objparam );
 
 //*************************************************************************************************
 //*************************************************************************************************
@@ -77,6 +77,7 @@ struct RigidBodyParameters {
    bool communicating_, infiniteMass_;
    id_t sid_, uid_;
    Vec3 gpos_, rpos_, v_, w_;
+   bool hasSuperBody_;
    Quat q_;
 };
 
@@ -98,7 +99,7 @@ void marshal( mpi::SendBuffer& buffer, const RigidBody& obj );
  * \param hasSuperBody False if body is not part of a union. Subordinate bodies in unions do not encode velocities but encode relative positions.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, RigidBodyParameters& objparam, bool hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, RigidBodyParameters& objparam );
 
 //*************************************************************************************************
 //*************************************************************************************************
@@ -130,7 +131,7 @@ void marshal( mpi::SendBuffer& buffer, const GeomPrimitive& obj );
  * \param hasSuperBody False if body is not part of a union. Passed on to rigid body unmarshalling.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, GeomPrimitiveParameters& objparam, bool hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, GeomPrimitiveParameters& objparam );
 
 }  // namespace communication
 }  // namespace pe

src/pe/communication/rigidbody/Box.cpp

@@ -49,8 +49,8 @@ void marshal( mpi::SendBuffer& buffer, const Box& obj ) {
  * \param hasSuperBody False if body is not part of a union. Passed on to rigid body unmarshalling.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, BoxParameters& objparam, bool hasSuperBody ) {
-   unmarshal( buffer, static_cast<GeomPrimitiveParameters&>( objparam ), hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, BoxParameters& objparam ) {
+   unmarshal( buffer, static_cast<GeomPrimitiveParameters&>( objparam ) );
    buffer >> objparam.lengths_;
 }
 //*************************************************************************************************

src/pe/communication/rigidbody/Box.h

@@ -58,21 +58,18 @@ void marshal( mpi::SendBuffer& buffer, const Box& obj );
  * \param hasSuperBody False if body is not part of a union. Passed on to rigid body unmarshalling.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, BoxParameters& objparam, bool hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, BoxParameters& objparam );
 //*************************************************************************************************
-
-
-template <>
-inline BodyID instantiate<Box>( mpi::RecvBuffer& buffer, const math::AABB& domain, const math::AABB& block )
+inline BoxID instantiate( mpi::RecvBuffer& buffer, const math::AABB& domain, const math::AABB& block, BoxID& newBody )
 {
    BoxParameters subobjparam;
-   unmarshal( buffer, subobjparam, false );
+   unmarshal( buffer, subobjparam );
    correctBodyPosition(domain, block.center(), subobjparam.gpos_);
-   BoxID obj = new Box( subobjparam.sid_, subobjparam.uid_, subobjparam.gpos_, subobjparam.rpos_, subobjparam.q_, subobjparam.lengths_, subobjparam.material_, false, subobjparam.communicating_, subobjparam.infiniteMass_ );
-   obj->setLinearVel( subobjparam.v_ );
-   obj->setAngularVel( subobjparam.w_ );
-   obj->MPITrait.setOwner( subobjparam.mpiTrait_.owner_ );
-   return obj;
+   newBody = new Box( subobjparam.sid_, subobjparam.uid_, subobjparam.gpos_, subobjparam.rpos_, subobjparam.q_, subobjparam.lengths_, subobjparam.material_, false, subobjparam.communicating_, subobjparam.infiniteMass_ );
+   newBody->setLinearVel( subobjparam.v_ );
+   newBody->setAngularVel( subobjparam.w_ );
+   newBody->MPITrait.setOwner( subobjparam.mpiTrait_.owner_ );
+   return newBody;
 }
 
 }  // namespace communication

src/pe/communication/rigidbody/Capsule.cpp

@@ -51,9 +51,9 @@ void marshal( mpi::SendBuffer& buffer, const Capsule& obj )
  * \param hasSuperBody False if body is not part of a union. Passed on to rigid body unmarshalling.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, CapsuleParameters& objparam, bool hasSuperBody )
+void unmarshal( mpi::RecvBuffer& buffer, CapsuleParameters& objparam )
 {
-   unmarshal( buffer, static_cast<GeomPrimitiveParameters&>( objparam ), hasSuperBody );
+   unmarshal( buffer, static_cast<GeomPrimitiveParameters&>( objparam ) );
    buffer >> objparam.radius_;
    buffer >> objparam.length_;
 }

src/pe/communication/rigidbody/Capsule.h

@@ -58,21 +58,20 @@ void marshal( mpi::SendBuffer& buffer, const Capsule& obj );
  * \param hasSuperBody False if body is not part of a union. Passed on to rigid body unmarshalling.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, CapsuleParameters& objparam, bool hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, CapsuleParameters& objparam );
 //*************************************************************************************************
 
 
-template <>
-inline BodyID instantiate<Capsule>( mpi::RecvBuffer& buffer, const math::AABB& domain, const math::AABB& block )
+inline CapsuleID instantiate( mpi::RecvBuffer& buffer, const math::AABB& domain, const math::AABB& block, CapsuleID& newBody )
 {
    CapsuleParameters subobjparam;
-   unmarshal( buffer, subobjparam, false );
+   unmarshal( buffer, subobjparam );
    correctBodyPosition(domain, block.center(), subobjparam.gpos_);
-   CapsuleID obj = new Capsule( subobjparam.sid_, subobjparam.uid_, subobjparam.gpos_, subobjparam.rpos_, subobjparam.q_, subobjparam.radius_, subobjparam.length_, subobjparam.material_, false, subobjparam.communicating_, subobjparam.infiniteMass_ );
-   obj->setLinearVel( subobjparam.v_ );
-   obj->setAngularVel( subobjparam.w_ );
-   obj->MPITrait.setOwner( subobjparam.mpiTrait_.owner_ );
-   return obj;
+   newBody = new Capsule( subobjparam.sid_, subobjparam.uid_, subobjparam.gpos_, subobjparam.rpos_, subobjparam.q_, subobjparam.radius_, subobjparam.length_, subobjparam.material_, false, subobjparam.communicating_, subobjparam.infiniteMass_ );
+   newBody->setLinearVel( subobjparam.v_ );
+   newBody->setAngularVel( subobjparam.w_ );
+   newBody->MPITrait.setOwner( subobjparam.mpiTrait_.owner_ );
+   return newBody;
 }
 
 }  // namespace communication

src/pe/communication/rigidbody/Sphere.cpp

@@ -49,8 +49,8 @@ void marshal( mpi::SendBuffer& buffer, const Sphere& obj ) {
  * \param hasSuperBody False if body is not part of a union. Passed on to rigid body unmarshalling.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, SphereParameters& objparam, bool hasSuperBody ) {
-   unmarshal( buffer, static_cast<GeomPrimitiveParameters&>( objparam ), hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, SphereParameters& objparam ) {
+   unmarshal( buffer, static_cast<GeomPrimitiveParameters&>( objparam ) );
    buffer >> objparam.radius_;
 }
 //*************************************************************************************************

src/pe/communication/rigidbody/Sphere.h

@@ -58,21 +58,20 @@ void marshal( mpi::SendBuffer& buffer, const Sphere& obj );
  * \param hasSuperBody False if body is not part of a union. Passed on to rigid body unmarshalling.
  * \return void
  */
-void unmarshal( mpi::RecvBuffer& buffer, SphereParameters& objparam, bool hasSuperBody );
+void unmarshal( mpi::RecvBuffer& buffer, SphereParameters& objparam );
 //*************************************************************************************************
 
 
-template <>
-inline BodyID instantiate<Sphere>( mpi::RecvBuffer& buffer, const math::AABB& domain, const math::AABB& block )
+inline SphereID instantiate( mpi::RecvBuffer& buffer, const math::AABB& domain, const math::AABB& block, SphereID& newBody )
 {
    SphereParameters subobjparam;
-   unmarshal( buffer, subobjparam, false );
+   unmarshal( buffer, subobjparam );
    correctBodyPosition(domain, block.center(), subobjparam.gpos_);
-   SphereID obj = new Sphere( subobjparam.sid_, subobjparam.uid_, subobjparam.gpos_, subobjparam.rpos_, subobjparam.q_, subobjparam.radius_, subobjparam.material_, false, subobjparam.communicating_, subobjparam.infiniteMass_ );
-   obj->setLinearVel( subobjparam.v_ );
-   obj->setAngularVel( subobjparam.w_ );
-   obj->MPITrait.setOwner( subobjparam.mpiTrait_.owner_ );
-   return obj;
+   newBody = new Sphere( subobjparam.sid_, subobjparam.uid_, subobjparam.gpos_, subobjparam.rpos_, subobjparam.q_, subobjparam.radius_, subobjparam.material_, false, subobjparam.communicating_, subobjparam.infiniteMass_ );
+   newBody->setLinearVel( subobjparam.v_ );
+   newBody->setAngularVel( subobjparam.w_ );
+   newBody->MPITrait.setOwner( subobjparam.mpiTrait_.owner_ );
+   return newBody;
 }
 
 }  // namespace communication

src/pe/communication/rigidbody/Union.h

+//======================================================================================================================
+//
+//  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 Union.h
+//! \ingroup pe
+//! \author Tobias Preclik
+//! \author Sebastian Eibl <sebastian.eibl@fau.de>
+//! \brief Marshalling of objects for data transmission or storage.
+//
+//======================================================================================================================
+
+#pragma once
+
+
+//*************************************************************************************************
+// Includes
+//*************************************************************************************************
+
+#include "pe/communication/Instantiate.h"
+#include "pe/communication/Marshalling.h"
+#include "pe/rigidbody/Union.h"
+
+namespace walberla {
+namespace pe {
+namespace communication {
+
+//forward declaration
+template < typename BodyTypeTuple >
+struct MarshalDynamically;
+
+//forward declaration
+template < typename BodyTypeTuple >
+struct UnmarshalDynamically;
+
+struct UnionParameters : public GeomPrimitiveParameters {
+   real_t         m_;
+   Mat3           I_;
+   math::AABB     aabb_;
+   size_t         size_;
+};
+
+//*************************************************************************************************
+/*!\brief Marshalling a box primitive.
+ *
+ * \param buffer The buffer to be filled.
+ * \param obj The object to be marshalled.
+ * \return void
+ */
+template < typename BodyTypeTuple >
+void marshal( mpi::SendBuffer& buffer, const Union<BodyTypeTuple>& obj )
+{
+   // Material of union is not relevant for reconstruction thus marshal RigidBody instead of GeomPrimitive.
+   marshal( buffer, static_cast<const RigidBody&>( obj ) );
+
+   buffer << obj.getMass();               // Encoding the total mass
+   buffer << obj.getBodyInertia();        // Encoding the moment of inertia
+
+   // Preparing the axis-aligned bounding box
+   buffer << obj.getAABB();               // Encoding the axis-aligned bounding box
+
+   buffer << static_cast<size_t> (obj.size());                  // Encoding the number of contained bodies
+
+   // Encoding the contained primitives
+   const typename Union<BodyTypeTuple>::ConstIterator begin( obj.begin() );
+   const typename Union<BodyTypeTuple>::ConstIterator end  ( obj.end()   );
+   for(  typename Union<BodyTypeTuple>::ConstIterator body=begin; body!=end; ++body )
+   {
+      buffer << body->getTypeID();
+      MarshalDynamically<BodyTypeTuple>::execute( buffer, **body );
+   }
+}
+
+//*************************************************************************************************
+
+//*************************************************************************************************
+/*!\brief Unmarshalling a box primitive.
+ *
+ * \param buffer The buffer from where to read.
+ * \param objparam The object to be reconstructed.
+ * \param hasSuperBody False if body is not part of a union. Passed on to rigid body unmarshalling.
+ * \return void
+ */
+inline
+void unmarshal( mpi::RecvBuffer& buffer, UnionParameters& objparam )
+{
+   // Material of union is not relevant for reconstruction thus marshal RigidBody instead of GeomPrimitive.
+   unmarshal( buffer, static_cast<RigidBodyParameters&>( objparam ) );
+
+   buffer >> objparam.m_;
+   buffer >> objparam.I_;
+   buffer >> objparam.aabb_;
+   //std::cout << "mass of union: " << objparam.m_ << std::endl;
+
+   // Decode the number of contained bodies
+   buffer >> objparam.size_;
+   //std::cout << "size of union: " << size << std::endl;
+}
+//*************************************************************************************************
+
+
+template <typename BodyTypeTuple>
+inline Union<BodyTypeTuple>* instantiate( mpi::RecvBuffer& buffer, const math::AABB& domain, const math::AABB& block, Union<BodyTypeTuple>*& newBody )
+{
+   UnionParameters subobjparam;
+   unmarshal( buffer, subobjparam );
+   correctBodyPosition(domain, block.center(), subobjparam.gpos_);
+   newBody = new Union<BodyTypeTuple>( subobjparam.sid_,
+                                       subobjparam.uid_,
+                                       subobjparam.gpos_,
+                                       subobjparam.rpos_,
+                                       subobjparam.q_,
+                                       false,
+                                       subobjparam.communicating_,
+                                       subobjparam.infiniteMass_ );
+   newBody->setLinearVel( subobjparam.v_ );
+   newBody->setAngularVel( subobjparam.w_ );
+   newBody->MPITrait.setOwner( subobjparam.mpiTrait_.owner_ );
+
+   // Decoding the contained primitives
+   for( size_t i = 0; i < subobjparam.size_; ++i )
+   {
+      decltype ( static_cast<BodyID>(nullptr)->getTypeID() ) type;
+      buffer >> type;
+      BodyID obj = UnmarshalDynamically<BodyTypeTuple>::execute(buffer, type, domain, block);
+      obj->setRemote( true );
+      newBody->add(obj);
+   }
+
+   return newBody;
+}
+
+}  // namespace communication
+}  // namespace pe
+}  // namespace walberla

src/pe/contact/Contact.cpp

@@ -64,8 +64,8 @@ Contact::Contact( GeomID g1, GeomID g2, const Vec3& gpos, const Vec3& normal, re
    WALBERLA_ASSERT( !( b1_->hasInfiniteMass() && b2_->hasInfiniteMass() ), "Invalid contact between two rigid bodies with infinite masses" );
 
    // Registering the contact with both attached rigid bodies
-   if( !b1_->hasInfiniteMass() ) b1_->registerContact( this );
-   if( !b2_->hasInfiniteMass() ) b2_->registerContact( this );
+//   if( !b1_->hasInfiniteMass() ) b1_->registerContact( this );
+//   if( !b2_->hasInfiniteMass() ) b2_->registerContact( this );
 
    // Merging the contact graph
 //   if( !b1_->hasInfiniteMass() && !b2_->hasInfiniteMass() )

src/pe/cr/HCSITS.impl.h

@@ -1774,12 +1774,8 @@ inline void HardContactSemiImplicitTimesteppingSolvers::integratePositions( Body
    body->resetNode();
    body->clearContacts();
 
-   if( body->isAwake() ) {
-      if (body->getSystemID() == 6917529027641081896)
-      {
-         //WALBERLA_LOG_DEVEL("pos: " << body->getPosition() << "\nvel: " << body->getLinearVel() << "\nv: " << v);
-      }
-
+   if( body->isAwake() )
+   {
       // Calculating the translational displacement
       body->setPosition( body->getPosition() + v * dt );
 
@@ -1794,14 +1790,6 @@ inline void HardContactSemiImplicitTimesteppingSolvers::integratePositions( Body
       body->setLinearVel( v );
       body->setAngularVel( w );
 
-      //      if( body->getType() == unionType ) {
-      //         // Updating the contained bodies
-      //         UnionID u( static_cast<UnionID>( body ) );
-
-      //         for( BodyIterator subbody=u->begin(); subbody!=u->end(); ++subbody )
-      //            u->updateBody( *subbody, dq );
-      //      }
-
       // Setting the axis-aligned bounding box
       body->calcBoundingBox();