]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Avoid MPI calls in ~TimerOutput. 5558/head
authorDavid Wells <wellsd2@rpi.edu>
Thu, 30 Nov 2017 17:52:20 +0000 (12:52 -0500)
committerDavid Wells <wellsd2@rpi.edu>
Wed, 6 Dec 2017 13:36:24 +0000 (08:36 -0500)
~TimerOutput prints timing information after communicating across an MPI
communicator. This MPI call can lead to deadlocks if some, but not all,
proceses throw an exception. Get around this by skipping inter-process
MPI calls when there is an uncaught exception.

include/deal.II/base/timer.h
source/base/timer.cc
tests/base/timer_07.cc [new file with mode: 0644]
tests/base/timer_07.with_mpi=true.output [new file with mode: 0644]

index ae28d5c665d39234547a32ce1dbd5a91593ee1d0..6ba90bf7fe41a702ba658ad2bc65ad897888e40b 100644 (file)
@@ -1003,16 +1003,7 @@ TimerOutput::Scope::Scope(dealii::TimerOutput &timer_,
   timer.enter_section(section_name);
 }
 
-inline
-TimerOutput::Scope::~Scope()
-{
-  try
-    {
-      stop();
-    }
-  catch (...)
-    {}
-}
+
 
 inline
 void
index bd382920c8958d2f77bcf38d2d3eec51fd0eaf77..3168c461e720f9270b0174d01ebbb7a0eae52250 100644 (file)
@@ -349,17 +349,48 @@ TimerOutput::TimerOutput (MPI_Comm      mpi_communicator,
 
 TimerOutput::~TimerOutput()
 {
-  try
+  auto do_exit = [this]()
+  {
+    try
+      {
+        while (active_sections.size() > 0)
+          leave_subsection();
+      }
+    catch (...)
+      {}
+
+    if ( (output_frequency == summary || output_frequency == every_call_and_summary)
+         && output_is_enabled == true)
+      print_summary();
+  };
+
+  // avoid communicating with other processes if there is an uncaught
+  // exception
+#ifdef DEAL_II_WITH_MPI
+  if (std::uncaught_exception() && mpi_communicator != MPI_COMM_SELF)
     {
-      while (active_sections.size() > 0)
-        leave_subsection();
+      std::cerr << "---------------------------------------------------------"
+                << std::endl
+                << "TimerOutput objects finalize timed values printed to the"
+                << std::endl
+                << "screen by communicating over MPI in their destructors."
+                << std::endl
+                << "Since an exception is currently uncaught, this"
+                << std::endl
+                << "synchronization (and subsequent output) will be skipped to"
+                << std::endl
+                << "avoid a possible deadlock."
+                << std::endl
+                << "---------------------------------------------------------"
+                << std::endl;
     }
-  catch (...)
-    {}
-
-  if ( (output_frequency == summary || output_frequency == every_call_and_summary)
-       && output_is_enabled == true)
-    print_summary();
+  else
+    {
+      do_exit();
+    }
+#else
+  do_exit();
+#endif
 }
 
 
@@ -686,5 +717,15 @@ TimerOutput::reset ()
   timer_all.restart();
 }
 
+TimerOutput::Scope::~Scope()
+{
+  try
+    {
+      stop();
+    }
+  catch (...)
+    {}
+}
+
 
 DEAL_II_NAMESPACE_CLOSE
