]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Also allow to print quantiles
authorMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Tue, 4 Feb 2020 11:39:56 +0000 (12:39 +0100)
committerMartin Kronbichler <kronbichler@lnm.mw.tum.de>
Fri, 7 Feb 2020 07:16:36 +0000 (08:16 +0100)
include/deal.II/base/timer.h
source/base/timer.cc
tests/base/timer_09.cc
tests/base/timer_09.with_mpi=true.mpirun=2.output
tests/base/timer_10.cc [new file with mode: 0644]
tests/base/timer_10.with_mpi=true.mpirun=2.output [new file with mode: 0644]

index 8d67048f19819b8833b4fbc3f542dc60f432faf5..c5a53a4a27bff4b07bae84cd51abf74dd87b139d 100644 (file)
@@ -584,9 +584,10 @@ private:
  * approach allows to identify on which ranks certain slowdowns occur. In case
  * some imbalance between the MPI ranks from one section to the next can be
  * tolerated, this strategy can hence be advantageous over the barrier variant
- * as it does not synchronize the program in places where it is not
- * necessary. The usage of this variant is to initialize the output object
- * without any native print settings and without communicator,
+ * as it does not synchronize the program in places where it is not necessary,
+ * and rather tries to display the imbalance observed in various phases. In
+ * order to use this variant initialize the output object without any native
+ * print settings and without communicator,
  * @code
  *   TimerOutput timer (pcout,
  *                      TimerOutput::never,
@@ -594,11 +595,16 @@ private:
  * @endcode
  * and then call
  * @code
- *   timer.print_summary_statistics(MPI_COMM_WORLD);
+ *   timer.print_wall_time_statistics(MPI_COMM_WORLD);
  * @endcode
- * Here, the output is written to the <code>pcout</code> object of type
- * ConditionalOStream passed to the constructor, making sure the information
- * is only printed once. See step-67 for an example usage of this variant.
+ * where appropriate. Here, the output is written to the <code>pcout</code>
+ * object of type ConditionalOStream passed to the constructor, making sure
+ * the information is only printed once. See step-67 for an example usage of
+ * this variant. Besides the basic minimum, average, and maximum of times over
+ * all MPI ranks, the TimerOutput::print_wall_time_statistics() function also
+ * takes a second argument to specify output of quantiles, e.g., the time
+ * taken by the 10\% of the slowest and fastest ranks, respectively, to get
+ * additional insight into the statistical distribution.
  *
  * @ingroup utilities
  * @author M. Kronbichler, 2009.
@@ -858,9 +864,18 @@ public:
    * useful information when the TimerOutput object is constructed without an
    * MPI_Comm argument, to let individual sections run without being disturbed
    * by barriers.
+   *
+   * The optional argument `quantile` allows to add two additional columns to
+   * the output in terms of the distribution of run times. If quantile = 0.1,
+   * the value and rank of the 10% lowest data is printed as well as the value
+   * and rank at 90% of the distribution function, in addition to the minimum
+   * and the maximum. The value of `quantile` needs to be between 0 (no
+   * quantiles are printed besides the minimum and maximum) and 0.5 (when the
+   * median is given).
    */
   void
-  print_summary_of_wall_time_statistics(const MPI_Comm mpi_comm) const;
+  print_wall_time_statistics(const MPI_Comm mpi_comm,
+                             const double   print_quantile = 0.) const;
 
   /**
    * By calling this function, all output can be disabled. This function
index 26943b36bbd2ef3b75a6e33b868a6494f454369b..cf358d8941a2463649f7fc6fae3baa40613d711e 100644 (file)
@@ -862,8 +862,8 @@ TimerOutput::print_summary() const
 
 
 void
-TimerOutput::print_summary_of_wall_time_statistics(
-  const MPI_Comm mpi_comm) const
+TimerOutput::print_wall_time_statistics(const MPI_Comm mpi_comm,
+                                        const double   quantile) const
 {
   // we are going to change the precision and width of output below. store the
   // old values so we can restore it later on
@@ -873,6 +873,8 @@ TimerOutput::print_summary_of_wall_time_statistics(
 
   AssertDimension(sections.size(),
                   Utilities::MPI::max(sections.size(), mpi_comm));
+  Assert(quantile >= 0. && quantile <= 0.5,
+         ExcMessage("The quantile must be between 0 and 0.5"));
 
   // get the maximum width among all sections
   unsigned int max_width = 0;
@@ -880,40 +882,138 @@ TimerOutput::print_summary_of_wall_time_statistics(
     max_width =
       std::max(max_width, static_cast<unsigned int>(i.first.length()));
 
-  // 20 is the default width until | character
-  max_width = std::max(max_width + 1, static_cast<unsigned int>(20));
-  const std::string extra_dash  = std::string(max_width - 20, '-');
-  const std::string extra_space = std::string(max_width - 20, ' ');
+  // 17 is the default width until | character
+  max_width = std::max(max_width + 1, static_cast<unsigned int>(17));
+  const std::string extra_dash  = std::string(max_width - 17, '-');
+  const std::string extra_space = std::string(max_width - 17, ' ');
+
+  // function to print data in a nice table
+  const auto print_statistics = [&](const double given_time) {
+    const unsigned int n_ranks = Utilities::MPI::n_mpi_processes(mpi_comm);
+    if (n_ranks == 1 || quantile == 0.)
+      {
+        Utilities::MPI::MinMaxAvg data =
+          Utilities::MPI::min_max_avg(given_time, mpi_comm);
+
+        out_stream << std::setw(10) << std::setprecision(4) << std::right;
+        out_stream << data.min << "s ";
+        out_stream << std::setw(5) << std::right;
+        out_stream << data.min_index << (n_ranks > 99999 ? "" : " ") << "|";
+        out_stream << std::setw(10) << std::setprecision(4) << std::right;
+        out_stream << data.avg << "s |";
+        out_stream << std::setw(10) << std::setprecision(4) << std::right;
+        out_stream << data.max << "s ";
+        out_stream << std::setw(5) << std::right;
+        out_stream << data.max_index << (n_ranks > 99999 ? "" : " ") << "|\n";
+      }
+    else
+      {
+        const unsigned int my_rank = Utilities::MPI::this_mpi_process(mpi_comm);
+        std::vector<double> receive_data(my_rank == 0 ? n_ranks : 0);
+        std::vector<double> result(9);
+#ifdef DEAL_II_WITH_MPI
+        MPI_Gather(&given_time,
+                   1,
+                   MPI_DOUBLE,
+                   receive_data.data(),
+                   1,
+                   MPI_DOUBLE,
+                   0,
+                   mpi_comm);
+        if (my_rank == 0)
+          {
+            // fill the received data in a pair and sort; on the way, also
+            // compute the average
+            std::vector<std::pair<double, unsigned int>> data_rank;
+            data_rank.reserve(n_ranks);
+            for (unsigned int i = 0; i < n_ranks; ++i)
+              {
+                data_rank.emplace_back(receive_data[i], i);
+                result[4] += receive_data[i];
+              }
+            result[4] /= n_ranks;
+            std::sort(data_rank.begin(), data_rank.end());
+
+            const unsigned int quantile_index =
+              static_cast<unsigned int>(std::round(quantile * n_ranks));
+            AssertIndexRange(quantile_index, data_rank.size());
+            result[0] = data_rank[0].first;
+            result[1] = data_rank[0].second;
+            result[2] = data_rank[quantile_index].first;
+            result[3] = data_rank[quantile_index].second;
+            result[5] = data_rank[n_ranks - 1 - quantile_index].first;
+            result[6] = data_rank[n_ranks - 1 - quantile_index].second;
+            result[7] = data_rank[n_ranks - 1].first;
+            result[8] = data_rank[n_ranks - 1].second;
+          }
+        MPI_Bcast(result.data(), 9, MPI_DOUBLE, 0, mpi_comm);
+#endif
+        out_stream << std::setw(10) << std::setprecision(4) << std::right;
+        out_stream << result[0] << "s ";
+        out_stream << std::setw(5) << std::right;
+        out_stream << static_cast<unsigned int>(result[1])
+                   << (n_ranks > 99999 ? "" : " ") << "|";
+        out_stream << std::setw(10) << std::setprecision(4) << std::right;
+        out_stream << result[2] << "s ";
+        out_stream << std::setw(5) << std::right;
+        out_stream << static_cast<unsigned int>(result[3])
+                   << (n_ranks > 99999 ? "" : " ") << "|";
+        out_stream << std::setw(10) << std::setprecision(4) << std::right;
+        out_stream << result[4] << "s |";
+        out_stream << std::setw(10) << std::setprecision(4) << std::right;
+        out_stream << result[5] << "s ";
+        out_stream << std::setw(5) << std::right;
+        out_stream << static_cast<unsigned int>(result[6])
+                   << (n_ranks > 99999 ? "" : " ") << "|";
+        out_stream << std::setw(10) << std::setprecision(4) << std::right;
+        out_stream << result[7] << "s ";
+        out_stream << std::setw(5) << std::right;
+        out_stream << static_cast<unsigned int>(result[8])
+                   << (n_ranks > 99999 ? "" : " ") << "|\n";
+      }
+  };
 
   // in case we want to write out wallclock times
   {
-    double total_wall_time = timer_all.wall_time();
-
-    Utilities::MPI::MinMaxAvg data =
-      Utilities::MPI::min_max_avg(total_wall_time, mpi_comm);
     const unsigned int n_ranks = Utilities::MPI::n_mpi_processes(mpi_comm);
 
+    const std::string time_rank_column = "------------------+";
+    const std::string time_rank_space  = "                  |";
+
     // now generate a nice table
     out_stream << "\n\n"
-               << "+---------------------------------" << extra_dash
-               << "+------------+------+------------+------------+------+\n"
-               << "| Total wallclock time elapsed    " << extra_space << "|";
-    out_stream << std::setw(10) << std::setprecision(4) << std::right;
-    out_stream << data.min << "s |";
-    out_stream << std::setw(5) << std::right;
-    out_stream << data.min_index << (n_ranks > 99999 ? "" : " ") << "|";
-    out_stream << std::setw(10) << std::setprecision(4) << std::right;
-    out_stream << data.avg << "s |";
-    out_stream << std::setw(10) << std::setprecision(4) << std::right;
-    out_stream << data.max << "s |";
-    out_stream << std::setw(5) << std::right;
-    out_stream << data.max_index << (n_ranks > 99999 ? "" : " ") << "|\n";
-    out_stream << "|                                 " << extra_space
-               << "|            |  min |            |            |  max |\n";
-    out_stream << "| Section             " << extra_space << "| no. calls "
-               << "|   min time | rank |   avg time |   max time | rank |\n";
-    out_stream << "+---------------------" << extra_dash << "+-----------+"
-               << "------------+------+------------+------------+------+\n";
+               << "+------------------------------" << extra_dash << "+"
+               << time_rank_column
+               << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
+               << "------------+"
+               << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
+               << time_rank_column << "\n"
+               << "| Total wallclock time elapsed " << extra_space << "|";
+
+    print_statistics(timer_all.wall_time());
+
+    out_stream << "|                              " << extra_space << "|"
+               << time_rank_space
+               << (n_ranks > 1 && quantile > 0. ? time_rank_space : "")
+               << "             "
+               << (n_ranks > 1 && quantile > 0. ? time_rank_space : "")
+               << time_rank_space << "\n";
+    out_stream << "| Section          " << extra_space << "| no. calls "
+               << "|   min time  rank |";
+    if (n_ranks > 1 && quantile > 0.)
+      out_stream << " " << std::setw(5) << std::setprecision(2) << std::right
+                 << quantile << "-tile  rank |";
+    out_stream << "   avg time |";
+    if (n_ranks > 1 && quantile > 0.)
+      out_stream << " " << std::setw(5) << std::setprecision(2) << std::right
+                 << 1. - quantile << "-tile  rank |";
+    out_stream << "   max time  rank |\n";
+    out_stream << "+------------------------------" << extra_dash << "+"
+               << time_rank_column
+               << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
+               << "------------+"
+               << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
+               << time_rank_column << "\n";
     for (const auto &i : sections)
       {
         std::string name_out = i.first;
@@ -926,21 +1026,15 @@ TimerOutput::print_summary_of_wall_time_statistics(
         out_stream << "| ";
         out_stream << std::setw(9);
         out_stream << i.second.n_calls << " |";
-        data = Utilities::MPI::min_max_avg(i.second.total_wall_time, mpi_comm);
-        out_stream << std::setw(10) << std::setprecision(4) << std::right;
-        out_stream << data.min << "s |";
-        out_stream << std::setw(5) << std::right;
-        out_stream << data.min_index << (n_ranks > 99999 ? "" : " ") << "|";
-        out_stream << std::setw(10) << std::setprecision(4) << std::right;
-        out_stream << data.avg << "s |";
-        out_stream << std::setw(10) << std::setprecision(4) << std::right;
-        out_stream << data.max << "s |";
-        out_stream << std::setw(5) << std::right;
-        out_stream << data.max_index << (n_ranks > 99999 ? "" : " ") << "|\n";
+
+        print_statistics(i.second.total_wall_time);
       }
-    out_stream << "+---------------------" << extra_dash << "+-----------+"
-               << "------------+------+------------+------------+------+\n"
-               << std::endl;
+    out_stream << "+------------------------------" << extra_dash << "+"
+               << time_rank_column
+               << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
+               << "------------+"
+               << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
+               << time_rank_column << "\n";
   }
 
   // restore previous precision and width
index a021657983001cbc520a280382d2daa6d3f424c6..27e02b04b9c0ab811f7963cefa67adea759b2182 100644 (file)
@@ -1,6 +1,6 @@
 // ---------------------------------------------------------------------
 //
-// Copyright (C) 2018 by the deal.II authors
+// Copyright (C) 2020 by the deal.II authors
 //
 // This file is part of the deal.II library.
 //
@@ -14,7 +14,7 @@
 // ---------------------------------------------------------------------
 
 
-// test TimerOutput::print_summary_of_wall_time_statistics()
+// test TimerOutput::print_wall_time_statistics()
 
 #include <deal.II/base/timer.h>
 
@@ -48,7 +48,7 @@ test()
   burn(50);
   t.leave_subsection("hi");
 
-  t.print_summary_of_wall_time_statistics(MPI_COMM_WORLD);
+  t.print_wall_time_statistics(MPI_COMM_WORLD);
 
   t.enter_subsection(
     "this is a very long section name that previously did not work");
@@ -56,7 +56,7 @@ test()
   t.leave_subsection(
     "this is a very long section name that previously did not work");
 
-  t.print_summary_of_wall_time_statistics(MPI_COMM_WORLD);
+  t.print_wall_time_statistics(MPI_COMM_WORLD);
 
   std::string s = ss.str();
   std::replace_if(s.begin(), s.end(), ::isdigit, ' ');
index b6815a7c8e2b2c616598a87f2e956d1131f50f14..e3178fe9282e443e8b3f56ebaae8b5426bf7a424 100644 (file)
@@ -1,24 +1,22 @@
 
 DEAL::
 
-+---------------------------------+------------+------+------------+------------+------+
-| Total wallclock time elapsed    |          s |      |          s |          s |      |
-|                                 |            |  min |            |            |  max |
-| Section             | no  calls |   min time | rank |   avg time |   max time | rank |
-+---------------------+-----------+------------+------+------------+------------+------+
-| hi                  |           |          s |      |          s |          s |      |
-+---------------------+-----------+------------+------+------------+------------+------+
++------------------------------+------------------+------------+------------------+
+| Total wallclock time elapsed |          s       |          s |          s       |
+|                              |                  |                               |
+| Section          | no  calls |   min time  rank |   avg time |   max time  rank |
++------------------------------+------------------+------------+------------------+
+| hi               |           |          s       |          s |          s       |
++------------------------------+------------------+------------+------------------+
 
 
-
-+---------------------------------------------------------------------------+------------+------+------------+------------+------+
-| Total wallclock time elapsed                                              |          s |      |          s |          s |      |
-|                                                                           |            |  min |            |            |  max |
-| Section                                                       | no  calls |   min time | rank |   avg time |   max time | rank |
-+---------------------------------------------------------------+-----------+------------+------+------------+------------+------+
-| hi                                                            |           |          s |      |          s |          s |      |
-| this is a very long section name that previously did not work |           |          s |      |          s |          s |      |
-+---------------------------------------------------------------+-----------+------------+------+------------+------------+------+
-
++---------------------------------------------------------------------------+------------------+------------+------------------+
+| Total wallclock time elapsed                                              |          s       |          s |          s       |
+|                                                                           |                  |                               |
+| Section                                                       | no  calls |   min time  rank |   avg time |   max time  rank |
++---------------------------------------------------------------------------+------------------+------------+------------------+
+| hi                                                            |           |          s       |          s |          s       |
+| this is a very long section name that previously did not work |           |          s       |          s |          s       |
++---------------------------------------------------------------------------+------------------+------------+------------------+
 
 DEAL::
diff --git a/tests/base/timer_10.cc b/tests/base/timer_10.cc
new file mode 100644 (file)
index 0000000..6f82e53
--- /dev/null
@@ -0,0 +1,77 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2020 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// test TimerOutput::print_wall_time_statistics() with optional
+// quantile argument
+
+#include <deal.II/base/timer.h>
+
+#include <algorithm>
+#include <sstream>
+
+#include "../tests.h"
+
+// burn computer time
+double s = 0.;
+void
+burn(unsigned int n)
+{
+  for (unsigned int i = 0; i < n; ++i)
+    {
+      for (unsigned int j = 1; j < 100000; ++j)
+        {
+          s += 1. / j * i;
+        }
+    }
+}
+
+void
+test()
+{
+  std::stringstream ss;
+
+  TimerOutput t(ss, TimerOutput::never, TimerOutput::wall_times);
+
+  t.enter_subsection("hi");
+  burn(50);
+  t.leave_subsection("hi");
+
+  t.print_wall_time_statistics(MPI_COMM_WORLD, 0.1);
+
+  t.enter_subsection(
+    "this is a very long section name that previously did not work");
+  burn(50);
+  t.leave_subsection(
+    "this is a very long section name that previously did not work");
+
+  t.print_wall_time_statistics(MPI_COMM_WORLD, 0.1);
+
+  std::string s = ss.str();
+  std::replace_if(s.begin(), s.end(), ::isdigit, ' ');
+  std::replace_if(s.begin(), s.end(), [](char x) { return x == '.'; }, ' ');
+
+  deallog << s << std::endl << std::endl;
+}
+
+int
+main(int argc, char **argv)
+{
+  Utilities::MPI::MPI_InitFinalize mpi(argc, argv);
+
+  mpi_initlog();
+
+  test();
+}
diff --git a/tests/base/timer_10.with_mpi=true.mpirun=2.output b/tests/base/timer_10.with_mpi=true.mpirun=2.output
new file mode 100644 (file)
index 0000000..7340772
--- /dev/null
@@ -0,0 +1,22 @@
+
+DEAL::
+
++------------------------------+------------------+------------------+------------+------------------+------------------+
+| Total wallclock time elapsed |          s       |          s       |          s |          s       |          s       |
+|                              |                  |                  |                               |                  |
+| Section          | no  calls |   min time  rank |      -tile  rank |   avg time |      -tile  rank |   max time  rank |
++------------------------------+------------------+------------------+------------+------------------+------------------+
+| hi               |           |          s       |          s       |          s |          s       |          s       |
++------------------------------+------------------+------------------+------------+------------------+------------------+
+
+
++---------------------------------------------------------------------------+------------------+------------------+------------+------------------+------------------+
+| Total wallclock time elapsed                                              |          s       |          s       |          s |          s       |          s       |
+|                                                                           |                  |                  |                               |                  |
+| Section                                                       | no  calls |   min time  rank |      -tile  rank |   avg time |      -tile  rank |   max time  rank |
++---------------------------------------------------------------------------+------------------+------------------+------------+------------------+------------------+
+| hi                                                            |           |          s       |          s       |          s |          s       |          s       |
+| this is a very long section name that previously did not work |           |          s       |          s       |          s |          s       |          s       |
++---------------------------------------------------------------------------+------------------+------------------+------------+------------------+------------------+
+
+DEAL::

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.