Reference documentation for deal.II version GIT ca9e7e8105 2023-03-24 14:15:03+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 - 2022 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<decltype(
213  wall_times)::duration_type>(last_lap_wall_time_data.max);
214  cpu_times.last_lap_time =
215  internal::TimerImplementation::from_seconds<decltype(
216  cpu_times)::duration_type>(
218  internal::TimerImplementation::to_seconds(
219  cpu_times.last_lap_time),
221  .max);
222  }
224  cpu_times.accumulated_time += cpu_times.last_lap_time;
226  Utilities::MPI::min_max_avg(internal::TimerImplementation::to_seconds(
229  }
230  return internal::TimerImplementation::to_seconds(cpu_times.accumulated_time);
231 }
232 
233 
234 
235 double
237 {
238  if (running)
239  {
240  const double running_time = internal::TimerImplementation::to_seconds(
241  cpu_clock_type::now() - cpu_times.current_lap_start_time +
242  cpu_times.accumulated_time);
243  return Utilities::MPI::sum(running_time, mpi_communicator);
244  }
245  else
246  {
247  return Utilities::MPI::sum(internal::TimerImplementation::to_seconds(
248  cpu_times.accumulated_time),
250  }
251 }
252 
253 
254 
255 double
257 {
258  return internal::TimerImplementation::to_seconds(cpu_times.last_lap_time);
259 }
260 
261 
262 
263 double
265 {
266  wall_clock_type::duration current_elapsed_wall_time;
267  if (running)
268  current_elapsed_wall_time = wall_clock_type::now() -
271  else
272  current_elapsed_wall_time = wall_times.accumulated_time;
273 
274  return internal::TimerImplementation::to_seconds(current_elapsed_wall_time);
275 }
276 
277 
278 
279 double
281 {
282  return internal::TimerImplementation::to_seconds(wall_times.last_lap_time);
283 }
284 
285 
286 
287 void
289 {
290  wall_times.reset();
291  cpu_times.reset();
292  running = false;
293  internal::TimerImplementation::clear_timing_data(last_lap_wall_time_data);
294  internal::TimerImplementation::clear_timing_data(accumulated_wall_time_data);
295 }
296 
297 
298 
299 /* ---------------------------- TimerOutput -------------------------- */
300 
301 TimerOutput::TimerOutput(std::ostream & stream,
302  const OutputFrequency output_frequency,
303  const OutputType output_type)
304  : output_frequency(output_frequency)
305  , output_type(output_type)
306  , out_stream(stream, true)
307  , output_is_enabled(true)
308  , mpi_communicator(MPI_COMM_SELF)
309 {}
310 
311 
312 
314  const OutputFrequency output_frequency,
315  const OutputType output_type)
316  : output_frequency(output_frequency)
317  , output_type(output_type)
318  , out_stream(stream)
319  , output_is_enabled(true)
320  , mpi_communicator(MPI_COMM_SELF)
321 {}
322 
323 
324 
325 TimerOutput::TimerOutput(const MPI_Comm & mpi_communicator,
326  std::ostream & stream,
327  const OutputFrequency output_frequency,
328  const OutputType output_type)
329  : output_frequency(output_frequency)
330  , output_type(output_type)
331  , out_stream(stream, true)
332  , output_is_enabled(true)
333  , mpi_communicator(mpi_communicator)
334 {}
335 
336 
337 
338 TimerOutput::TimerOutput(const MPI_Comm & mpi_communicator,
339  ConditionalOStream & stream,
340  const OutputFrequency output_frequency,
341  const OutputType output_type)
342  : output_frequency(output_frequency)
343  , output_type(output_type)
344  , out_stream(stream)
345  , output_is_enabled(true)
346  , mpi_communicator(mpi_communicator)
347 {}
348 
349 
350 
352 {
353  auto do_exit = [this]() {
354  try
355  {
356  while (active_sections.size() > 0)
358  // don't print unless we leave all subsections
359  if ((output_frequency == summary ||
361  output_is_enabled == true)
362  print_summary();
363  }
364  catch (...)
365  {}
366  };
367 
368  // avoid communicating with other processes if there is an uncaught
369  // exception
370 #ifdef DEAL_II_WITH_MPI
371 # if __cpp_lib_uncaught_exceptions >= 201411
372  // std::uncaught_exception() is deprecated in c++17
373  if (std::uncaught_exceptions() > 0 && mpi_communicator != MPI_COMM_SELF)
374 # else
375  if (std::uncaught_exception() == true && mpi_communicator != MPI_COMM_SELF)
376 # endif
377  {
378  const unsigned int myid =
380  if (myid == 0)
381  std::cerr
382  << "---------------------------------------------------------\n"
383  << "TimerOutput objects finalize timed values printed to the\n"
384  << "screen by communicating over MPI in their destructors.\n"
385  << "Since an exception is currently uncaught, this\n"
386  << "synchronization (and subsequent output) will be skipped\n"
387  << "to avoid a possible deadlock.\n"
388  << "---------------------------------------------------------"
389  << std::endl;
390  }
391  else
392  {
393  do_exit();
394  }
395 #else
396  do_exit();
397 #endif
398 }
399 
400 
401 
402 void
403 TimerOutput::enter_subsection(const std::string &section_name)
404 {
405  std::lock_guard<std::mutex> lock(mutex);
406 
407  Assert(section_name.empty() == false, ExcMessage("Section string is empty."));
408 
409  Assert(std::find(active_sections.begin(),
410  active_sections.end(),
411  section_name) == active_sections.end(),
412  ExcMessage(std::string("Cannot enter the already active section <") +
413  section_name + ">."));
414 
415  if (sections.find(section_name) == sections.end())
416  {
417  if (mpi_communicator != MPI_COMM_SELF)
418  {
419  // create a new timer for this section. the second argument
420  // will ensure that we have an MPI barrier before starting
421  // and stopping a timer, and this ensures that we get the
422  // maximum run time for this section over all processors.
423  // The mpi_communicator from TimerOutput is passed to the
424  // Timer here, so this Timer will collect timing information
425  // among all processes inside mpi_communicator.
426  sections[section_name].timer = Timer(mpi_communicator, true);
427  }
428 
429 
430  sections[section_name].total_cpu_time = 0;
431  sections[section_name].total_wall_time = 0;
432  sections[section_name].n_calls = 0;
433  }
434 
435  sections[section_name].timer.reset();
436  sections[section_name].timer.start();
437  ++sections[section_name].n_calls;
438 
439  active_sections.push_back(section_name);
440 }
441 
442 
443 
444 void
445 TimerOutput::leave_subsection(const std::string &section_name)
446 {
447  Assert(!active_sections.empty(),
448  ExcMessage("Cannot exit any section because none has been entered!"));
449 
450  std::lock_guard<std::mutex> lock(mutex);
451 
452  if (!section_name.empty())
453  {
454  Assert(sections.find(section_name) != sections.end(),
455  ExcMessage("Cannot delete a section that was never created."));
456  Assert(std::find(active_sections.begin(),
457  active_sections.end(),
458  section_name) != active_sections.end(),
459  ExcMessage("Cannot delete a section that has not been entered."));
460  }
461 
462  // if no string is given, exit the last
463  // active section.
464  const std::string actual_section_name =
465  (section_name.empty() ? active_sections.back() : section_name);
466 
467  sections[actual_section_name].timer.stop();
468  sections[actual_section_name].total_wall_time +=
469  sections[actual_section_name].timer.last_wall_time();
470 
471  // Get cpu time. On MPI systems, if constructed with an mpi_communicator
472  // like MPI_COMM_WORLD, then the Timer will sum up the CPU time between
473  // processors among the provided mpi_communicator. Therefore, no
474  // communication is needed here.
475  const double cpu_time = sections[actual_section_name].timer.last_cpu_time();
476  sections[actual_section_name].total_cpu_time += cpu_time;
477 
478  // in case we have to print out something, do that here...
479  if ((output_frequency == every_call ||
481  output_is_enabled == true)
482  {
483  std::string output_time;
484  std::ostringstream cpu;
485  cpu << cpu_time << "s";
486  std::ostringstream wall;
487  wall << sections[actual_section_name].timer.last_wall_time() << "s";
488  if (output_type == cpu_times)
489  output_time = ", CPU time: " + cpu.str();
490  else if (output_type == wall_times)
491  output_time = ", wall time: " + wall.str() + ".";
492  else
493  output_time =
494  ", CPU/wall time: " + cpu.str() + " / " + wall.str() + ".";
495 
496  out_stream << actual_section_name << output_time << std::endl;
497  }
498 
499  // delete the index from the list of
500  // active ones
501  active_sections.erase(std::find(active_sections.begin(),
502  active_sections.end(),
503  actual_section_name));
504 }
505 
506 
507 
508 std::map<std::string, double>
510 {
511  std::map<std::string, double> output;
512  for (const auto &section : sections)
513  {
514  switch (kind)
515  {
516  case TimerOutput::OutputData::total_cpu_time:
517  output[section.first] = section.second.total_cpu_time;
518  break;
519  case TimerOutput::OutputData::total_wall_time:
520  output[section.first] = section.second.total_wall_time;
521  break;
522  case TimerOutput::OutputData::n_calls:
523  output[section.first] = section.second.n_calls;
524  break;
525  default:
526  Assert(false, ExcNotImplemented());
527  }
528  }
529  return output;
530 }
531 
532 
533 
534 void
536 {
537  // we are going to change the precision and width of output below. store the
538  // old values so the get restored when exiting this function
539  const boost::io::ios_base_all_saver restore_stream(out_stream.get_stream());
540 
541  // get the maximum width among all sections
542  unsigned int max_width = 0;
543  for (const auto &i : sections)
544  max_width = std::max(max_width, static_cast<unsigned int>(i.first.size()));
545 
546  // 32 is the default width until | character
547  max_width = std::max(max_width + 1, static_cast<unsigned int>(32));
548  const std::string extra_dash = std::string(max_width - 32, '-');
549  const std::string extra_space = std::string(max_width - 32, ' ');
550 
552  {
553  // in case we want to write CPU times
554  if (output_type != wall_times)
555  {
556  double total_cpu_time =
558 
559  // check that the sum of all times is less or equal than the total
560  // time. otherwise, we might have generated a lot of overhead in this
561  // function.
562  double check_time = 0.;
563  for (const auto &i : sections)
564  check_time += i.second.total_cpu_time;
565 
566  const double time_gap = check_time - total_cpu_time;
567  if (time_gap > 0.0)
568  total_cpu_time = check_time;
569 
570  // generate a nice table
571  out_stream << "\n\n"
572  << "+---------------------------------------------"
573  << extra_dash << "+------------"
574  << "+------------+\n"
575  << "| Total CPU time elapsed since start "
576  << extra_space << "|";
577  out_stream << std::setw(10) << std::setprecision(3) << std::right;
578  out_stream << total_cpu_time << "s | |\n";
579  out_stream << "| "
580  << extra_space << "| "
581  << "| |\n";
582  out_stream << "| Section " << extra_space
583  << "| no. calls |";
584  out_stream << std::setw(10);
585  out_stream << std::setprecision(3);
586  out_stream << " CPU time "
587  << " | % of total |\n";
588  out_stream << "+---------------------------------" << extra_dash
589  << "+-----------+------------"
590  << "+------------+";
591  for (const auto &i : sections)
592  {
593  std::string name_out = i.first;
594 
595  // resize the array so that it is always of the same size
596  unsigned int pos_non_space = name_out.find_first_not_of(' ');
597  name_out.erase(0, pos_non_space);
598  name_out.resize(max_width, ' ');
599  out_stream << std::endl;
600  out_stream << "| " << name_out;
601  out_stream << "| ";
602  out_stream << std::setw(9);
603  out_stream << i.second.n_calls << " |";
604  out_stream << std::setw(10);
605  out_stream << std::setprecision(3);
606  out_stream << i.second.total_cpu_time << "s |";
607  out_stream << std::setw(10);
608  if (total_cpu_time != 0)
609  {
610  // if run time was less than 0.1%, just print a zero to avoid
611  // printing silly things such as "2.45e-6%". otherwise print
612  // the actual percentage
613  const double fraction =
614  i.second.total_cpu_time / total_cpu_time;
615  if (fraction > 0.001)
616  {
617  out_stream << std::setprecision(2);
618  out_stream << fraction * 100;
619  }
620  else
621  out_stream << 0.0;
622 
623  out_stream << "% |";
624  }
625  else
626  out_stream << 0.0 << "% |";
627  }
628  out_stream << std::endl
629  << "+---------------------------------" << extra_dash
630  << "+-----------+"
631  << "------------+------------+\n"
632  << std::endl;
633 
634  if (time_gap > 0.0)
635  out_stream
636  << std::endl
637  << "Note: The sum of counted times is " << time_gap
638  << " seconds larger than the total time.\n"
639  << "(Timer function may have introduced too much overhead, or different\n"
640  << "section timers may have run at the same time.)" << std::endl;
641  }
642 
643  // in case we want to write out wallclock times
644  if (output_type != cpu_times)
645  {
647 
648  // now generate a nice table
649  out_stream << "\n\n"
650  << "+---------------------------------------------"
651  << extra_dash << "+------------"
652  << "+------------+\n"
653  << "| Total wallclock time elapsed since start "
654  << extra_space << "|";
655  out_stream << std::setw(10) << std::setprecision(3) << std::right;
656  out_stream << total_wall_time << "s | |\n";
657  out_stream << "| "
658  << extra_space << "| "
659  << "| |\n";
660  out_stream << "| Section " << extra_space
661  << "| no. calls |";
662  out_stream << std::setw(10);
663  out_stream << std::setprecision(3);
664  out_stream << " wall time | % of total |\n";
665  out_stream << "+---------------------------------" << extra_dash
666  << "+-----------+------------"
667  << "+------------+";
668  for (const auto &i : sections)
669  {
670  std::string name_out = i.first;
671 
672  // resize the array so that it is always of the same size
673  unsigned int pos_non_space = name_out.find_first_not_of(' ');
674  name_out.erase(0, pos_non_space);
675  name_out.resize(max_width, ' ');
676  out_stream << std::endl;
677  out_stream << "| " << name_out;
678  out_stream << "| ";
679  out_stream << std::setw(9);
680  out_stream << i.second.n_calls << " |";
681  out_stream << std::setw(10);
682  out_stream << std::setprecision(3);
683  out_stream << i.second.total_wall_time << "s |";
684  out_stream << std::setw(10);
685 
686  if (total_wall_time != 0)
687  {
688  // if run time was less than 0.1%, just print a zero to avoid
689  // printing silly things such as "2.45e-6%". otherwise print
690  // the actual percentage
691  const double fraction =
692  i.second.total_wall_time / total_wall_time;
693  if (fraction > 0.001)
694  {
695  out_stream << std::setprecision(2);
696  out_stream << fraction * 100;
697  }
698  else
699  out_stream << 0.0;
700 
701  out_stream << "% |";
702  }
703  else
704  out_stream << 0.0 << "% |";
705  }
706  out_stream << std::endl
707  << "+---------------------------------" << extra_dash
708  << "+-----------+"
709  << "------------+------------+\n"
710  << std::endl;
711  }
712  }
713  else
714  // output_type == cpu_and_wall_times_grouped
715  {
716  const double total_wall_time = timer_all.wall_time();
717  double total_cpu_time =
719 
720  // check that the sum of all times is less or equal than the total time.
721  // otherwise, we might have generated a lot of overhead in this function.
722  double check_time = 0.;
723 
724  for (const auto &i : sections)
725  check_time += i.second.total_cpu_time;
726 
727  const double time_gap = check_time - total_cpu_time;
728  if (time_gap > 0.0)
729  total_cpu_time = check_time;
730 
731  // generate a nice table
732  out_stream << "\n\n+---------------------------------------------"
733  << extra_dash << "+"
734  << "------------+------------+"
735  << "------------+------------+" << '\n'
736  << "| Total CPU/wall time elapsed since start "
737  << extra_space << "|" << std::setw(10) << std::setprecision(3)
738  << std::right << total_cpu_time << "s | |"
739  << total_wall_time << "s | |"
740  << "\n| "
741  << extra_space << "|"
742  << " | |"
743  << " | |"
744  << "\n| Section " << extra_space
745  << "| no. calls |"
746  << " CPU time | % of total |"
747  << " wall time | % of total |"
748  << "\n+---------------------------------" << extra_dash
749  << "+-----------+"
750  << "------------+------------+"
751  << "------------+------------+" << std::endl;
752 
753  for (const auto &i : sections)
754  {
755  std::string name_out = i.first;
756 
757  // resize the array so that it is always of the same size
758  unsigned int pos_non_space = name_out.find_first_not_of(' ');
759  name_out.erase(0, pos_non_space);
760  name_out.resize(max_width, ' ');
761  out_stream << "| " << name_out << "| ";
762 
763  out_stream << std::setw(9);
764  out_stream << i.second.n_calls << " |";
765 
766  if (output_type != wall_times)
767  {
768  out_stream << std::setw(10);
769  out_stream << std::setprecision(3);
770  out_stream << i.second.total_cpu_time << "s |";
771  out_stream << std::setw(10);
772  if (total_cpu_time != 0)
773  {
774  // if run time was less than 0.1%, just print a zero to avoid
775  // printing silly things such as "2.45e-6%". otherwise print
776  // the actual percentage
777  const double fraction =
778  i.second.total_cpu_time / total_cpu_time;
779  if (fraction > 0.001)
780  {
781  out_stream << std::setprecision(2);
782  out_stream << fraction * 100;
783  }
784  else
785  out_stream << 0.0;
786 
787  out_stream << "% |";
788  }
789  else
790  out_stream << 0.0 << "% |";
791  }
792 
793  if (output_type != cpu_times)
794  {
795  out_stream << std::setw(10);
796  out_stream << std::setprecision(3);
797  out_stream << i.second.total_wall_time << "s |";
798  out_stream << std::setw(10);
799 
800  if (total_wall_time != 0)
801  {
802  // if run time was less than 0.1%, just print a zero to avoid
803  // printing silly things such as "2.45e-6%". otherwise print
804  // the actual percentage
805  const double fraction =
806  i.second.total_wall_time / total_wall_time;
807  if (fraction > 0.001)
808  {
809  out_stream << std::setprecision(2);
810  out_stream << fraction * 100;
811  }
812  else
813  out_stream << 0.0;
814 
815  out_stream << "% |";
816  }
817  else
818  out_stream << 0.0 << "% |";
819  }
820  out_stream << std::endl;
821  }
822 
823  out_stream << "+---------------------------------" << extra_dash
824  << "+-----------+"
825  << "------------+------------+"
826  << "------------+------------+" << std::endl
827  << std::endl;
828 
829  if (output_type != wall_times && time_gap > 0.0)
830  out_stream
831  << std::endl
832  << "Note: The sum of counted times is " << time_gap
833  << " seconds larger than the total time.\n"
834  << "(Timer function may have introduced too much overhead, or different\n"
835  << "section timers may have run at the same time.)" << std::endl;
836  }
837 }
838 
839 
840 
841 void
843  const double quantile) const
844 {
845  // we are going to change the precision and width of output below. store the
846  // old values so the get restored when exiting this function
847  const boost::io::ios_base_all_saver restore_stream(out_stream.get_stream());
848 
849  AssertDimension(sections.size(),
850  Utilities::MPI::max(sections.size(), mpi_comm));
851  Assert(quantile >= 0. && quantile <= 0.5,
852  ExcMessage("The quantile must be between 0 and 0.5"));
853 
854  // get the maximum width among all sections
855  unsigned int max_width = 0;
856  for (const auto &i : sections)
857  max_width = std::max(max_width, static_cast<unsigned int>(i.first.size()));
858 
859  // 17 is the default width until | character
860  max_width = std::max(max_width + 1, static_cast<unsigned int>(17));
861  const std::string extra_dash = std::string(max_width - 17, '-');
862  const std::string extra_space = std::string(max_width - 17, ' ');
863 
864  // function to print data in a nice table
865  const auto print_statistics = [&](const double given_time) {
866  const unsigned int n_ranks = Utilities::MPI::n_mpi_processes(mpi_comm);
867  if (n_ranks == 1 || quantile == 0.)
868  {
870  Utilities::MPI::min_max_avg(given_time, mpi_comm);
871 
872  out_stream << std::setw(10) << std::setprecision(4) << std::right;
873  out_stream << data.min << "s ";
874  out_stream << std::setw(5) << std::right;
875  out_stream << data.min_index << (n_ranks > 99999 ? "" : " ") << "|";
876  out_stream << std::setw(10) << std::setprecision(4) << std::right;
877  out_stream << data.avg << "s |";
878  out_stream << std::setw(10) << std::setprecision(4) << std::right;
879  out_stream << data.max << "s ";
880  out_stream << std::setw(5) << std::right;
881  out_stream << data.max_index << (n_ranks > 99999 ? "" : " ") << "|\n";
882  }
883  else
884  {
885  const unsigned int my_rank = Utilities::MPI::this_mpi_process(mpi_comm);
886  std::vector<double> receive_data(my_rank == 0 ? n_ranks : 0);
887  std::vector<double> result(9);
888 #ifdef DEAL_II_WITH_MPI
889  int ierr = MPI_Gather(&given_time,
890  1,
891  MPI_DOUBLE,
892  receive_data.data(),
893  1,
894  MPI_DOUBLE,
895  0,
896  mpi_comm);
897  AssertThrowMPI(ierr);
898  if (my_rank == 0)
899  {
900  // fill the received data in a pair and sort; on the way, also
901  // compute the average
902  std::vector<std::pair<double, unsigned int>> data_rank;
903  data_rank.reserve(n_ranks);
904  for (unsigned int i = 0; i < n_ranks; ++i)
905  {
906  data_rank.emplace_back(receive_data[i], i);
907  result[4] += receive_data[i];
908  }
909  result[4] /= n_ranks;
910  std::sort(data_rank.begin(), data_rank.end());
911 
912  const unsigned int quantile_index =
913  static_cast<unsigned int>(std::round(quantile * n_ranks));
914  AssertIndexRange(quantile_index, data_rank.size());
915  result[0] = data_rank[0].first;
916  result[1] = data_rank[0].second;
917  result[2] = data_rank[quantile_index].first;
918  result[3] = data_rank[quantile_index].second;
919  result[5] = data_rank[n_ranks - 1 - quantile_index].first;
920  result[6] = data_rank[n_ranks - 1 - quantile_index].second;
921  result[7] = data_rank[n_ranks - 1].first;
922  result[8] = data_rank[n_ranks - 1].second;
923  }
924  ierr = MPI_Bcast(result.data(), 9, MPI_DOUBLE, 0, mpi_comm);
925  AssertThrowMPI(ierr);
926 #endif
927  out_stream << std::setw(10) << std::setprecision(4) << std::right;
928  out_stream << result[0] << "s ";
929  out_stream << std::setw(5) << std::right;
930  out_stream << static_cast<unsigned int>(result[1])
931  << (n_ranks > 99999 ? "" : " ") << "|";
932  out_stream << std::setw(10) << std::setprecision(4) << std::right;
933  out_stream << result[2] << "s ";
934  out_stream << std::setw(5) << std::right;
935  out_stream << static_cast<unsigned int>(result[3])
936  << (n_ranks > 99999 ? "" : " ") << "|";
937  out_stream << std::setw(10) << std::setprecision(4) << std::right;
938  out_stream << result[4] << "s |";
939  out_stream << std::setw(10) << std::setprecision(4) << std::right;
940  out_stream << result[5] << "s ";
941  out_stream << std::setw(5) << std::right;
942  out_stream << static_cast<unsigned int>(result[6])
943  << (n_ranks > 99999 ? "" : " ") << "|";
944  out_stream << std::setw(10) << std::setprecision(4) << std::right;
945  out_stream << result[7] << "s ";
946  out_stream << std::setw(5) << std::right;
947  out_stream << static_cast<unsigned int>(result[8])
948  << (n_ranks > 99999 ? "" : " ") << "|\n";
949  }
950  };
951 
952  // in case we want to write out wallclock times
953  {
954  const unsigned int n_ranks = Utilities::MPI::n_mpi_processes(mpi_comm);
955 
956  const std::string time_rank_column = "------------------+";
957  const std::string time_rank_space = " |";
958 
959  // now generate a nice table
960  out_stream << '\n'
961  << "+------------------------------" << extra_dash << "+"
962  << time_rank_column
963  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
964  << "------------+"
965  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
966  << time_rank_column << '\n'
967  << "| Total wallclock time elapsed " << extra_space << "|";
968 
969  print_statistics(timer_all.wall_time());
970 
971  out_stream << "| " << extra_space << "|"
972  << time_rank_space
973  << (n_ranks > 1 && quantile > 0. ? time_rank_space : "")
974  << " "
975  << (n_ranks > 1 && quantile > 0. ? time_rank_space : "")
976  << time_rank_space << '\n';
977  out_stream << "| Section " << extra_space << "| no. calls "
978  << "| min time rank |";
979  if (n_ranks > 1 && quantile > 0.)
980  out_stream << " " << std::setw(5) << std::setprecision(2) << std::right
981  << quantile << "-tile rank |";
982  out_stream << " avg time |";
983  if (n_ranks > 1 && quantile > 0.)
984  out_stream << " " << std::setw(5) << std::setprecision(2) << std::right
985  << 1. - quantile << "-tile rank |";
986  out_stream << " max time rank |\n";
987  out_stream << "+------------------------------" << extra_dash << "+"
988  << time_rank_column
989  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
990  << "------------+"
991  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
992  << time_rank_column << '\n';
993  for (const auto &i : sections)
994  {
995  std::string name_out = i.first;
996 
997  // resize the array so that it is always of the same size
998  unsigned int pos_non_space = name_out.find_first_not_of(' ');
999  name_out.erase(0, pos_non_space);
1000  name_out.resize(max_width, ' ');
1001  out_stream << "| " << name_out;
1002  out_stream << "| ";
1003  out_stream << std::setw(9);
1004  out_stream << i.second.n_calls << " |";
1005 
1006  print_statistics(i.second.total_wall_time);
1007  }
1008  out_stream << "+------------------------------" << extra_dash << "+"
1009  << time_rank_column
1010  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
1011  << "------------+"
1012  << (n_ranks > 1 && quantile > 0. ? time_rank_column : "")
1013  << time_rank_column << '\n';
1014  }
1015 }
1016 
1017 
1018 
1019 void
1021 {
1022  output_is_enabled = false;
1023 }
1024 
1025 
1026 
1027 void
1029 {
1030  output_is_enabled = true;
1031 }
1032 
1033 void
1035 {
1036  std::lock_guard<std::mutex> lock(mutex);
1037  sections.clear();
1038  active_sections.clear();
1039  timer_all.restart();
1040 }
1041 
1043 {
1044  try
1045  {
1046  stop();
1047  }
1048  catch (...)
1049  {}
1050 }
1051 
1052 
std::ostream & get_stream() const
void reset()
Definition: timer.cc:1034
MPI_Comm mpi_communicator
Definition: timer.h:882
void print_summary() const
Definition: timer.cc:535
@ 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:1020
std::map< std::string, Section > sections
Definition: timer.h:858
void enable_output()
Definition: timer.cc:1028
void leave_subsection(const std::string &section_name="")
Definition: timer.cc:445
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
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:509
bool output_is_enabled
Definition: timer.h:869
~TimerOutput()
Definition: timer.cc:351
Timer timer_all
Definition: timer.h:841
void enter_subsection(const std::string &section_name)
Definition: timer.cc:403
void print_wall_time_statistics(const MPI_Comm &mpi_comm, const double print_quantile=0.) const
Definition: timer.cc:842
TimerOutput(std::ostream &stream, const OutputFrequency output_frequency, const OutputType output_type)
Definition: timer.cc:301
Definition: timer.h:118
double last_cpu_time() const
Definition: timer.cc:256
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:236
Timer()
Definition: timer.cc:159
double wall_time() const
Definition: timer.cc:264
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:288
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:280
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:474
#define DEAL_II_WARNING(desc)
Definition: config.h:563
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:475
#define Assert(cond, exc)
Definition: exceptions.h:1586
static ::ExceptionBase & ExcNotImplemented()
#define AssertDimension(dim1, dim2)
Definition: exceptions.h:1759
#define AssertThrowMPI(error_code)
Definition: exceptions.h:1886
#define AssertIndexRange(index, range)
Definition: exceptions.h:1827
static ::ExceptionBase & ExcMessage(std::string arg1)
static const types::blas_int zero
static const char T
unsigned int this_mpi_process(const MPI_Comm &mpi_communicator)
Definition: mpi.cc:161
T sum(const T &t, const MPI_Comm &mpi_communicator)
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
Definition: mpi.cc:150
MinMaxAvg min_max_avg(const double my_value, const MPI_Comm &mpi_communicator)
Definition: mpi.cc:124
T max(const T &t, const MPI_Comm &mpi_communicator)
static const unsigned int invalid_unsigned_int
Definition: types.h:213
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