diff --git a/tests/base/timer_07.cc b/tests/base/timer_07.cc
new file mode 100644 (file)
index 0000000..7fa11dc
--- /dev/null
@@ -0,0 +1,108 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2017 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+#include "../tests.h"
+#include <deal.II/base/timer.h>
+
+#include <sstream>
+
+DeclException0(Timer07Exception);
+
+int main (int argc, char **argv)
+{
+  initlog();
+
+  // capture cerr for testing purposes
+  std::stringstream captured_cerr;
+  std::streambuf *old_cerr = std::cerr.rdbuf(captured_cerr.rdbuf());
+
+  Utilities::MPI::MPI_InitFinalize mpi_initialization(argc, argv);
+  // Test printing from TimerOutput
+  try
+    {
+      std::cerr << "TimerOutput\n";
+      TimerOutput timer_out(MPI_COMM_WORLD,
+                            std::cerr,
+                            TimerOutput::summary,
+                            TimerOutput::cpu_times);
+      timer_out.enter_section("Section1");
+
+      throw Timer07Exception();
+    }
+  catch (const Timer07Exception &exc)
+    {
+    }
+
+  // The Scope should still exit correctly
+  try
+    {
+      std::cerr << "TimerOutput::Scope\n";
+      TimerOutput timer_out(MPI_COMM_WORLD,
+                            std::cerr,
+                            TimerOutput::summary,
+                            TimerOutput::cpu_times);
+      TimerOutput::Scope timer_scope(timer_out, "Section1");
+
+      throw Timer07Exception();
+    }
+  catch (const Timer07Exception &exc)
+    {
+    }
+
+  // Test that no errors are printed for MPI_COMM_SELF since no communication occurs
+  try
+    {
+      std::cerr << "TimerOutput::Scope with MPI_COMM_SELF\n";
+      TimerOutput timer_out(MPI_COMM_SELF,
+                            std::cerr,
+                            TimerOutput::summary,
+                            TimerOutput::cpu_times);
+      TimerOutput::Scope timer_scope(timer_out, "Section1");
+
+      throw Timer07Exception();
+    }
+  catch (const Timer07Exception &exc)
+    {
+    }
+
+  // convert numbers to xs to avoid printing time data
+  auto is_digit = [](const char c) -> bool
+  {
+    return std::isdigit(c);
+  };
+  std::string output = captured_cerr.str();
+  std::string::iterator next_number = std::find_if(output.begin(), output.end(), is_digit);
+  while (next_number != output.end())
+    {
+      // convert everything between the |s to xs so that we have consistent output.
+      const std::string::iterator start_pipe
+        = std::find(std::string::reverse_iterator(next_number),
+                    output.rend(),
+                    '|').base();
+      Assert(start_pipe != output.end(), ExcInternalError());
+      const std::string::iterator end_pipe = std::find(next_number, output.end(), '|');
+      Assert(end_pipe != output.end(), ExcInternalError());
+      Assert(end_pipe - start_pipe > 1, ExcInternalError());
+
+      std::fill(start_pipe + 1, end_pipe - 1, 'x');
+      next_number = std::find_if(next_number, output.end(), is_digit);
+    }
+
+
+  deallog << output << std::endl;
+
+  // restore stream
+  std::cerr.rdbuf(old_cerr);
+}
diff --git a/tests/base/timer_07.with_mpi=true.output b/tests/base/timer_07.with_mpi=true.output
new file mode 100644 (file)
index 0000000..5822b12
--- /dev/null
@@ -0,0 +1,29 @@
+
+DEAL::TimerOutput
+---------------------------------------------------------
+TimerOutput objects finalize timed values printed to the
+screen by communicating over MPI in their destructors.
+Since an exception is currently uncaught, this
+synchronization (and subsequent output) will be skipped to
+avoid a possible deadlock.
+---------------------------------------------------------
+TimerOutput::Scope
+---------------------------------------------------------
+TimerOutput objects finalize timed values printed to the
+screen by communicating over MPI in their destructors.
+Since an exception is currently uncaught, this
+synchronization (and subsequent output) will be skipped to
+avoid a possible deadlock.
+---------------------------------------------------------
+TimerOutput::Scope with MPI_COMM_SELF
+
+
++---------------------------------------------+------------+------------+
+| Total CPU time elapsed since start          | xxxxxxxxxx |            |
+|                                             |            |            |
+| Section                         | no. calls |  CPU time  | % of total |
++---------------------------------+-----------+------------+------------+
+| xxxxxxxxxxxxxxxxxxxxxxxxxxxxxxx | xxxxxxxxx | xxxxxxxxxx | xxxxxxxxxx |
++---------------------------------+-----------+------------+------------+
+
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.