From: kronbichler Date: Fri, 23 Jan 2009 09:32:46 +0000 (+0000) Subject: Apply the changes that I accidentally submitted to the old server yesterday. Now... X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=c7c14c2fe4f4d64d8086226a3b2fb29716f01be8;p=dealii-svn.git Apply the changes that I accidentally submitted to the old server yesterday. Now the timer class makes sense in MPI, too, and the FEValues::get_function_xxx have a shortcut when the vector contains zeros. git-svn-id: https://svn.dealii.org/trunk@18277 0785d39b-7218-0410-832d-ea1e28bc413d --- diff --git a/deal.II/base/include/base/timer.h b/deal.II/base/include/base/timer.h index 698c66438d..ada33f4056 100644 --- a/deal.II/base/include/base/timer.h +++ b/deal.II/base/include/base/timer.h @@ -16,6 +16,10 @@ #include #include +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI +#include +#endif + #include #include #include @@ -73,6 +77,24 @@ class Timer */ Timer (); +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + /** + * Constructor that takes an MPI + * communicator as input. A timer + * constructed this way will sum up the + * CPU times over all processors in the + * MPI network when requested by the + * operator (). + * + * Starts the timer at 0 sec. + * + * This constructor is only available + * if the deal.II compiler is an MPI + * compiler. + */ + Timer (MPI_Comm mpi_communicator); +#endif + /** * Re-start the timer at the point where * it was stopped. This way a cumulative @@ -170,6 +192,14 @@ class Timer * running. */ bool running; + +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + /** + * Store whether the timer is presently + * running. + */ + MPI_Comm mpi_communicator; +#endif }; @@ -220,9 +250,43 @@ class TimerOutput * ConditionalOStream to write output to. */ TimerOutput (ConditionalOStream &stream, + const enum OutputFrequency output_frequency, + const enum OutputType output_type); + +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + /** + * Constructor that takes an MPI + * communicator as input. A timer + * constructed this way will sum up the + * CPU times over all processors in the + * MPI network for calculating the CPU + * time. + * + * Meant for using std::cout as output + * stream. + */ + TimerOutput (MPI_Comm mpi_comm, + std::ostream &stream, const enum OutputFrequency output_frequency, const enum OutputType output_type); + /** + * Constructor that takes an MPI + * communicator as input. A timer + * constructed this way will sum up the + * CPU times over all processors in the + * MPI network for calculating the CPU + * time. + * + * Constructor that takes a + * ConditionalOStream to write output to. + */ + TimerOutput (MPI_Comm mpi_comm, + ConditionalOStream &stream, + const enum OutputFrequency output_frequency, + const enum OutputType output_type); +#endif + /** * Destructor. Calls print_summary() in * case the option for writing the @@ -298,6 +362,14 @@ class TimerOutput * function. */ std::list active_sections; + +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + /** + * Store whether the timer is presently + * running. + */ + MPI_Comm mpi_communicator; +#endif }; DEAL_II_NAMESPACE_CLOSE diff --git a/deal.II/base/source/timer.cc b/deal.II/base/source/timer.cc index d0688febb1..81f4eaeae3 100644 --- a/deal.II/base/source/timer.cc +++ b/deal.II/base/source/timer.cc @@ -22,9 +22,7 @@ #include #include -#ifdef DEAL_II_COMPILER_SUPPORTS_MPI -#include -#endif +#include // on SunOS 4.x, getrusage is stated in the man pages and exists, but @@ -39,16 +37,37 @@ DEAL_II_NAMESPACE_OPEN + // in case we use an MPI compiler, need + // to create a communicator just for the + // current process Timer::Timer() : cumulative_time (0.), cumulative_wall_time (0.) +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + , mpi_communicator (MPI_COMM_SELF) +#endif { start(); } + // in case we use an MPI compiler, use + // the communicator given from input +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI +Timer::Timer(MPI_Comm mpi_communicator) + : + cumulative_time (0.), + cumulative_wall_time (0.), + mpi_communicator (mpi_communicator) +{ + start(); +} +#endif + + + void Timer::start () { running = true; @@ -107,10 +126,52 @@ double Timer::operator() () const const double dtime_children = usage_children.ru_utime.tv_sec + 1.e-6 * usage_children.ru_utime.tv_usec; + // in case of MPI, need to get the time + // passed by summing the time over all + // processes in the network. works also + // in case we just want to have the time + // of a single thread, since then the + // communicator is MPI_COMM_SELF +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + double local_time = dtime - start_time + dtime_children + - start_time_children + cumulative_time; + + int mpiInitialized; + MPI_Initialized(&mpiInitialized); + + if ( mpiInitialized ) + { + double global_time = 0.; + MPI_Allreduce (&local_time, &global_time, 1, MPI_DOUBLE, MPI_SUM, + mpi_communicator); + return global_time; + } + else + return local_time; +#else return dtime - start_time + dtime_children - start_time_children + cumulative_time; +#endif } else +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + { + int mpiInitialized; + MPI_Initialized(&mpiInitialized); + + if ( mpiInitialized ) + { + double local_time = cumulative_time; + double global_time = 0.; + MPI_Allreduce (&local_time, &global_time, 1, MPI_DOUBLE, MPI_SUM, + mpi_communicator); + return global_time; + } + else + return cumulative_time; + } +#else return cumulative_time; +#endif } @@ -148,6 +209,9 @@ TimerOutput::TimerOutput (std::ostream &stream, output_frequency (output_frequency), output_type (output_type), out_stream (stream, true) +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + , mpi_communicator (MPI_COMM_SELF) +#endif {} @@ -159,8 +223,39 @@ TimerOutput::TimerOutput (ConditionalOStream &stream, output_frequency (output_frequency), output_type (output_type), out_stream (stream) +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + , mpi_communicator (MPI_COMM_SELF) +#endif +{} + + +#ifdef DEAL_II_COMPILER_SUPPORTS_MPI + +TimerOutput::TimerOutput (MPI_Comm mpi_communicator, + std::ostream &stream, + const enum OutputFrequency output_frequency, + const enum OutputType output_type) + : + output_frequency (output_frequency), + output_type (output_type), + out_stream (stream, true), + mpi_communicator (mpi_communicator) +{} + + + +TimerOutput::TimerOutput (MPI_Comm mpi_communicator, + ConditionalOStream &stream, + const enum OutputFrequency output_frequency, + const enum OutputType output_type) + : + output_frequency (output_frequency), + output_type (output_type), + out_stream (stream), + mpi_communicator (mpi_communicator) {} +#endif TimerOutput::~TimerOutput() @@ -177,6 +272,9 @@ TimerOutput::~TimerOutput() void TimerOutput::enter_section (const std::string §ion_name) { + Threads::Mutex mutex; + Threads::Mutex::ScopedLock lock (mutex); + Assert (section_name.empty() == false, ExcMessage ("Section string is empty.")); @@ -203,6 +301,9 @@ TimerOutput::enter_section (const std::string §ion_name) void TimerOutput::exit_section (const std::string §ion_name) { + Threads::Mutex mutex; + Threads::Mutex::ScopedLock lock (mutex); + if (section_name != "") { Assert (sections.find (section_name) != sections.end(), @@ -222,13 +323,11 @@ TimerOutput::exit_section (const std::string §ion_name) sections[actual_section_name].total_wall_time += sections[actual_section_name].timer.wall_time(); - // get cpu time. on MPI - // systems, add the local - // contributions. - // - // TODO: this should rather be - // in the Timer class itself, - // shouldn't it? + // get cpu time. on MPI systems, add + // the local contributions. we could + // do that also in the Timer class + // itself, but we didn't initialize + // the Timers here according to that double cpu_time = sections[actual_section_name].timer(); { @@ -242,7 +341,7 @@ TimerOutput::exit_section (const std::string §ion_name) if( mpiInitialized ) { MPI_Allreduce (&cpu_time, &total_cpu_time, 1, MPI_DOUBLE, MPI_SUM, - MPI_COMM_WORLD); + mpi_communicator); cpu_time = total_cpu_time; } #endif @@ -294,8 +393,8 @@ TimerOutput::print_summary () if( mpiInitialized ) { - MPI_Allreduce (&my_cpu_time, &total_cpu_time, 1, MPI_DOUBLE, MPI_SUM, - MPI_COMM_WORLD); + MPI_Allreduce (&my_cpu_time, &total_cpu_time, 1, MPI_DOUBLE, + MPI_SUM, mpi_communicator); } else total_cpu_time = my_cpu_time; diff --git a/deal.II/deal.II/source/fe/fe_values.cc b/deal.II/deal.II/source/fe/fe_values.cc index 566ccff831..0f06613d7a 100644 --- a/deal.II/deal.II/source/fe/fe_values.cc +++ b/deal.II/deal.II/source/fe/fe_values.cc @@ -255,6 +255,9 @@ namespace FEValuesViews if (shape_function_data[shape_function].is_nonzero_shape_function_component) { const double value = dof_values(shape_function); + if (value == 0.) + continue; + const double * shape_value_ptr = &fe_values.shape_values(shape_function_data[shape_function].row_index, 0); for (unsigned int q_point=0; q_point * shape_gradient_ptr = &fe_values.shape_gradients[shape_function_data[shape_function]. row_index][0]; @@ -334,6 +340,9 @@ namespace FEValuesViews if (shape_function_data[shape_function].is_nonzero_shape_function_component) { const double value = dof_values(shape_function); + if (value == 0.) + continue; + const Tensor<2,spacedim> * shape_hessian_ptr = &fe_values.shape_hessians[shape_function_data[shape_function]. row_index][0]; @@ -374,6 +383,9 @@ namespace FEValuesViews if (shape_function_data[shape_function].is_nonzero_shape_function_component) { const double value = dof_values(shape_function); + if (value == 0.) + continue; + const unsigned int row_index = shape_function_data[shape_function].row_index; for (unsigned int q_point=0; q_point::get_function_values ( for (unsigned int shape_func=0; shape_funcshape_values(shape_func, 0); for (unsigned int point=0; point::get_function_values ( for (unsigned int shape_func=0; shape_funcshape_values(shape_func, 0); for (unsigned int point=0; point::get_function_values ( for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { @@ -1987,6 +2025,8 @@ void FEValuesBase::get_function_values ( for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { @@ -2089,6 +2129,8 @@ void FEValuesBase::get_function_values ( for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { @@ -2187,9 +2229,11 @@ FEValuesBase::get_function_gradients ( for (unsigned int shape_func=0; shape_func *shape_gradient_ptr = &this->shape_gradients[shape_func][0]; - for (unsigned int point=0; point::get_function_gradients ( for (unsigned int shape_func=0; shape_func *shape_gradient_ptr = &this->shape_gradients[shape_func][0]; - for (unsigned int point=0; point::get_function_gradients ( for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { @@ -2397,6 +2445,8 @@ void FEValuesBase::get_function_gradients ( for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { @@ -2479,9 +2529,11 @@ get_function_hessians (const InputVector &fe_function, for (unsigned int shape_func=0; shape_func *shape_hessians_ptr = &this->shape_hessians[shape_func][0]; - for (unsigned int point=0; point::get_function_hessians ( for (unsigned int shape_func=0; shape_func *shape_hessians_ptr = &this->shape_hessians[shape_func][0]; - for (unsigned int point=0; pointis_primitive(shape_func)) { @@ -2680,6 +2736,8 @@ void FEValuesBase::get_function_hessians ( for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { @@ -2765,9 +2823,11 @@ void FEValuesBase::get_function_laplacians ( for (unsigned int shape_func=0; shape_func *shape_hessian_ptr = &this->shape_hessians[shape_func][0]; - for (unsigned int point=0; point::get_function_laplacians ( for (unsigned int shape_func=0; shape_func *shape_hessian_ptr = &this->shape_hessians[shape_func][0]; - for (unsigned int point=0; point::get_function_laplacians ( for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { @@ -2958,6 +3022,8 @@ void FEValuesBase::get_function_laplacians ( for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) { @@ -3062,6 +3128,8 @@ void FEValuesBase::get_function_laplacians ( for (unsigned int shape_func=0; shape_funcis_primitive(shape_func)) {