Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
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>
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
45namespace 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
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
139template <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
148template <typename clock_type_>
149void
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
165Timer::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
176void
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
193double
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
235double
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
255double
257{
258 return internal::TimerImplementation::to_seconds(cpu_times.last_lap_time);
259}
260
261
262
263double
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
279double
281{
282 return internal::TimerImplementation::to_seconds(wall_times.last_lap_time);
283}
284
285
286
287void
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}
297
298
299/* ---------------------------- TimerOutput -------------------------- */
300
301TimerOutput::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
325TimerOutput::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
338TimerOutput::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)
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
402void
403TimerOutput::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
444void
445TimerOutput::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...
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
508std::map<std::string, double>
510{
511 std::map<std::string, double> output;
512 for (const auto &section : sections)
513 {
514 switch (kind)
515 {
517 output[section.first] = section.second.total_cpu_time;
518 break;
520 output[section.first] = section.second.total_wall_time;
521 break;
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
534void
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)
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)
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
841void
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
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
1019void
1021{
1022 output_is_enabled = false;
1023}
1024
1025
1026
1027void
1029{
1030 output_is_enabled = true;
1031}
1032
1033void
1035{
1036 std::lock_guard<std::mutex> lock(mutex);
1037 sections.clear();
1038 active_sections.clear();
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
@ 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:842
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
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:472
#define DEAL_II_WARNING(desc)
Definition config.h:561
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
static ::ExceptionBase & ExcNotImplemented()
#define Assert(cond, exc)
#define AssertDimension(dim1, dim2)
#define AssertThrowMPI(error_code)
#define AssertIndexRange(index, range)
static ::ExceptionBase & ExcMessage(std::string arg1)
T sum(const T &t, const MPI_Comm mpi_communicator)
unsigned int n_mpi_processes(const MPI_Comm mpi_communicator)
Definition mpi.cc:150
T max(const T &t, const MPI_Comm mpi_communicator)
unsigned int this_mpi_process(const MPI_Comm mpi_communicator)
Definition mpi.cc:161
MinMaxAvg min_max_avg(const double my_value, const MPI_Comm mpi_communicator)
Definition mpi.cc:124
static const unsigned int invalid_unsigned_int
Definition types.h:213
STL namespace.
::VectorizedArray< Number, width > max(const ::VectorizedArray< Number, width > &, const ::VectorizedArray< Number, width > &)
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