Teuchos Package Browser (Single Doxygen Collection)  Version of the Day
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
Teuchos_PerformanceMonitorBase.cpp
Go to the documentation of this file.
1 // @HEADER
2 // ***********************************************************************
3 //
4 // Teuchos: Common Tools Package
5 // Copyright (2004) Sandia Corporation
6 //
7 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
8 // license for use of this work by or on behalf of the U.S. Government.
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 
43 #include "Teuchos_CommHelpers.hpp"
44 #include "Teuchos_DefaultComm.hpp"
45 #include <iterator> // std::back_inserter
46 
47 namespace Teuchos {
48 
49  namespace {
50 
51  // Pack the given array of strings into a single string with an
52  // offsets array. This is a helper routine of \c sendStrings().
53  // For strings[k], offsets[k] gives its start position in
54  // packedString, and offsets[k+1]-ofsets[k] gives its length.
55  // Thus, packedString.substr (offsets[k], offsets[k+1]-offsets[k])
56  // == strings[k].
57  //
58  // \param packedString [out] The packed string. It will be
59  // resized on output as necessary.
60  //
61  // \param offsets [out] Array of offsets, of length one more than
62  // the number of strings in the \c strings array. Thus, the
63  // offsets array always has positive length.
64  //
65  // \param strings [in] The array of strings to pack. It may have
66  // zero length. In that case, on output, the offsets array will
67  // have length one, offsets[0] = 0, and packedString will have
68  // length zero.
69  //
70  // Although std::string could be considered an array, it does not
71  // have a size_type typedef. Instead, it uses size_t for that
72  // purpose.
73  void
74  packStringsForSend (std::string& packedString,
75  Array<size_t>& offsets,
76  const Array<std::string>& strings)
77  {
78  using std::cerr;
79  using std::endl;
80  using std::string;
81 
82  const bool debug = false;
83 
84  // Compute index offsets in the packed string.
85  offsets.resize (strings.size() + 1);
86  size_t totalLength = 0;
87  Array<size_t>::size_type offsetsIndex = 0;
88  for (Array<string>::const_iterator it = strings.begin();
89  it != strings.end(); ++it, ++offsetsIndex)
90  {
91  offsets[offsetsIndex] = totalLength;
92  totalLength += it->size();
93  }
94  offsets[offsetsIndex] = totalLength;
95 
96  // Pack the array of strings into the packed string.
97  packedString.resize (totalLength);
98  string::iterator packedStringIter = packedString.begin();
99  for (Array<string>::const_iterator it = strings.begin();
100  it != strings.end(); ++it)
101  packedStringIter = std::copy (it->begin(), it->end(), packedStringIter);
102 
103  if (debug)
104  {
105  std::ostringstream out;
106  RCP<const Comm<int> > pComm = DefaultComm<int>::getComm ();
107  out << "Proc " << pComm->getRank() << ": in pack: offsets = [";
108  for (Array<size_t>::const_iterator it = offsets.begin();
109  it != offsets.end(); ++it)
110  {
111  out << *it;
112  if (it + 1 != offsets.end())
113  out << ", ";
114  }
115  out << "], packedString = " << packedString << endl;
116  cerr << out.str();
117  }
118  }
119 
120  // \brief Send an array of strings.
121  //
122  // Teuchos::send() (or rather, Teuchos::SerializationTraits)
123  // doesn't know how to send an array of strings. This function
124  // packs an array of strings into a single string with an offsets
125  // array, and sends the offsets array (and the packed string, if
126  // it is not empty).
127  void
128  sendStrings (const Comm<int>& comm, // in
129  const Array<std::string>& strings, // in
130  const int destRank) // in
131  {
132  // Pack the string array into the packed string, and compute
133  // offsets.
134  std::string packedString;
135  Array<size_t> offsets;
136  packStringsForSend (packedString, offsets, strings);
137  TEUCHOS_TEST_FOR_EXCEPTION(offsets.size() == 0, std::logic_error,
138  "packStringsForSend() returned a zero-length offsets "
139  "array on MPI Proc " << comm.getRank() << ", to be "
140  "sent to Proc " << destRank << ". The offsets array "
141  "should always have positive length. Please report "
142  "this bug to the Teuchos developers.");
143 
144  // Send the count of offsets.
145  send (comm, offsets.size(), destRank);
146 
147  // Send the array of offsets. There is always at least one
148  // element in the offsets array, so we can always take the
149  // address of the first element.
150  const int offsetsSendCount = static_cast<int> (offsets.size());
151  send (comm, offsetsSendCount, &offsets[0], destRank);
152 
153  // Now send the packed string. It may be empty if the strings
154  // array has zero elements or if all the strings in the array
155  // are empty. If the packed string is empty, we don't send
156  // anything, since the receiving process already knows (from the
157  // offsets array) not to expect anything.
158  const int stringSendCount = static_cast<int> (packedString.size());
159  if (stringSendCount > 0)
160  send (comm, stringSendCount, &packedString[0], destRank);
161  }
162 
163  void
164  unpackStringsAfterReceive (Array<std::string>& strings,
165  const std::string& packedString,
166  const Array<size_t> offsets)
167  {
168  const bool debug = false;
169  if (debug)
170  {
171  using std::cerr;
172  using std::endl;
173 
174  std::ostringstream out;
175  RCP<const Comm<int> > pComm = DefaultComm<int>::getComm ();
176  out << "Proc " << pComm->getRank() << ": in unpack: offsets = [";
177  for (Array<size_t>::const_iterator it = offsets.begin();
178  it != offsets.end(); ++it)
179  {
180  out << *it;
181  if (it + 1 != offsets.end())
182  out << ", ";
183  }
184  out << "], packedString = " << packedString << endl;
185  cerr << out.str();
186  }
187  TEUCHOS_TEST_FOR_EXCEPTION(offsets.size() == 0, std::logic_error,
188  "The offsets array has length zero, which does not "
189  "make sense. Even when sending / receiving zero "
190  "strings, the offsets array should have one entry "
191  "(namely, zero).");
192  const Array<size_t>::size_type numStrings = offsets.size() - 1;
193  strings.resize (numStrings);
194  for (Array<size_t>::size_type k = 0; k < numStrings; ++k)
195  { // Exclusive index range in the packed string in which to
196  // find the current string.
197  const size_t start = offsets[k];
198  const size_t end = offsets[k+1];
199  strings[k] = packedString.substr (start, end - start);
200  }
201  }
202 
203  // Function corresponding to \c sendStrings() that receives an
204  // array of strings (Array<std::string>) in packed form.
205  void
206  receiveStrings (const Comm<int>& comm,
207  const int sourceRank,
208  Array<std::string>& strings)
209  {
210  // Receive the number of offsets. There should always be at
211  // least 1 offset.
212  Array<size_t>::size_type numOffsets = 0;
213  receive (comm, sourceRank, &numOffsets);
214  TEUCHOS_TEST_FOR_EXCEPTION(numOffsets == 0, std::logic_error,
215  "Invalid number of offsets numOffsets=" << numOffsets
216  << " received on MPI Rank " << comm.getRank()
217  << " from Rank " << sourceRank << ". Please report "
218  "this bug to the Teuchos developers.");
219 
220  // Receive the array of offsets.
221  Array<size_t> offsets (numOffsets);
222  const int offsetsRecvCount = static_cast<int> (numOffsets);
223  receive (comm, sourceRank, offsetsRecvCount, &offsets[0]);
224 
225  // If the packed string is nonempty, receive the packed string,
226  // and unpack it. The last entry of offsets is the length of
227  // the packed string.
228  std::string packedString (offsets.back(), ' ');
229  const int stringRecvCount = static_cast<int> (offsets.back());
230  if (stringRecvCount > 0)
231  {
232  receive (comm, sourceRank, stringRecvCount, &packedString[0]);
233  unpackStringsAfterReceive (strings, packedString, offsets);
234  }
235  }
236 
237 
238  void
239  broadcastStringsHelper (const Comm<int>& comm,
240  const int myRank,
241  const int left,
242  const int right,
243  Array<std::string>& globalNames)
244  {
245  // If left >= right, there is only one process, so we don't need
246  // to do anything.
247  //
248  // If left < right, then split the inclusive interval [left,
249  // right] into [left, mid-1] and [mid, right]. Send from left
250  // to mid, and recurse on the two subintervals.
251  if (left >= right)
252  return;
253  else
254  {
255  const int mid = left + (right - left + 1) / 2;
256 
257  // This could be optimized further on the sending rank, by
258  // prepacking the strings so that they don't have to be
259  // packed more than once.
260  if (myRank == left)
261  sendStrings (comm, globalNames, mid);
262  else if (myRank == mid)
263  receiveStrings (comm, left, globalNames);
264 
265  // Recurse on [left, mid-1] or [mid, right], depending on myRank.
266  if (myRank >= left && myRank <= mid-1)
267  broadcastStringsHelper (comm, myRank, left, mid-1, globalNames);
268  else if (myRank >= mid && myRank <= right)
269  broadcastStringsHelper (comm, myRank, mid, right, globalNames);
270  else
271  return; // Don't recurse if not participating.
272  }
273  }
274 
275 
276  void
277  broadcastStrings (const Comm<int>& comm,
278  Array<std::string>& globalNames)
279  {
280  const int myRank = comm.getRank();
281  const int left = 0;
282  const int right = comm.getSize() - 1;
283 
284  broadcastStringsHelper (comm, myRank, left, right, globalNames);
285  }
286 
287  // \brief Helper function for \c mergeCounterNamesHelper().
288  //
289  // The \c mergeCounterNamesHelper() function implements (using a
290  // parallel reduction) the set union resp. intersection (depending
291  // on the \c setOp argument) of the MPI process' sets of counter
292  // names. This function implements the binary associative
293  // operator which computes the set union resp. intersection of two
294  // sets: the "left" process' intermediate reduction result (global
295  // counter names), and the "mid" process' local counter names.
296  //
297  // \param comm [in] Communicator for which \c
298  // mergeCounterNamesHelper() was called.
299  //
300  // \param myRank [in] Rank of the calling MPI process; must be
301  // either == left or == mid.
302  //
303  // \param left [in] The "left" input argument of
304  // \c mergeCounterNamesHelper().
305  //
306  // \param mid [in] The value of "mid" in the implementation of \c
307  // mergeCounterNamesHelper().
308  //
309  // \param globalNames [in/out] Only accessed if myRank == left.
310  // If so, on input: the intermediate reduction result of the
311  // union resp. intersection (depending on \c setOp). On output:
312  // the union resp. intersection of the input value of the "left"
313  // MPI process' globalNames with the "mid" MPI process'
314  // localNames.
315  //
316  // \param setOp [in] If Intersection, compute the set intersection
317  // of counter names, else if Union, compute the set union of
318  // counter names.
319  void
320  mergeCounterNamesPair (const Comm<int>& comm,
321  const int myRank,
322  const int left,
323  const int mid,
324  Array<std::string>& globalNames,
325  const ECounterSetOp setOp)
326  {
327  using std::cerr;
328  using std::endl;
329  using std::string;
330 
331  const bool debug = false;
332 
333  if (myRank == left)
334  { // Receive names from the other process, and merge its names
335  // with the names on this process.
336  Array<string> otherNames;
337  receiveStrings (comm, mid, otherNames);
338 
339  if (debug)
340  {
341  // Buffering locally in an ostringstream before writing to
342  // the shared stderr sometimes helps avoid interleaved
343  // debugging output.
344  std::ostringstream out;
345  out << "Proc " << myRank << ": in mergePair: otherNames = [";
346  for (Array<std::string>::const_iterator it = otherNames.begin();
347  it != otherNames.end(); ++it)
348  {
349  out << "\"" << *it << "\"";
350  if (it + 1 != otherNames.end())
351  out << ", ";
352  }
353  out << "]" << endl;
354  cerr << out.str();
355  }
356 
357  // Assume that both globalNames and otherNames are sorted.
358  // Compute the set intersection / union as specified by the
359  // enum.
360  Array<string> newNames;
361  if ( std::is_sorted(globalNames.begin(), globalNames.end()) &&
362  std::is_sorted(otherNames.begin(), otherNames.end())) {
363  if (setOp == Intersection)
364  std::set_intersection (globalNames.begin(), globalNames.end(),
365  otherNames.begin(), otherNames.end(),
366  std::back_inserter (newNames));
367  else if (setOp == Union)
368  std::set_union (globalNames.begin(), globalNames.end(),
369  otherNames.begin(), otherNames.end(),
370  std::back_inserter (newNames));
371  else
372  TEUCHOS_TEST_FOR_EXCEPTION(setOp != Intersection && setOp != Union,
373  std::logic_error,
374  "Invalid set operation enum value. Please "
375  "report this bug to the Teuchos developers.");
376  globalNames.swap (newNames);
377  } else { // Need a brute force merge
378  unsortedMergePair(otherNames, globalNames, setOp);
379  }
380  }
381  else if (myRank == mid)
382  sendStrings (comm, globalNames, left);
383  else
384  TEUCHOS_TEST_FOR_EXCEPTION(myRank != left && myRank != mid,
385  std::logic_error,
386  "myRank=" << myRank << " is neither left=" << left
387  << " nor mid=" << mid << ". Please report this "
388  "bug to the Teuchos developers.");
389  }
390 
391 
392  // Recursive helper function for \c mergeCounterNames().
393  //
394  // This function implements the set union resp. intersection
395  // (depending on the \c setOp argument) of the MPI process' sets
396  // of counter names, using a parallel reduction. (Since the
397  // Teuchos comm wrappers as of 11 July 2011 lack a wrapper for
398  // MPI_Reduce(), we hand-roll the reduction using a binary tree
399  // via recursion. We don't need an all-reduce in this case.)
400  void
401  mergeCounterNamesHelper (const Comm<int>& comm,
402  const int myRank,
403  const int left,
404  const int right, // inclusive range [left, right]
405  const Array<std::string>& localNames,
406  Array<std::string>& globalNames,
407  const ECounterSetOp setOp)
408  {
409  // Correctness proof:
410  //
411  // 1. Both set intersection and set union are associative (and
412  // indeed even commutative) operations.
413  // 2. mergeCounterNamesHelper() is just a reduction by binary tree.
414  // 3. Reductions may use any tree shape as long as the binary
415  // operation is associative.
416  //
417  // Recursive "reduction" algorithm:
418  //
419  // Let mid be the midpoint of the inclusive interval [left,
420  // right]. If the (intersection, union) of [left, mid-1] and
421  // the (intersection, union) of [mid, right] are both computed
422  // correctly, then the (intersection, union) of these two sets
423  // is the (intersection, union) of [left, right].
424  //
425  // The base case is left == right: the (intersection, union) of
426  // one set is simply that set, so copy localNames into
427  // globalNames.
428  //
429  // We include another base case for safety: left > right,
430  // meaning that the set of processes is empty, so we do nothing
431  // (the (intersection, union) of an empty set of sets is the
432  // empty set).
433  if (left > right)
434  return;
435  else if (left == right)
436  {
437  Array<string> newNames;
438  newNames.reserve (localNames.size());
439  std::copy (localNames.begin(), localNames.end(),
440  std::back_inserter (newNames));
441  globalNames.swap (newNames);
442  }
443  else
444  { // You're sending messages across the network, so don't
445  // bother to optimize away a few branches here.
446  //
447  // Recurse on [left, mid-1] or [mid, right], depending on myRank.
448  const int mid = left + (right - left + 1) / 2;
449  if (myRank >= left && myRank <= mid-1)
450  mergeCounterNamesHelper (comm, myRank, left, mid-1,
451  localNames, globalNames, setOp);
452  else if (myRank >= mid && myRank <= right)
453  mergeCounterNamesHelper (comm, myRank, mid, right,
454  localNames, globalNames, setOp);
455  else
456  return; // Don't call mergeCounterNamesPair() if not participating.
457 
458  // Combine the results of the recursive step.
459  if (myRank == left || myRank == mid)
460  mergeCounterNamesPair (comm, myRank, left, mid,
461  globalNames, setOp);
462  }
463  }
464 
465  } // namespace (anonymous)
466 
477  void unsortedMergePair(const Array<std::string>& localNames,
478  Array<std::string>& globalNames,
479  const ECounterSetOp setOp){
480  if (setOp == Union) {
481  for (int i=0; i<localNames.size();++i) {
482  bool found=false;
483  // If the name is not in globalNames add it
484  for (int j=0;j<globalNames.size() && !found; ++j)
485  if (localNames[i] == globalNames[j])
486  found=true;
487  if (!found)
488  globalNames.push_back(localNames[i]);
489  }
490  } else if (setOp == Intersection) {
491  for (int i=0; i<globalNames.size();++i) {
492  bool found=false;
493  // If the name is not in localNames remove it
494  for (int j=0;j<localNames.size() && !found; ++j)
495  if (localNames[j] == globalNames[i])
496  found=true;
497  if (!found) {
498  globalNames.remove(i);
499  --i;
500  }
501  }
502  } else
503  TEUCHOS_TEST_FOR_EXCEPTION(setOp != Intersection && setOp != Union,
504  std::logic_error,
505  "Invalid set operation enum value. Please "
506  "report this bug to the Teuchos developers.");
507  }
508 
509 
510  void
512  const Array<std::string>& localNames,
513  Array<std::string>& globalNames,
514  const ECounterSetOp setOp)
515  {
516  const int myRank = comm.getRank();
517  const int left = 0;
518  const int right = comm.getSize() - 1;
519  Array<std::string> theGlobalNames;
520  mergeCounterNamesHelper (comm, myRank, left, right,
521  localNames, theGlobalNames, setOp);
522 
523  // Proc 0 has the list of counter names. Now broadcast it back to
524  // all the procs.
525  broadcastStrings (comm, theGlobalNames);
526 
527  // "Transactional" semantics ensure strong exception safety for
528  // output.
529  globalNames.swap (theGlobalNames);
530  }
531 
532 } // namespace Teuchos
void remove(int i)
Remove the i-th element from the array, with optional boundschecking.
virtual int getSize() const =0
Returns the number of processes that make up this communicator.
virtual int getRank() const =0
Returns the rank of this process.
#define TEUCHOS_TEST_FOR_EXCEPTION(throw_exception_test, Exception, msg)
Macro for throwing an exception with breakpointing to ease debugging.
static Teuchos::RCP< const Comm< OrdinalType > > getComm()
Return the default global communicator.
void send(const Packet sendBuffer[], const Ordinal count, const int destRank, const int tag, const Comm< Ordinal > &comm)
Variant of send() that takes a tag (and restores the correct order of arguments). ...
void mergeCounterNames(const Comm< int > &comm, const Array< std::string > &localNames, Array< std::string > &globalNames, const ECounterSetOp setOp)
Merge counter names over all processors.
friend void swap(Array< T2 > &a1, Array< T2 > &a2)
void unsortedMergePair(const Array< std::string > &localNames, Array< std::string > &globalNames, const ECounterSetOp setOp)
std::vector< std::string >::const_iterator const_iterator
The type of a const forward iterator.
void push_back(const value_type &x)
size_type size() const
Common capabilities for collecting and reporting performance data collectively across MPI processes...
ECounterSetOp
Set operation type for mergeCounterNames() to perform.