EllipsoidVtkOutput.h 3.95 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36
//======================================================================================================================
//
//  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>
37
#include <array>
38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62

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_;
63
   std::vector< std::array<real_t,6> > tensorGlyphs_;
64 65 66 67 68 69
};


void EllipsoidVtkOutput::push( std::ostream& os, const uint_t data, const uint_t point, const uint_t component )
{
   WALBERLA_ASSERT_LESS( point, bodies_.size() );
70
   WALBERLA_ASSERT_LESS( component, 6u );
71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100

   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() );
101
   WALBERLA_ASSERT_LESS( component, 6u );
102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129

   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