22 #include <boost/io/ios_state.hpp> 31 #include <type_traits> 33 #ifdef DEAL_II_HAVE_SYS_RESOURCE_H 34 # include <sys/resource.h> 47 namespace TimerImplementation
56 struct is_duration : std::false_type
62 template <
typename Rep,
typename Period>
63 struct is_duration<std::chrono::duration<Rep, Period>> : std::true_type
73 from_seconds(
const double time)
76 "The template type should be a duration type.");
77 return T(std::lround(T::period::den * (time / T::period::num)));
84 template <
typename Rep,
typename Period>
86 to_seconds(
const std::chrono::duration<Rep, Period> duration)
88 return Period::num *
double(duration.count()) / Period::den;
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>();
113 double system_cpu_duration = 0.0;
115 FILETIME cpuTime, sysTime, createTime, exitTime;
116 const auto succeeded = GetProcessTimes(
117 GetCurrentProcess(), &createTime, &exitTime, &sysTime, &cpuTime);
120 system_cpu_duration =
121 (
double)(((
unsigned long long)cpuTime.dwHighDateTime << 32) |
122 cpuTime.dwLowDateTime) /
126 #elif defined(DEAL_II_HAVE_SYS_RESOURCE_H) 128 getrusage(RUSAGE_SELF, &usage);
129 system_cpu_duration = usage.ru_utime.tv_sec + 1.e-6 * usage.ru_utime.tv_usec;
134 internal::TimerImplementation::from_seconds<duration>(system_cpu_duration));
139 template <
typename clock_type_>
148 template <
typename clock_type_>
152 current_lap_start_time = clock_type::now();
160 :
Timer(MPI_COMM_SELF, false)
167 , mpi_communicator(mpi_communicator)
180 #ifdef DEAL_II_WITH_MPI 212 internal::TimerImplementation::from_seconds<decltype(
215 internal::TimerImplementation::from_seconds<decltype(
218 internal::TimerImplementation::to_seconds(
230 return internal::TimerImplementation::to_seconds(
cpu_times.accumulated_time);
240 const double running_time = internal::TimerImplementation::to_seconds(
258 return internal::TimerImplementation::to_seconds(
cpu_times.last_lap_time);
266 wall_clock_type::duration current_elapsed_wall_time;
268 current_elapsed_wall_time = wall_clock_type::now() -
274 return internal::TimerImplementation::to_seconds(current_elapsed_wall_time);
304 : output_frequency(output_frequency)
305 , output_type(output_type)
306 , out_stream(stream, true)
307 , output_is_enabled(true)
316 : output_frequency(output_frequency)
317 , output_type(output_type)
326 std::ostream & stream,
329 : output_frequency(output_frequency)
330 , output_type(output_type)
333 , mpi_communicator(mpi_communicator)
342 : output_frequency(output_frequency)
343 , output_type(output_type)
346 , mpi_communicator(mpi_communicator)
353 auto do_exit = [
this]() {
370 #ifdef DEAL_II_WITH_MPI 371 # if __cpp_lib_uncaught_exceptions >= 201411 378 const unsigned int myid =
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 <<
"---------------------------------------------------------" 405 std::lock_guard<std::mutex> lock(
mutex);
407 Assert(section_name.empty() ==
false,
ExcMessage(
"Section string is empty."));
412 ExcMessage(std::string(
"Cannot enter the already active section <") +
413 section_name +
">."));
430 sections[section_name].total_cpu_time = 0;
431 sections[section_name].total_wall_time = 0;
435 sections[section_name].timer.reset();
436 sections[section_name].timer.start();
448 ExcMessage(
"Cannot exit any section because none has been entered!"));
450 std::lock_guard<std::mutex> lock(
mutex);
452 if (!section_name.empty())
455 ExcMessage(
"Cannot delete a section that was never created."));
459 ExcMessage(
"Cannot delete a section that has not been entered."));
464 const std::string actual_section_name =
467 sections[actual_section_name].timer.stop();
468 sections[actual_section_name].total_wall_time +=
469 sections[actual_section_name].timer.last_wall_time();
475 const double cpu_time =
sections[actual_section_name].timer.last_cpu_time();
476 sections[actual_section_name].total_cpu_time += cpu_time;
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";
489 output_time =
", CPU time: " + cpu.str();
491 output_time =
", wall time: " + wall.str() +
".";
494 ", CPU/wall time: " + cpu.str() +
" / " + wall.str() +
".";
496 out_stream << actual_section_name << output_time << std::endl;
503 actual_section_name));
508 std::map<std::string, double>
511 std::map<std::string, double> output;
512 for (
const auto §ion :
sections)
516 case TimerOutput::OutputData::total_cpu_time:
517 output[section.first] = section.second.total_cpu_time;
519 case TimerOutput::OutputData::total_wall_time:
520 output[section.first] = section.second.total_wall_time;
522 case TimerOutput::OutputData::n_calls:
523 output[section.first] = section.second.n_calls;
542 unsigned int max_width = 0;
545 std::max(max_width, static_cast<unsigned int>(i.first.length()));
548 max_width =
std::max(max_width + 1, static_cast<unsigned int>(32));
549 const std::string extra_dash = std::string(max_width - 32,
'-');
550 const std::string extra_space = std::string(max_width - 32,
' ');
563 double check_time = 0.;
564 for (
const auto &i : sections)
565 check_time += i.second.total_cpu_time;
569 total_cpu_time = check_time;
573 <<
"+---------------------------------------------" 574 << extra_dash <<
"+------------" 575 <<
"+------------+\n" 576 <<
"| Total CPU time elapsed since start " 577 << extra_space <<
"|";
578 out_stream << std::setw(10) << std::setprecision(3) << std::right;
581 << extra_space <<
"| " 588 <<
" | % of total |\n";
589 out_stream <<
"+---------------------------------" << extra_dash
590 <<
"+-----------+------------" 592 for (
const auto &i : sections)
594 std::string name_out = i.first;
597 unsigned int pos_non_space = name_out.find_first_not_of(
' ');
598 name_out.erase(0, pos_non_space);
599 name_out.resize(max_width,
' ');
607 out_stream << i.second.total_cpu_time <<
"s |";
609 if (total_cpu_time != 0)
614 const double fraction =
616 if (fraction > 0.001)
630 <<
"+---------------------------------" << extra_dash
632 <<
"------------+------------+\n" 638 <<
"Note: The sum of counted times is " << time_gap
639 <<
" seconds larger than the total time.\n" 640 <<
"(Timer function may have introduced too much overhead, or different\n" 641 <<
"section timers may have run at the same time.)" << std::endl;
651 <<
"+---------------------------------------------" 652 << extra_dash <<
"+------------" 653 <<
"+------------+\n" 654 <<
"| Total wallclock time elapsed since start " 655 << extra_space <<
"|";
656 out_stream << std::setw(10) << std::setprecision(3) << std::right;
659 << extra_space <<
"| " 666 out_stream <<
"+---------------------------------" << extra_dash
667 <<
"+-----------+------------" 669 for (
const auto &i : sections)
671 std::string name_out = i.first;
674 unsigned int pos_non_space = name_out.find_first_not_of(
' ');
675 name_out.erase(0, pos_non_space);
676 name_out.resize(max_width,
' ');
684 out_stream << i.second.total_wall_time <<
"s |";
687 if (total_wall_time != 0)
692 const double fraction =
694 if (fraction > 0.001)
708 <<
"+---------------------------------" << extra_dash
710 <<
"------------+------------+\n" 723 double check_time = 0.;
725 for (
const auto &i : sections)
726 check_time += i.second.total_cpu_time;
730 total_cpu_time = check_time;
733 out_stream <<
"\n\n+---------------------------------------------" 735 <<
"------------+------------+" 736 <<
"------------+------------+" 738 <<
"| Total CPU/wall time elapsed since start " 739 << extra_space <<
"|" << std::setw(10) << std::setprecision(3)
740 << std::right << total_cpu_time <<
"s | |" 741 << total_wall_time <<
"s | |" 743 << extra_space <<
"|" 746 <<
"\n| Section " << extra_space
748 <<
" CPU time | % of total |" 749 <<
" wall time | % of total |" 750 <<
"\n+---------------------------------" << extra_dash
752 <<
"------------+------------+" 753 <<
"------------+------------+" << std::endl;
755 for (
const auto &i : sections)
757 std::string name_out = i.first;
760 unsigned int pos_non_space = name_out.find_first_not_of(
' ');
761 name_out.erase(0, pos_non_space);
762 name_out.resize(max_width,
' ');
772 out_stream << i.second.total_cpu_time <<
"s |";
774 if (total_cpu_time != 0)
779 const double fraction =
781 if (fraction > 0.001)
799 out_stream << i.second.total_wall_time <<
"s |";
802 if (total_wall_time != 0)
807 const double fraction =
809 if (fraction > 0.001)
825 out_stream <<
"+---------------------------------" << extra_dash
827 <<
"------------+------------+" 828 <<
"------------+------------+" << std::endl
834 <<
"Note: The sum of counted times is " << time_gap
835 <<
" seconds larger than the total time.\n" 836 <<
"(Timer function may have introduced too much overhead, or different\n" 837 <<
"section timers may have run at the same time.)" << std::endl;
845 const double quantile)
const 853 Assert(quantile >= 0. && quantile <= 0.5,
854 ExcMessage(
"The quantile must be between 0 and 0.5"));
857 unsigned int max_width = 0;
860 std::max(max_width, static_cast<unsigned int>(i.first.length()));
863 max_width =
std::max(max_width + 1, static_cast<unsigned int>(17));
864 const std::string extra_dash = std::string(max_width - 17,
'-');
865 const std::string extra_space = std::string(max_width - 17,
' ');
868 const auto print_statistics = [&](
const double given_time) {
870 if (n_ranks == 1 || quantile == 0.)
875 out_stream << std::setw(10) << std::setprecision(4) << std::right;
879 out_stream << std::setw(10) << std::setprecision(4) << std::right;
881 out_stream << std::setw(10) << std::setprecision(4) << std::right;
889 std::vector<double> receive_data(my_rank == 0 ? n_ranks : 0);
890 std::vector<double> result(9);
891 #ifdef DEAL_II_WITH_MPI 892 int ierr = MPI_Gather(&given_time,
905 std::vector<std::pair<double, unsigned int>> data_rank;
906 data_rank.reserve(n_ranks);
907 for (
unsigned int i = 0; i < n_ranks; ++i)
909 data_rank.emplace_back(receive_data[i], i);
910 result[4] += receive_data[i];
912 result[4] /= n_ranks;
913 std::sort(data_rank.begin(), data_rank.end());
915 const unsigned int quantile_index =
916 static_cast<unsigned int>(std::round(quantile * n_ranks));
918 result[0] = data_rank[0].first;
919 result[1] = data_rank[0].second;
920 result[2] = data_rank[quantile_index].first;
921 result[3] = data_rank[quantile_index].second;
922 result[5] = data_rank[n_ranks - 1 - quantile_index].first;
923 result[6] = data_rank[n_ranks - 1 - quantile_index].second;
924 result[7] = data_rank[n_ranks - 1].first;
925 result[8] = data_rank[n_ranks - 1].second;
927 ierr = MPI_Bcast(result.data(), 9, MPI_DOUBLE, 0, mpi_comm);
930 out_stream << std::setw(10) << std::setprecision(4) << std::right;
933 out_stream << static_cast<unsigned int>(result[1])
934 << (n_ranks > 99999 ?
"" :
" ") <<
"|";
935 out_stream << std::setw(10) << std::setprecision(4) << std::right;
938 out_stream << static_cast<unsigned int>(result[3])
939 << (n_ranks > 99999 ?
"" :
" ") <<
"|";
940 out_stream << std::setw(10) << std::setprecision(4) << std::right;
942 out_stream << std::setw(10) << std::setprecision(4) << std::right;
945 out_stream << static_cast<unsigned int>(result[6])
946 << (n_ranks > 99999 ?
"" :
" ") <<
"|";
947 out_stream << std::setw(10) << std::setprecision(4) << std::right;
950 out_stream << static_cast<unsigned int>(result[8])
951 << (n_ranks > 99999 ?
"" :
" ") <<
"|\n";
959 const std::string time_rank_column =
"------------------+";
960 const std::string time_rank_space =
" |";
964 <<
"+------------------------------" << extra_dash <<
"+" 966 << (n_ranks > 1 && quantile > 0. ? time_rank_column :
"")
968 << (n_ranks > 1 && quantile > 0. ? time_rank_column :
"")
969 << time_rank_column <<
"\n" 970 <<
"| Total wallclock time elapsed " << extra_space <<
"|";
976 << (n_ranks > 1 && quantile > 0. ? time_rank_space :
"")
978 << (n_ranks > 1 && quantile > 0. ? time_rank_space :
"")
979 << time_rank_space <<
"\n";
980 out_stream <<
"| Section " << extra_space <<
"| no. calls " 981 <<
"| min time rank |";
982 if (n_ranks > 1 && quantile > 0.)
983 out_stream <<
" " << std::setw(5) << std::setprecision(2) << std::right
984 << quantile <<
"-tile rank |";
986 if (n_ranks > 1 && quantile > 0.)
987 out_stream <<
" " << std::setw(5) << std::setprecision(2) << std::right
988 << 1. - quantile <<
"-tile rank |";
990 out_stream <<
"+------------------------------" << extra_dash <<
"+" 992 << (n_ranks > 1 && quantile > 0. ? time_rank_column :
"")
994 << (n_ranks > 1 && quantile > 0. ? time_rank_column :
"")
995 << time_rank_column <<
"\n";
996 for (
const auto &i : sections)
998 std::string name_out = i.first;
1001 unsigned int pos_non_space = name_out.find_first_not_of(
' ');
1002 name_out.erase(0, pos_non_space);
1003 name_out.resize(max_width,
' ');
1009 print_statistics(i.second.total_wall_time);
1011 out_stream <<
"+------------------------------" << extra_dash <<
"+" 1013 << (n_ranks > 1 && quantile > 0. ? time_rank_column :
"")
1015 << (n_ranks > 1 && quantile > 0. ? time_rank_column :
"")
1016 << time_rank_column <<
"\n";
1039 std::lock_guard<std::mutex> lock(
mutex);
void print_summary() const
static const unsigned int invalid_unsigned_int
double last_wall_time() const
#define AssertDimension(dim1, dim2)
MPI_Comm mpi_communicator
MPI_Comm mpi_communicator
#define AssertIndexRange(index, range)
static time_point now() noexcept
Utilities::MPI::MinMaxAvg last_lap_wall_time_data
static ::ExceptionBase & ExcMessage(std::string arg1)
T sum(const T &t, const MPI_Comm &mpi_communicator)
#define Assert(cond, exc)
duration_type accumulated_time
#define DEAL_II_NAMESPACE_CLOSE
typename clock_type::duration duration_type
ClockMeasurements< wall_clock_type > wall_times
duration_type last_lap_time
std::list< std::string > active_sections
void leave_subsection(const std::string §ion_name="")
unsigned int n_mpi_processes(const MPI_Comm &mpi_communicator)
std::map< std::string, double > get_summary_data(const OutputData kind) const
#define AssertThrowMPI(error_code)
double last_cpu_time() const
TimerOutput(std::ostream &stream, const OutputFrequency output_frequency, const OutputType output_type)
#define DEAL_II_NAMESPACE_OPEN
std::map< std::string, Section > sections
ConditionalOStream out_stream
std::ostream & get_stream() const
unsigned int this_mpi_process(const MPI_Comm &mpi_communicator)
std::chrono::time_point< CPUClock, duration > time_point
static ::ExceptionBase & ExcNotImplemented()
static const types::blas_int zero
void enter_subsection(const std::string §ion_name)
MinMaxAvg min_max_avg(const double my_value, const MPI_Comm &mpi_communicator)
OutputFrequency output_frequency
time_point_type current_lap_start_time
void print_wall_time_statistics(const MPI_Comm &mpi_comm, const double print_quantile=0.) const
T max(const T &t, const MPI_Comm &mpi_communicator)
#define DEAL_II_WARNING(desc)
Utilities::MPI::MinMaxAvg accumulated_wall_time_data
ClockMeasurements< cpu_clock_type > cpu_times