1 #ifndef DUNE_ALU3DGRIDGRID_HH
2 #define DUNE_ALU3DGRIDGRID_HH
11 #include <dune/common/bigunsignedint.hh>
12 #include <dune/common/deprecated.hh>
13 #include <dune/common/static_assert.hh>
15 #include <dune/geometry/referenceelements.hh>
37 #include <dune/common/mpihelper.hh>
39 #if ALU3DGRID_PARALLEL
40 #include <dune/common/mpicollectivecommunication.hh>
42 #include <dune/common/collectivecommunication.hh>
49 template<
int cd,
int dim,
class Gr
idImp>
50 class ALU3dGridEntity;
51 template<
int cd, PartitionIteratorType pitype,
class Gr
idImp >
52 class ALU3dGridLevelIterator;
53 template<
int cd,
class Gr
idImp >
54 class ALU3dGridEntityPointerBase;
55 template<
int cd,
class Gr
idImp >
56 class ALU3dGridEntitySeed;
57 template<
int cd,
class Gr
idImp >
58 class ALU3dGridEntityPointer;
59 template<
int mydim,
int coorddim,
class Gr
idImp>
60 class ALU3dGridGeometry;
61 template<
class Gr
idImp>
62 class ALU3dGridHierarchicIterator;
63 template<
class Gr
idImp>
64 class ALU3dGridIntersectionIterator;
65 template<
class Gr
idImp>
67 template<
int codim, PartitionIteratorType pitype,
class Gr
idImp>
69 template <
int mydim,
int coorddim,
class Gr
idImp>
71 template <
class Gr
idImp>
73 template< ALU3dGr
idElementType,
class >
75 template< ALU3dGr
idElementType,
class >
77 template< ALU3dGr
idElementType,
class >
79 template <
class EntityImp>
83 template <
class Gr
idImp,
class GeometryImp,
int nChild>
85 template< ALU3dGr
idElementType elType,
class Comm >
93 #if ALU3DGRID_PARALLEL
94 template< ALU3dGr
idElementType elType,
class Comm = MPI_Comm >
96 #else // #if ALU3DGRID_PARALLEL
97 template< ALU3dGr
idElementType elType,
class Comm = No_Comm >
99 #endif // #else // #if ALU3DGRID_PARALLEL
101 template <
class Comm >
107 template <
class Comm>
113 namespace DefaultIndexSetHelper
116 template< ALU3dGr
idElementType elType,
class Comm,
class Index >
125 return (container.getData( index ).index() >= 0);
136 template< ALU3dGr
idElementType elType,
class Comm >
139 template< ALU3dGr
idElementType elType >
155 if( macroName.empty() )
163 static int getRank ( No_Comm comm ) {
return 0; }
167 ALU3DSPACE Gitter::Geometric::BuilderIF* builder =
168 dynamic_cast< ALU3DSPACE Gitter::Geometric::BuilderIF*
>( &grid.container() );
170 DUNE_THROW(InvalidStateException,
"dynamic_cast of ALUGrid builder failed");
179 #if ALU3DGRID_PARALLEL
180 template< ALU3dGr
idElementType elType >
186 typedef ALU3DSPACE GitterDunePll GitterImplType;
188 typedef Dune::CollectiveCommunication< MPI_Comm > CollectiveCommunication;
191 : ccobj_( comm ), mpAccess_( comm )
194 int nlinks ()
const {
return mpAccess_.nlinks(); }
196 GitterImplType *createALUGrid (
const std::string ¯oName,
ALU3DSPACE ProjectVertex *projection )
198 return new GitterImplType( macroName.c_str(), mpAccess_, projection );
201 static MPI_Comm defaultComm () {
return MPI_COMM_WORLD; }
203 static int getRank ( MPI_Comm comm )
206 MPI_Comm_rank( comm, &rank );
210 static typename ALU3DSPACE Gitter::Geometric::BuilderIF &getBuilder ( GitterImplType &grid )
212 ALU3DSPACE Gitter::Geometric::BuilderIF* builder =
213 dynamic_cast< ALU3DSPACE Gitter::Geometric::BuilderIF*
>( &grid.containerPll() );
215 DUNE_THROW(InvalidStateException,
"dynamic_cast of ALUGrid builder failed");
219 static void duneNotifyMacroGridChanges ( GitterImplType &grid )
221 grid.duneNotifyMacroGridChanges();
224 CollectiveCommunication ccobj_;
227 #endif // #if ALU3DGRID_PARALLEL
234 template< ALU3dGr
idElementType elType,
class Comm >
295 template< PartitionIteratorType pitype >
306 template< PartitionIteratorType pitype >
364 template< ALU3dGr
idElementType elType,
class Comm >
367 ALU3dGridFamily< elType, Comm > >,
375 typedef ThisType MyType;
530 ALU3dGrid (
const std::string ¯oTriangFilename,
541 static inline std::string
name ();
549 template<
int cd, PartitionIteratorType pitype>
554 template<
int cd, PartitionIteratorType pitype>
556 lend (
int level)
const;
568 lend (
int level)
const;
572 template <
int codim, PartitionIteratorType pitype>
574 leafbegin(
int level)
const;
577 template <
int codim, PartitionIteratorType pitype>
579 leafend(
int level)
const;
584 leafbegin(
int level)
const;
589 leafend(
int level)
const;
605 template <
int codim, PartitionIteratorType pitype>
610 template <
int codim, PartitionIteratorType pitype>
626 template <
int codim, PartitionIteratorType pitype>
628 createLeafIteratorBegin (
int level)
const;
631 template <
int codim, PartitionIteratorType pitype>
633 createLeafIteratorEnd(
int level)
const;
637 int size (
int level,
int cd)
const;
640 int size (
int codim)
const;
706 template <
class DataHandle>
709 template<
class DataHandleImpl,
class Data >
713 LBHandle lbHandle( *
this, dataHandle );
718 int ghostSize (
int level,
int codim)
const;
730 template<
class DataHandleImp,
class DataTypeImp>
737 template<
class DataHandleImp,
class DataTypeImp>
761 template<
class Gr
idImp,
class DataHandle >
767 template<
class Gr
idImp,
class DataHandle >
775 template <GrapeIOFileFormatType ftype>
784 bool writeMacroGrid(
const std::string path,
const std::string filename )
const ;
788 template <GrapeIOFileFormatType ftype>
802 bool mark(
int refCount ,
const typename Traits::template Codim<0>::Entity & e);
805 int getMark(
const typename Traits::template Codim<0>::Entity & e)
const;
810 return Communications::defaultComm();
815 template<
class IntersectionType >
816 static const typename BaseType
817 :: template ReturnImplementationType< IntersectionType >
818 :: ImplementationType &
843 return Communications::getBuilder(
myGrid() );
849 Communications::duneNotifyMacroGridChanges(
myGrid() );
855 template <
class EntitySeed >
856 typename Traits :: template Codim< EntitySeed :: codimension >
:: EntityPointer
859 enum { codim = EntitySeed :: codimension };
862 return ALUPointer(
factory(), seed ) ;
882 assert( level >= 0 );
891 assert( codim >= 1 );
892 assert( codim <= 3 );
898 assert( codim >= 1 );
899 assert( codim <= 3 );
901 assert( level >= 0 );
908 assert( level >= 0 );
918 const ThisType &
operator= (
const ThisType & );
946 assert( segmentIndex < (
int)
bndVec_->size() );
947 return (*
bndVec_)[ segmentIndex ];
958 #ifdef USE_SMP_PARALLEL
1022 #ifdef USE_SMP_PARALLEL
1023 std::vector< GridObjectFactoryType > factoryVec_;
1049 const std::string filename );
1052 namespace Capabilities
1055 template< ALU3dGr
idElementType elType,
class Comm,
int cdim >
1058 static const bool v =
true;
1061 template< ALU3dGr
idElementType elType,
class Comm >
1064 static const bool v =
true;
1067 template< ALU3dGr
idElementType elType,
class Comm >
1070 static const bool v =
true;
1073 template< ALU3dGr
idElementType elType,
class Comm >
1076 static const bool v =
true;
1084 #if COMPILE_ALUGRID_INLINE
1085 #include "grid_imp.cc"