Reference documentation for deal.II version GIT 29f9da0a34 2023-12-07 10:00:01+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
timer.cc
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 1998 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
17 #include <deal.II/base/mpi.h>
19 #include <deal.II/base/timer.h>
20 #include <deal.II/base/utilities.h>
21 
22 #include <boost/io/ios_state.hpp>
23 
24 #include <algorithm>
25 #include <chrono>
26 #include <iomanip>
27 #include <iostream>
28 #include <map>
29 #include <sstream>
30 #include <string>
31 #include <type_traits>
32 
33 #ifdef DEAL_II_HAVE_SYS_RESOURCE_H
34 # include <sys/resource.h>
35 #endif
36 
37 #ifdef DEAL_II_MSVC
38 # include <windows.h>
39 #endif
40 
41 
42 
44 
45 namespace internal
46 {
47  namespace TimerImplementation
48  {
49  namespace
50  {
55  template <typename T>
56  struct is_duration : std::false_type
57  {};
58 
62  template <typename Rep, typename Period>
63  struct is_duration<std::chrono::duration<Rep, Period>> : std::true_type
64  {};
65 
71  template <typename T>
72  T
73  from_seconds(const double time)
74  {
75  static_assert(is_duration<T>::value,
76  "The template type should be a duration type.");
77  return T(std::lround(T::period::den * (time / T::period::num)));
78  }
79 
84  template <typename Rep, typename Period>
85  double
86  to_seconds(const std::chrono::duration<Rep, Period> duration)
87  {
88  return Period::num * double(duration.count()) / Period::den;
89  }
90 
94  void
95  clear_timing_data(Utilities::MPI::MinMaxAvg &data)
96  {
97  data.sum = numbers::signaling_nan<double>();
98  data.min = numbers::signaling_nan<double>();
99  data.max = numbers::signaling_nan<double>();
100  data.avg = numbers::signaling_nan<double>();
103  }
104  } // namespace
105  } // namespace TimerImplementation
106 } // namespace internal
107 
108 
109 
111 CPUClock::now() noexcept
112 {
113  double system_cpu_duration = 0.0;
114 #ifdef DEAL_II_MSVC
115  FILETIME cpuTime, sysTime, createTime, exitTime;
116  const auto succeeded = GetProcessTimes(
117  GetCurrentProcess(), &createTime, &exitTime, &sysTime, &cpuTime);
118  if (succeeded)
119  {
120  system_cpu_duration =
121  (double)(((unsigned long long)cpuTime.dwHighDateTime << 32) |
122  cpuTime.dwLowDateTime) /
123  1e7;
124  }
125  // keep the zero value if GetProcessTimes didn't work
126 #elif defined(DEAL_II_HAVE_SYS_RESOURCE_H)
127  rusage usage;
128  getrusage(RUSAGE_SELF, &usage);
129  system_cpu_duration = usage.ru_utime.tv_sec + 1.e-6 * usage.ru_utime.tv_usec;
130 #else
131  DEAL_II_WARNING("Unsupported platform. Porting not finished.")
132 #endif
133  return time_point(
134  internal::TimerImplementation::from_seconds<duration>(system_cpu_duration));
135 }
136 
137 
138 
139 template <typename clock_type_>
141  : current_lap_start_time(clock_type::now())
142  , accumulated_time(duration_type::zero())
143  , last_lap_time(duration_type::zero())
144 {}
145 
146 
147 
148 template <typename clock_type_>
149 void
151 {
152  current_lap_start_time = clock_type::now();
153  accumulated_time = duration_type::zero();
154  last_lap_time = duration_type::zero();
155 }
156 
157 
158 
160  : Timer(MPI_COMM_SELF, /*sync_lap_times=*/false)
161 {}
162 
163 
164 
165 Timer::Timer(const MPI_Comm mpi_communicator, const bool sync_lap_times_)
166  : running(false)
167  , mpi_communicator(mpi_communicator)
168  , sync_lap_times(sync_lap_times_)
169 {
170  reset();
171  start();
172 }
173 
174 
175 
176 void
178 {
179  running = true;
180 #ifdef DEAL_II_WITH_MPI
181  if (sync_lap_times)
182  {
183  const int ierr = MPI_Barrier(mpi_communicator);
184  AssertThrowMPI(ierr);
185  }
186 #endif
187  wall_times.current_lap_start_time = wall_clock_type::now();
188  cpu_times.current_lap_start_time = cpu_clock_type::now();
189 }
190 
191 
192 
193 double
195 {
196  if (running)
197  {
198  running = false;
199 
201  wall_clock_type::now() - wall_times.current_lap_start_time;
202  cpu_times.last_lap_time =
203  cpu_clock_type::now() - cpu_times.current_lap_start_time;
204 
206  Utilities::MPI::min_max_avg(internal::TimerImplementation::to_seconds(
209  if (sync_lap_times)
210  {
212  internal::TimerImplementation::from_seconds<
213  decltype(wall_times)::duration_type>(last_lap_wall_time_data.max);
214  cpu_times.last_lap_time = internal::TimerImplementation::from_seconds<
215  decltype(cpu_times)::duration_type>(
217  internal::TimerImplementation::to_seconds(
218  cpu_times.last_lap_time),
220  .max);
221  }
223  cpu_times.accumulated_time += cpu_times.last_lap_time;
225  Utilities::MPI::min_max_avg(internal::TimerImplementation::to_seconds(
228  }
229  return internal::TimerImplementation::to_seconds(cpu_times.accumulated_time);
230 }
231 
232 
233 
234 double
236 {
237  if (running)
238  {
239  const double running_time = internal::TimerImplementation::to_seconds(
240  cpu_clock_type::now() - cpu_times.current_lap_start_time +
241  cpu_times.accumulated_time);
242  return Utilities::MPI::sum(running_time, mpi_communicator);
243  }
244  else
245  {
246  return Utilities::MPI::sum(internal::TimerImplementation::to_seconds(
247  cpu_times.accumulated_time),
249  }
250 }
251 
252 
253 
254 double
256 {
257  return internal::TimerImplementation::to_seconds(cpu_times.last_lap_time);
258 }
259 
260 
261 
262 double
264 {
265  wall_clock_type::duration current_elapsed_wall_time;
266  if (running)
267  current_elapsed_wall_time = wall_clock_type::now() -
270  else
271  current_elapsed_wall_time = wall_times.accumulated_time;
272 
273  return internal::TimerImplementation::to_seconds(current_elapsed_wall_time);
274 }
275 
276 
277 
278 double
280 {
281  return internal::TimerImplementation::to_seconds(wall_times.last_lap_time);
282 }
283 
284 
285 
286 void
288 {
290  cpu_times.reset();
291  running = false;
292  internal::TimerImplementation::clear_timing_data(last_lap_wall_time_data);
293  internal::TimerImplementation::clear_timing_data(accumulated_wall_time_data);
294 }
295 
296 
297 
298 /* ---------------------------- TimerOutput -------------------------- */
299 
300 TimerOutput::TimerOutput(std::ostream &stream,
301  const OutputFrequency output_frequency,
302  const OutputType output_type)
303  : output_frequency(output_frequency)
304  , output_type(output_type)
305  , out_stream(stream, true)
306  , output_is_enabled(true)
307  , mpi_communicator(MPI_COMM_SELF)
308 {}
309 
310 
311 
313  const OutputFrequency output_frequency,
314  const OutputType output_type)
315  : output_frequency(output_frequency)
316  , output_type(output_type)
317  , out_stream(stream)
318  , output_is_enabled(true)
319  , mpi_communicator(MPI_COMM_SELF)
320 {}
321 
322 
323 
324 TimerOutput::TimerOutput(const MPI_Comm mpi_communicator,
325  std::ostream &stream,
326  const OutputFrequency output_frequency,
327  const OutputType output_type)
328  : output_frequency(output_frequency)
329  , output_type(output_type)
330  , out_stream(stream, true)
331  , output_is_enabled(true)
332  , mpi_communicator(mpi_communicator)
333 {}
334 
335 
336 
337 TimerOutput::TimerOutput(const MPI_Comm mpi_communicator,
338  ConditionalOStream &stream,
339  const OutputFrequency output_frequency,
340  const OutputType output_type)
341  : output_frequency(output_frequency)
342  , output_type(output_type)
343  , out_stream(stream)
344  , output_is_enabled(true)
345  , mpi_communicator(mpi_communicator)
346 {}
347 
348 
349 
351 {
352  auto do_exit = [this]() {
353  try
354  {
355  while (active_sections.size() > 0)
357  // don't print unless we leave all subsections
358  if ((output_frequency == summary ||
360  output_is_enabled == true)
361  print_summary();
362  }
363  catch (...)
364  {}
365  };
366 
367  // avoid communicating with other processes if there is an uncaught
368  // exception
369 #ifdef DEAL_II_WITH_MPI
370 # if __cpp_lib_uncaught_exceptions >= 201411
371  // std::uncaught_exception() is deprecated in c++17
372  if (std::uncaught_exceptions() > 0 && mpi_communicator != MPI_COMM_SELF)
373 # else
374  if (std::uncaught_exception() == true && mpi_communicator != MPI_COMM_SELF)
375 # endif
376  {
377  const unsigned int myid =
379  if (myid == 0)
380  std::cerr
381  << "---------------------------------------------------------\n"
382  << "TimerOutput objects finalize timed values printed to the\n"
383  << "screen by communicating over MPI in their destructors.\n"
384  << "Since an exception is currently uncaught, this\n"
385  << "synchronization (and subsequent output) will be skipped\n"
386  << "to avoid a possible deadlock.\n"
387  << "---------------------------------------------------------"
388  << std::endl;
389  }
390  else
391  {
392  do_exit();
393  }
394 #else
395  do_exit();
396 #endif
397 }
398 
399 
400 
401 void
402 TimerOutput::enter_subsection(const std::string &section_name)
403 {
404  std::lock_guard<std::mutex> lock(mutex);
405 
406  Assert(section_name.empty() == false, ExcMessage("Section string is empty."));
407 
408  Assert(std::find(active_sections.begin(),
409  active_sections.end(),
410  section_name) == active_sections.end(),
411  ExcMessage(std::string("Cannot enter the already active section <") +
412  section_name + ">."));
413 
414  if (sections.find(section_name) == sections.end())
415  {
416  if (mpi_communicator != MPI_COMM_SELF)
417  {
418  // create a new timer for this section. the second argument
419  // will ensure that we have an MPI barrier before starting
420  // and stopping a timer, and this ensures that we get the
421  // maximum run time for this section over all processors.
422  // The mpi_communicator from TimerOutput is passed to the
423  // Timer here, so this Timer will collect timing information
424  // among all processes inside mpi_communicator.
425  sections[section_name].timer = Timer(mpi_communicator, true);
426  }
427 
428 
429  sections[section_name].total_cpu_time = 0;
430  sections[section_name].total_wall_time = 0;
431  sections[section_name].n_calls = 0;
432  }
433 
434  sections[section_name].timer.reset();
435  sections[section_name].timer.start();
436  ++sections[section_name].n_calls;
437 
438  active_sections.push_back(section_name);
439 }
440 
441 
442 
443 void
444 TimerOutput::leave_subsection(const std::string &section_name)
445 {
446  Assert(!active_sections.empty(),
447  ExcMessage("Cannot exit any section because none has been entered!"));
448 
449  std::lock_guard<std::mutex> lock(mutex);
450 
451  if (!section_name.empty())
452  {
453  Assert(sections.find(section_name) != sections.end(),
454  ExcMessage("Cannot delete a section that was never created."));
455  Assert(std::find(active_sections.begin(),
456  active_sections.end(),
457  section_name) != active_sections.end(),
458  ExcMessage("Cannot delete a section that has not been entered."));
459  }
460 
461  // if no string is given, exit the last
462  // active section.
463  const std::string actual_section_name =
464  (section_name.empty() ? active_sections.back() : section_name);
465 
466  sections[actual_section_name].timer.stop();
467  sections[actual_section_name].total_wall_time +=
468  sections[actual_section_name].timer.last_wall_time();
469 
470  // Get cpu time. On MPI systems, if constructed with an mpi_communicator
471  // like MPI_COMM_WORLD, then the Timer will sum up the CPU time between
472  // processors among the provided mpi_communicator. Therefore, no
473  // communication is needed here.
474  const double cpu_time = sections[actual_section_name].timer.last_cpu_time();
475  sections[actual_section_name].total_cpu_time += cpu_time;
476 
477  // in case we have to print out something, do that here...
478  if ((output_frequency == every_call ||
480  output_is_enabled == true)
481  {
482  std::string output_time;
483  std::ostringstream cpu;
484  cpu << cpu_time << "s";
485  std::ostringstream wall;
486  wall << sections[actual_section_name].timer.last_wall_time() << "s";
487  if (output_type == cpu_times)
488  output_time = ", CPU time: " + cpu.str();
489  else if (output_type == wall_times)
490  output_time = ", wall time: " + wall.str() + ".";
491  else
492  output_time =
493  ", CPU/wall time: " + cpu.str() + " / " + wall.str() + ".";
494 
495  out_stream << actual_section_name << output_time << std::endl;
496  }
497 
498  // delete the index from the list of
499  // active ones
500  active_sections.erase(std::find(active_sections.begin(),
501  active_sections.end(),
502  actual_section_name));
503 }
504 
505 
506 
507 std::map<std::string, double>
509 {
510  std::map<std::string, double> output;
511  for (const auto &section : sections)
512  {
513  switch (kind)
514  {
515  case TimerOutput::OutputData::total_cpu_time:
516  output[section.first] = section.second.total_cpu_time;
517  break;
518  case TimerOutput::OutputData::total_wall_time:
519  output[section.first] = section.second.total_wall_time;
520  break;
521  case TimerOutput::OutputData::n_calls:
522  output[section.first] = section.second.n_calls;
523  break;
524  default:
525  Assert(false, ExcNotImplemented());
526  }
527  }
528  return output;
529 }
530 
531 
532 
533 void
535 {
536  // we are going to change the precision and width of output below. store the
537  // old values so the get restored when exiting this function
538  const boost::io::ios_base_all_saver restore_stream(out_stream.get_stream());
539 
540  // get the maximum width among all sections
541  unsigned int max_width = 0;
542  for (const auto &i : sections)
543  max_width = std::max(max_width, static_cast<unsigned int>(i.first.size()));
544 
545  // 32 is the default width until | character
546  max_width = std::max(max_width + 1, static_cast<unsigned int>(32));
547  const std::string extra_dash = std::string(max_width - 32, '-');
548  const std::string extra_space = std::string(max_width - 32, ' ');
549 
551  {
552  // in case we want to write CPU times
553  if (output_type != wall_times)
554  {
555  double total_cpu_time =
557 
558  // check that the sum of all times is less or equal than the total
559  // time. otherwise, we might have generated a lot of overhead in this
560  // function.
561  double check_time = 0.;
562  for (const auto &i : sections)
563  check_time += i.second.total_cpu_time;
564 
565  const double time_gap = check_time - total_cpu_time;
566  if (time_gap > 0.0)
567  total_cpu_time = check_time;
568 
569  // generate a nice table
570  out_stream << "\n\n"
571  << "+---------------------------------------------"
572  << extra_dash << "+------------"
573  << "+------------+\n"
574  << "| Total CPU time elapsed since start "
575  << extra_space << "|";
576  out_stream << std::setw(10) << std::setprecision(3) << std::right;
577  out_stream << total_cpu_time << "s | |\n";
578  out_stream << "| "
579  << extra_space << "| "
580  << "| |\n";
581  out_stream << "| Section " << extra_space
582  << "| no. calls |";
583  out_stream << std::setw(10);
584  out_stream << std::setprecision(3);
585  out_stream << " CPU time "
586  << " | % of total |\n";
587  out_stream << "+---------------------------------" << extra_dash
588  << "+-----------+------------"
589  << "+------------+";
590  for (const auto &i : sections)
591  {
592  std::string name_out = i.first;
593 
594  // resize the array so that it is always of the same size
595  unsigned int pos_non_space = name_out.find_first_not_of(' ');
596  name_out.erase(0, pos_non_space);
597  name_out.resize(max_width, ' ');
598  out_stream << std::endl;
599  out_stream << "| " << name_out;
600  out_stream << "| ";
601  out_stream << std::setw(9);
602  out_stream << i.second.n_calls << " |";
603  out_stream << std::setw(10);
604  out_stream << std::setprecision(3);
605  out_stream << i.second.total_cpu_time << "s |";
606  out_stream << std::setw(10);
607  if (total_cpu_time != 0)
608  {
609  // if run time was less than 0.1%, just print a zero to avoid
610  // printing silly things such as "2.45e-6%". otherwise print
611  // the actual percentage
612  const double fraction =
613  i.second.total_cpu_time / total_cpu_time;
614  if (fraction > 0.001)
615  {
616  out_stream << std::setprecision(2);
617  out_stream << fraction * 100;
618  }
619  else
620  out_stream << 0.0;
621 
622  out_stream << "% |";
623  }
624  else
625  out_stream << 0.0 << "% |";
626  }
627  out_stream << std::endl
628  << "+---------------------------------" << extra_dash
629  << "+-----------+"
630  << "------------+------------+\n"
631  << std::endl;
632 
633  if (time_gap > 0.0)
634  out_stream
635  << std::endl
636  << "Note: The sum of counted times is " << time_gap
637  << " seconds larger than the total time.\n"
638  << "(Timer function may have introduced too much overhead, or different\n"
639  << "section timers may have run at the same time.)" << std::endl;
640  }
641 
642  // in case we want to write out wallclock times
643  if (output_type != cpu_times)
644  {
646 
647  // now generate a nice table
648  out_stream << "\n\n"
649  << "+---------------------------------------------"
650  << extra_dash << "+------------"
651  << "+------------+\n"
652  << "| Total wallclock time elapsed since start "
653  << extra_space << "|";
654  out_stream << std::setw(10) << std::setprecision(3) << std::right;
655  out_stream << total_wall_time << "s | |\n";
656  out_stream << "| "
657  << extra_space << "| "
658  << "| |\n";
659  out_stream << "| Section " << extra_space
660  << "| no. calls |";
661  out_stream << std::setw(10);
662  out_stream << std::setprecision(3);
663  out_stream << " wall time | % of total |\n";
664  out_stream << "+---------------------------------" << extra_dash
665  << "+-----------+------------"
666  << "+------------+";
667  for (const auto &i : sections)
668  {
669  std::string name_out = i.first;
670 
671  // resize the array so that it is always of the same size
672  unsigned int pos_non_space = name_out.find_first_not_of(' ');
673  name_out.erase(0, pos_non_space);
674  name_out.resize(max_width, ' ');
675  out_stream << std::endl;
676  out_stream << "| " << name_out;
677  out_stream << "| ";
678  out_stream << std::setw(9);
679  out_stream << i.second.n_calls << " |";
680  out_stream << std::setw(10);
681  out_stream << std::setprecision(3);
682  out_stream << i.second.total_wall_time << "s |";
683  out_stream << std::setw(10);
684 
685  if (total_wall_time != 0)
686  {
687  // if run time was less than 0.1%, just print a zero to avoid
688  // printing silly things such as "2.45e-6%". otherwise print
689  // the actual percentage
690  const double fraction =
691  i.second.total_wall_time / total_wall_time;
692  if (fraction > 0.001)
693  {
694  out_stream << std::setprecision(2);
695  out_stream << fraction * 100;
696  }
697  else
698  out_stream << 0.0;
699 
700  out_stream << "% |";
701  }
702  else
703  out_stream << 0.0 << "% |";
704  }
705  out_stream << std::endl
706  << "+---------------------------------" << extra_dash
707  << "+-----------+"
708  << "------------+------------+\n"
709  << std::endl;
710  }
711  }
712  else
713  // output_type == cpu_and_wall_times_grouped
714  {
715  const double total_wall_time = timer_all.wall_time();
716  double total_cpu_time =
718 
719  // check that the sum of all times is less or equal than the total time.
720  // otherwise, we might have generated a lot of overhead in this function.
721  double check_time = 0.;
722 
723  for (const auto &i : sections)
724  check_time += i.second.total_cpu_time;
725 
726  const double time_gap = check_time - total_cpu_time;
727  if (time_gap > 0.0)
728  total_cpu_time = check_time;
729 
730  // generate a nice table
731  out_stream << "\n\n+---------------------------------------------"
732  << extra_dash << "+"
733  << "------------+------------+"
734  << "------------+------------+" << '\n'
735  << "| Total CPU/wall time elapsed since start "
736  << extra_space << "|" << std::setw(10) << std::setprecision(3)
737  << std::right << total_cpu_time << "s | |"
738  << total_wall_time << "s | |"
739  << "\n| "
740  << extra_space << "|"
741  << " | |"
742  << " | |"
743  << "\n| Section " << extra_space
744  << "| no. calls |"
745  << " CPU time | % of total |"
746  << " wall time | % of total |"
747  << "\n+---------------------------------" << extra_dash
748  << "+-----------+"
749  << "------------+------------+"
750  << "------------+------------+" << std::endl;
751 
752  for (const auto &i : sections)
753  {
754  std::string name_out = i.first;
755 
756  // resize the array so that it is always of the same size
757  unsigned int pos_non_space = name_out.find_first_not_of(' ');
758  name_out.erase(0, pos_non_space);
759  name_out.resize(max_width, ' ');
760  out_stream << "| " << name_out << "| ";
761 
762  out_stream << std::setw(9);
763  out_stream << i.second.n_calls << " |";
764 
765  if (output_type != wall_times)
766  {
767  out_stream << std::setw(10);
768  out_stream << std::setprecision(3);
769  out_stream << i.second.total_cpu_time << "s |";
770  out_stream << std::setw(10);
771  if (total_cpu_time != 0)
772  {
773  // if run time was less than 0.1%, just print a zero to avoid
774  // printing silly things such as "2.45e-6%". otherwise print
775  // the actual percentage
776  const double fraction =
777  i.second.total_cpu_time / total_cpu_time;
778  if (fraction > 0.001)
779  {
780  out_stream << std::setprecision(2);
781  out_stream << fraction * 100;
782  }
783  else
784  out_stream << 0.0;
785 
786  out_stream << "% |";
787  }
788  else
789  out_stream << 0.0 << "% |";
790  }
791 
792  if (output_type != cpu_times)
793  {
794  out_stream << std::setw(10);
795  out_stream << std::setprecision(3);
796  out_stream << i.second.total_wall_time << "s |";
797  out_stream << std::setw(10);
798 
799  if (total_wall_time != 0)
800  {
801  // if run time was less than 0.1%, just print a zero to avoid
802  // printing silly things such as "2.45e-6%". otherwise print
803  // the actual percentage
804  const double fraction =
805  i.second.total_wall_time / total_wall_time;
806  if (fraction > 0.001)
807  {
808  out_stream << std::setprecision(2);
809  out_stream << fraction * 100;
810  }
811  else
812  out_stream << 0.0;
813 
814  out_stream << "% |";
815  }
816  else
817  out_stream << 0.0 << "% |";
818  }
819  out_stream << std::endl;
820  }
821 
822  out_stream << "+---------------------------------" << extra_dash
823  << "+-----------+"
824  << "------------+------------+"
825  << "------------+------------+" << std::endl
826  << std::endl;
827 
828  if (output_type != wall_times && time_gap > 0.0)
829  out_stream
830  << std::endl
831  << "Note: The sum of counted times is " << time_gap
832  << " seconds larger than the total time.\n"
833  << "(Timer function may have introduced too much overhead, or different\n"
834  << "section timers may have run at the same time.)" << std::endl;
835  }
836 }
837 
838 
839 
840 void
842  const double quantile) const
843 {
844  // we are going to change the precision and width of output below. store the
845  // old values so the get restored when exiting this function
846  const boost::io::ios_base_all_saver restore_stream(out_stream.get_stream());
847 
848  AssertDimension(sections.size(),
849  Utilities::MPI::max(sections.size(), mpi_comm));
850  Assert(quantile >= 0. && quantile <= 0.5,
851  ExcMessage("The quantile must be between 0 and 0.5"));
852 
853  // get the maximum width among all sections
854  unsigned int max_width = 0;
855  for (const auto &i : sections)
856  max_width = std::max(max_width, static_cast<unsigned int>(i.first.size()));
857 
858  // 17 is the default width until | character
859  max_width = std::max(max_width + 1, static_cast<unsigned int>(17));
860  const std::string extra_dash = std::string(max_width - 17, '-');
861  const std::string extra_space = std::string(max_width - 17, ' ');
862 
863  // function to print data in a nice table
864  const auto print_statistics = [&](const double given_time) {
865  const unsigned int n_ranks = Utilities::MPI::n_mpi_processes(mpi_comm);
866  if (n_ranks == 1 || quantile == 0.)
867  {
869  Utilities::MPI::min_max_avg(given_time, mpi_comm);
870 
871  out_stream << std::setw(10) << std::setprecision(4) << std::right;
872  out_stream << data.min << "s ";
873  out_stream << std::setw(5) << std::right;
874  out_stream << data.min_index << (n_ranks > 99999 ? "" : " ") << "|";
875  out_stream << std::setw(10) << std::setprecision(4) << std::right;
876  out_stream << data.avg << "s |";
877  out_stream << std::setw(10) << std::setprecision(4) << std::right;
878  out_stream << data.max << "s ";
879  out_stream << std::setw(5) << std::right;
880  out_stream << data.max_index << (n_ranks > 99999 ? "" : " ") << "|\n";
881  }
882  else
883  {
884  const unsigned int my_rank = Utilities::MPI::this_mpi_process(mpi_comm);
885  std::vector<double> receive_data(my_rank == 0 ? n_ranks : 0);
886  std::vector<double> result(9);
887 #ifdef DEAL_II_WITH_MPI
888  int ierr = MPI_Gather(&given_time,
889  1,
890  MPI_DOUBLE,
891  receive_data.data(),
892  1,
893  MPI_DOUBLE,
894  0,
895  mpi_comm);
896  AssertThrowMPI(ierr);
897  if (my_rank == 0)
898  {
899  // fill the received data in a pair and sort; on the way, also
900  // compute the average
901  std::vector<std::pair<double, unsigned int>> data_rank;
902  data_rank.reserve(n_ranks);
903  for (unsigned int i = 0; i < n_ranks; ++i)
904  {
905  data_rank.emplace_back(receive_data[i], i);
906  result[4] += receive_data[i];
907  }
908  result[4] /= n_ranks;
909  std::sort(data_rank.begin(), data_rank.end());
910 
911  const unsigned int quantile_index =
912  static_cast<unsigned int>(std::round(quantile * n_ranks));
913  AssertIndexRange(quantile_index, data_rank.size());
914  result[0] = data_rank[0].first;
915  result[1] = data_rank[0].second;
916  result[2] = data_rank[quantile_index].first;
917  result[3] = data_rank[quantile_index].second;
918  result[5] = data_rank[n_ranks - 1 - quantile_index].first;
919  result[6] = data_rank[n_ranks - 1 - quantile_index].second;
920  result[7] = data_rank[n_ranks - 1].first;
921  result[8] = data_rank[n_ranks - 1].second;
922  }
923  ierr = MPI_Bcast(result.data(), 9, MPI_DOUBLE, 0, mpi_comm);
924  AssertThrowMPI(ierr);
925 #endif
926  out_stream << std::setw(10) << std::setprecision(4) << std::right;
927  out_stream << result[0] << "s ";
928  out_stream << std::setw(5) << std::right;
929  out_stream << static_cast<unsigned int>(result[1])
930  << (n_ranks > 99999 ? "" : " ") << "|";
931  out_stream << std::setw(10) << std::setprecision(4) << std::right;
932  out_stream << result[2] << "s ";
933  out_stream << std::setw(5) << std::right;
934  out_stream << static_cast<unsigned int>(result[3])
935  << (n_ranks > 99999 ? "" : " ") << "|";
936  out_stream << std::setw(10) << std::setprecision(4) << std::right;
937  out_stream << result[4] << "s |";
938  out_stream << std::setw(10) << std::setprecision(4) << std::right;
939  out_stream << result[5] << "s ";
940  out_stream << std::setw(5) << std::right;
941  out_stream << static_cast<unsigned int>(result[6])
942  << (n_ranks > 99999 ? "" : " ") << "|";
943  out_stream << std::setw(10) << std::setprecision(4) << std::right;
944  out_stream << result[7] << "s ";
945  out_stream << std::setw(5) << std::right;
946  out_stream << static_cast<unsigned int>(result[8])
947  << (n_ranks > 99999 ? "" : " ") << "|\n";
948  }
949  };
950 
951  // in case we want to write out wallclock times
952  {
953  const unsigned int n_ranks = Utilities::MPI::n_mpi_processes(mpi_comm);
954 
955  const std::string time_rank_column = "------------------+";
956  const std::string time_rank_space = " |";
957 
958  // now generate a nice table
959  out_stream << '\n'
960  << "+------------------------------" << extra_dash << "+"
961  << time_rank_column
962  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
963  << "------------+"
964  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
965  << time_rank_column << '\n'
966  << "| Total wallclock time elapsed " << extra_space << "|";
967 
968  print_statistics(timer_all.wall_time());
969 
970  out_stream << "| " << extra_space << "|"
971  << time_rank_space
972  << (n_ranks > 1 && quantile > 0. ? time_rank_space : "")
973  << " "
974  << (n_ranks > 1 && quantile > 0. ? time_rank_space : "")
975  << time_rank_space << '\n';
976  out_stream << "| Section " << extra_space << "| no. calls "
977  << "| min time rank |";
978  if (n_ranks > 1 && quantile > 0.)
979  out_stream << " " << std::setw(5) << std::setprecision(2) << std::right
980  << quantile << "-tile rank |";
981  out_stream << " avg time |";
982  if (n_ranks > 1 && quantile > 0.)
983  out_stream << " " << std::setw(5) << std::setprecision(2) << std::right
984  << 1. - quantile << "-tile rank |";
985  out_stream << " max time rank |\n";
986  out_stream << "+------------------------------" << extra_dash << "+"
987  << time_rank_column
988  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
989  << "------------+"
990  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
991  << time_rank_column << '\n';
992  for (const auto &i : sections)
993  {
994  std::string name_out = i.first;
995 
996  // resize the array so that it is always of the same size
997  unsigned int pos_non_space = name_out.find_first_not_of(' ');
998  name_out.erase(0, pos_non_space);
999  name_out.resize(max_width, ' ');
1000  out_stream << "| " << name_out;
1001  out_stream << "| ";
1002  out_stream << std::setw(9);
1003  out_stream << i.second.n_calls << " |";
1004 
1005  print_statistics(i.second.total_wall_time);
1006  }
1007  out_stream << "+------------------------------" << extra_dash << "+"
1008  << time_rank_column
1009  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
1010  << "------------+"
1011  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
1012  << time_rank_column << '\n';
1013  }
1014 }
1015 
1016 
1017 
1018 void
1020 {
1021  output_is_enabled = false;
1022 }
1023 
1024 
1025 
1026 void
1028 {
1029  output_is_enabled = true;
1030 }
1031 
1032 void
1034 {
1035  std::lock_guard<std::mutex> lock(mutex);
1036  sections.clear();
1037  active_sections.clear();
1038  timer_all.restart();
1039 }
1040 
1042 {
1043  try
1044  {
1045  stop();
1046  }
1047  catch (...)
1048  {}
1049 }
1050 
1051 
std::ostream & get_stream() const
void reset()
Definition: timer.cc:1033
MPI_Comm mpi_communicator
Definition: timer.h:882
void print_summary() const
Definition: timer.cc:534
@ cpu_and_wall_times_grouped
Definition: timer.h:660
@ cpu_times
Definition: timer.h:648
@ wall_times
Definition: timer.h:652
OutputFrequency output_frequency
Definition: timer.h:829
void disable_output()
Definition: timer.cc:1019
std::map< std::string, Section > sections
Definition: timer.h:858
void enable_output()
Definition: timer.cc:1027
void leave_subsection(const std::string &section_name="")
Definition: timer.cc:444
OutputType output_type
Definition: timer.h:834
OutputFrequency
Definition: timer.h:600
@ every_call
Definition: timer.h:604
@ summary
Definition: timer.h:608
@ every_call_and_summary
Definition: timer.h:612
void print_wall_time_statistics(const MPI_Comm mpi_comm, const double print_quantile=0.) const
Definition: timer.cc:841
std::list< std::string > active_sections
Definition: timer.h:877
Threads::Mutex mutex
Definition: timer.h:888
ConditionalOStream out_stream
Definition: timer.h:863
@ total_wall_time
Definition: timer.h:632
@ total_cpu_time
Definition: timer.h:628
std::map< std::string, double > get_summary_data(const OutputData kind) const
Definition: timer.cc:508
bool output_is_enabled
Definition: timer.h:869
~TimerOutput()
Definition: timer.cc:350
Timer timer_all
Definition: timer.h:841
void enter_subsection(const std::string &section_name)
Definition: timer.cc:402
TimerOutput(std::ostream &stream, const OutputFrequency output_frequency, const OutputType output_type)
Definition: timer.cc:300
Definition: timer.h:118
double last_cpu_time() const
Definition: timer.cc:255
bool sync_lap_times
Definition: timer.h:335
void start()
Definition: timer.cc:177
bool running
Definition: timer.h:322
Utilities::MPI::MinMaxAvg accumulated_wall_time_data
Definition: timer.h:350
double cpu_time() const
Definition: timer.cc:235
Timer()
Definition: timer.cc:159
double wall_time() const
Definition: timer.cc:263
MPI_Comm mpi_communicator
Definition: timer.h:329
ClockMeasurements< wall_clock_type > wall_times
Definition: timer.h:312
ClockMeasurements< cpu_clock_type > cpu_times
Definition: timer.h:317
void reset()
Definition: timer.cc:287
double stop()
Definition: timer.cc:194
Utilities::MPI::MinMaxAvg last_lap_wall_time_data
Definition: timer.h:342
void restart()
Definition: timer.h:897
double last_wall_time() const
Definition: timer.cc:279
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_WARNING(desc)
Definition: config.h:566
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
#define Assert(cond, exc)
Definition: exceptions.h:1631
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1820
#define AssertThrowMPI(error_code)
Definition: exceptions.h:1947
#define AssertIndexRange(index, range)
Definition: exceptions.h:1888
static ::ExceptionBase & ExcMessage(std::string arg1)
static const types::blas_int zero
static const char T
T sum(const T &t, const MPI_Comm mpi_communicator)
unsigned int n_mpi_processes(const MPI_Comm mpi_communicator)
Definition: mpi.cc:149
T max(const T &t, const MPI_Comm mpi_communicator)
unsigned int this_mpi_process(const MPI_Comm mpi_communicator)
Definition: mpi.cc:164
MinMaxAvg min_max_avg(const double my_value, const MPI_Comm mpi_communicator)
Definition: mpi.cc:123
static const unsigned int invalid_unsigned_int
Definition: types.h:221
static time_point now() noexcept
Definition: timer.cc:111
std::chrono::time_point< CPUClock, duration > time_point
Definition: timer.h:59
clock_type_ clock_type
Definition: timer.h:257
time_point_type current_lap_start_time
Definition: timer.h:273
typename clock_type::duration duration_type
Definition: timer.h:267
duration_type accumulated_time
Definition: timer.h:278
duration_type last_lap_time
Definition: timer.h:283
unsigned int max_index
Definition: mpi.h:980
unsigned int min_index
Definition: mpi.h:970