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_TPETRA_MPI
44 # include <Teuchos_DefaultMpiComm.hpp>// this includes mpi.h too
45 #else
46 # include <Teuchos_DefaultSerialComm.hpp>
47 #endif // HAVE_TPETRA_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  // Tpetra's default communicator, wrapped in a Teuchos wrapper.
56  // After Tpetra::finalize() is called, this GOES AWAY (is set to null).
57  Teuchos::RCP<const Teuchos::Comm<int> > wrappedDefaultComm_;
58 
59 #ifdef HAVE_TPETRA_MPI
60  // Initialize MPI, if needed, and check for errors.
61  void initMpi (int* argc, char*** argv)
62  {
63  int isInitialized = 0;
64  int err = MPI_Initialized (&isInitialized);
65 
66  TEUCHOS_TEST_FOR_EXCEPTION(
67  err != MPI_SUCCESS, std::runtime_error, "MPI_Initialized failed with "
68  "error code " << err << " != MPI_SUCCESS. This probably indicates "
69  "that your MPI library is corrupted or that it is incorrectly linked "
70  "to your program, since this function should otherwise always "
71  "succeed.");
72 
73  if (isInitialized == 0) { // MPI not yet initialized
74  // Tpetra doesn't currently need to call MPI_Init_thread, since
75  // with Tpetra, only one thread ever calls MPI functions. If we
76  // ever want to explore MPI_THREAD_MULTIPLE, here would be the
77  // place to call MPI_Init_thread.
78  err = MPI_Init (argc, argv);
79  }
80 
81  TEUCHOS_TEST_FOR_EXCEPTION(
82  err != MPI_SUCCESS, std::runtime_error, "MPI_Init failed with error "
83  "code " << err << " != MPI_SUCCESS. If MPI was set up correctly, this "
84  "should not happen, since we have already checked that MPI_Init (or "
85  "MPI_Init_thread) has not yet been called. This may indicate that "
86  "your MPI library is corrupted or that it is incorrectly linked to "
87  "your program.");
88  }
89 #endif // HAVE_TPETRA_MPI
90 
91  } // namespace (anonymous)
92 
93  bool isInitialized () {
94  return tpetraIsInitialized_;
95  }
96 
97  Teuchos::RCP<const Teuchos::Comm<int> > getDefaultComm ()
98  {
99  TEUCHOS_TEST_FOR_EXCEPTION(
100  tpetraIsInitialized_, std::runtime_error, "Tpetra::getDefaultComm: "
101  "You must call Tpetra::initialize before you may get a default "
102  "communicator.");
103  TEUCHOS_TEST_FOR_EXCEPTION(
104  wrappedDefaultComm_.is_null (), std::logic_error,
105  "Tpetra::getDefaultComm: wrappedDefaultComm_ is null! "
106  "This should never happen, because Tpetra claims that it has been "
107  "initialized. Please report this bug to the Tpetra developers.");
108  return wrappedDefaultComm_;
109  }
110 
111 
112  void initialize (int* argc, char*** argv)
113  {
114 #ifdef HAVE_TPETRA_MPI
115  initMpi (argc, argv);
116 #endif // HAVE_TPETRA_MPI
117 
118  if (! tpetraIsInitialized_) {
119  // Unlike MPI_Init, Kokkos promises not to modify argc and argv.
120  Kokkos::initialize (*argc, *argv);
121 
122 #ifdef HAVE_TPETRA_MPI
123  wrappedDefaultComm_ =
124  Teuchos::rcp (new Teuchos::MpiComm<int> (MPI_COMM_WORLD));
125 #else
126  wrappedDefaultComm_ = Teuchos::rcp (new Teuchos::SerialComm<int> ());
127 #endif // HAVE_TPETRA_MPI
128  tpetraIsInitialized_ = true;
129  }
130  }
131 
132 
133 #ifdef HAVE_TPETRA_MPI
134  void initialize (int* argc, char*** argv, MPI_Comm comm)
135  {
136 #ifdef HAVE_TPETRA_MPI
137  initMpi (argc, argv);
138 #endif // HAVE_TPETRA_MPI
139 
140  // Set the default communicator. What if users have already
141  // called initialize() before, but with a different default
142  // communicator? There are two possible things we could do here:
143  //
144  // 1. Test via MPI_Comm_compare whether comm differs from the
145  // raw MPI communicator in wrappedDefaultComm_ (if indeed it
146  // is an MpiComm).
147  // 2. Accept that the user might want to change the default
148  // communicator, and let them do it.
149  //
150  // I prefer #2. Perhaps it would be sensible to print a warning
151  // here, but on which process? Would we use the old or the new
152  // communicator to find that process' rank? We don't want to
153  // use MPI_COMM_WORLD's Process 0, since neither communicator
154  // might include that process. Furthermore, in some
155  // environments, only Process 0 in MPI_COMM_WORLD is allowed to
156  // do I/O. Thus, we just let the change go without a warning.
157  wrappedDefaultComm_ = Teuchos::rcp (new Teuchos::MpiComm<int> (comm));
158 
159  if (! tpetraIsInitialized_) {
160  // Unlike MPI_Init, Kokkos promises not to modify argc and argv.
161  Kokkos::initialize (*argc, *argv);
162  tpetraIsInitialized_ = true;
163  }
164  }
165 #endif // HAVE_TPETRA_MPI
166 
167 
168  void
169  initialize (int* argc, char*** argv,
170  const Teuchos::RCP<const Teuchos::Comm<int> >& comm)
171  {
172 #ifdef HAVE_TPETRA_MPI
173  initMpi (argc, argv);
174 #endif // HAVE_TPETRA_MPI
175 
176  // Set the default communicator. What if users have already
177  // called initialize() before, but with a different default
178  // communicator? There are two possible things we could do here:
179  //
180  // 1. Test via MPI_Comm_compare whether comm differs from the
181  // raw MPI communicator in wrappedDefaultComm_ (if indeed it
182  // is an MpiComm).
183  // 2. Accept that the user might want to change the default
184  // communicator, and let them do it.
185  //
186  // I prefer #2. Perhaps it would be sensible to print a warning
187  // here, but on which process? Would we use the old or the new
188  // communicator to find that process' rank? We don't want to
189  // use MPI_COMM_WORLD's Process 0, since neither communicator
190  // might include that process. Furthermore, in some
191  // environments, only Process 0 in MPI_COMM_WORLD is allowed to
192  // do I/O. Thus, we just let the change go without a warning.
193  wrappedDefaultComm_ = comm;
194 
195  if (! tpetraIsInitialized_) {
196  // Unlike MPI_Init, Kokkos promises not to modify argc and argv.
197  Kokkos::initialize (*argc, *argv);
198  tpetraIsInitialized_ = true;
199  }
200  }
201 
202 
203  void finalize ()
204  {
205  using std::cerr;
206  using std::endl;
207 
208  if (! tpetraIsInitialized_) {
209  // If the user didn't call initialize(), then we shouldn't call
210  // MPI_Finalize() or otherwise shut anything else off. For
211  // example, the user might want to use Kokkos without Tpetra, so
212  // they might have called Kokkos::initialize() without calling
213  // Tpetra::initialize(). In that case, we shouldn't call
214  // Kokkos::initialize() here.
215  return;
216  }
217 
218  Kokkos::finalize ();
219 
220  // Make sure that no outstanding references to the communicator
221  // remain. If users gave initialize() an MPI_Comm, _they_ are
222  // responsible for freeing it before calling finalize().
223  wrappedDefaultComm_ = Teuchos::null;
224 
225 #ifdef HAVE_TPETRA_MPI
226  int isInitialized = 0;
227  int err = MPI_Initialized (&isInitialized);
228 
229  // finalize() is a kind of destructor, so it's a bad idea to throw
230  // an exception. Instead, just print out an error message and
231  // hope for the best.
232  if (err != MPI_SUCCESS) {
233  cerr << "MPI_Initialized failed with error code " << err << " != "
234  "MPI_SUCCESS. This probably indicates that your MPI library is "
235  "corrupted or that it is incorrectly linked to your program, since "
236  "this function should otherwise always succeed." << endl;
237  }
238  else if (isInitialized != 0) {
239  // This must be called by the same thread that called MPI_Init
240  // (possibly, but not necessarily, in Tpetra::initialize()).
241  err = MPI_Finalize ();
242 
243  // finalize() is a kind of destructor, so it's a bad idea to
244  // throw an exception. Instead, just print out an error message
245  // and hope for the best.
246  if (err != MPI_SUCCESS) {
247  cerr << "MPI_Finalize failed with error code " << err << " != "
248  "MPI_SUCCESS. I'm not sure how to recover from this safely, "
249  "so I will keep going and hope for the best." << endl;
250  }
251  }
252 #endif // HAVE_TPETRA_MPI
253 
254  tpetraIsInitialized_ = false; // it's not anymore.
255  }
256 } // namespace Tpetra
bool isInitialized()
Whether Tpetra is in an initialized state.
Definition: Tpetra_Core.cpp:93
void initialize(int *argc, char ***argv)
Initialize Tpetra.
Namespace Tpetra contains the class and methods constituting the Tpetra library.
Functions for initializing and finalizing Tpetra.
Teuchos::RCP< const Teuchos::Comm< int > > getDefaultComm()
Get Tpetra&#39;s default communicator.
Definition: Tpetra_Core.cpp:97