dune-fem  2.4.1-rc
mpimanager.hh
Go to the documentation of this file.
1 #ifndef DUNE_FEM_MPIMANAGER_HH
2 #define DUNE_FEM_MPIMANAGER_HH
3 
4 #include <dune/common/parallel/mpicollectivecommunication.hh>
5 #include <dune/common/parallel/mpihelper.hh>
6 
8 
9 #if HAVE_PETSC
11 #endif
12 
13 namespace Dune
14 {
15 
16  namespace Fem
17  {
18 
19  struct MPIManager
20  {
21  typedef Dune::CollectiveCommunication< MPIHelper::MPICommunicator >
23 
24  private:
25  MPIManager ()
26  : helper_( 0 ),
27  comm_( 0 )
28  {}
29 
30  ~MPIManager ()
31  {
32  delete comm_;
33  }
34 
35  // prohibit copying and assignment
36  MPIManager ( const MPIManager & );
37  MPIManager &operator= ( const MPIManager & );
38 
39  static MPIManager &instance ()
40  {
41  static MPIManager instance;
42  return instance;
43  }
44 
45 #if HAVE_PETSC
46  class PETSc
47  {
48  ~PETSc() { ::Dune::Petsc::finalize(); }
49  public:
50  static void initialize( const bool verbose, int &argc, char **&argv )
51  {
52  // needed for later calling Petsc::finalize to the right time
53  static PETSc petsc;
54  ::Dune::Petsc::initialize( verbose, argc, argv );
55  }
56  };
57 #endif
58 
59  public:
60  static void initialize ( int &argc, char **&argv )
61  {
62  MPIHelper *&helper = instance().helper_;
63  CollectiveCommunication *&comm = instance().comm_;
64 
65  // the following initalization is only enabled for
66  // MPI-thread parallel programs
67 #if HAVE_MPI && MPI_2
68 #ifdef USE_SMP_PARALLEL
69  int provided;
70  // use MPI_Init_thread for hybrid parallel programs
71  int is_initialized = MPI_Init_thread(&argc, &argv, MPI_THREAD_FUNNELED, &provided );
72 
73  if( is_initialized != MPI_SUCCESS )
74  DUNE_THROW(InvalidStateException,"MPI_Init_thread failed!");
75 
76 #if not defined NDEBUG && defined DUNE_DEVEL_MODE
77  // for OpenMPI provided seems to be MPI_THREAD_SINGLE
78  // but the bybrid version still works. On BlueGene systems
79  // the MPI_THREAD_FUNNELED is really needed
80  if( provided != MPI_THREAD_FUNNELED )
81  {
82  if( provided == MPI_THREAD_SINGLE )
83  dwarn << "MPI thread support = single (instead of funneled)!" << std::endl;
84  else
85  dwarn << "WARNING: MPI thread support = " << provided << " != MPI_THREAD_FUNNELED " << MPI_THREAD_FUNNELED << std::endl;
86  }
87 #endif // end NDEBUG
88 
89 #endif // end USE_SMP_PARALLEL
90 #endif // end HAVE_MPI && MPI_2
91 
92  if( (helper != 0) || (comm != 0) )
93  DUNE_THROW( InvalidStateException, "MPIManager has already been initialized." );
94 
95  // if not already called, this will call MPI_Init
96  helper = &MPIHelper::instance( argc, argv );
97  comm = new CollectiveCommunication( helper->getCommunicator() );
98 
99 #if HAVE_PETSC
100  // initialize PETSc if pressent
101  PETSc::initialize( rank() == 0, argc, argv );
102 #endif
103 
104  // initialize static variables of QuadratureStorageRegistry
106  }
107 
108  static const CollectiveCommunication &comm ()
109  {
110  const CollectiveCommunication *const comm = instance().comm_;
111  if( comm == 0 )
112  DUNE_THROW( InvalidStateException, "MPIManager has not been initialized." );
113  return *comm;
114  }
115 
116  static int rank ()
117  {
118  return comm().rank();
119  }
120 
121  static int size ()
122  {
123  return comm().size();
124  }
125 
126  private:
127  MPIHelper *helper_;
129  };
130 
131  } // namespace Fem
132 
133 } // namespace Dune
134 
135 #endif // #ifndef DUNE_FEM_MPIMANAGER_HH
static int rank()
Definition: mpimanager.hh:116
static int size()
Definition: mpimanager.hh:121
static const CollectiveCommunication & comm()
Definition: mpimanager.hh:108
static void initialize()
initialize static variables
Definition: registry.hh:63
static void initialize(int &argc, char **&argv)
Definition: mpimanager.hh:60
Definition: coordinate.hh:4
Definition: mpimanager.hh:19
Dune::CollectiveCommunication< MPIHelper::MPICommunicator > CollectiveCommunication
Definition: mpimanager.hh:22