diff --git a/src/core/math/Matrix3.h b/src/core/math/Matrix3.h
index 42a36dc9e2e95727dc0bb66117124c33162c31f8..27a29b9450d6c9d95630440b7ffc346d6e0fa630 100644
--- a/src/core/math/Matrix3.h
+++ b/src/core/math/Matrix3.h
@@ -1611,6 +1611,23 @@ inline const Matrix3<Type> skewSymCrossProduct( const Vector3<Type>& vec, const
//*************************************************************************************************
+//*************************************************************************************************
+/*!\brief Dyadic product of two vectors (\f$ M = u \otimes v \f$).
+ *
+ * \param vec1 The first vector argument.
+ * \param vec2 The second vector argument.
+ * \return The matrix \f$ u \otimes v \f$.
+ */
+template< typename Type >
+inline const Matrix3<Type> dyadicProduct( const Vector3<Type>& vec1, const Vector3<Type>& vec2 )
+{
+ return Matrix3<Type>( vec1[0] * vec2[0], vec1[0] * vec2[1], vec1[0] * vec2[2],
+ vec1[1] * vec2[0], vec1[1] * vec2[1], vec1[1] * vec2[2],
+ vec1[2] * vec2[0], vec1[2] * vec2[1], vec1[2] * vec2[2]);
+}
+//*************************************************************************************************
+
+
//**********************************************************************************************************************
/*!\fn std::ostream& operator<<( std::ostream& os, const Matrix3<Type>& m )
// \brief Global output operator for 3x3 matrices.
diff --git a/src/pe/vtk/EllipsoidVtkOutput.cpp b/src/pe/vtk/EllipsoidVtkOutput.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..ec43aa33d6c643517b730e09d4118bc2fc289e6f
--- /dev/null
+++ b/src/pe/vtk/EllipsoidVtkOutput.cpp
@@ -0,0 +1,91 @@
+//======================================================================================================================
+//
+// 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 EllipsoidVtkOutput.cpp
+//! \author Christoph Rettinger <christoph.rettinger@fau.de>
+//
+//======================================================================================================================
+
+#include "core/DataTypes.h"
+#include "EllipsoidVtkOutput.h"
+#include "pe/rigidbody/BodyStorage.h"
+
+#include "core/selectable/IsSetSelected.h"
+#include "core/uid/GlobalState.h"
+
+
+namespace walberla {
+namespace pe {
+
+std::vector< EllipsoidVtkOutput::Attributes > EllipsoidVtkOutput::getAttributes() const
+{
+ std::vector< Attributes > attributes;
+ attributes.push_back( Attributes( vtk::typeToString< float >(), "mass", uint_c(1) ) );
+ attributes.push_back( Attributes( vtk::typeToString< float >(), "tensorGlyph", uint_c(9) ) );
+ attributes.push_back( Attributes( vtk::typeToString< float >(), "velocity", uint_c(3) ) );
+ attributes.push_back( Attributes( vtk::typeToString< int >(), "rank", uint_c(1) ) );
+ attributes.push_back( Attributes( vtk::typeToString< id_t >(), "id", uint_c(1) ) );
+ attributes.push_back( Attributes( vtk::typeToString< id_t >(), "uid", uint_c(1) ) );
+
+ return attributes;
+}
+
+void EllipsoidVtkOutput::configure()
+{
+ bodies_.clear();
+ tensorGlyphs_.clear();
+
+ for( auto blockIt = blockStorage_.begin(); blockIt != blockStorage_.end(); ++blockIt )
+ {
+
+ const Storage& bs = *(blockIt->getData<const Storage>( storageID_ ));
+
+ for( auto it = bs[0].begin(); it != bs[0].end(); ++it )
+ {
+ if (it->getTypeID() == Ellipsoid::getStaticTypeID())
+ {
+ auto ellipsoid = static_cast<ConstEllipsoidID> (*it);
+ bodies_.push_back(ellipsoid);
+
+ // compute tensor glyph for visualization with ParaView (tensorGlyph)
+ Mat3 rotMat = ellipsoid->getRotation();
+ Vector3<real_t> directionVectorX(rotMat[0], rotMat[3], rotMat[6]);
+ Vector3<real_t> directionVectorY(rotMat[1], rotMat[4], rotMat[7]);
+ Vector3<real_t> directionVectorZ(rotMat[2], rotMat[5], rotMat[8]);
+ Vector3<real_t> semiAxes = ellipsoid->getSemiAxes();
+ Mat3 axa = math::dyadicProduct(directionVectorX, directionVectorX);
+ Mat3 bxb = math::dyadicProduct(directionVectorY, directionVectorY);
+ Mat3 cxc = math::dyadicProduct(directionVectorZ, directionVectorZ);
+ Mat3 tensor = axa * semiAxes[0] + bxb * semiAxes[1] + cxc * semiAxes[2];
+ tensorGlyphs_.push_back(tensor);
+ }
+ }
+ }
+}
+
+std::vector< math::Vector3< real_t > > EllipsoidVtkOutput::getPoints()
+{
+ std::vector< math::Vector3< real_t > > result( bodies_.size() );
+ auto resultIt = result.begin();
+ for( auto it = bodies_.begin(); it != bodies_.end(); ++it, ++resultIt )
+ {
+ *resultIt = (*it)->getPosition();
+ }
+ return result;
+}
+
+} // namespace pe
+} // namespace walberla
+
diff --git a/src/pe/vtk/EllipsoidVtkOutput.h b/src/pe/vtk/EllipsoidVtkOutput.h
new file mode 100644
index 0000000000000000000000000000000000000000..43bcec3e34d1f98fb0199acd73bf96e4520d1630
--- /dev/null
+++ b/src/pe/vtk/EllipsoidVtkOutput.h
@@ -0,0 +1,130 @@
+//======================================================================================================================
+//
+// 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 EllipsoidVtkOutput.h
+//! \author Christoph Rettinger <christoph.rettinger@fau.de>
+//
+//======================================================================================================================
+
+#pragma once
+
+#include "pe/rigidbody/RigidBody.h"
+#include "pe/rigidbody/Ellipsoid.h"
+
+#include "core/DataTypes.h"
+#include "core/Set.h"
+#include "core/uid/SUID.h"
+
+#include "domain_decomposition/BlockStorage.h"
+
+#include "vtk/Base64Writer.h"
+#include "vtk/PointDataSource.h"
+#include "vtk/UtilityFunctions.h"
+
+#include <vector>
+
+#include <boost/array.hpp>
+
+namespace walberla {
+namespace pe {
+
+class EllipsoidVtkOutput : public vtk::PointDataSource
+{
+public:
+ EllipsoidVtkOutput( ConstBlockDataID storageID, const BlockStorage & blockStorage)
+ : storageID_( storageID )
+ , blockStorage_( blockStorage ) { }
+
+ std::vector< Attributes > getAttributes() const;
+
+ void configure();
+
+ std::vector< Vector3< real_t > > getPoints();
+
+ inline void push( std::ostream& os , const uint_t /*data*/, const uint_t point, const uint_t component );
+ inline void push( vtk::Base64Writer& b64, const uint_t /*data*/, const uint_t point, const uint_t component );
+
+private:
+
+ ConstBlockDataID storageID_;
+ const BlockStorage & blockStorage_;
+ std::vector< ConstEllipsoidID > bodies_;
+ std::vector< Mat3 > tensorGlyphs_;
+};
+
+
+void EllipsoidVtkOutput::push( std::ostream& os, const uint_t data, const uint_t point, const uint_t component )
+{
+ WALBERLA_ASSERT_LESS( point, bodies_.size() );
+ WALBERLA_ASSERT_LESS( component, 9u );
+
+ switch( data )
+ {
+ case 0:
+ vtk::toStream( os, numeric_cast< float >(bodies_.at( point )->getMass()) );
+ break;
+ case 1:
+ vtk::toStream( os, numeric_cast< float >(tensorGlyphs_.at(point)[component]) );
+ break;
+ case 2:
+ vtk::toStream( os, numeric_cast< float >(bodies_.at( point )->getLinearVel()[component]) );
+ break;
+ case 3:
+ vtk::toStream( os, numeric_cast< int >(bodies_.at( point )->MPITrait.getOwner().rank_) );
+ break;
+ case 4:
+ vtk::toStream( os, numeric_cast< id_t >(bodies_.at( point )->getTopSuperBody()->getSystemID()) );
+ break;
+ case 5:
+ vtk::toStream( os, numeric_cast< id_t >(bodies_.at( point )->getTopSuperBody()->getID()) );
+ break;
+ }
+
+}
+
+
+
+void EllipsoidVtkOutput::push( vtk::Base64Writer& b64, const uint_t data, const uint_t point, const uint_t component )
+{
+ WALBERLA_ASSERT_LESS( point, bodies_.size() );
+ WALBERLA_ASSERT_LESS( component, 9u );
+
+ switch( data )
+ {
+ case 0:
+ b64 << numeric_cast< float >(bodies_.at( point )->getMass());
+ break;
+ case 1:
+ b64 << numeric_cast< float >(tensorGlyphs_.at(point)[component]);
+ break;
+ case 2:
+ b64 << numeric_cast< float >(bodies_.at( point )->getLinearVel()[component]);
+ break;
+ case 3:
+ b64 << numeric_cast< int >(bodies_.at( point )->MPITrait.getOwner().rank_);
+ break;
+ case 4:
+ b64 << numeric_cast< id_t >(bodies_.at( point )->getTopSuperBody()->getSystemID());
+ break;
+ case 5:
+ b64 << numeric_cast< id_t >(bodies_.at( point )->getTopSuperBody()->getID());
+ break;
+ }
+}
+
+
+} // namespace pe
+} // namespace walberla
+