Tpetra parallel linear algebra  Version of the Day
Tpetra_Core.cpp
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Tpetra: Templated Linear Algebra Services Package
5 // Copyright (2008) Sandia Corporation
6 //
7 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
8 // the U.S. Government retains certain rights in this software.
9 //
10 // Redistribution and use in source and binary forms, with or without
11 // modification, are permitted provided that the following conditions are
12 // met:
13 //
14 // 1. Redistributions of source code must retain the above copyright
15 // notice, this list of conditions and the following disclaimer.
16 //
17 // 2. Redistributions in binary form must reproduce the above copyright
18 // notice, this list of conditions and the following disclaimer in the
19 // documentation and/or other materials provided with the distribution.
20 //
21 // 3. Neither the name of the Corporation nor the names of the
22 // contributors may be used to endorse or promote products derived from
23 // this software without specific prior written permission.
24 //
25 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
26 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
27 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
28 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
29 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
30 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
31 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
32 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
33 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
34 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
35 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
36 //
37 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
38 //
39 // ************************************************************************
40 // @HEADER
41 
42 #include <Tpetra_Core.hpp>
43 #ifdef HAVE_TPETRACORE_MPI
44 # include <Teuchos_DefaultMpiComm.hpp>// this includes mpi.h too
45 #else
46 # include <Teuchos_DefaultSerialComm.hpp>
47 #endif // HAVE_TPETRACORE_MPI
48 #include <Kokkos_Core.hpp>
49 
50 namespace Tpetra {
51  namespace { // (anonymous)
52  // Whether one of the Tpetra::initialize() functions has been called before.
53  bool tpetraIsInitialized_ = false;
54 
55  // Whether Tpetra initialized Kokkos. Tpetra::finalize only
56  // finalizes Kokkos if it initialized Kokkos. Otherwise,
57  // something else initialized Kokkos and is responsible for
58  // finalizing it.
59  bool tpetraInitializedKokkos_ = false;
60 
61 #ifdef HAVE_TPETRACORE_MPI
62  // Whether Tpetra initialized MPI. Tpetra::finalize only
63  // finalizes MPI if it initialized MPI. Otherwise, something else
64  // initialized MPI and is responsible for finalizing it.
65  bool tpetraInitializedMpi_ = false;
66 #endif // HAVE_TPETRACORE_MPI
67 
68  // Tpetra's default communicator, wrapped in a Teuchos wrapper.
69  // After Tpetra::finalize() is called, this GOES AWAY (is set to null).
70  Teuchos::RCP<const Teuchos::Comm<int> > wrappedDefaultComm_;
71 
72  // Initialize Kokkos, if it needs initialization.
73  // This takes the same arguments as (the first two of) initialize().
74  void initKokkos (int* argc, char*** argv)
75  {
76  if (! tpetraInitializedKokkos_) {
77  // Kokkos doesn't have a global is_initialized(). However,
78  // Kokkos::initialize() always initializes the default execution
79  // space, so it suffices to check whether that was initialized.
80  const bool kokkosIsInitialized =
81  Kokkos::DefaultExecutionSpace::is_initialized ();
82  if (! kokkosIsInitialized) {
83  // Unlike MPI_Init, Kokkos promises not to modify argc and argv.
84  Kokkos::initialize (*argc, *argv);
85  tpetraInitializedKokkos_ = true;
86  }
87  }
88 
89  const bool kokkosIsInitialized =
90  Kokkos::DefaultExecutionSpace::is_initialized ();
91  TEUCHOS_TEST_FOR_EXCEPTION
92  (! kokkosIsInitialized, std::logic_error, "At the end of initKokkos, "
93  "Kokkos is not initialized. Please report this bug to the Tpetra "
94  "developers.");
95  }
96 
97 #ifdef HAVE_TPETRACORE_MPI
98  bool mpiIsInitialized (const bool throwExceptionOnFailure = true)
99  {
100  int isInitialized = 0;
101  int err = MPI_Initialized (&isInitialized);
102 
103  if (throwExceptionOnFailure) {
104  TEUCHOS_TEST_FOR_EXCEPTION
105  (err != MPI_SUCCESS, std::runtime_error, "MPI_Initialized failed with "
106  "error code " << err << " != MPI_SUCCESS. This probably indicates "
107  "that your MPI library is corrupted or that it is incorrectly linked "
108  "to your program, since this function should otherwise always "
109  "succeed.");
110  }
111  else if (err != MPI_SUCCESS) {
112  using std::cerr;
113  using std::endl;
114  cerr << "MPI_Initialized failed with error code " << err << " != "
115  "MPI_SUCCESS. This probably indicates that your MPI library is "
116  "corrupted or that it is incorrectly linked to your program, since "
117  "this function should otherwise always succeed." << endl;
118  return false; // the best we can do in this case is nothing
119  }
120  return isInitialized != 0;
121  }
122 
123  // Initialize MPI, if needed, and check for errors. This takes
124  // the same arguments as MPI_Init and the first two arguments of
125  // initialize().
126  void initMpi (int* argc, char*** argv)
127  {
128  if (! tpetraInitializedMpi_) {
129  int err = MPI_SUCCESS;
130  const bool initialized = mpiIsInitialized ();
131  if (! initialized) { // MPI not yet initialized
132  // Tpetra doesn't currently need to call MPI_Init_thread, since
133  // with Tpetra, only one thread ever calls MPI functions. If we
134  // ever want to explore MPI_THREAD_MULTIPLE, here would be the
135  // place to call MPI_Init_thread.
136  err = MPI_Init (argc, argv);
137  }
138  TEUCHOS_TEST_FOR_EXCEPTION
139  (err != MPI_SUCCESS, std::runtime_error, "MPI_Init failed with error "
140  "code " << err << " != MPI_SUCCESS. If MPI was set up correctly, this "
141  "should not happen, since we have already checked that MPI_Init (or "
142  "MPI_Init_thread) has not yet been called. This may indicate that "
143  "your MPI library is corrupted or that it is incorrectly linked to "
144  "your program.");
145  tpetraInitializedMpi_ = true;
146  }
147 
148  const bool initialized = mpiIsInitialized ();
149  TEUCHOS_TEST_FOR_EXCEPTION
150  (! initialized, std::logic_error, "At the end of initMpi, MPI is not "
151  "initialized. Please report this bug to the Tpetra developers.");
152  }
153 #endif // HAVE_TPETRACORE_MPI
154 
155  } // namespace (anonymous)
156 
157  bool isInitialized () {
158  return tpetraIsInitialized_;
159  }
160 
161  Teuchos::RCP<const Teuchos::Comm<int> > getDefaultComm ()
162  {
163  // It's OK to call this function if Tpetra::initialize hasn't been
164  // called. Users aren't obligated to call Tpetra::initialize, as
165  // long as they take responsibility for initializing Kokkos and
166  // MPI. However, MPI must be initialized. (Kokkos need not be
167  // initialized for this method to be called.)
168 
169 #ifdef HAVE_TPETRACORE_MPI
170  const bool mpiInitd = mpiIsInitialized (false); // don't throw on failure
171  TEUCHOS_TEST_FOR_EXCEPTION
172  (! mpiInitd, std::runtime_error,
173  "Tpetra::getDefaultComm: MPI has not been initialized. Before "
174  "calling this method, you must either initialize MPI (by calling "
175  "MPI_Init), or you must call Tpetra::initialize (which initializes "
176  "MPI, if it has not yet been initialized).");
177 #endif // HAVE_TPETRACORE_MPI
178 
179  // This serves as lazy initialization. Thus, the various
180  // Tpetra::initialize functions do not need to set up
181  // wrappedDefaultComm_.
182  if (wrappedDefaultComm_.is_null ()) {
183  Teuchos::RCP<const Teuchos::Comm<int> > comm;
184 #ifdef HAVE_TPETRACORE_MPI
185  comm = Teuchos::rcp (new Teuchos::MpiComm<int> (MPI_COMM_WORLD));
186 #else
187  comm = Teuchos::rcp (new Teuchos::SerialComm<int> ());
188 #endif // HAVE_TPETRACORE_MPI
189  wrappedDefaultComm_ = comm;
190  }
191  return wrappedDefaultComm_;
192  }
193 
194  void initialize (int* argc, char*** argv)
195  {
196  if (! tpetraIsInitialized_) {
197 #ifdef HAVE_TPETRACORE_MPI
198  initMpi (argc, argv); // initialize MPI, if needed
199 #endif // HAVE_TPETRACORE_MPI
200  initKokkos (argc, argv); // initialize Kokkos, if needed
201  }
202  tpetraIsInitialized_ = true;
203  }
204 
205 #ifdef HAVE_TPETRACORE_MPI
206  void initialize (int* argc, char*** argv, MPI_Comm comm)
207  {
208  initialize (argc, argv);
209  // Set the default communicator. What if users have already
210  // called initialize() before, but with a different default
211  // communicator? There are two possible things we could do here:
212  //
213  // 1. Test via MPI_Comm_compare whether comm differs from the
214  // raw MPI communicator in wrappedDefaultComm_ (if indeed it
215  // is an MpiComm).
216  // 2. Accept that the user might want to change the default
217  // communicator, and let them do it.
218  //
219  // I prefer #2. Perhaps it would be sensible to print a warning
220  // here, but on which process? Would we use the old or the new
221  // communicator to find that process' rank? We don't want to use
222  // MPI_COMM_WORLD's Process 0, since neither communicator might
223  // include that process. Furthermore, in some environments, only
224  // Process 0 in MPI_COMM_WORLD is allowed to do I/O. Thus, we
225  // just let the change go without a warning.
226  wrappedDefaultComm_ = Teuchos::rcp (new Teuchos::MpiComm<int> (comm));
227  }
228 #endif // HAVE_TPETRACORE_MPI
229 
230  void
231  initialize (int* argc, char*** argv,
232  const Teuchos::RCP<const Teuchos::Comm<int> >& comm)
233  {
234  initialize (argc, argv);
235  // Set the default communicator. What if users have already
236  // called initialize() before, but with a different default
237  // communicator? There are two possible things we could do here:
238  //
239  // 1. Test via MPI_Comm_compare whether comm differs from the
240  // raw MPI communicator in wrappedDefaultComm_ (if indeed it
241  // is an MpiComm).
242  // 2. Accept that the user might want to change the default
243  // communicator, and let them do it.
244  //
245  // I prefer #2. Perhaps it would be sensible to print a warning
246  // here, but on which process? Would we use the old or the new
247  // communicator to find that process' rank? We don't want to use
248  // MPI_COMM_WORLD's Process 0, since neither communicator might
249  // include that process. Furthermore, in some environments, only
250  // Process 0 in MPI_COMM_WORLD is allowed to do I/O. Thus, we
251  // just let the change go without a warning.
252  wrappedDefaultComm_ = comm;
253  }
254 
255  void finalize ()
256  {
257  if (! tpetraIsInitialized_) {
258  return; // user didn't call initialize(), so do nothing at all
259  }
260 
261  // Tpetra should only finalize Kokkos if it initialized Kokkos.
262  // See Github Issue #434.
263  if (tpetraInitializedKokkos_) {
264  Kokkos::finalize ();
265  }
266 
267  // Make sure that no outstanding references to the communicator
268  // remain. If users gave initialize() an MPI_Comm, _they_ are
269  // responsible for freeing it before calling finalize().
270  wrappedDefaultComm_ = Teuchos::null;
271 
272 #ifdef HAVE_TPETRACORE_MPI
273  // Tpetra should only finalize MPI if it initialized MPI.
274  // See Github Issue #434.
275  if (tpetraInitializedMpi_) {
276  // finalize() is a kind of destructor, so it's a bad idea to throw
277  // an exception on error. Instead, we print an error message and
278  // do nothing.
279  const bool throwOnFail = false;
280  const bool mpiInitd = mpiIsInitialized (throwOnFail);
281 
282  if (mpiInitd) {
283  // This must be called by the same thread that called MPI_Init
284  // (possibly, but not necessarily, in Tpetra::initialize()).
285  const int err = MPI_Finalize ();
286 
287  // finalize() is a kind of destructor, so it's a bad idea to
288  // throw an exception. Instead, just print out an error message
289  // and hope for the best.
290  if (err != MPI_SUCCESS) {
291  std::cerr << "MPI_Finalize failed with error code " << err << " != "
292  "MPI_SUCCESS. I'm not sure how to recover from this safely, "
293  "so I will keep going and hope for the best." << std::endl;
294  }
295  }
296  }
297 #endif // HAVE_TPETRACORE_MPI
298 
299  tpetraIsInitialized_ = false; // it's not anymore.
300  }
301 } // namespace Tpetra
bool isInitialized()
Whether Tpetra is in an initialized state.
void initialize(int *argc, char ***argv)
Initialize Tpetra.
Namespace Tpetra contains the class and methods constituting the Tpetra library.
void finalize()
Finalize Tpetra.
Functions for initializing and finalizing Tpetra.
Teuchos::RCP< const Teuchos::Comm< int > > getDefaultComm()
Get Tpetra&#39;s default communicator.