From: Mark Abraham Date: Sat, 28 Sep 2013 13:45:30 +0000 (+0200) Subject: Encapsulate gmx_wallclock_accounting_t into new timing module X-Git-Url: http://biod.pnpi.spb.ru/gitweb/?a=commitdiff_plain;h=690c59fdfd0f8cff4207fc7b61c2c284ab402566;p=alexxy%2Fgromacs.git Encapsulate gmx_wallclock_accounting_t into new timing module Started a timing module where things like gmx_wallcycle.c and gmx_cyclecounter.c will eventually also live. Renamed s/runtime/wallclock_accounting/ to avoid ambiguity with "compile time," "simulation time," etc. This patch preserves the release-4-6 behaviour of the reported times with all supported combinations of MPI, thread-MPI, GPUs, OpenMP. I think the old behaviour of -ntomp_pme != -ntomp was buggy with respect to timing, and this patch fixes that. Change-Id: I22cba17cdf866e17d9d0732de204b85b0ab31994 --- diff --git a/CMakeLists.txt b/CMakeLists.txt index be9a1f04f4..0a7ff76667 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -364,7 +364,7 @@ check_function_exists(sqrtf HAVE_SQRTF) include(CheckLibraryExists) check_library_exists(m sqrt "" HAVE_LIBM) -check_library_exists(rt clock_gettime "" HAVE_LIBRT) +check_library_exists(rt clock_gettime "" HAVE_CLOCK_GETTIME) include(CheckTypeSize) @@ -1060,8 +1060,9 @@ if(GMX_LOAD_PLUGINS) endif(GMX_LOAD_PLUGINS) set(VMD_QUIETLY TRUE CACHE INTERNAL "") -# Link real-time library for POSIX timers -if(HAVE_TIME_H AND HAVE_UNISTD_H AND HAVE_LIBRT) +# Link real-time library for POSIX timers. The check for clock_gettime +# confirms the linkability of rt. +if(HAVE_TIME_H AND HAVE_UNISTD_H AND HAVE_CLOCK_GETTIME) list(APPEND GMX_EXTRA_LIBRARIES rt) endif() diff --git a/src/config.h.cmakein b/src/config.h.cmakein index 3ce464188c..abc42377a2 100644 --- a/src/config.h.cmakein +++ b/src/config.h.cmakein @@ -235,6 +235,9 @@ /* Define to 1 if you have the MSVC _aligned_malloc() function. */ #cmakedefine HAVE__ALIGNED_MALLOC +/* Define to 1 if you have the clock_gettime() function. */ +#cmakedefine HAVE_CLOCK_GETTIME + /* Define to 1 if you have the gettimeofday() function. */ #cmakedefine HAVE_GETTIMEOFDAY diff --git a/src/gromacs/CMakeLists.txt b/src/gromacs/CMakeLists.txt index cfd23d4cf8..32a82ff0f2 100644 --- a/src/gromacs/CMakeLists.txt +++ b/src/gromacs/CMakeLists.txt @@ -53,6 +53,7 @@ add_subdirectory(fft) add_subdirectory(linearalgebra) add_subdirectory(onlinehelp) add_subdirectory(options) +add_subdirectory(timing) add_subdirectory(utility) if (NOT GMX_BUILD_MDRUN_ONLY) add_subdirectory(legacyheaders) diff --git a/src/gromacs/gmxana/gmx_membed.c b/src/gromacs/gmxana/gmx_membed.c index 3d4ae1829b..753f548797 100644 --- a/src/gromacs/gmxana/gmx_membed.c +++ b/src/gromacs/gmxana/gmx_membed.c @@ -142,7 +142,7 @@ int gmx_membed(int argc, char *argv[]) { "-start", FALSE, etBOOL, {&bStart}, "Call mdrun with membed options" }, { "-stepout", FALSE, etINT, {&nstepout}, - "HIDDENFrequency of writing the remaining runtime" }, + "HIDDENFrequency of writing the remaining wall clock time for the run" }, { "-v", FALSE, etBOOL, {&bVerbose}, "Be loud and noisy" }, { "-mdrun_path", FALSE, etSTR, {&mdrun_path}, diff --git a/src/gromacs/gmxana/gmx_tune_pme.c b/src/gromacs/gmxana/gmx_tune_pme.c index b3683936ab..04b9119a3c 100644 --- a/src/gromacs/gmxana/gmx_tune_pme.c +++ b/src/gromacs/gmxana/gmx_tune_pme.c @@ -58,7 +58,7 @@ #include "gmx_ana.h" #include "names.h" #include "perf_est.h" -#include "sim_util.h" +#include "gromacs/timing/walltime_accounting.h" /* Enum for situations that can occur during log file parsing, the diff --git a/src/gromacs/gmxlib/nrnb.c b/src/gromacs/gmxlib/nrnb.c index 9c4917e821..748785e5dd 100644 --- a/src/gromacs/gmxlib/nrnb.c +++ b/src/gromacs/gmxlib/nrnb.c @@ -387,7 +387,7 @@ void print_perf(FILE *out, double time_per_thread, double time_per_node, gmx_large_int_t nsteps, real delta_t, double nbfs, double mflop) { - real runtime; + real wallclocktime; fprintf(out, "\n"); @@ -404,15 +404,15 @@ void print_perf(FILE *out, double time_per_thread, double time_per_node, } if (delta_t > 0) { - mflop = mflop/time_per_node; - runtime = nsteps*delta_t; + mflop = mflop/time_per_node; + wallclocktime = nsteps*delta_t; if (getenv("GMX_DETAILED_PERF_STATS") == NULL) { fprintf(out, "%12s %12s %12s\n", "", "(ns/day)", "(hour/ns)"); fprintf(out, "%12s %12.3f %12.3f\n", "Performance:", - runtime*24*3.6/time_per_node, 1000*time_per_node/(3600*runtime)); + wallclocktime*24*3.6/time_per_node, 1000*time_per_node/(3600*wallclocktime)); } else { @@ -421,7 +421,7 @@ void print_perf(FILE *out, double time_per_thread, double time_per_node, "(ns/day)", "(hour/ns)"); fprintf(out, "%12s %12.3f %12.3f %12.3f %12.3f\n", "Performance:", nbfs/time_per_node, (mflop > 1000) ? (mflop/1000) : mflop, - runtime*24*3.6/time_per_node, 1000*time_per_node/(3600*runtime)); + wallclocktime*24*3.6/time_per_node, 1000*time_per_node/(3600*wallclocktime)); } } else diff --git a/src/gromacs/legacyheaders/mdrun.h b/src/gromacs/legacyheaders/mdrun.h index b6a3707cab..ed921d3440 100644 --- a/src/gromacs/legacyheaders/mdrun.h +++ b/src/gromacs/legacyheaders/mdrun.h @@ -137,7 +137,7 @@ typedef double gmx_integrator_t (FILE *log, t_commrec *cr, real cpt_period, real max_hours, const char *deviceOptions, unsigned long Flags, - gmx_runtime_t *runtime); + gmx_walltime_accounting_t walltime_accounting); /* ROUTINES from md.c */ diff --git a/src/gromacs/legacyheaders/pme.h b/src/gromacs/legacyheaders/pme.h index e405926981..510178bb56 100644 --- a/src/gromacs/legacyheaders/pme.h +++ b/src/gromacs/legacyheaders/pme.h @@ -98,7 +98,7 @@ int gmx_pme_do(gmx_pme_t pme, int gmx_pmeonly(gmx_pme_t pme, t_commrec *cr, t_nrnb *mynrnb, gmx_wallcycle_t wcycle, - gmx_runtime_t *runtime, + gmx_walltime_accounting_t walltime_accounting, real ewaldcoeff, t_inputrec *ir); /* Called on the nodes that do PME exclusively (as slaves) diff --git a/src/gromacs/legacyheaders/sim_util.h b/src/gromacs/legacyheaders/sim_util.h index c38c2f183a..7c943fe41f 100644 --- a/src/gromacs/legacyheaders/sim_util.h +++ b/src/gromacs/legacyheaders/sim_util.h @@ -36,12 +36,12 @@ #ifndef _sim_util_h #define _sim_util_h -#include #include "typedefs.h" #include "enxio.h" #include "mdebin.h" #include "update.h" #include "vcm.h" +#include "gromacs/timing/walltime_accounting.h" #ifdef __cplusplus extern "C" { @@ -64,16 +64,6 @@ typedef struct { typedef struct gmx_global_stat *gmx_global_stat_t; -/*! /brief Manages measuring wall clock times for simulations */ -typedef struct { - double start_time_stamp; //!< Seconds since the epoch recorded at the start of the simulation - double start_time_stamp_per_thread; //!< Seconds since the epoch recorded at the start of the simulation for this thread - double elapsed_run_time; //!< Total seconds elapsed over the simulation - double elapsed_run_time_per_thread; //!< Total seconds elapsed over the simulation running this thread - gmx_large_int_t nsteps_done; //!< Used by integrators to report the amount of work they did -} gmx_runtime_t; - - void do_pbc_first(FILE *log, matrix box, t_forcerec *fr, t_graph *graph, rvec x[]); @@ -136,24 +126,16 @@ int do_per_step(gmx_large_int_t step, gmx_large_int_t nstep); /* ROUTINES from sim_util.c */ -double gmx_gettime(); - -void print_time(FILE *out, gmx_runtime_t *runtime, +void print_time(FILE *out, gmx_walltime_accounting_t walltime_accounting, gmx_large_int_t step, t_inputrec *ir, t_commrec *cr); -void runtime_start(gmx_runtime_t *runtime); - -void runtime_end(gmx_runtime_t *runtime); - -double runtime_get_elapsed_time(gmx_runtime_t *runtime); - void print_date_and_time(FILE *log, int pid, const char *title, - const gmx_runtime_t *runtime); + const gmx_walltime_accounting_t walltime_accounting); void finish_run(FILE *log, t_commrec *cr, t_inputrec *inputrec, t_nrnb nrnb[], gmx_wallcycle_t wcycle, - gmx_runtime_t *runtime, + gmx_walltime_accounting_t walltime_accounting, wallclock_gpu_t *gputimes, gmx_bool bWriteStat); diff --git a/src/gromacs/mdlib/minimize.c b/src/gromacs/mdlib/minimize.c index c15fe93bd5..cf67098ff4 100644 --- a/src/gromacs/mdlib/minimize.c +++ b/src/gromacs/mdlib/minimize.c @@ -80,6 +80,7 @@ #include "gromacs/linearalgebra/mtxio.h" #include "gromacs/linearalgebra/sparsematrix.h" +#include "gromacs/timing/walltime_accounting.h" typedef struct { t_state s; @@ -102,25 +103,27 @@ static em_state_t *init_em_state() return ems; } -static void print_em_start(FILE *fplog, t_commrec *cr, gmx_runtime_t *runtime, - gmx_wallcycle_t wcycle, - const char *name) +static void print_em_start(FILE *fplog, + t_commrec *cr, + gmx_walltime_accounting_t walltime_accounting, + gmx_wallcycle_t wcycle, + const char *name) { char buf[STRLEN]; - runtime_start(runtime); + walltime_accounting_start(walltime_accounting); sprintf(buf, "Started %s", name); print_date_and_time(fplog, cr->nodeid, buf, NULL); wallcycle_start(wcycle, ewcRUN); } -static void em_time_end(gmx_runtime_t *runtime, - gmx_wallcycle_t wcycle) +static void em_time_end(gmx_walltime_accounting_t walltime_accounting, + gmx_wallcycle_t wcycle) { wallcycle_stop(wcycle, ewcRUN); - runtime_end(runtime); + walltime_accounting_end(walltime_accounting); } static void sp_header(FILE *out, const char *minimizer, real ftol, int nsteps) @@ -463,7 +466,8 @@ void init_em(FILE *fplog, const char *title, } static void finish_em(t_commrec *cr, gmx_mdoutf_t *outf, - gmx_runtime_t *runtime, gmx_wallcycle_t wcycle) + gmx_walltime_accounting_t walltime_accounting, + gmx_wallcycle_t wcycle) { if (!(cr->duty & DUTY_PME)) { @@ -473,7 +477,7 @@ static void finish_em(t_commrec *cr, gmx_mdoutf_t *outf, done_mdoutf(outf); - em_time_end(runtime, wcycle); + em_time_end(walltime_accounting, wcycle); } static void swap_em_state(em_state_t *ems1, em_state_t *ems2) @@ -967,7 +971,7 @@ double do_cg(FILE *fplog, t_commrec *cr, real gmx_unused cpt_period, real gmx_unused max_hours, const char gmx_unused *deviceOptions, unsigned long gmx_unused Flags, - gmx_runtime_t *runtime) + gmx_walltime_accounting_t walltime_accounting) { const char *CG = "Polak-Ribiere Conjugate Gradients"; @@ -1008,7 +1012,7 @@ double do_cg(FILE *fplog, t_commrec *cr, nfile, fnm, &outf, &mdebin); /* Print to log file */ - print_em_start(fplog, cr, runtime, wcycle, CG); + print_em_start(fplog, cr, walltime_accounting, wcycle, CG); /* Max number of steps */ number_steps = inputrec->nsteps; @@ -1554,10 +1558,10 @@ double do_cg(FILE *fplog, t_commrec *cr, fprintf(fplog, "\nPerformed %d energy evaluations in total.\n", neval); } - finish_em(cr, outf, runtime, wcycle); + finish_em(cr, outf, walltime_accounting, wcycle); /* To print the actual number of steps we needed somewhere */ - runtime->nsteps_done = step; + walltime_accounting_set_nsteps_done(walltime_accounting, step); return 0; } /* That's all folks */ @@ -1581,7 +1585,7 @@ double do_lbfgs(FILE *fplog, t_commrec *cr, real gmx_unused cpt_period, real gmx_unused max_hours, const char gmx_unused *deviceOptions, unsigned long gmx_unused Flags, - gmx_runtime_t *runtime) + gmx_walltime_accounting_t walltime_accounting) { static const char *LBFGS = "Low-Memory BFGS Minimizer"; em_state_t ems; @@ -1677,7 +1681,7 @@ double do_lbfgs(FILE *fplog, t_commrec *cr, end = mdatoms->homenr + start; /* Print to log file */ - print_em_start(fplog, cr, runtime, wcycle, LBFGS); + print_em_start(fplog, cr, walltime_accounting, wcycle, LBFGS); do_log = do_ene = do_x = do_f = TRUE; @@ -2336,10 +2340,10 @@ double do_lbfgs(FILE *fplog, t_commrec *cr, fprintf(fplog, "\nPerformed %d energy evaluations in total.\n", neval); } - finish_em(cr, outf, runtime, wcycle); + finish_em(cr, outf, walltime_accounting, wcycle); /* To print the actual number of steps we needed somewhere */ - runtime->nsteps_done = step; + walltime_accounting_set_nsteps_done(walltime_accounting, step); return 0; } /* That's all folks */ @@ -2363,7 +2367,7 @@ double do_steep(FILE *fplog, t_commrec *cr, real gmx_unused cpt_period, real gmx_unused max_hours, const char gmx_unused *deviceOptions, unsigned long gmx_unused Flags, - gmx_runtime_t *runtime) + gmx_walltime_accounting_t walltime_accounting) { const char *SD = "Steepest Descents"; em_state_t *s_min, *s_try; @@ -2396,7 +2400,7 @@ double do_steep(FILE *fplog, t_commrec *cr, nfile, fnm, &outf, &mdebin); /* Print to log file */ - print_em_start(fplog, cr, runtime, wcycle, SD); + print_em_start(fplog, cr, walltime_accounting, wcycle, SD); /* Set variables for stepsize (in nm). This is the largest * step that we are going to make in any direction. @@ -2561,12 +2565,12 @@ double do_steep(FILE *fplog, t_commrec *cr, s_min->epot, s_min->fmax, s_min->a_fmax, fnormn); } - finish_em(cr, outf, runtime, wcycle); + finish_em(cr, outf, walltime_accounting, wcycle); /* To print the actual number of steps we needed somewhere */ inputrec->nsteps = count; - runtime->nsteps_done = count; + walltime_accounting_set_nsteps_done(walltime_accounting, count); return 0; } /* That's all folks */ @@ -2590,7 +2594,7 @@ double do_nm(FILE *fplog, t_commrec *cr, real gmx_unused cpt_period, real gmx_unused max_hours, const char gmx_unused *deviceOptions, unsigned long gmx_unused Flags, - gmx_runtime_t *runtime) + gmx_walltime_accounting_t walltime_accounting) { const char *NM = "Normal Mode Analysis"; gmx_mdoutf_t *outf; @@ -2698,7 +2702,7 @@ double do_nm(FILE *fplog, t_commrec *cr, where(); /* Write start time and temperature */ - print_em_start(fplog, cr, runtime, wcycle, NM); + print_em_start(fplog, cr, walltime_accounting, wcycle, NM); /* fudge nr of steps to nr of atoms */ inputrec->nsteps = natoms*2; @@ -2857,9 +2861,9 @@ double do_nm(FILE *fplog, t_commrec *cr, gmx_mtxio_write(ftp2fn(efMTX, nfile, fnm), sz, sz, full_matrix, sparse_matrix); } - finish_em(cr, outf, runtime, wcycle); + finish_em(cr, outf, walltime_accounting, wcycle); - runtime->nsteps_done = natoms*2; + walltime_accounting_set_nsteps_done(walltime_accounting, natoms*2); return 0; } diff --git a/src/gromacs/mdlib/pme.c b/src/gromacs/mdlib/pme.c index 0558c901f8..20ed264458 100644 --- a/src/gromacs/mdlib/pme.c +++ b/src/gromacs/mdlib/pme.c @@ -4143,7 +4143,7 @@ void gmx_pme_calc_energy(gmx_pme_t pme, int n, rvec *x, real *q, real *V) static void reset_pmeonly_counters(gmx_wallcycle_t wcycle, - gmx_runtime_t *runtime, + gmx_walltime_accounting_t walltime_accounting, t_nrnb *nrnb, t_inputrec *ir, gmx_large_int_t step) { @@ -4158,7 +4158,7 @@ static void reset_pmeonly_counters(gmx_wallcycle_t wcycle, } ir->init_step = step; wallcycle_start(wcycle, ewcRUN); - runtime_start(runtime); + walltime_accounting_start(walltime_accounting); } @@ -4199,7 +4199,7 @@ static void gmx_pmeonly_switch(int *npmedata, gmx_pme_t **pmedata, int gmx_pmeonly(gmx_pme_t pme, t_commrec *cr, t_nrnb *nrnb, gmx_wallcycle_t wcycle, - gmx_runtime_t *runtime, + gmx_walltime_accounting_t walltime_accounting, real ewaldcoeff, t_inputrec *ir) { @@ -4255,7 +4255,7 @@ int gmx_pmeonly(gmx_pme_t pme, if (ret == pmerecvqxRESETCOUNTERS) { /* Reset the cycle and flop counters */ - reset_pmeonly_counters(wcycle, runtime, nrnb, ir, step); + reset_pmeonly_counters(wcycle, walltime_accounting, nrnb, ir, step); } } while (ret == pmerecvqxSWITCHGRID || ret == pmerecvqxRESETCOUNTERS); @@ -4271,7 +4271,7 @@ int gmx_pmeonly(gmx_pme_t pme, if (count == 0) { wallcycle_start(wcycle, ewcRUN); - runtime_start(runtime); + walltime_accounting_start(walltime_accounting); } wallcycle_start(wcycle, ewcPMEMESH); @@ -4293,7 +4293,7 @@ int gmx_pmeonly(gmx_pme_t pme, } /***** end of quasi-loop, we stop with the break above */ while (TRUE); - runtime_end(runtime); + walltime_accounting_end(walltime_accounting); return 0; } diff --git a/src/gromacs/mdlib/sim_util.c b/src/gromacs/mdlib/sim_util.c index 03347e3831..3d77d2c314 100644 --- a/src/gromacs/mdlib/sim_util.c +++ b/src/gromacs/mdlib/sim_util.c @@ -38,10 +38,6 @@ #endif #include -#include -#ifdef HAVE_UNISTD_H -#include -#endif #ifdef HAVE_SYS_TIME_H #include #endif @@ -91,6 +87,7 @@ #include "nbnxn_kernels/nbnxn_kernel_gpu_ref.h" #include "gromacs/utility/gmxmpi.h" +#include "gromacs/timing/walltime_accounting.h" #include "adress.h" #include "qmmm.h" @@ -98,55 +95,15 @@ #include "nbnxn_cuda_data_mgmt.h" #include "nbnxn_cuda/nbnxn_cuda.h" -double -gmx_gettime() -{ -#if _POSIX_TIMERS > 0 - /* Mac and Windows do not support this */ - struct timespec t; - double seconds; - - clock_gettime(CLOCK_REALTIME, &t); - seconds = (double) t.tv_sec + 1e-9*(double)t.tv_nsec; - return seconds; -#elif defined HAVE_GETTIMEOFDAY - struct timeval t; - double seconds; - - gettimeofday(&t, NULL); - seconds = (double) t.tv_sec + 1e-6*(double)t.tv_usec; - - return seconds; -#else - double seconds; - - seconds = time(NULL); - - return seconds; -#endif -} - -double -gmx_gettime_per_thread() -{ -#if _POSIX_THREAD_CPUTIME > 0 - struct timespec t; - double seconds; - - clock_gettime(CLOCK_THREAD_CPUTIME_ID, &t); - seconds = (double) t.tv_sec + 1e-9*(double)t.tv_nsec; - return seconds; -#else - return gmx_gettime(); -#endif -} - -void print_time(FILE *out, gmx_runtime_t *runtime, gmx_large_int_t step, - t_inputrec *ir, t_commrec gmx_unused *cr) +void print_time(FILE *out, + gmx_walltime_accounting_t walltime_accounting, + gmx_large_int_t step, + t_inputrec *ir, + t_commrec gmx_unused *cr) { time_t finish; char timebuf[STRLEN]; - double dt, time_per_step; + double dt, elapsed_seconds, time_per_step; char buf[48]; #ifndef GMX_THREAD_MPI @@ -159,10 +116,9 @@ void print_time(FILE *out, gmx_runtime_t *runtime, gmx_large_int_t step, if ((step >= ir->nstlist)) { double seconds_since_epoch = gmx_gettime(); - dt = seconds_since_epoch - runtime->start_time_stamp; - time_per_step = dt/(step - ir->init_step + 1); - - dt = (ir->nsteps + ir->init_step - step) * time_per_step; + elapsed_seconds = seconds_since_epoch - walltime_accounting_get_start_time_stamp(walltime_accounting); + time_per_step = elapsed_seconds/(step - ir->init_step + 1); + dt = (ir->nsteps + ir->init_step - step) * time_per_step; if (ir->nsteps >= 0) { @@ -176,7 +132,7 @@ void print_time(FILE *out, gmx_runtime_t *runtime, gmx_large_int_t step, } else { - fprintf(out, ", remaining runtime: %5d s ", (int)dt); + fprintf(out, ", remaining wall clock time: %5d s ", (int)dt); } } else @@ -195,31 +151,8 @@ void print_time(FILE *out, gmx_runtime_t *runtime, gmx_large_int_t step, fflush(out); } -void runtime_start(gmx_runtime_t *runtime) -{ - runtime->start_time_stamp = gmx_gettime(); - runtime->start_time_stamp_per_thread = gmx_gettime_per_thread(); - runtime->elapsed_run_time = 0; -} - -void runtime_end(gmx_runtime_t *runtime) -{ - double now, now_per_thread; - - now = gmx_gettime(); - now_per_thread = gmx_gettime_per_thread(); - - runtime->elapsed_run_time = now - runtime->start_time_stamp; - runtime->elapsed_run_time_per_thread = now_per_thread - runtime->start_time_stamp_per_thread; -} - -double runtime_get_elapsed_time(gmx_runtime_t *runtime) -{ - return gmx_gettime() - runtime->start_time_stamp; -} - void print_date_and_time(FILE *fplog, int nodeid, const char *title, - const gmx_runtime_t *runtime) + const gmx_walltime_accounting_t walltime_accounting) { int i; char timebuf[STRLEN]; @@ -228,9 +161,9 @@ void print_date_and_time(FILE *fplog, int nodeid, const char *title, if (fplog) { - if (runtime != NULL) + if (walltime_accounting != NULL) { - tmptime = (time_t) runtime->start_time_stamp; + tmptime = (time_t) walltime_accounting_get_start_time_stamp(walltime_accounting); gmx_ctime_r(&tmptime, timebuf, STRLEN); } else @@ -2443,7 +2376,7 @@ void do_pbc_mtop(FILE *fplog, int ePBC, matrix box, void finish_run(FILE *fplog, t_commrec *cr, t_inputrec *inputrec, t_nrnb nrnb[], gmx_wallcycle_t wcycle, - gmx_runtime_t *runtime, + gmx_walltime_accounting_t walltime_accounting, wallclock_gpu_t *gputimes, gmx_bool bWriteStat) { @@ -2451,8 +2384,10 @@ void finish_run(FILE *fplog, t_commrec *cr, t_nrnb *nrnb_tot = NULL; real delta_t; double nbfs, mflop; - double elapsed_run_time_over_all_ranks = 0; - double elapsed_run_time_per_thread_over_all_ranks = 0; + double elapsed_time, + elapsed_time_over_all_ranks, + elapsed_time_over_all_threads, + elapsed_time_over_all_threads_over_all_ranks; wallcycle_sum(cr, wcycle); if (cr->nnodes > 1) @@ -2468,27 +2403,27 @@ void finish_run(FILE *fplog, t_commrec *cr, nrnb_tot = nrnb; } + elapsed_time = walltime_accounting_get_elapsed_time(walltime_accounting); + elapsed_time_over_all_ranks = elapsed_time; + elapsed_time_over_all_threads = walltime_accounting_get_elapsed_time_over_all_threads(walltime_accounting); + elapsed_time_over_all_threads_over_all_ranks = elapsed_time_over_all_threads; #ifdef GMX_MPI if (cr->nnodes > 1) { - /* reduce elapsed_run_time over all MPI ranks in the current simulation */ - MPI_Allreduce(&runtime->elapsed_run_time, - &elapsed_run_time_over_all_ranks, + /* reduce elapsed_time over all MPI ranks in the current simulation */ + MPI_Allreduce(&elapsed_time, + &elapsed_time_over_all_ranks, 1, MPI_DOUBLE, MPI_SUM, cr->mpi_comm_mysim); - elapsed_run_time_over_all_ranks /= cr->nnodes; - /* reduce elapsed_run_time_per_thread over all MPI ranks in the current simulation */ - MPI_Allreduce(&runtime->elapsed_run_time_per_thread, - &elapsed_run_time_per_thread_over_all_ranks, + elapsed_time_over_all_ranks /= cr->nnodes; + /* Reduce elapsed_time_over_all_threads over all MPI ranks in the + * current simulation. */ + MPI_Allreduce(&elapsed_time_over_all_threads, + &elapsed_time_over_all_threads_over_all_ranks, 1, MPI_DOUBLE, MPI_SUM, cr->mpi_comm_mysim); } - else #endif - { - elapsed_run_time_over_all_ranks = runtime->elapsed_run_time; - elapsed_run_time_per_thread_over_all_ranks = runtime->elapsed_run_time_per_thread; - } if (SIMMASTER(cr)) { @@ -2533,7 +2468,8 @@ void finish_run(FILE *fplog, t_commrec *cr, if (SIMMASTER(cr)) { - wallcycle_print(fplog, cr->nnodes, cr->npmenodes, runtime->elapsed_run_time, + wallcycle_print(fplog, cr->nnodes, cr->npmenodes, + elapsed_time_over_all_ranks, wcycle, gputimes); if (EI_DYNAMICS(inputrec->eI)) @@ -2547,15 +2483,17 @@ void finish_run(FILE *fplog, t_commrec *cr, if (fplog) { - print_perf(fplog, elapsed_run_time_per_thread_over_all_ranks, - elapsed_run_time_over_all_ranks, - runtime->nsteps_done, delta_t, nbfs, mflop); + print_perf(fplog, elapsed_time_over_all_threads_over_all_ranks, + elapsed_time_over_all_ranks, + walltime_accounting_get_nsteps_done(walltime_accounting), + delta_t, nbfs, mflop); } if (bWriteStat) { - print_perf(stderr, elapsed_run_time_per_thread_over_all_ranks, - elapsed_run_time_over_all_ranks, - runtime->nsteps_done, delta_t, nbfs, mflop); + print_perf(stderr, elapsed_time_over_all_threads_over_all_ranks, + elapsed_time_over_all_ranks, + walltime_accounting_get_nsteps_done(walltime_accounting), + delta_t, nbfs, mflop); } } } diff --git a/src/gromacs/mdlib/tpi.c b/src/gromacs/mdlib/tpi.c index ec5da04f54..7666a24b7f 100644 --- a/src/gromacs/mdlib/tpi.c +++ b/src/gromacs/mdlib/tpi.c @@ -77,6 +77,7 @@ #include "gmxfio.h" #include "pme.h" #include "gbutil.h" +#include "gromacs/timing/walltime_accounting.h" #ifdef GMX_X86_SSE2 #include "gmx_x86_sse2.h" @@ -131,7 +132,7 @@ double do_tpi(FILE *fplog, t_commrec *cr, real gmx_unused cpt_period, real gmx_unused max_hours, const char gmx_unused *deviceOptions, unsigned long gmx_unused Flags, - gmx_runtime_t *runtime) + gmx_walltime_accounting_t walltime_accounting) { const char *TPI = "Test Particle Insertion"; gmx_localtop_t *top; @@ -254,9 +255,10 @@ double do_tpi(FILE *fplog, t_commrec *cr, snew(f, top_global->natoms); /* Print to log file */ - runtime_start(runtime); + walltime_accounting_start(walltime_accounting); print_date_and_time(fplog, cr->nodeid, - "Started Test Particle Insertion", runtime); + "Started Test Particle Insertion", + walltime_accounting); wallcycle_start(wcycle, ewcRUN); /* The last charge group is the group to be inserted */ @@ -786,7 +788,7 @@ double do_tpi(FILE *fplog, t_commrec *cr, bNotLastFrame = read_next_frame(oenv, status, &rerun_fr); } /* End of the loop */ - runtime_end(runtime); + walltime_accounting_end(walltime_accounting); close_trj(status); @@ -833,7 +835,7 @@ double do_tpi(FILE *fplog, t_commrec *cr, sfree(sum_UgembU); - runtime->nsteps_done = frame*inputrec->nsteps; + walltime_accounting_set_nsteps_done(walltime_accounting, frame*inputrec->nsteps); return 0; } diff --git a/src/gromacs/timing/CMakeLists.txt b/src/gromacs/timing/CMakeLists.txt new file mode 100644 index 0000000000..4d3d87fc3a --- /dev/null +++ b/src/gromacs/timing/CMakeLists.txt @@ -0,0 +1,42 @@ +# +# This file is part of the GROMACS molecular simulation package. +# +# Copyright (c) 2013, by the GROMACS development team, led by +# David van der Spoel, Berk Hess, Erik Lindahl, and including many +# others, as listed in the AUTHORS file in the top-level source +# directory and at http://www.gromacs.org. +# +# GROMACS is free software; you can redistribute it and/or +# modify it under the terms of the GNU Lesser General Public License +# as published by the Free Software Foundation; either version 2.1 +# of the License, or (at your option) any later version. +# +# GROMACS is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +# Lesser General Public License for more details. +# +# You should have received a copy of the GNU Lesser General Public +# License along with GROMACS; if not, see +# http://www.gnu.org/licenses, or write to the Free Software Foundation, +# Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA. +# +# If you want to redistribute modifications to GROMACS, please +# consider that scientific software is very special. Version +# control is crucial - bugs must be traceable. We will be happy to +# consider code for inclusion in the official distribution, but +# derived work must not be called official GROMACS. Details are found +# in the README & COPYING files - if they are missing, get the +# official version at http://www.gromacs.org. +# +# To help us fund GROMACS development, we humbly ask that you cite +# the research papers on the package. Check out http://www.gromacs.org. + +file(GLOB TIMING_SOURCES *.cpp *.c) +set(LIBGROMACS_SOURCES ${LIBGROMACS_SOURCES} ${TIMING_SOURCES} PARENT_SCOPE) + +# No installed headers for this module + +if (BUILD_TESTING) +# add_subdirectory(tests) +endif (BUILD_TESTING) diff --git a/src/gromacs/timing/walltime_accounting.c b/src/gromacs/timing/walltime_accounting.c new file mode 100644 index 0000000000..20e976c0d8 --- /dev/null +++ b/src/gromacs/timing/walltime_accounting.c @@ -0,0 +1,237 @@ +/* + * + * This source code is part of + * + * G R O M A C S + * + * GROningen MAchine for Chemical Simulations + * + * VERSION 3.2.0 + * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others. + * Copyright (c) 2013, The GROMACS development team, + * check out http://www.gromacs.org for more information. + + * This program is free software; you can redistribute it and/or + * modify it under the terms of the GNU General Public License + * as published by the Free Software Foundation; either version 2 + * of the License, or (at your option) any later version. + * + * If you want to redistribute modifications, please consider that + * scientific software is very special. Version control is crucial - + * bugs must be traceable. We will be happy to consider code for + * inclusion in the official distribution, but derived work must not + * be called official GROMACS. Details are found in the README & COPYING + * files - if they are missing, get the official version at www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the papers on the package - you can find them in the top README file. + * + * For more info, check our website at http://www.gromacs.org + * + * And Hey: + * GROwing Monsters And Cloning Shrimps + */ +#include "gromacs/timing/walltime_accounting.h" + +#ifdef HAVE_CONFIG_H +#include +#endif + +#include "gromacs/legacyheaders/smalloc.h" +#include "gromacs/legacyheaders/types/simple.h" + +#include +#ifdef HAVE_UNISTD_H +#include +#endif +#ifdef HAVE_SYS_TIME_H +#include +#endif + +/* TODO in future: convert gmx_walltime_accounting to a class, + * resolve who should have responsibility for recording the number of + * steps done, consider whether parts of finish_time, print_perf, + * wallcycle_print belong in this module. + * + * If/when any kind of task parallelism is implemented (even OpenMP + * regions simultaneously assigned to different tasks), consider + * whether this data structure (and/or cycle counters) should be + * maintained on a per-OpenMP-thread basis. */ + +/*! \brief Manages caching wall-clock time measurements for + * simulations */ +typedef struct gmx_walltime_accounting { + //! Seconds since the epoch recorded at the start of the simulation + double start_time_stamp; + //! Seconds since the epoch recorded at the start of the simulation for this thread + double start_time_stamp_per_thread; + //! Total seconds elapsed over the simulation + double elapsed_time; + //! Total seconds elapsed over the simulation running this thread + double elapsed_time_over_all_threads; + /*! \brief Number of OpenMP threads that will be launched by this + * MPI rank. + * + * This is used to scale elapsed_time_over_all_threads so + * that any combination of real MPI, thread MPI and OpenMP (even + * mdrun -ntomp_pme) processes/threads would (when run at maximum + * efficiency) return values such that the sum of + * elapsed_time_over_all_threads over all threads was constant + * with respect to parallelism implementation. */ + int numOpenMPThreads; + //! Set by integrators to report the amount of work they did + gmx_large_int_t nsteps_done; +} t_gmx_walltime_accounting; + +/*! \brief Calls system timing routines (e.g. clock_gettime) to get + * the (fractional) number of seconds elapsed since the epoch when + * this thread was executing. + * + * This can be used to measure system load. This can be unreliable if + * threads migrate between sockets. If thread-specific timers are not + * supported by the OS (e.g. if the OS is not POSIX-compliant), this + * function is implemented by gmx_gettime. */ +static double gmx_gettime_per_thread(); + +// TODO In principle, all this should get protected by checks that +// walltime_accounting is not null. In practice, that NULL condition +// does not happen, and future refactoring will likely enforce it by +// having the gmx_walltime_accounting_t object be owned by the runner +// object. When these become member functions, existence will be +// guaranteed. + +gmx_walltime_accounting_t +walltime_accounting_init(int numOpenMPThreads) +{ + gmx_walltime_accounting_t walltime_accounting; + + snew(walltime_accounting, 1); + walltime_accounting->start_time_stamp = 0; + walltime_accounting->start_time_stamp_per_thread = 0; + walltime_accounting->elapsed_time = 0; + walltime_accounting->nsteps_done = 0; + walltime_accounting->numOpenMPThreads = numOpenMPThreads; + + return walltime_accounting; +} + +void +walltime_accounting_destroy(gmx_walltime_accounting_t walltime_accounting) +{ + sfree(walltime_accounting); +} + +void +walltime_accounting_start(gmx_walltime_accounting_t walltime_accounting) +{ + walltime_accounting->start_time_stamp = gmx_gettime(); + walltime_accounting->start_time_stamp_per_thread = gmx_gettime_per_thread(); + walltime_accounting->elapsed_time = 0; + walltime_accounting->nsteps_done = 0; +} + +void +walltime_accounting_end(gmx_walltime_accounting_t walltime_accounting) +{ + double now, now_per_thread; + + now = gmx_gettime(); + now_per_thread = gmx_gettime_per_thread(); + + walltime_accounting->elapsed_time = now - walltime_accounting->start_time_stamp; + walltime_accounting->elapsed_time_over_all_threads = now_per_thread - walltime_accounting->start_time_stamp_per_thread; + /* For thread-MPI, the per-thread CPU timer makes this just + * work. For OpenMP threads, the per-thread CPU timer measurement + * needs to be multiplied by the number of OpenMP threads used, + * under the current assumption that all regions ever opened + * within a process are of the same size, and each thread should + * keep one core busy. + */ + walltime_accounting->elapsed_time_over_all_threads *= walltime_accounting->numOpenMPThreads; +} + +double +walltime_accounting_get_current_elapsed_time(gmx_walltime_accounting_t walltime_accounting) +{ + return gmx_gettime() - walltime_accounting->start_time_stamp; +} + +double +walltime_accounting_get_elapsed_time(gmx_walltime_accounting_t walltime_accounting) +{ + return walltime_accounting->elapsed_time; +} + +double +walltime_accounting_get_elapsed_time_over_all_threads(gmx_walltime_accounting_t walltime_accounting) +{ + return walltime_accounting->elapsed_time_over_all_threads; +} + +double +walltime_accounting_get_start_time_stamp(gmx_walltime_accounting_t walltime_accounting) +{ + return walltime_accounting->start_time_stamp; +} + +double +walltime_accounting_get_nsteps_done(gmx_walltime_accounting_t walltime_accounting) +{ + return walltime_accounting->nsteps_done; +} + +void +walltime_accounting_set_nsteps_done(gmx_walltime_accounting_t walltime_accounting, + gmx_large_int_t nsteps_done) +{ + walltime_accounting->nsteps_done = nsteps_done; +} + +double +gmx_gettime() +{ +#if defined HAVE_CLOCK_GETTIME && _POSIX_TIMERS >= 0 + /* Mac and Windows do not support this. For added fun, Windows + * defines _POSIX_TIMERS without actually providing the + * implementation. */ + struct timespec t; + double seconds; + + clock_gettime(CLOCK_REALTIME, &t); + seconds = (double) t.tv_sec + 1e-9*(double)t.tv_nsec; + + return seconds; +#elif defined HAVE_GETTIMEOFDAY + // Note that gettimeofday() is deprecated by POSIX, but since Mac + // and Windows do not yet support POSIX, we are still stuck. + struct timeval t; + double seconds; + + gettimeofday(&t, NULL); + seconds = (double) t.tv_sec + 1e-6*(double)t.tv_usec; + + return seconds; +#else + double seconds; + + seconds = time(NULL); + + return seconds; +#endif +} + +static double +gmx_gettime_per_thread() +{ +#if defined HAVE_CLOCK_GETTIME && _POSIX_THREAD_CPUTIME >= 0 + struct timespec t; + double seconds; + + clock_gettime(CLOCK_THREAD_CPUTIME_ID, &t); + seconds = (double) t.tv_sec + 1e-9*(double)t.tv_nsec; + + return seconds; +#else + return gmx_gettime(); +#endif +} diff --git a/src/gromacs/timing/walltime_accounting.h b/src/gromacs/timing/walltime_accounting.h new file mode 100644 index 0000000000..c2dde1bc55 --- /dev/null +++ b/src/gromacs/timing/walltime_accounting.h @@ -0,0 +1,110 @@ +/* + * + * This source code is part of + * + * G R O M A C S + * + * GROningen MAchine for Chemical Simulations + * + * VERSION 3.2.0 + * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others. + * Copyright (c) 1991-2000, University of Groningen, The Netherlands. + * Copyright (c) 2001-2004, The GROMACS development team, + * check out http://www.gromacs.org for more information. + + * This program is free software; you can redistribute it and/or + * modify it under the terms of the GNU General Public License + * as published by the Free Software Foundation; either version 2 + * of the License, or (at your option) any later version. + * + * If you want to redistribute modifications, please consider that + * scientific software is very special. Version control is crucial - + * bugs must be traceable. We will be happy to consider code for + * inclusion in the official distribution, but derived work must not + * be called official GROMACS. Details are found in the README & COPYING + * files - if they are missing, get the official version at www.gromacs.org. + * + * To help us fund GROMACS development, we humbly ask that you cite + * the papers on the package - you can find them in the top README file. + * + * For more info, check our website at http://www.gromacs.org + * + * And Hey: + * Gromacs Runs On Most of All Computer Systems + */ + +#ifndef GMX_TIMING_WALLTIME_ACCOUNTING_H +#define GMX_TIMING_WALLTIME_ACCOUNTING_H + +#include "gromacs/legacyheaders/types/simple.h" + +#ifdef __cplusplus +extern "C" { +#endif + +#if 0 +} +#endif + +/*! Contains per-process and per-thread data about elapsed wall-clock + * times and integration steps performed. */ +typedef struct gmx_walltime_accounting *gmx_walltime_accounting_t; + +//! Constructor +gmx_walltime_accounting_t +walltime_accounting_init(int numOpenMPThreads); + +//! Destructor +void +walltime_accounting_destroy(gmx_walltime_accounting_t walltime_accounting); + +/*! Record initial time stamps, e.g. at run end or counter + * re-initalization time */ +void +walltime_accounting_start(gmx_walltime_accounting_t walltime_accounting); + +/*! Measure and cache the elapsed wall-clock time since + * walltime_accounting_start */ +void +walltime_accounting_end(gmx_walltime_accounting_t walltime_accounting); + +/*! Measure and return the elapsed wall-clock time since + * walltime_accounting_start */ +double +walltime_accounting_get_current_elapsed_time(gmx_walltime_accounting_t walltime_accounting); + +//! Get the cached wall-clock time for this node +double +walltime_accounting_get_elapsed_time(gmx_walltime_accounting_t walltime_accounting); + +//! Get the cached wall-clock time, multiplied by the number of OpenMP threads +double +walltime_accounting_get_elapsed_time_over_all_threads(gmx_walltime_accounting_t walltime_accounting); + +//! Get the cached initial time stamp for this node +double +walltime_accounting_get_start_time_stamp(gmx_walltime_accounting_t walltime_accounting); + +//! Get the number of integration steps done +double +walltime_accounting_get_nsteps_done(gmx_walltime_accounting_t walltime_accounting); + +/*! Set the number of integration steps done + * + * TODO consider whether this should get done in walltime_accounting_end */ +void +walltime_accounting_set_nsteps_done(gmx_walltime_accounting_t walltime_accounting, + gmx_large_int_t nsteps_done); + +/*! \brief Calls system timing routines (e.g. clock_gettime) to get the + * (fractional) number of seconds elapsed since the epoch. + * + * Resolution is implementation-dependent, but typically nanoseconds + * or microseconds. */ +double gmx_gettime(); + +#ifdef __cplusplus +} +#endif + +#endif /* GMX_TIMING_WALLTIME_ACCOUNTING_H */ diff --git a/src/programs/mdrun/md.c b/src/programs/mdrun/md.c index 86a0dc5176..504a4ba514 100644 --- a/src/programs/mdrun/md.c +++ b/src/programs/mdrun/md.c @@ -89,6 +89,7 @@ #include "nbnxn_cuda_data_mgmt.h" #include "gromacs/utility/gmxmpi.h" +#include "gromacs/timing/walltime_accounting.h" #ifdef GMX_FAHCORE #include "corewrap.h" @@ -98,7 +99,7 @@ static void reset_all_counters(FILE *fplog, t_commrec *cr, gmx_large_int_t step, gmx_large_int_t *step_rel, t_inputrec *ir, gmx_wallcycle_t wcycle, t_nrnb *nrnb, - gmx_runtime_t *runtime, + gmx_walltime_accounting_t walltime_accounting, nbnxn_cuda_ptr_t cu_nbv) { char sbuf[STEPSTRSIZE]; @@ -123,8 +124,8 @@ static void reset_all_counters(FILE *fplog, t_commrec *cr, ir->nsteps -= *step_rel; *step_rel = 0; wallcycle_start(wcycle, ewcRUN); - runtime_start(runtime); - print_date_and_time(fplog, cr->nodeid, "Restarted time", runtime); + walltime_accounting_start(walltime_accounting); + print_date_and_time(fplog, cr->nodeid, "Restarted time", walltime_accounting); } double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], @@ -142,11 +143,11 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], real cpt_period, real max_hours, const char gmx_unused *deviceOptions, unsigned long Flags, - gmx_runtime_t *runtime) + gmx_walltime_accounting_t walltime_accounting) { gmx_mdoutf_t *outf; gmx_large_int_t step, step_rel; - double elapsed_run_time; + double elapsed_time; double t, t0, lam0[efptNR]; gmx_bool bGStatEveryStep, bGStat, bCalcVir, bCalcEner; gmx_bool bNS, bNStList, bSimAnn, bStopCM, bRerunMD, bNotLastFrame = FALSE, @@ -668,8 +669,8 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], } /* Set and write start time */ - runtime_start(runtime); - print_date_and_time(fplog, cr->nodeid, "Started mdrun", runtime); + walltime_accounting_start(walltime_accounting); + print_date_and_time(fplog, cr->nodeid, "Started mdrun", walltime_accounting); wallcycle_start(wcycle, ewcRUN); if (fplog) { @@ -1337,7 +1338,7 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], copy_mat(state->fvir_prev, force_vir); } - elapsed_run_time = runtime_get_elapsed_time(runtime); + elapsed_time = walltime_accounting_get_current_elapsed_time(walltime_accounting); /* Check whether everything is still allright */ if (((int)gmx_get_stop_condition() > handled_stop_condition) @@ -1374,7 +1375,7 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], handled_stop_condition = (int)gmx_get_stop_condition(); } else if (MASTER(cr) && (bNS || ir->nstlist <= 0) && - (max_hours > 0 && elapsed_run_time > max_hours*60.0*60.0*0.99) && + (max_hours > 0 && elapsed_time > max_hours*60.0*60.0*0.99) && gs.sig[eglsSTOPCOND] == 0 && gs.set[eglsSTOPCOND] == 0) { /* Signal to terminate the run */ @@ -1387,7 +1388,7 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], } if (bResetCountersHalfMaxH && MASTER(cr) && - elapsed_run_time > max_hours*60.0*60.0*0.495) + elapsed_time > max_hours*60.0*60.0*0.495) { gs.sig[eglsRESETCOUNTERS] = 1; } @@ -1424,7 +1425,7 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], if (MASTER(cr) && ((bGStat || !PAR(cr)) && cpt_period >= 0 && (cpt_period == 0 || - elapsed_run_time >= nchkpt*cpt_period*60.0)) && + elapsed_time >= nchkpt*cpt_period*60.0)) && gs.set[eglsCHKPT] == 0) { gs.sig[eglsCHKPT] = 1; @@ -1782,15 +1783,14 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], /* Gets written into the state at the beginning of next loop*/ state->fep_state = lamnew; } - - /* Remaining runtime */ + /* Print the remaining wall clock time for the run */ if (MULTIMASTER(cr) && (do_verbose || gmx_got_usr_signal()) && !bPMETuneRunning) { if (shellfc) { fprintf(stderr, "\n"); } - print_time(stderr, runtime, step, ir, cr); + print_time(stderr, walltime_accounting, step, ir, cr); } /* Replica exchange */ @@ -1917,7 +1917,7 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], gs.set[eglsRESETCOUNTERS] != 0) { /* Reset all the counters related to performance over the run */ - reset_all_counters(fplog, cr, step, &step_rel, ir, wcycle, nrnb, runtime, + reset_all_counters(fplog, cr, step, &step_rel, ir, wcycle, nrnb, walltime_accounting, fr->nbv != NULL && fr->nbv->bUseGPU ? fr->nbv->cu_nbv : NULL); wcycle_set_reset_counters(wcycle, -1); if (!(cr->duty & DUTY_PME)) @@ -1926,7 +1926,7 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], gmx_pme_send_resetcounters(cr, step); } /* Correct max_hours for the elapsed time */ - max_hours -= elapsed_run_time/(60.0*60.0); + max_hours -= elapsed_time/(60.0*60.0); bResetCountersHalfMaxH = FALSE; gs.set[eglsRESETCOUNTERS] = 0; } @@ -1935,8 +1935,8 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], /* End of main MD loop */ debug_gmx(); - /* Stop the time */ - runtime_end(runtime); + /* Stop measuring walltime */ + walltime_accounting_end(walltime_accounting); if (bRerunMD && MASTER(cr)) { @@ -1987,7 +1987,7 @@ double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[], print_replica_exchange_statistics(fplog, repl_ex); } - runtime->nsteps_done = step_rel; + walltime_accounting_set_nsteps_done(walltime_accounting, step_rel); return 0; } diff --git a/src/programs/mdrun/mdrun.cpp b/src/programs/mdrun/mdrun.cpp index cac5894f19..f98dea1622 100644 --- a/src/programs/mdrun/mdrun.cpp +++ b/src/programs/mdrun/mdrun.cpp @@ -527,7 +527,7 @@ int gmx_mdrun(int argc, char *argv[]) { "-confout", FALSE, etBOOL, {&bConfout}, "HIDDENWrite the last configuration with [TT]-c[tt] and force checkpointing at the last step" }, { "-stepout", FALSE, etINT, {&nstepout}, - "HIDDENFrequency of writing the remaining runtime" }, + "HIDDENFrequency of writing the remaining wall clock time for the run" }, { "-resetstep", FALSE, etINT, {&resetstep}, "HIDDENReset cycle counters after these many time steps" }, { "-resethway", FALSE, etBOOL, {&bResetCountersHalfWay}, diff --git a/src/programs/mdrun/runner.c b/src/programs/mdrun/runner.c index 266f018fcd..5bac378eb4 100644 --- a/src/programs/mdrun/runner.c +++ b/src/programs/mdrun/runner.c @@ -921,38 +921,38 @@ int mdrunner(gmx_hw_opt_t *hw_opt, int repl_ex_seed, real pforce, real cpt_period, real max_hours, const char *deviceOptions, unsigned long Flags) { - gmx_bool bForceUseGPU, bTryUseGPU; - double nodetime = 0, realtime; - t_inputrec *inputrec; - t_state *state = NULL; - matrix box; - gmx_ddbox_t ddbox = {0}; - int npme_major, npme_minor; - real tmpr1, tmpr2; - t_nrnb *nrnb; - gmx_mtop_t *mtop = NULL; - t_mdatoms *mdatoms = NULL; - t_forcerec *fr = NULL; - t_fcdata *fcd = NULL; - real ewaldcoeff = 0; - gmx_pme_t *pmedata = NULL; - gmx_vsite_t *vsite = NULL; - gmx_constr_t constr; - int i, m, nChargePerturbed = -1, status, nalloc; - char *gro; - gmx_wallcycle_t wcycle; - gmx_bool bReadRNG, bReadEkin; - int list; - gmx_runtime_t runtime; - int rc; - gmx_large_int_t reset_counters; - gmx_edsam_t ed = NULL; - t_commrec *cr_old = cr; - int nthreads_pme = 1; - int nthreads_pp = 1; - gmx_membed_t membed = NULL; - gmx_hw_info_t *hwinfo = NULL; - master_inf_t minf = {-1, FALSE}; + gmx_bool bForceUseGPU, bTryUseGPU; + double nodetime = 0, realtime; + t_inputrec *inputrec; + t_state *state = NULL; + matrix box; + gmx_ddbox_t ddbox = {0}; + int npme_major, npme_minor; + real tmpr1, tmpr2; + t_nrnb *nrnb; + gmx_mtop_t *mtop = NULL; + t_mdatoms *mdatoms = NULL; + t_forcerec *fr = NULL; + t_fcdata *fcd = NULL; + real ewaldcoeff = 0; + gmx_pme_t *pmedata = NULL; + gmx_vsite_t *vsite = NULL; + gmx_constr_t constr; + int i, m, nChargePerturbed = -1, status, nalloc; + char *gro; + gmx_wallcycle_t wcycle; + gmx_bool bReadRNG, bReadEkin; + int list; + gmx_walltime_accounting_t walltime_accounting = NULL; + int rc; + gmx_large_int_t reset_counters; + gmx_edsam_t ed = NULL; + t_commrec *cr_old = cr; + int nthreads_pme = 1; + int nthreads_pp = 1; + gmx_membed_t membed = NULL; + gmx_hw_info_t *hwinfo = NULL; + master_inf_t minf = {-1, FALSE}; /* CAUTION: threads may be started later on in this function, so cr doesn't reflect the final parallel state right now */ @@ -1545,6 +1545,9 @@ int mdrunner(gmx_hw_opt_t *hw_opt, if (cr->duty & DUTY_PP) { + /* Assumes uniform use of the number of OpenMP threads */ + walltime_accounting = walltime_accounting_init(gmx_omp_nthreads_get(emntDefault)); + if (inputrec->ePull != epullNO) { /* Initialize pull code */ @@ -1584,7 +1587,7 @@ int mdrunner(gmx_hw_opt_t *hw_opt, cpt_period, max_hours, deviceOptions, Flags, - &runtime); + walltime_accounting); if (inputrec->ePull != epullNO) { @@ -1600,7 +1603,8 @@ int mdrunner(gmx_hw_opt_t *hw_opt, else { /* do PME only */ - gmx_pmeonly(*pmedata, cr, nrnb, wcycle, &runtime, ewaldcoeff, inputrec); + walltime_accounting = walltime_accounting_init(gmx_omp_nthreads_get(emntPME)); + gmx_pmeonly(*pmedata, cr, nrnb, wcycle, walltime_accounting, ewaldcoeff, inputrec); } wallcycle_stop(wcycle, ewcRUN); @@ -1609,7 +1613,7 @@ int mdrunner(gmx_hw_opt_t *hw_opt, * if rerunMD, don't write last frame again */ finish_run(fplog, cr, - inputrec, nrnb, wcycle, &runtime, + inputrec, nrnb, wcycle, walltime_accounting, fr != NULL && fr->nbv != NULL && fr->nbv->bUseGPU ? nbnxn_cuda_get_timings(fr->nbv->cu_nbv) : NULL, EI_DYNAMICS(inputrec->eI) && !MULTISIM(cr)); @@ -1636,7 +1640,8 @@ int mdrunner(gmx_hw_opt_t *hw_opt, gmx_hardware_info_free(hwinfo); /* Does what it says */ - print_date_and_time(fplog, cr->nodeid, "Finished mdrun", &runtime); + print_date_and_time(fplog, cr->nodeid, "Finished mdrun", walltime_accounting); + walltime_accounting_destroy(walltime_accounting); /* Close logfile already here if we were appending to it */ if (MASTER(cr) && (Flags & MD_APPENDFILES))