1 #ifndef DUNE_ALU3DGRID_HSFC_HH
2 #define DUNE_ALU3DGRID_HSFC_HH
8 #include <dune/common/exceptions.hh>
10 #include <dune/common/parallel/mpihelper.hh>
11 #include <dune/common/parallel/collectivecommunication.hh>
12 #include <dune/common/parallel/mpicollectivecommunication.hh>
14 #include <dune/alugrid/impl/parallel/zcurve.hh>
16 #include <dune/alugrid/impl/parallel/aluzoltan.hh>
19 extern double Zoltan_HSFC_InvHilbert3d (Zoltan_Struct *zz,
double *coord);
20 extern double Zoltan_HSFC_InvHilbert2d (Zoltan_Struct *zz,
double *coord);
26 template <
class Coordinate>
30 typedef Dune :: CollectiveCommunication< typename Dune :: MPIHelper :: MPICommunicator >
31 CollectiveCommunication ;
34 typedef void Zoltan_Struct;
35 typedef CollectiveCommunication Zoltan;
39 typedef double zoltan_hsfc_inv_t(Zoltan_Struct *zz,
double *coord);
41 static const int dimension = Coordinate::dimension;
46 zoltan_hsfc_inv_t* hsfcInv_;
52 const Coordinate& upper,
53 const CollectiveCommunication& comm =
54 CollectiveCommunication( Dune::MPIHelper::getCommunicator() ) )
58 hsfcInv_( dimension == 3 ? Zoltan_HSFC_InvHilbert3d : Zoltan_HSFC_InvHilbert2d ),
72 assert( point.size() == (
unsigned int)dimension );
74 Coordinate center( 0 ) ;
76 for(
int d=0; d<dimension; ++d )
77 center[ d ] = (point[ d ] - lower_[ d ]) / length_[ d ];
80 return hsfcInv_( zz_.Get_C_Handle(), ¢er[ 0 ] );
82 DUNE_THROW(Dune::SystemError,
"Zoltan not found, cannot use Zoltan's Hilbert curve");
88 double index(
const Coordinate& point )
const
94 template<
class Gr
idView >
97 typedef typename GridView :: template Codim< 0 > :: Iterator Iterator ;
98 typedef typename Iterator :: Entity :: Geometry :: GlobalCoordinate GlobalCoordinate ;
100 std::vector< GlobalCoordinate > vertices;
101 vertices.reserve( view.indexSet().size( 0 ) );
103 const Iterator endit = view.template end< 0 > ();
104 for(Iterator it = view.template begin< 0 > (); it != endit; ++ it )
106 GlobalCoordinate center = (*it).geometry().center();
107 vertices.push_back( center );
110 std::stringstream gnufilename;
111 gnufilename << name <<
".gnu";
112 if( view.grid().comm().size() > 1 )
113 gnufilename <<
"." << view.grid().comm().rank();
115 std::ofstream gnuFile ( gnufilename.str() );
118 for(
size_t i=0; i<vertices.size(); ++i )
120 gnuFile << vertices[ i ] << std::endl;
127 std::stringstream vtkfilename;
128 vtkfilename << name <<
".vtk";
129 if( view.grid().comm().size() > 1 )
130 vtkfilename <<
"." << view.grid().comm().rank();
131 std::ofstream vtkFile ( vtkfilename.str() );
135 vtkFile <<
"# vtk DataFile Version 1.0" << std::endl;
136 vtkFile <<
"Line representation of vtk" << std::endl;
137 vtkFile <<
"ASCII" << std::endl;
138 vtkFile <<
"DATASET POLYDATA" << std::endl;
139 vtkFile <<
"POINTS "<< vertices.size() <<
" FLOAT" << std::endl;
141 for(
size_t i=0; i<vertices.size(); ++i )
143 vtkFile << vertices[ i ];
144 for(
int d=GlobalCoordinate::dimension; d<3; ++d )
146 vtkFile << std::endl;
150 vtkFile <<
"LINES " << vertices.size()-1 <<
" " << (vertices.size()-1)*3 << std::endl;
152 for(
size_t i=0; i<vertices.size()-1; ++i )
153 vtkFile <<
"2 " << i <<
" " << i+1 << std::endl;
164 template <
class Coordinate>
167 typedef ::ALUGridSFC::ZCurve< long int, Coordinate::dimension> ZCurveOrderingType;
171 typedef Dune :: CollectiveCommunication< typename MPIHelper :: MPICommunicator >
172 CollectiveCommunication ;
190 const Coordinate& lower,
191 const Coordinate& upper,
192 const CollectiveCommunication& comm =
193 CollectiveCommunication( Dune::MPIHelper::getCommunicator() ) )
194 :
zCurve_ ( lower, upper, comm )
200 template <
class OtherComm>
202 const Coordinate& lower,
203 const Coordinate& upper,
204 const OtherComm& otherComm )
206 otherComm.size() > 1 ? CollectiveCommunication(
Dune::MPIHelper::getCommunicator() ) :
207 CollectiveCommunication(
Dune::MPIHelper::getLocalCommunicator() ) )
209 otherComm.size() > 1 ? CollectiveCommunication(
Dune::MPIHelper::getCommunicator() ) :
210 CollectiveCommunication(
Dune::MPIHelper::getLocalCommunicator() ) )
216 double index(
const Coordinate& point )
const
220 return double(
zCurve_.index( point ) );
228 DUNE_THROW(NotImplemented,
"Wrong space filling curve ordering selected");
236 #endif // #ifndef DUNE_ALU3DGRID_HSFC_HH