2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2009,2010,2011,2012,2013,2014,2015,2016,2017,2018,2019, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
49 #ifdef HAVE_SYS_TIME_H
53 #include "gromacs/commandline/pargs.h"
54 #include "gromacs/ewald/pme.h"
55 #include "gromacs/fft/calcgrid.h"
56 #include "gromacs/fileio/checkpoint.h"
57 #include "gromacs/fileio/tpxio.h"
58 #include "gromacs/math/utilities.h"
59 #include "gromacs/math/vec.h"
60 #include "gromacs/mdlib/perf_est.h"
61 #include "gromacs/mdtypes/commrec.h"
62 #include "gromacs/mdtypes/inputrec.h"
63 #include "gromacs/mdtypes/md_enums.h"
64 #include "gromacs/mdtypes/state.h"
65 #include "gromacs/pbcutil/pbc.h"
66 #include "gromacs/taskassignment/usergpuids.h"
67 #include "gromacs/timing/walltime_accounting.h"
68 #include "gromacs/topology/topology.h"
69 #include "gromacs/utility/arraysize.h"
70 #include "gromacs/utility/baseversion.h"
71 #include "gromacs/utility/cstringutil.h"
72 #include "gromacs/utility/fatalerror.h"
73 #include "gromacs/utility/futil.h"
74 #include "gromacs/utility/gmxassert.h"
75 #include "gromacs/utility/path.h"
76 #include "gromacs/utility/smalloc.h"
77 #include "gromacs/utility/stringutil.h"
79 /* Enum for situations that can occur during log file parsing, the
80 * corresponding string entries can be found in do_the_tests() in
81 * const char* ParseLog[] */
82 /* TODO clean up CamelCasing of these enum names */
88 eParselogResetProblem,
92 eParselogLargePrimeFactor,
93 eParselogMismatchOfNumberOfPPRanksAndAvailableGPUs,
102 int nPMEnodes; /* number of PME-only nodes used in this test */
103 int nx, ny, nz; /* DD grid */
104 int guessPME; /* if nPMEnodes == -1, this is the guessed number of PME nodes */
105 double *Gcycles; /* This can contain more than one value if doing multiple tests */
109 float *PME_f_load; /* PME mesh/force load average*/
110 float PME_f_load_Av; /* Average average ;) ... */
111 char *mdrun_cmd_line; /* Mdrun command line used for this test */
117 int nr_inputfiles; /* The number of tpr and mdp input files */
118 int64_t orig_sim_steps; /* Number of steps to be done in the real simulation */
119 int64_t orig_init_step; /* Init step for the real simulation */
120 real *rcoulomb; /* The coulomb radii [0...nr_inputfiles] */
121 real *rvdw; /* The vdW radii */
122 real *rlist; /* Neighbourlist cutoff radius */
123 int *nkx, *nky, *nkz;
124 real *fsx, *fsy, *fsz; /* Fourierspacing in x,y,z dimension */
128 static void sep_line(FILE *fp)
130 fprintf(fp, "\n------------------------------------------------------------\n");
134 /* Wrapper for system calls */
135 static int gmx_system_call(char *command)
137 return ( system(command) );
141 /* Check if string starts with substring */
142 static gmx_bool str_starts(const char *string, const char *substring)
144 return ( std::strncmp(string, substring, std::strlen(substring)) == 0);
148 static void cleandata(t_perf *perfdata, int test_nr)
150 perfdata->Gcycles[test_nr] = 0.0;
151 perfdata->ns_per_day[test_nr] = 0.0;
152 perfdata->PME_f_load[test_nr] = 0.0;
156 static void remove_if_exists(const char *fn)
160 fprintf(stdout, "Deleting %s\n", fn);
166 static void finalize(const char *fn_out)
172 fp = fopen(fn_out, "r");
173 fprintf(stdout, "\n\n");
175 while (fgets(buf, STRLEN-1, fp) != nullptr)
177 fprintf(stdout, "%s", buf);
180 fprintf(stdout, "\n\n");
185 eFoundNothing, eFoundDDStr, eFoundAccountingStr, eFoundCycleStr
188 static int parse_logfile(const char *logfile, const char *errfile,
189 t_perf *perfdata, int test_nr, int presteps, int64_t cpt_steps,
193 char line[STRLEN], dumstring[STRLEN], dumstring2[STRLEN];
194 const char matchstrdd[] = "Domain decomposition grid";
195 const char matchstrcr[] = "resetting all time and cycle counters";
196 const char matchstrbal[] = "Average PME mesh/force load:";
197 const char matchstring[] = "R E A L C Y C L E A N D T I M E A C C O U N T I N G";
198 const char errSIG[] = "signal, stopping at the next";
200 float dum1, dum2, dum3, dum4;
203 int64_t resetsteps = -1;
204 gmx_bool bFoundResetStr = FALSE;
205 gmx_bool bResetChecked = FALSE;
208 if (!gmx_fexist(logfile))
210 fprintf(stderr, "WARNING: Could not find logfile %s.\n", logfile);
211 cleandata(perfdata, test_nr);
212 return eParselogNotFound;
215 fp = fopen(logfile, "r");
216 perfdata->PME_f_load[test_nr] = -1.0;
217 perfdata->guessPME = -1;
219 iFound = eFoundNothing;
222 iFound = eFoundDDStr; /* Skip some case statements */
225 while (fgets(line, STRLEN, fp) != nullptr)
227 /* Remove leading spaces */
230 /* Check for TERM and INT signals from user: */
231 if (std::strstr(line, errSIG) != nullptr)
234 cleandata(perfdata, test_nr);
235 return eParselogTerm;
238 /* Check whether cycle resetting worked */
239 if (presteps > 0 && !bFoundResetStr)
241 if (std::strstr(line, matchstrcr) != nullptr)
243 sprintf(dumstring, "step %s", "%" SCNd64);
244 sscanf(line, dumstring, &resetsteps);
245 bFoundResetStr = TRUE;
246 if (resetsteps == presteps+cpt_steps)
248 bResetChecked = TRUE;
252 sprintf(dumstring, "%" PRId64, resetsteps);
253 sprintf(dumstring2, "%" PRId64, presteps+cpt_steps);
254 fprintf(stderr, "WARNING: Time step counters were reset at step %s,\n"
255 " though they were supposed to be reset at step %s!\n",
256 dumstring, dumstring2);
261 /* Look for strings that appear in a certain order in the log file: */
265 /* Look for domain decomp grid and separate PME nodes: */
266 if (str_starts(line, matchstrdd))
268 sscanf(line, "Domain decomposition grid %d x %d x %d, separate PME ranks %d",
269 &(perfdata->nx), &(perfdata->ny), &(perfdata->nz), &npme);
270 if (perfdata->nPMEnodes == -1)
272 perfdata->guessPME = npme;
274 else if (perfdata->nPMEnodes != npme)
276 gmx_fatal(FARGS, "PME ranks from command line and output file are not identical");
278 iFound = eFoundDDStr;
280 /* Catch a few errors that might have occurred: */
281 else if (str_starts(line, "There is no domain decomposition for"))
284 return eParselogNoDDGrid;
286 else if (str_starts(line, "The number of ranks you selected"))
289 return eParselogLargePrimeFactor;
291 else if (str_starts(line, "reading tpx file"))
294 return eParselogTPXVersion;
296 else if (str_starts(line, "The -dd or -npme option request a parallel simulation"))
299 return eParselogNotParallel;
303 /* Even after the "Domain decomposition grid" string was found,
304 * it could be that mdrun had to quit due to some error. */
305 if (str_starts(line, "Incorrect launch configuration: mismatching number of"))
308 return eParselogMismatchOfNumberOfPPRanksAndAvailableGPUs;
310 else if (str_starts(line, "Some of the requested GPUs do not exist"))
313 return eParselogGpuProblem;
315 /* Look for PME mesh/force balance (not necessarily present, though) */
316 else if (str_starts(line, matchstrbal))
318 sscanf(&line[std::strlen(matchstrbal)], "%f", &(perfdata->PME_f_load[test_nr]));
320 /* Look for matchstring */
321 else if (str_starts(line, matchstring))
323 iFound = eFoundAccountingStr;
326 case eFoundAccountingStr:
327 /* Already found matchstring - look for cycle data */
328 if (str_starts(line, "Total "))
330 sscanf(line, "Total %*f %lf", &(perfdata->Gcycles[test_nr]));
331 iFound = eFoundCycleStr;
335 /* Already found cycle data - look for remaining performance info and return */
336 if (str_starts(line, "Performance:"))
338 ndum = sscanf(line, "%s %f %f %f %f", dumstring, &dum1, &dum2, &dum3, &dum4);
339 /* (ns/day) is the second last entry, depending on whether GMX_DETAILED_PERF_STATS was set in print_perf(), nrnb.c */
340 perfdata->ns_per_day[test_nr] = (ndum == 5) ? dum3 : dum1;
342 if (bResetChecked || presteps == 0)
348 return eParselogResetProblem;
355 /* Close the log file */
358 /* Check why there is no performance data in the log file.
359 * Did a fatal errors occur? */
360 if (gmx_fexist(errfile))
362 fp = fopen(errfile, "r");
363 while (fgets(line, STRLEN, fp) != nullptr)
365 if (str_starts(line, "Fatal error:") )
367 if (fgets(line, STRLEN, fp) != nullptr)
369 fprintf(stderr, "\nWARNING: An error occurred during this benchmark:\n"
373 cleandata(perfdata, test_nr);
374 return eParselogFatal;
381 fprintf(stderr, "WARNING: Could not find stderr file %s.\n", errfile);
384 /* Giving up ... we could not find out why there is no performance data in
386 fprintf(stdout, "No performance data in log file.\n");
387 cleandata(perfdata, test_nr);
389 return eParselogNoPerfData;
393 static gmx_bool analyze_data(
402 int *index_tpr, /* OUT: Nr of mdp file with best settings */
403 int *npme_optimal) /* OUT: Optimal number of PME nodes */
406 int line = 0, line_win = -1;
407 int k_win = -1, i_win = -1, winPME;
408 double s = 0.0; /* standard deviation */
411 char str_PME_f_load[13];
412 gmx_bool bCanUseOrigTPR;
413 gmx_bool bRefinedCoul, bRefinedVdW, bRefinedGrid;
419 fprintf(fp, "Summary of successful runs:\n");
420 fprintf(fp, "Line tpr PME ranks Gcycles Av. Std.dev. ns/day PME/f");
423 fprintf(fp, " DD grid");
429 for (k = 0; k < ntprs; k++)
431 for (i = 0; i < ntests; i++)
433 /* Select the right dataset: */
434 pd = &(perfdata[k][i]);
436 pd->Gcycles_Av = 0.0;
437 pd->PME_f_load_Av = 0.0;
438 pd->ns_per_day_Av = 0.0;
440 if (pd->nPMEnodes == -1)
442 sprintf(strbuf, "(%3d)", pd->guessPME);
446 sprintf(strbuf, " ");
449 /* Get the average run time of a setting */
450 for (j = 0; j < nrepeats; j++)
452 pd->Gcycles_Av += pd->Gcycles[j];
453 pd->PME_f_load_Av += pd->PME_f_load[j];
455 pd->Gcycles_Av /= nrepeats;
456 pd->PME_f_load_Av /= nrepeats;
458 for (j = 0; j < nrepeats; j++)
460 if (pd->ns_per_day[j] > 0.0)
462 pd->ns_per_day_Av += pd->ns_per_day[j];
466 /* Somehow the performance number was not aquired for this run,
467 * therefor set the average to some negative value: */
468 pd->ns_per_day_Av = -1.0F*nrepeats;
472 pd->ns_per_day_Av /= nrepeats;
475 if (pd->PME_f_load_Av > 0.0)
477 sprintf(str_PME_f_load, "%12.3f", pd->PME_f_load_Av);
481 sprintf(str_PME_f_load, "%s", " - ");
485 /* We assume we had a successful run if both averages are positive */
486 if (pd->Gcycles_Av > 0.0 && pd->ns_per_day_Av > 0.0)
488 /* Output statistics if repeats were done */
491 /* Calculate the standard deviation */
493 for (j = 0; j < nrepeats; j++)
495 s += gmx::square( pd->Gcycles[j] - pd->Gcycles_Av );
500 fprintf(fp, "%4d %3d %4d%s %12.3f %12.3f %12.3f %s",
501 line, k, pd->nPMEnodes, strbuf, pd->Gcycles_Av, s,
502 pd->ns_per_day_Av, str_PME_f_load);
505 fprintf(fp, " %3d %3d %3d", pd->nx, pd->ny, pd->nz);
509 /* Store the index of the best run found so far in 'winner': */
510 if ( (k_win == -1) || (pd->Gcycles_Av < perfdata[k_win][i_win].Gcycles_Av) )
523 gmx_fatal(FARGS, "None of the runs was successful! Check %s for problems.", fn);
528 winPME = perfdata[k_win][i_win].nPMEnodes;
532 /* We stuck to a fixed number of PME-only nodes */
533 sprintf(strbuf, "settings No. %d", k_win);
537 /* We have optimized the number of PME-only nodes */
540 sprintf(strbuf, "%s", "the automatic number of PME ranks");
544 sprintf(strbuf, "%d PME ranks", winPME);
547 fprintf(fp, "Best performance was achieved with %s", strbuf);
548 if ((nrepeats > 1) && (ntests > 1))
550 fprintf(fp, " (see line %d)", line_win);
554 /* Only mention settings if they were modified: */
555 bRefinedCoul = !gmx_within_tol(info->rcoulomb[k_win], info->rcoulomb[0], GMX_REAL_EPS);
556 bRefinedVdW = !gmx_within_tol(info->rvdw[k_win], info->rvdw[0], GMX_REAL_EPS);
557 bRefinedGrid = !(info->nkx[k_win] == info->nkx[0] &&
558 info->nky[k_win] == info->nky[0] &&
559 info->nkz[k_win] == info->nkz[0]);
561 if (bRefinedCoul || bRefinedVdW || bRefinedGrid)
563 fprintf(fp, "Optimized PME settings:\n");
564 bCanUseOrigTPR = FALSE;
568 bCanUseOrigTPR = TRUE;
573 fprintf(fp, " New Coulomb radius: %f nm (was %f nm)\n", info->rcoulomb[k_win], info->rcoulomb[0]);
578 fprintf(fp, " New Van der Waals radius: %f nm (was %f nm)\n", info->rvdw[k_win], info->rvdw[0]);
583 fprintf(fp, " New Fourier grid xyz: %d %d %d (was %d %d %d)\n", info->nkx[k_win], info->nky[k_win], info->nkz[k_win],
584 info->nkx[0], info->nky[0], info->nkz[0]);
587 if (bCanUseOrigTPR && ntprs > 1)
589 fprintf(fp, "and original PME settings.\n");
594 /* Return the index of the mdp file that showed the highest performance
595 * and the optimal number of PME nodes */
597 *npme_optimal = winPME;
599 return bCanUseOrigTPR;
603 /* Get the commands we need to set up the runs from environment variables */
604 static void get_program_paths(gmx_bool bThreads, char *cmd_mpirun[], char *cmd_mdrun[])
607 const char def_mpirun[] = "mpirun";
609 const char empty_mpirun[] = "";
611 /* Get the commands we need to set up the runs from environment variables */
614 if ( (cp = getenv("MPIRUN")) != nullptr)
616 *cmd_mpirun = gmx_strdup(cp);
620 *cmd_mpirun = gmx_strdup(def_mpirun);
625 *cmd_mpirun = gmx_strdup(empty_mpirun);
628 if (*cmd_mdrun == nullptr)
630 /* The use of MDRUN is deprecated, but made available in 5.1
631 for backward compatibility. It may be removed in a future
633 if ( (cp = getenv("MDRUN" )) != nullptr)
635 *cmd_mdrun = gmx_strdup(cp);
639 gmx_fatal(FARGS, "The way to call mdrun must be set in the -mdrun command-line flag.");
644 /* Check that the commands will run mdrun (perhaps via mpirun) by
645 * running a very quick test simulation. Requires MPI environment or
646 * GPU support to be available if applicable. */
647 /* TODO implement feature to parse the log file to get the list of
648 compatible GPUs from mdrun, if the user of gmx tune-pme has not
650 static void check_mdrun_works(gmx_bool bThreads,
651 const char *cmd_mpirun,
653 const char *cmd_mdrun,
654 gmx_bool bNeedGpuSupport)
656 char *command = nullptr;
660 const char filename[] = "benchtest.log";
662 /* This string should always be identical to the one in copyrite.c,
663 * gmx_print_version_info() in the GMX_MPI section */
664 const char match_mpi[] = "MPI library: MPI";
665 const char match_mdrun[] = "Executable: ";
666 const char match_nogpu[] = "GPU support: disabled";
667 gmx_bool bMdrun = FALSE;
668 gmx_bool bMPI = FALSE;
669 gmx_bool bHaveGpuSupport = TRUE;
671 /* Run a small test to see whether mpirun + mdrun work */
672 fprintf(stdout, "Making sure that mdrun can be executed. ");
675 snew(command, std::strlen(cmd_mdrun) + std::strlen(cmd_np) + std::strlen(filename) + 50);
676 sprintf(command, "%s%s -version -maxh 0.001 1> %s 2>&1", cmd_mdrun, cmd_np, filename);
680 snew(command, std::strlen(cmd_mpirun) + std::strlen(cmd_np) + std::strlen(cmd_mdrun) + std::strlen(filename) + 50);
681 sprintf(command, "%s%s%s -version -maxh 0.001 1> %s 2>&1", cmd_mpirun, cmd_np, cmd_mdrun, filename);
683 fprintf(stdout, "Trying '%s' ... ", command);
684 make_backup(filename);
685 gmx_system_call(command);
687 /* Check if we find the characteristic string in the output: */
688 if (!gmx_fexist(filename))
690 gmx_fatal(FARGS, "Output from test run could not be found.");
693 fp = fopen(filename, "r");
694 /* We need to scan the whole output file, since sometimes the queuing system
695 * also writes stuff to stdout/err */
698 cp = fgets(line, STRLEN, fp);
701 if (str_starts(line, match_mdrun) )
705 if (str_starts(line, match_mpi) )
709 if (str_starts(line, match_nogpu) )
711 bHaveGpuSupport = FALSE;
721 gmx_fatal(FARGS, "Need a threaded version of mdrun. This one\n"
723 "seems to have been compiled with MPI instead.",
731 gmx_fatal(FARGS, "Need an MPI-enabled version of mdrun. This one\n"
733 "seems to have been compiled without MPI support.",
740 gmx_fatal(FARGS, "Cannot execute mdrun. Please check %s for problems!",
744 if (bNeedGpuSupport && !bHaveGpuSupport)
746 gmx_fatal(FARGS, "The mdrun executable did not have the expected GPU support.");
749 fprintf(stdout, "passed.\n");
756 /* Handles the no-GPU case by emitting an empty string. */
757 static std::string make_gpu_id_command_line(const char *eligible_gpu_ids)
759 /* If the user has given no eligible GPU IDs, or we're trying the
760 * default behaviour, then there is nothing for tune_pme to give
761 * to mdrun -gpu_id */
762 if (eligible_gpu_ids != nullptr)
764 return gmx::formatString("-gpu_id %s", eligible_gpu_ids);
768 return std::string();
771 static void launch_simulation(
772 gmx_bool bLaunch, /* Should the simulation be launched? */
773 FILE *fp, /* General log file */
774 gmx_bool bThreads, /* whether to use threads */
775 char *cmd_mpirun, /* Command for mpirun */
776 char *cmd_np, /* Switch for -np or -ntmpi or empty */
777 char *cmd_mdrun, /* Command for mdrun */
778 char *args_for_mdrun, /* Arguments for mdrun */
779 const char *simulation_tpr, /* This tpr will be simulated */
780 int nPMEnodes, /* Number of PME ranks to use */
781 const char *eligible_gpu_ids) /* Available GPU IDs for
782 * constructing mdrun command lines */
787 /* Make enough space for the system call command,
788 * (200 extra chars for -npme ... etc. options should suffice): */
789 snew(command, std::strlen(cmd_mpirun)+std::strlen(cmd_mdrun)+std::strlen(cmd_np)+std::strlen(args_for_mdrun)+std::strlen(simulation_tpr)+200);
791 auto cmd_gpu_ids = make_gpu_id_command_line(eligible_gpu_ids);
793 /* Note that the -passall options requires args_for_mdrun to be at the end
794 * of the command line string */
797 sprintf(command, "%s%s-npme %d -s %s %s %s",
798 cmd_mdrun, cmd_np, nPMEnodes, simulation_tpr, args_for_mdrun, cmd_gpu_ids.c_str());
802 sprintf(command, "%s%s%s -npme %d -s %s %s %s",
803 cmd_mpirun, cmd_np, cmd_mdrun, nPMEnodes, simulation_tpr, args_for_mdrun, cmd_gpu_ids.c_str());
806 fprintf(fp, "%s this command line to launch the simulation:\n\n%s", bLaunch ? "Using" : "Please use", command);
810 /* Now the real thing! */
813 fprintf(stdout, "\nLaunching simulation with best parameters now.\nExecuting '%s'", command);
816 gmx_system_call(command);
821 static void modify_PMEsettings(
822 int64_t simsteps, /* Set this value as number of time steps */
823 int64_t init_step, /* Set this value as init_step */
824 const char *fn_best_tpr, /* tpr file with the best performance */
825 const char *fn_sim_tpr) /* name of tpr file to be launched */
831 t_inputrec irInstance;
832 t_inputrec *ir = &irInstance;
833 read_tpx_state(fn_best_tpr, ir, &state, &mtop);
835 /* Reset nsteps and init_step to the value of the input .tpr file */
836 ir->nsteps = simsteps;
837 ir->init_step = init_step;
839 /* Write the tpr file which will be launched */
840 sprintf(buf, "Writing optimized simulation file %s with nsteps=%s.\n", fn_sim_tpr, "%" PRId64);
841 fprintf(stdout, buf, ir->nsteps);
843 write_tpx_state(fn_sim_tpr, ir, &state, &mtop);
846 static gmx_bool can_scale_rvdw(int vdwtype)
848 return (evdwCUT == vdwtype ||
852 #define EPME_SWITCHED(e) ((e) == eelPMESWITCH || (e) == eelPMEUSERSWITCH)
854 /* Make additional TPR files with more computational load for the
855 * direct space processors: */
856 static void make_benchmark_tprs(
857 const char *fn_sim_tpr, /* READ : User-provided tpr file */
858 char *fn_bench_tprs[], /* WRITE: Names of benchmark tpr files */
859 int64_t benchsteps, /* Number of time steps for benchmark runs */
860 int64_t statesteps, /* Step counter in checkpoint file */
861 real rmin, /* Minimal Coulomb radius */
862 real rmax, /* Maximal Coulomb radius */
863 bool bScaleRvdw, /* Scale rvdw along with rcoulomb */
864 const int *ntprs, /* No. of TPRs to write, each with a different
865 rcoulomb and fourierspacing */
866 t_inputinfo *info, /* Contains information about mdp file options */
867 FILE *fp) /* Write the output here */
872 real nlist_buffer; /* Thickness of the buffer regions for PME-switch potentials */
875 gmx_bool bNote = FALSE;
876 real add; /* Add this to rcoul for the next test */
877 real fac = 1.0; /* Scaling factor for Coulomb radius */
878 real fourierspacing; /* Basic fourierspacing from tpr */
881 sprintf(buf, "Making benchmark tpr file%s with %s time step%s",
882 *ntprs > 1 ? "s" : "", "%" PRId64, benchsteps > 1 ? "s" : "");
883 fprintf(stdout, buf, benchsteps);
886 sprintf(buf, " (adding %s steps from checkpoint file)", "%" PRId64);
887 fprintf(stdout, buf, statesteps);
888 benchsteps += statesteps;
890 fprintf(stdout, ".\n");
892 t_inputrec irInstance;
893 t_inputrec *ir = &irInstance;
894 read_tpx_state(fn_sim_tpr, ir, &state, &mtop);
896 /* Check if some kind of PME was chosen */
897 if (EEL_PME(ir->coulombtype) == FALSE)
899 gmx_fatal(FARGS, "Can only do optimizations for simulations with %s electrostatics.",
903 /* Check if rcoulomb == rlist, which is necessary for plain PME. */
904 if ( (ir->cutoff_scheme != ecutsVERLET) &&
905 (eelPME == ir->coulombtype) && !(ir->rcoulomb == ir->rlist))
907 gmx_fatal(FARGS, "%s requires rcoulomb (%f) to be equal to rlist (%f).",
908 EELTYPE(eelPME), ir->rcoulomb, ir->rlist);
910 /* For other PME types, rcoulomb is allowed to be smaller than rlist */
911 else if (ir->rcoulomb > ir->rlist)
913 gmx_fatal(FARGS, "%s requires rcoulomb (%f) to be equal to or smaller than rlist (%f)",
914 EELTYPE(ir->coulombtype), ir->rcoulomb, ir->rlist);
917 if (bScaleRvdw && ir->rvdw != ir->rcoulomb)
919 fprintf(stdout, "NOTE: input rvdw != rcoulomb, will not scale rvdw\n");
923 /* Reduce the number of steps for the benchmarks */
924 info->orig_sim_steps = ir->nsteps;
925 ir->nsteps = benchsteps;
926 /* We must not use init_step from the input tpr file for the benchmarks */
927 info->orig_init_step = ir->init_step;
930 /* For PME-switch potentials, keep the radial distance of the buffer region */
931 nlist_buffer = ir->rlist - ir->rcoulomb;
933 /* Determine length of triclinic box vectors */
934 for (d = 0; d < DIM; d++)
937 for (i = 0; i < DIM; i++)
939 box_size[d] += state.box[d][i]*state.box[d][i];
941 box_size[d] = std::sqrt(box_size[d]);
944 if (ir->fourier_spacing > 0)
946 info->fsx[0] = ir->fourier_spacing;
947 info->fsy[0] = ir->fourier_spacing;
948 info->fsz[0] = ir->fourier_spacing;
952 /* Reconstruct fourierspacing per dimension from the number of grid points and box size */
953 info->fsx[0] = box_size[XX]/ir->nkx;
954 info->fsy[0] = box_size[YY]/ir->nky;
955 info->fsz[0] = box_size[ZZ]/ir->nkz;
958 /* If no value for the fourierspacing was provided on the command line, we
959 * use the reconstruction from the tpr file */
960 if (ir->fourier_spacing > 0)
962 /* Use the spacing from the tpr */
963 fourierspacing = ir->fourier_spacing;
967 /* Use the maximum observed spacing */
968 fourierspacing = std::max(std::max(info->fsx[0], info->fsy[0]), info->fsz[0]);
971 fprintf(stdout, "Calculating PME grid points on the basis of a fourierspacing of %f nm\n", fourierspacing);
973 /* For performance comparisons the number of particles is useful to have */
974 fprintf(fp, " Number of particles : %d\n", mtop.natoms);
976 /* Print information about settings of which some are potentially modified: */
977 fprintf(fp, " Coulomb type : %s\n", EELTYPE(ir->coulombtype));
978 fprintf(fp, " Grid spacing x y z : %f %f %f\n",
979 box_size[XX]/ir->nkx, box_size[YY]/ir->nky, box_size[ZZ]/ir->nkz);
980 fprintf(fp, " Van der Waals type : %s\n", EVDWTYPE(ir->vdwtype));
981 if (ir_vdw_switched(ir))
983 fprintf(fp, " rvdw_switch : %f nm\n", ir->rvdw_switch);
985 if (EPME_SWITCHED(ir->coulombtype))
987 fprintf(fp, " rlist : %f nm\n", ir->rlist);
990 /* Print a descriptive line about the tpr settings tested */
991 fprintf(fp, "\nWill try these real/reciprocal workload settings:\n");
992 fprintf(fp, " No. scaling rcoulomb");
993 fprintf(fp, " nkx nky nkz");
994 fprintf(fp, " spacing");
995 if (can_scale_rvdw(ir->vdwtype))
997 fprintf(fp, " rvdw");
999 if (EPME_SWITCHED(ir->coulombtype))
1001 fprintf(fp, " rlist");
1003 fprintf(fp, " tpr file\n");
1005 /* Loop to create the requested number of tpr input files */
1006 for (j = 0; j < *ntprs; j++)
1008 /* The first .tpr is the provided one, just need to modify nsteps,
1009 * so skip the following block */
1012 /* Determine which Coulomb radii rc to use in the benchmarks */
1013 add = (rmax-rmin)/(*ntprs-1);
1014 if (gmx_within_tol(rmin, info->rcoulomb[0], GMX_REAL_EPS))
1016 ir->rcoulomb = rmin + j*add;
1018 else if (gmx_within_tol(rmax, info->rcoulomb[0], GMX_REAL_EPS))
1020 ir->rcoulomb = rmin + (j-1)*add;
1024 /* rmin != rcoul != rmax, ergo test between rmin and rmax */
1025 add = (rmax-rmin)/(*ntprs-2);
1026 ir->rcoulomb = rmin + (j-1)*add;
1029 /* Determine the scaling factor fac */
1030 fac = ir->rcoulomb/info->rcoulomb[0];
1032 /* Scale the Fourier grid spacing */
1033 ir->nkx = ir->nky = ir->nkz = 0;
1034 calcFftGrid(nullptr, state.box, fourierspacing*fac, minimalPmeGridSize(ir->pme_order),
1035 &ir->nkx, &ir->nky, &ir->nkz);
1037 /* Adjust other radii since various conditions need to be fulfilled */
1038 if (eelPME == ir->coulombtype)
1040 /* plain PME, rcoulomb must be equal to rlist TODO only in the group scheme? */
1041 ir->rlist = ir->rcoulomb;
1045 /* rlist must be >= rcoulomb, we keep the size of the buffer region */
1046 ir->rlist = ir->rcoulomb + nlist_buffer;
1049 if (bScaleRvdw && can_scale_rvdw(ir->vdwtype))
1051 if (ecutsVERLET == ir->cutoff_scheme ||
1052 evdwPME == ir->vdwtype)
1054 /* With either the Verlet cutoff-scheme or LJ-PME,
1055 the van der Waals radius must always equal the
1057 ir->rvdw = ir->rcoulomb;
1061 /* For vdw cutoff, rvdw >= rlist */
1062 ir->rvdw = std::max(info->rvdw[0], ir->rlist);
1065 } /* end of "if (j != 0)" */
1067 /* for j==0: Save the original settings
1068 * for j >0: Save modified radii and Fourier grids */
1069 info->rcoulomb[j] = ir->rcoulomb;
1070 info->rvdw[j] = ir->rvdw;
1071 info->nkx[j] = ir->nkx;
1072 info->nky[j] = ir->nky;
1073 info->nkz[j] = ir->nkz;
1074 info->rlist[j] = ir->rlist;
1075 info->fsx[j] = fac*fourierspacing;
1076 info->fsy[j] = fac*fourierspacing;
1077 info->fsz[j] = fac*fourierspacing;
1079 /* Write the benchmark tpr file */
1080 fn_bench_tprs[j] = gmx_strdup(gmx::Path::concatenateBeforeExtension(fn_sim_tpr, gmx::formatString("_bench%.2d", j)).c_str());
1082 fprintf(stdout, "Writing benchmark tpr %s with nsteps=", fn_bench_tprs[j]);
1083 fprintf(stdout, "%" PRId64, ir->nsteps);
1086 fprintf(stdout, ", scaling factor %f\n", fac);
1090 fprintf(stdout, ", unmodified settings\n");
1093 write_tpx_state(fn_bench_tprs[j], ir, &state, &mtop);
1095 /* Write information about modified tpr settings to log file */
1096 fprintf(fp, "%4d%10f%10f", j, fac, ir->rcoulomb);
1097 fprintf(fp, "%5d%5d%5d", ir->nkx, ir->nky, ir->nkz);
1098 fprintf(fp, " %9f ", info->fsx[j]);
1099 if (can_scale_rvdw(ir->vdwtype))
1101 fprintf(fp, "%10f", ir->rvdw);
1103 if (EPME_SWITCHED(ir->coulombtype))
1105 fprintf(fp, "%10f", ir->rlist);
1107 fprintf(fp, " %-14s\n", fn_bench_tprs[j]);
1109 /* Make it clear to the user that some additional settings were modified */
1110 if (!gmx_within_tol(ir->rvdw, info->rvdw[0], GMX_REAL_EPS)
1111 || !gmx_within_tol(ir->rlist, info->rlist[0], GMX_REAL_EPS) )
1118 fprintf(fp, "\nNote that in addition to the Coulomb radius and the Fourier grid\n"
1119 "other input settings were also changed (see table above).\n"
1120 "Please check if the modified settings are appropriate.\n");
1127 /* Rename the files we want to keep to some meaningful filename and
1128 * delete the rest */
1129 static void cleanup(const t_filenm *fnm, int nfile, int k, int nnodes,
1130 int nPMEnodes, int nr, gmx_bool bKeepStderr)
1132 char numstring[STRLEN];
1133 const char *fn = nullptr;
1138 fprintf(stdout, "Cleaning up, deleting benchmark temp files ...\n");
1140 for (i = 0; i < nfile; i++)
1142 opt = const_cast<char *>(fnm[i].opt);
1143 if (std::strcmp(opt, "-p") == 0)
1145 /* do nothing; keep this file */
1148 else if (std::strcmp(opt, "-bg") == 0)
1150 /* Give the log file a nice name so one can later see which parameters were used */
1151 numstring[0] = '\0';
1154 sprintf(numstring, "_%d", nr);
1156 std::string newfilename = gmx::formatString("%s_no%d_np%d_npme%d%s", opt2fn("-bg", nfile, fnm), k, nnodes, nPMEnodes, numstring);
1157 if (gmx_fexist(opt2fn("-bg", nfile, fnm)))
1159 fprintf(stdout, "renaming log file to %s\n", newfilename.c_str());
1160 make_backup(newfilename);
1161 rename(opt2fn("-bg", nfile, fnm), newfilename.c_str());
1164 else if (std::strcmp(opt, "-err") == 0)
1166 /* This file contains the output of stderr. We want to keep it in
1167 * cases where there have been problems. */
1168 fn = opt2fn(opt, nfile, fnm);
1169 numstring[0] = '\0';
1172 sprintf(numstring, "_%d", nr);
1174 std::string newfilename = gmx::formatString("%s_no%d_np%d_npme%d%s", fn, k, nnodes, nPMEnodes, numstring);
1179 fprintf(stdout, "Saving stderr output in %s\n", newfilename.c_str());
1180 make_backup(newfilename);
1181 rename(fn, newfilename.c_str());
1185 fprintf(stdout, "Deleting %s\n", fn);
1190 /* Delete the files which are created for each benchmark run: (options -b*) */
1191 else if ( (0 == std::strncmp(opt, "-b", 2)) && (opt2bSet(opt, nfile, fnm) || !is_optional(&fnm[i])) )
1193 remove_if_exists(opt2fn(opt, nfile, fnm));
1200 eNpmeAuto, eNpmeAll, eNpmeReduced, eNpmeSubset, eNpmeNr
1203 /* Create a list of numbers of PME nodes to test */
1204 static void make_npme_list(
1205 const char *npmevalues_opt, /* Make a complete list with all
1206 * possibilities or a short list that keeps only
1207 * reasonable numbers of PME nodes */
1208 int *nentries, /* Number of entries we put in the nPMEnodes list */
1209 int *nPMEnodes[], /* Each entry contains the value for -npme */
1210 int nnodes, /* Total number of nodes to do the tests on */
1211 int minPMEnodes, /* Minimum number of PME nodes */
1212 int maxPMEnodes) /* Maximum number of PME nodes */
1215 int min_factor = 1; /* We request that npp and npme have this minimal
1216 * largest common factor (depends on npp) */
1217 int nlistmax; /* Max. list size */
1218 int nlist; /* Actual number of entries in list */
1222 /* Do we need to check all possible values for -npme or is a reduced list enough? */
1223 if (!std::strcmp(npmevalues_opt, "all") )
1227 else if (!std::strcmp(npmevalues_opt, "subset") )
1229 eNPME = eNpmeSubset;
1231 else /* "auto" or "range" */
1237 else if (nnodes < 128)
1239 eNPME = eNpmeReduced;
1243 eNPME = eNpmeSubset;
1247 /* Calculate how many entries we could possibly have (in case of -npme all) */
1250 nlistmax = maxPMEnodes - minPMEnodes + 3;
1251 if (0 == minPMEnodes)
1261 /* Now make the actual list which is at most of size nlist */
1262 snew(*nPMEnodes, nlistmax);
1263 nlist = 0; /* start counting again, now the real entries in the list */
1264 for (i = 0; i < nlistmax - 2; i++)
1266 npme = maxPMEnodes - i;
1277 /* For 2d PME we want a common largest factor of at least the cube
1278 * root of the number of PP nodes */
1279 min_factor = static_cast<int>(std::cbrt(npp));
1282 gmx_fatal(FARGS, "Unknown option for eNPME in make_npme_list");
1284 if (gmx_greatest_common_divisor(npp, npme) >= min_factor)
1286 (*nPMEnodes)[nlist] = npme;
1290 /* We always test 0 PME nodes and the automatic number */
1291 *nentries = nlist + 2;
1292 (*nPMEnodes)[nlist ] = 0;
1293 (*nPMEnodes)[nlist+1] = -1;
1295 fprintf(stderr, "Will try the following %d different values for -npme:\n", *nentries);
1296 for (i = 0; i < *nentries-1; i++)
1298 fprintf(stderr, "%d, ", (*nPMEnodes)[i]);
1300 fprintf(stderr, "and %d (auto).\n", (*nPMEnodes)[*nentries-1]);
1304 /* Allocate memory to store the performance data */
1305 static void init_perfdata(t_perf *perfdata[], int ntprs, int datasets, int repeats)
1310 for (k = 0; k < ntprs; k++)
1312 snew(perfdata[k], datasets);
1313 for (i = 0; i < datasets; i++)
1315 for (j = 0; j < repeats; j++)
1317 snew(perfdata[k][i].Gcycles, repeats);
1318 snew(perfdata[k][i].ns_per_day, repeats);
1319 snew(perfdata[k][i].PME_f_load, repeats);
1326 /* Check for errors on mdrun -h */
1327 static void make_sure_it_runs(char *mdrun_cmd_line, int length, FILE *fp,
1328 const t_filenm *fnm, int nfile)
1330 char *command, *msg;
1333 snew(command, length + 15);
1334 snew(msg, length + 500);
1336 fprintf(stdout, "Making sure the benchmarks can be executed by running just 1 step...\n");
1337 sprintf(command, "%s -nsteps 1 -quiet", mdrun_cmd_line);
1338 fprintf(stdout, "Executing '%s' ...\n", command);
1339 ret = gmx_system_call(command);
1343 /* To prevent confusion, do not again issue a gmx_fatal here since we already
1344 * get the error message from mdrun itself */
1346 "Cannot run the first benchmark simulation! Please check the error message of\n"
1347 "mdrun for the source of the problem. Did you provide a command line\n"
1348 "argument that neither gmx tune_pme nor mdrun understands? If you're\n"
1349 "sure your command line should work, you can bypass this check with \n"
1350 "gmx tune_pme -nocheck. The failing command was:\n"
1351 "\n%s\n\n", command);
1353 fprintf(stderr, "%s", msg);
1355 fprintf(fp, "%s", msg);
1359 fprintf(stdout, "Benchmarks can be executed!\n");
1361 /* Clean up the benchmark output files we just created */
1362 fprintf(stdout, "Cleaning up ...\n");
1363 remove_if_exists(opt2fn("-bc", nfile, fnm));
1364 remove_if_exists(opt2fn("-be", nfile, fnm));
1365 remove_if_exists(opt2fn("-bcpo", nfile, fnm));
1366 remove_if_exists(opt2fn("-bg", nfile, fnm));
1367 remove_if_exists(opt2fn("-bo", nfile, fnm));
1368 remove_if_exists(opt2fn("-bx", nfile, fnm));
1374 static void do_the_tests(
1375 FILE *fp, /* General tune_pme output file */
1376 char **tpr_names, /* Filenames of the input files to test */
1377 int maxPMEnodes, /* Max fraction of nodes to use for PME */
1378 int minPMEnodes, /* Min fraction of nodes to use for PME */
1379 int npme_fixed, /* If >= -1, test fixed number of PME
1381 const char *npmevalues_opt, /* Which -npme values should be tested */
1382 t_perf **perfdata, /* Here the performace data is stored */
1383 int *pmeentries, /* Entries in the nPMEnodes list */
1384 int repeats, /* Repeat each test this often */
1385 int nnodes, /* Total number of nodes = nPP + nPME */
1386 int nr_tprs, /* Total number of tpr files to test */
1387 gmx_bool bThreads, /* Threads or MPI? */
1388 char *cmd_mpirun, /* mpirun command string */
1389 char *cmd_np, /* "-np", "-n", whatever mpirun needs */
1390 char *cmd_mdrun, /* mdrun command string */
1391 char *cmd_args_bench, /* arguments for mdrun in a string */
1392 const t_filenm *fnm, /* List of filenames from command line */
1393 int nfile, /* Number of files specified on the cmdl. */
1394 int presteps, /* DLB equilibration steps, is checked */
1395 int64_t cpt_steps, /* Time step counter in the checkpoint */
1396 gmx_bool bCheck, /* Check whether benchmark mdrun works */
1397 const char *eligible_gpu_ids) /* GPU IDs for
1398 * constructing mdrun command lines */
1400 int i, nr, k, ret, count = 0, totaltests;
1401 int *nPMEnodes = nullptr;
1402 t_perf *pd = nullptr;
1404 char *command, *cmd_stub;
1406 gmx_bool bResetProblem = FALSE;
1407 gmx_bool bFirst = TRUE;
1409 /* This string array corresponds to the eParselog enum type at the start
1411 const char* ParseLog[] = {
1413 "Logfile not found!",
1414 "No timings, logfile truncated?",
1415 "Run was terminated.",
1416 "Counters were not reset properly.",
1417 "No DD grid found for these settings.",
1418 "TPX version conflict!",
1419 "mdrun was not started in parallel!",
1420 "Number of PP ranks has a prime factor that is too large.",
1421 "The number of PP ranks did not suit the number of GPUs.",
1422 "Some GPUs were not detected or are incompatible.",
1423 "An error occurred."
1425 char str_PME_f_load[13];
1428 /* Allocate space for the mdrun command line. 100 extra characters should
1429 be more than enough for the -npme etcetera arguments */
1430 cmdline_length = std::strlen(cmd_mpirun)
1431 + std::strlen(cmd_np)
1432 + std::strlen(cmd_mdrun)
1433 + std::strlen(cmd_args_bench)
1434 + std::strlen(tpr_names[0]) + 100;
1435 snew(command, cmdline_length);
1436 snew(cmd_stub, cmdline_length);
1438 /* Construct the part of the command line that stays the same for all tests: */
1441 sprintf(cmd_stub, "%s%s", cmd_mdrun, cmd_np);
1445 sprintf(cmd_stub, "%s%s%s ", cmd_mpirun, cmd_np, cmd_mdrun);
1448 /* Create a list of numbers of PME nodes to test */
1449 if (npme_fixed < -1)
1451 make_npme_list(npmevalues_opt, pmeentries, &nPMEnodes,
1452 nnodes, minPMEnodes, maxPMEnodes);
1458 nPMEnodes[0] = npme_fixed;
1459 fprintf(stderr, "Will use a fixed number of %d PME-only ranks.\n", nPMEnodes[0]);
1464 fprintf(fp, "\nNo benchmarks done since number of repeats (-r) is 0.\n");
1466 finalize(opt2fn("-p", nfile, fnm));
1470 /* Allocate one dataset for each tpr input file: */
1471 init_perfdata(perfdata, nr_tprs, *pmeentries, repeats);
1473 /*****************************************/
1474 /* Main loop over all tpr files to test: */
1475 /*****************************************/
1476 totaltests = nr_tprs*(*pmeentries)*repeats;
1477 for (k = 0; k < nr_tprs; k++)
1479 fprintf(fp, "\nIndividual timings for input file %d (%s):\n", k, tpr_names[k]);
1480 fprintf(fp, "PME ranks Gcycles ns/day PME/f Remark\n");
1481 /* Loop over various numbers of PME nodes: */
1482 for (i = 0; i < *pmeentries; i++)
1484 pd = &perfdata[k][i];
1486 auto cmd_gpu_ids = make_gpu_id_command_line(eligible_gpu_ids);
1488 /* Loop over the repeats for each scenario: */
1489 for (nr = 0; nr < repeats; nr++)
1491 pd->nPMEnodes = nPMEnodes[i];
1493 /* Add -npme and -s to the command line and save it. Note that
1494 * the -passall (if set) options requires cmd_args_bench to be
1495 * at the end of the command line string */
1496 snew(pd->mdrun_cmd_line, cmdline_length);
1497 sprintf(pd->mdrun_cmd_line, "%s-npme %d -s %s %s %s",
1498 cmd_stub, pd->nPMEnodes, tpr_names[k], cmd_args_bench, cmd_gpu_ids.c_str());
1500 /* To prevent that all benchmarks fail due to a show-stopper argument
1501 * on the mdrun command line, we make a quick check first.
1502 * This check can be turned off in cases where the automatically chosen
1503 * number of PME-only ranks leads to a number of PP ranks for which no
1504 * decomposition can be found (e.g. for large prime numbers) */
1505 if (bFirst && bCheck)
1507 /* TODO This check is really for a functional
1508 * .tpr, and if we need it, it should take place
1509 * for every .tpr, and the logic for it should be
1510 * immediately inside the loop over k, not in
1511 * this inner loop. */
1512 char *temporary_cmd_line;
1514 snew(temporary_cmd_line, cmdline_length);
1515 /* TODO -npme 0 is more likely to succeed at low
1516 parallelism than the default of -npme -1, but
1517 is more likely to fail at the scaling limit
1518 when the PP domains may be too small. "mpirun
1519 -np 1 mdrun" is probably a reasonable thing to
1520 do for this check, but it'll be easier to
1521 implement that after some refactoring of how
1522 the number of MPI ranks is managed. */
1523 sprintf(temporary_cmd_line, "%s-npme 0 -nb cpu -s %s %s",
1524 cmd_stub, tpr_names[k], cmd_args_bench);
1525 make_sure_it_runs(temporary_cmd_line, cmdline_length, fp, fnm, nfile);
1529 /* Do a benchmark simulation: */
1532 sprintf(buf, ", pass %d/%d", nr+1, repeats);
1538 fprintf(stdout, "\n=== Progress %2.0f%%, tpr %d/%d, run %d/%d%s:\n",
1539 (100.0*count)/totaltests,
1540 k+1, nr_tprs, i+1, *pmeentries, buf);
1541 make_backup(opt2fn("-err", nfile, fnm));
1542 sprintf(command, "%s 1> /dev/null 2>%s", pd->mdrun_cmd_line, opt2fn("-err", nfile, fnm));
1543 fprintf(stdout, "%s\n", pd->mdrun_cmd_line);
1544 gmx_system_call(command);
1546 /* Collect the performance data from the log file; also check stderr
1547 * for fatal errors */
1548 ret = parse_logfile(opt2fn("-bg", nfile, fnm), opt2fn("-err", nfile, fnm),
1549 pd, nr, presteps, cpt_steps, nnodes);
1550 if ((presteps > 0) && (ret == eParselogResetProblem))
1552 bResetProblem = TRUE;
1555 if (-1 == pd->nPMEnodes)
1557 sprintf(buf, "(%3d)", pd->guessPME);
1565 if (pd->PME_f_load[nr] > 0.0)
1567 sprintf(str_PME_f_load, "%12.3f", pd->PME_f_load[nr]);
1571 sprintf(str_PME_f_load, "%s", " - ");
1574 /* Write the data we got to disk */
1575 fprintf(fp, "%4d%s %12.3f %12.3f %s %s", pd->nPMEnodes,
1576 buf, pd->Gcycles[nr], pd->ns_per_day[nr], str_PME_f_load, ParseLog[ret]);
1577 if (!(ret == eParselogOK || ret == eParselogNoDDGrid || ret == eParselogNotFound) )
1579 fprintf(fp, " Check %s file for problems.", ret == eParselogFatal ? "err" : "log");
1585 /* Do some cleaning up and delete the files we do not need any more */
1586 cleanup(fnm, nfile, k, nnodes, pd->nPMEnodes, nr, ret == eParselogFatal);
1588 /* If the first run with this number of processors already failed, do not try again: */
1589 if (pd->Gcycles[0] <= 0.0 && repeats > 1)
1591 fprintf(stdout, "Skipping remaining passes of unsuccessful setting, see log file for details.\n");
1592 count += repeats-(nr+1);
1595 } /* end of repeats loop */
1596 } /* end of -npme loop */
1597 } /* end of tpr file loop */
1602 fprintf(fp, "WARNING: The cycle and time step counters could not be reset properly. ");
1610 static void check_input(
1617 real maxPMEfraction,
1618 real minPMEfraction,
1620 int64_t bench_nsteps,
1621 const t_filenm *fnm,
1631 /* Make sure the input file exists */
1632 if (!gmx_fexist(opt2fn("-s", nfile, fnm)))
1634 gmx_fatal(FARGS, "File %s not found.", opt2fn("-s", nfile, fnm));
1637 /* Make sure that the checkpoint file is not overwritten during benchmarking */
1638 if ( (0 == std::strcmp(opt2fn("-cpi", nfile, fnm), opt2fn("-bcpo", nfile, fnm)) ) && (sim_part > 1) )
1640 gmx_fatal(FARGS, "Checkpoint input (-cpi) and benchmark checkpoint output (-bcpo) files must not be identical.\n"
1641 "The checkpoint input file must not be overwritten during the benchmarks.\n");
1644 /* Make sure that repeats is >= 0 (if == 0, only write tpr files) */
1647 gmx_fatal(FARGS, "Number of repeats < 0!");
1650 /* Check number of nodes */
1653 gmx_fatal(FARGS, "Number of ranks/threads must be a positive integer.");
1656 /* Automatically choose -ntpr if not set */
1666 /* Set a reasonable scaling factor for rcoulomb */
1669 *rmax = rcoulomb * 1.2;
1672 fprintf(stderr, "Will test %d tpr file%s.\n", *ntprs, *ntprs == 1 ? "" : "s");
1678 fprintf(stderr, "Note: Choose ntpr>1 to shift PME load between real and reciprocal space.\n");
1682 /* Make shure that rmin <= rcoulomb <= rmax */
1691 if (!(*rmin <= *rmax) )
1693 gmx_fatal(FARGS, "Please choose the Coulomb radii such that rmin <= rmax.\n"
1694 "rmin = %g, rmax = %g, actual rcoul from .tpr file = %g\n", *rmin, *rmax, rcoulomb);
1696 /* Add test scenarios if rmin or rmax were set */
1699 if (!gmx_within_tol(*rmin, rcoulomb, GMX_REAL_EPS) && (*ntprs == 1) )
1702 fprintf(stderr, "NOTE: Setting -rmin to %g changed -ntpr to %d\n",
1705 if (!gmx_within_tol(*rmax, rcoulomb, GMX_REAL_EPS) && (*ntprs == 1) )
1708 fprintf(stderr, "NOTE: Setting -rmax to %g changed -ntpr to %d\n",
1713 /* If one of rmin, rmax is set, we need 2 tpr files at minimum */
1714 if (!gmx_within_tol(*rmax, rcoulomb, GMX_REAL_EPS) || !gmx_within_tol(*rmin, rcoulomb, GMX_REAL_EPS) )
1716 *ntprs = std::max(*ntprs, 2);
1719 /* If both rmin, rmax are set, we need 3 tpr files at minimum */
1720 if (!gmx_within_tol(*rmax, rcoulomb, GMX_REAL_EPS) && !gmx_within_tol(*rmin, rcoulomb, GMX_REAL_EPS) )
1722 *ntprs = std::max(*ntprs, 3);
1727 fprintf(stderr, "NOTE: Your rmin, rmax setting changed -ntpr to %d\n", *ntprs);
1732 if (gmx_within_tol(*rmin, rcoulomb, GMX_REAL_EPS) && gmx_within_tol(rcoulomb, *rmax, GMX_REAL_EPS)) /* We have just a single rc */
1734 fprintf(stderr, "WARNING: Resetting -ntpr to 1 since no Coulomb radius scaling is requested.\n"
1735 "Please set rmin < rmax to test Coulomb radii in the [rmin, rmax] interval\n"
1736 "with correspondingly adjusted PME grid settings\n");
1741 /* Check whether max and min fraction are within required values */
1742 if (maxPMEfraction > 0.5 || maxPMEfraction < 0)
1744 gmx_fatal(FARGS, "-max must be between 0 and 0.5");
1746 if (minPMEfraction > 0.5 || minPMEfraction < 0)
1748 gmx_fatal(FARGS, "-min must be between 0 and 0.5");
1750 if (maxPMEfraction < minPMEfraction)
1752 gmx_fatal(FARGS, "-max must be larger or equal to -min");
1755 /* Check whether the number of steps - if it was set - has a reasonable value */
1756 if (bench_nsteps < 0)
1758 gmx_fatal(FARGS, "Number of steps must be positive.");
1761 if (bench_nsteps > 10000 || bench_nsteps < 100)
1763 fprintf(stderr, "WARNING: steps=");
1764 fprintf(stderr, "%" PRId64, bench_nsteps);
1765 fprintf(stderr, ". Are you sure you want to perform so %s steps for each benchmark?\n", (bench_nsteps < 100) ? "few" : "many");
1770 gmx_fatal(FARGS, "Cannot have a negative number of presteps.\n");
1773 /* Check for rcoulomb scaling if more than one .tpr file is tested */
1776 if (*rmin/rcoulomb < 0.75 || *rmax/rcoulomb > 1.25)
1778 fprintf(stderr, "WARNING: Applying extreme scaling factor. I hope you know what you are doing.\n");
1782 /* If a fixed number of PME nodes is set we do rcoulomb and PME gird tuning
1783 * only. We need to check whether the requested number of PME-only nodes
1785 if (npme_fixed > -1)
1787 /* No more than 50% of all nodes can be assigned as PME-only nodes. */
1788 if (2*npme_fixed > nnodes)
1790 gmx_fatal(FARGS, "Cannot have more than %d PME-only ranks for a total of %d ranks (you chose %d).\n",
1791 nnodes/2, nnodes, npme_fixed);
1793 if ((npme_fixed > 0) && (5*npme_fixed < nnodes))
1795 fprintf(stderr, "WARNING: Only %g percent of the ranks are assigned as PME-only ranks.\n",
1796 (100.0*npme_fixed)/nnodes);
1798 if (opt2parg_bSet("-min", npargs, pa) || opt2parg_bSet("-max", npargs, pa))
1800 fprintf(stderr, "NOTE: The -min, -max, and -npme options have no effect when a\n"
1801 " fixed number of PME-only ranks is requested with -fix.\n");
1807 /* Returns TRUE when "opt" is needed at launch time */
1808 static gmx_bool is_launch_file(char *opt, gmx_bool bSet)
1810 if (0 == std::strncmp(opt, "-swap", 5))
1815 /* Apart from the input .tpr and the output log files we need all options that
1816 * were set on the command line and that do not start with -b */
1817 if (0 == std::strncmp(opt, "-b", 2) || 0 == std::strncmp(opt, "-s", 2)
1818 || 0 == std::strncmp(opt, "-err", 4) || 0 == std::strncmp(opt, "-p", 2) )
1827 /* Returns TRUE when "opt" defines a file which is needed for the benchmarks runs */
1828 static gmx_bool is_bench_file(char *opt, gmx_bool bSet, gmx_bool bOptional, gmx_bool bIsOutput)
1830 /* Apart from the input .tpr, all files starting with "-b" are for
1831 * _b_enchmark files exclusively */
1832 if (0 == std::strncmp(opt, "-s", 2))
1837 if (0 == std::strncmp(opt, "-b", 2) || 0 == std::strncmp(opt, "-s", 2))
1839 return !bOptional || bSet;
1849 return bSet; /* These are additional input files like -cpi -ei */
1855 /* Adds 'buf' to 'str' */
1856 static void add_to_string(char **str, const char *buf)
1861 len = std::strlen(*str) + std::strlen(buf) + 1;
1863 std::strcat(*str, buf);
1867 /* Create the command line for the benchmark as well as for the real run */
1868 static void create_command_line_snippets(
1869 gmx_bool bAppendFiles,
1870 gmx_bool bKeepAndNumCPT,
1871 gmx_bool bResetHWay,
1875 char *cmd_args_bench[], /* command line arguments for benchmark runs */
1876 char *cmd_args_launch[], /* command line arguments for simulation run */
1877 char extra_args[], /* Add this to the end of the command line */
1878 char *deffnm) /* Default file names, or NULL if not set */
1883 char strbuf[STRLEN];
1886 /* strlen needs at least '\0' as a string: */
1887 snew(*cmd_args_bench, 1);
1888 snew(*cmd_args_launch, 1);
1889 *cmd_args_launch[0] = '\0';
1890 *cmd_args_bench[0] = '\0';
1893 /*******************************************/
1894 /* 1. Process other command line arguments */
1895 /*******************************************/
1898 /* Add equilibration steps to benchmark options */
1899 sprintf(strbuf, "-resetstep %d ", presteps);
1900 add_to_string(cmd_args_bench, strbuf);
1902 /* These switches take effect only at launch time */
1905 sprintf(strbuf, "-deffnm %s ", deffnm);
1906 add_to_string(cmd_args_launch, strbuf);
1910 add_to_string(cmd_args_launch, "-noappend ");
1914 add_to_string(cmd_args_launch, "-cpnum ");
1918 add_to_string(cmd_args_launch, "-resethway ");
1921 /********************/
1922 /* 2. Process files */
1923 /********************/
1924 for (i = 0; i < nfile; i++)
1926 opt = const_cast<char *>(fnm[i].opt);
1927 name = opt2fn(opt, nfile, fnm);
1929 /* Strbuf contains the options, now let's sort out where we need that */
1930 sprintf(strbuf, "%s %s ", opt, name);
1932 if (is_bench_file(opt, opt2bSet(opt, nfile, fnm), is_optional(&fnm[i]), is_output(&fnm[i])) )
1934 /* All options starting with -b* need the 'b' removed,
1935 * therefore overwrite strbuf */
1936 if (0 == std::strncmp(opt, "-b", 2))
1938 sprintf(strbuf, "-%s %s ", &opt[2], name);
1941 add_to_string(cmd_args_bench, strbuf);
1944 if (is_launch_file(opt, opt2bSet(opt, nfile, fnm)) )
1946 add_to_string(cmd_args_launch, strbuf);
1950 add_to_string(cmd_args_bench, extra_args);
1951 add_to_string(cmd_args_launch, extra_args);
1955 /* Set option opt */
1956 static void setopt(const char *opt, int nfile, t_filenm fnm[])
1960 for (i = 0; (i < nfile); i++)
1962 if (std::strcmp(opt, fnm[i].opt) == 0)
1964 fnm[i].flag |= ffSET;
1970 /* This routine inspects the tpr file and ...
1971 * 1. checks for output files that get triggered by a tpr option. These output
1972 * files are marked as 'set' to allow for a proper cleanup after each
1974 * 2. returns the PME:PP load ratio
1975 * 3. returns rcoulomb from the tpr */
1976 static float inspect_tpr(int nfile, t_filenm fnm[], real *rcoulomb)
1978 gmx_bool bTpi; /* Is test particle insertion requested? */
1979 gmx_bool bFree; /* Is a free energy simulation requested? */
1980 gmx_bool bNM; /* Is a normal mode analysis requested? */
1981 gmx_bool bSwap; /* Is water/ion position swapping requested? */
1986 /* Check tpr file for options that trigger extra output files */
1987 t_inputrec irInstance;
1988 t_inputrec *ir = &irInstance;
1989 read_tpx_state(opt2fn("-s", nfile, fnm), ir, &state, &mtop);
1990 bFree = (efepNO != ir->efep );
1991 bNM = (eiNM == ir->eI );
1992 bSwap = (eswapNO != ir->eSwapCoords);
1993 bTpi = EI_TPI(ir->eI);
1995 /* Set these output files on the tuning command-line */
1998 setopt("-pf", nfile, fnm);
1999 setopt("-px", nfile, fnm);
2003 setopt("-dhdl", nfile, fnm);
2007 setopt("-tpi", nfile, fnm);
2008 setopt("-tpid", nfile, fnm);
2012 setopt("-mtx", nfile, fnm);
2016 setopt("-swap", nfile, fnm);
2019 *rcoulomb = ir->rcoulomb;
2021 /* Return the estimate for the number of PME nodes */
2022 float npme = pme_load_estimate(mtop, *ir, state.box);
2027 static void couple_files_options(int nfile, t_filenm fnm[])
2030 gmx_bool bSet, bBench;
2035 for (i = 0; i < nfile; i++)
2037 opt = const_cast<char *>(fnm[i].opt);
2038 bSet = ((fnm[i].flag & ffSET) != 0);
2039 bBench = (0 == std::strncmp(opt, "-b", 2));
2041 /* Check optional files */
2042 /* If e.g. -eo is set, then -beo also needs to be set */
2043 if (is_optional(&fnm[i]) && bSet && !bBench)
2045 sprintf(buf, "-b%s", &opt[1]);
2046 setopt(buf, nfile, fnm);
2048 /* If -beo is set, then -eo also needs to be! */
2049 if (is_optional(&fnm[i]) && bSet && bBench)
2051 sprintf(buf, "-%s", &opt[2]);
2052 setopt(buf, nfile, fnm);
2058 #define BENCHSTEPS (1000)
2060 int gmx_tune_pme(int argc, char *argv[])
2062 const char *desc[] = {
2063 "For a given number [TT]-np[tt] or [TT]-ntmpi[tt] of ranks, [THISMODULE] systematically",
2064 "times [gmx-mdrun] with various numbers of PME-only ranks and determines",
2065 "which setting is fastest. It will also test whether performance can",
2066 "be enhanced by shifting load from the reciprocal to the real space",
2067 "part of the Ewald sum. ",
2068 "Simply pass your [REF].tpr[ref] file to [THISMODULE] together with other options",
2069 "for [gmx-mdrun] as needed.[PAR]",
2070 "[THISMODULE] needs to call [gmx-mdrun] and so requires that you",
2071 "specify how to call mdrun with the argument to the [TT]-mdrun[tt]",
2072 "parameter. Depending how you have built GROMACS, values such as",
2073 "'gmx mdrun', 'gmx_d mdrun', or 'mdrun_mpi' might be needed.[PAR]",
2074 "The program that runs MPI programs can be set in the environment variable",
2075 "MPIRUN (defaults to 'mpirun'). Note that for certain MPI frameworks,",
2076 "you need to provide a machine- or hostfile. This can also be passed",
2077 "via the MPIRUN variable, e.g.[PAR]",
2078 "[TT]export MPIRUN=\"/usr/local/mpirun -machinefile hosts\"[tt]",
2079 "Note that in such cases it is normally necessary to compile",
2080 "and/or run [THISMODULE] without MPI support, so that it can call",
2081 "the MPIRUN program.[PAR]",
2082 "Before doing the actual benchmark runs, [THISMODULE] will do a quick",
2083 "check whether [gmx-mdrun] works as expected with the provided parallel settings",
2084 "if the [TT]-check[tt] option is activated (the default).",
2085 "Please call [THISMODULE] with the normal options you would pass to",
2086 "[gmx-mdrun] and add [TT]-np[tt] for the number of ranks to perform the",
2087 "tests on, or [TT]-ntmpi[tt] for the number of threads. You can also add [TT]-r[tt]",
2088 "to repeat each test several times to get better statistics. [PAR]",
2089 "[THISMODULE] can test various real space / reciprocal space workloads",
2090 "for you. With [TT]-ntpr[tt] you control how many extra [REF].tpr[ref] files will be",
2091 "written with enlarged cutoffs and smaller Fourier grids respectively.",
2092 "Typically, the first test (number 0) will be with the settings from the input",
2093 "[REF].tpr[ref] file; the last test (number [TT]ntpr[tt]) will have the Coulomb cutoff",
2094 "specified by [TT]-rmax[tt] with a somewhat smaller PME grid at the same time. ",
2095 "In this last test, the Fourier spacing is multiplied with [TT]rmax[tt]/rcoulomb. ",
2096 "The remaining [REF].tpr[ref] files will have equally-spaced Coulomb radii (and Fourier ",
2097 "spacings) between these extremes. [BB]Note[bb] that you can set [TT]-ntpr[tt] to 1",
2098 "if you just seek the optimal number of PME-only ranks; in that case",
2099 "your input [REF].tpr[ref] file will remain unchanged.[PAR]",
2100 "For the benchmark runs, the default of 1000 time steps should suffice for most",
2101 "MD systems. The dynamic load balancing needs about 100 time steps",
2102 "to adapt to local load imbalances, therefore the time step counters",
2103 "are by default reset after 100 steps. For large systems (>1M atoms), as well as ",
2104 "for a higher accuracy of the measurements, you should set [TT]-resetstep[tt] to a higher value.",
2105 "From the 'DD' load imbalance entries in the md.log output file you",
2106 "can tell after how many steps the load is sufficiently balanced. Example call:[PAR]",
2107 "[TT]gmx tune_pme -np 64 -s protein.tpr -launch[tt][PAR]",
2108 "After calling [gmx-mdrun] several times, detailed performance information",
2109 "is available in the output file [TT]perf.out[tt].",
2110 "[BB]Note[bb] that during the benchmarks, a couple of temporary files are written",
2111 "(options [TT]-b*[tt]), these will be automatically deleted after each test.[PAR]",
2112 "If you want the simulation to be started automatically with the",
2113 "optimized parameters, use the command line option [TT]-launch[tt].[PAR]",
2114 "Basic support for GPU-enabled [TT]mdrun[tt] exists. Give a string containing the IDs",
2115 "of the GPUs that you wish to use in the optimization in the [TT]-gpu_id[tt]",
2116 "command-line argument. This works exactly like [TT]mdrun -gpu_id[tt], does not imply a mapping,",
2117 "and merely declares the eligible set of GPU devices. [TT]gmx-tune_pme[tt] will construct calls to",
2118 "mdrun that use this set appropriately. [TT]gmx-tune_pme[tt] does not support",
2119 "[TT]-gputasks[tt].[PAR]",
2124 int pmeentries = 0; /* How many values for -npme do we actually test for each tpr file */
2125 real maxPMEfraction = 0.50;
2126 real minPMEfraction = 0.25;
2127 int maxPMEnodes, minPMEnodes;
2128 float guessPMEratio; /* guessed PME:PP ratio based on the tpr file */
2129 float guessPMEnodes;
2130 int npme_fixed = -2; /* If >= -1, use only this number
2131 * of PME-only nodes */
2133 real rmin = 0.0, rmax = 0.0; /* min and max value for rcoulomb if scaling is requested */
2134 real rcoulomb = -1.0; /* Coulomb radius as set in .tpr file */
2135 gmx_bool bScaleRvdw = TRUE;
2136 int64_t bench_nsteps = BENCHSTEPS;
2137 int64_t new_sim_nsteps = -1; /* -1 indicates: not set by the user */
2138 int64_t cpt_steps = 0; /* Step counter in .cpt input file */
2139 int presteps = 1500; /* Do a full cycle reset after presteps steps */
2140 gmx_bool bOverwrite = FALSE, bKeepTPR;
2141 gmx_bool bLaunch = FALSE;
2142 char *ExtraArgs = nullptr;
2143 char **tpr_names = nullptr;
2144 const char *simulation_tpr = nullptr;
2145 char *deffnm = nullptr;
2146 int best_npme, best_tpr;
2147 int sim_part = 1; /* For benchmarks with checkpoint files */
2151 /* Default program names if nothing else is found */
2152 char *cmd_mpirun = nullptr, *cmd_mdrun = nullptr;
2153 char *cmd_args_bench, *cmd_args_launch;
2154 char *cmd_np = nullptr;
2156 /* IDs of GPUs that are eligible for computation */
2157 char *eligible_gpu_ids = nullptr;
2159 t_perf **perfdata = nullptr;
2164 /* Print out how long the tuning took */
2167 static t_filenm fnm[] = {
2169 { efOUT, "-p", "perf", ffWRITE },
2170 { efLOG, "-err", "bencherr", ffWRITE },
2171 { efTPR, "-so", "tuned", ffWRITE },
2173 { efTPR, "-s", nullptr, ffREAD },
2174 { efTRN, "-o", nullptr, ffWRITE },
2175 { efCOMPRESSED, "-x", nullptr, ffOPTWR },
2176 { efCPT, "-cpi", nullptr, ffOPTRD },
2177 { efCPT, "-cpo", nullptr, ffOPTWR },
2178 { efSTO, "-c", "confout", ffWRITE },
2179 { efEDR, "-e", "ener", ffWRITE },
2180 { efLOG, "-g", "md", ffWRITE },
2181 { efXVG, "-dhdl", "dhdl", ffOPTWR },
2182 { efXVG, "-field", "field", ffOPTWR },
2183 { efXVG, "-table", "table", ffOPTRD },
2184 { efXVG, "-tablep", "tablep", ffOPTRD },
2185 { efXVG, "-tableb", "table", ffOPTRD },
2186 { efTRX, "-rerun", "rerun", ffOPTRD },
2187 { efXVG, "-tpi", "tpi", ffOPTWR },
2188 { efXVG, "-tpid", "tpidist", ffOPTWR },
2189 { efEDI, "-ei", "sam", ffOPTRD },
2190 { efXVG, "-eo", "edsam", ffOPTWR },
2191 { efXVG, "-px", "pullx", ffOPTWR },
2192 { efXVG, "-pf", "pullf", ffOPTWR },
2193 { efXVG, "-ro", "rotation", ffOPTWR },
2194 { efLOG, "-ra", "rotangles", ffOPTWR },
2195 { efLOG, "-rs", "rotslabs", ffOPTWR },
2196 { efLOG, "-rt", "rottorque", ffOPTWR },
2197 { efMTX, "-mtx", "nm", ffOPTWR },
2198 { efXVG, "-swap", "swapions", ffOPTWR },
2199 /* Output files that are deleted after each benchmark run */
2200 { efTRN, "-bo", "bench", ffWRITE },
2201 { efXTC, "-bx", "bench", ffWRITE },
2202 { efCPT, "-bcpo", "bench", ffWRITE },
2203 { efSTO, "-bc", "bench", ffWRITE },
2204 { efEDR, "-be", "bench", ffWRITE },
2205 { efLOG, "-bg", "bench", ffWRITE },
2206 { efXVG, "-beo", "benchedo", ffOPTWR },
2207 { efXVG, "-bdhdl", "benchdhdl", ffOPTWR },
2208 { efXVG, "-bfield", "benchfld", ffOPTWR },
2209 { efXVG, "-btpi", "benchtpi", ffOPTWR },
2210 { efXVG, "-btpid", "benchtpid", ffOPTWR },
2211 { efXVG, "-bdevout", "benchdev", ffOPTWR },
2212 { efXVG, "-brunav", "benchrnav", ffOPTWR },
2213 { efXVG, "-bpx", "benchpx", ffOPTWR },
2214 { efXVG, "-bpf", "benchpf", ffOPTWR },
2215 { efXVG, "-bro", "benchrot", ffOPTWR },
2216 { efLOG, "-bra", "benchrota", ffOPTWR },
2217 { efLOG, "-brs", "benchrots", ffOPTWR },
2218 { efLOG, "-brt", "benchrott", ffOPTWR },
2219 { efMTX, "-bmtx", "benchn", ffOPTWR },
2220 { efNDX, "-bdn", "bench", ffOPTWR },
2221 { efXVG, "-bswap", "benchswp", ffOPTWR }
2224 gmx_bool bThreads = FALSE;
2228 const char *procstring[] =
2229 { nullptr, "np", "n", "none", nullptr };
2230 const char *npmevalues_opt[] =
2231 { nullptr, "auto", "all", "subset", nullptr };
2233 gmx_bool bAppendFiles = TRUE;
2234 gmx_bool bKeepAndNumCPT = FALSE;
2235 gmx_bool bResetCountersHalfWay = FALSE;
2236 gmx_bool bBenchmark = TRUE;
2237 gmx_bool bCheck = TRUE;
2239 gmx_output_env_t *oenv = nullptr;
2242 /***********************/
2243 /* tune_pme options: */
2244 /***********************/
2245 { "-mdrun", FALSE, etSTR, {&cmd_mdrun},
2246 "Command line to run a simulation, e.g. 'gmx mdrun' or 'mdrun_mpi'" },
2247 { "-np", FALSE, etINT, {&nnodes},
2248 "Number of ranks to run the tests on (must be > 2 for separate PME ranks)" },
2249 { "-npstring", FALSE, etENUM, {procstring},
2250 "Name of the [TT]$MPIRUN[tt] option that specifies the number of ranks to use ('np', or 'n'; use 'none' if there is no such option)" },
2251 { "-ntmpi", FALSE, etINT, {&nthreads},
2252 "Number of MPI-threads to run the tests on (turns MPI & mpirun off)"},
2253 { "-r", FALSE, etINT, {&repeats},
2254 "Repeat each test this often" },
2255 { "-max", FALSE, etREAL, {&maxPMEfraction},
2256 "Max fraction of PME ranks to test with" },
2257 { "-min", FALSE, etREAL, {&minPMEfraction},
2258 "Min fraction of PME ranks to test with" },
2259 { "-npme", FALSE, etENUM, {npmevalues_opt},
2260 "Within -min and -max, benchmark all possible values for [TT]-npme[tt], or just a reasonable subset. "
2261 "Auto neglects -min and -max and chooses reasonable values around a guess for npme derived from the .tpr"},
2262 { "-fix", FALSE, etINT, {&npme_fixed},
2263 "If >= -1, do not vary the number of PME-only ranks, instead use this fixed value and only vary rcoulomb and the PME grid spacing."},
2264 { "-rmax", FALSE, etREAL, {&rmax},
2265 "If >0, maximal rcoulomb for -ntpr>1 (rcoulomb upscaling results in fourier grid downscaling)" },
2266 { "-rmin", FALSE, etREAL, {&rmin},
2267 "If >0, minimal rcoulomb for -ntpr>1" },
2268 { "-scalevdw", FALSE, etBOOL, {&bScaleRvdw},
2269 "Scale rvdw along with rcoulomb"},
2270 { "-ntpr", FALSE, etINT, {&ntprs},
2271 "Number of [REF].tpr[ref] files to benchmark. Create this many files with different rcoulomb scaling factors depending on -rmin and -rmax. "
2272 "If < 1, automatically choose the number of [REF].tpr[ref] files to test" },
2273 { "-steps", FALSE, etINT64, {&bench_nsteps},
2274 "Take timings for this many steps in the benchmark runs" },
2275 { "-resetstep", FALSE, etINT, {&presteps},
2276 "Let dlb equilibrate this many steps before timings are taken (reset cycle counters after this many steps)" },
2277 { "-nsteps", FALSE, etINT64, {&new_sim_nsteps},
2278 "If non-negative, perform this many steps in the real run (overwrites nsteps from [REF].tpr[ref], add [REF].cpt[ref] steps)" },
2279 { "-launch", FALSE, etBOOL, {&bLaunch},
2280 "Launch the real simulation after optimization" },
2281 { "-bench", FALSE, etBOOL, {&bBenchmark},
2282 "Run the benchmarks or just create the input [REF].tpr[ref] files?" },
2283 { "-check", FALSE, etBOOL, {&bCheck},
2284 "Before the benchmark runs, check whether mdrun works in parallel" },
2285 { "-gpu_id", FALSE, etSTR, {&eligible_gpu_ids},
2286 "List of unique GPU device IDs that are eligible for use" },
2287 /******************/
2288 /* mdrun options: */
2289 /******************/
2290 /* We let tune_pme parse and understand these options, because we need to
2291 * prevent that they appear on the mdrun command line for the benchmarks */
2292 { "-append", FALSE, etBOOL, {&bAppendFiles},
2293 "Append to previous output files when continuing from checkpoint instead of adding the simulation part number to all file names (for launch only)" },
2294 { "-cpnum", FALSE, etBOOL, {&bKeepAndNumCPT},
2295 "Keep and number checkpoint files (launch only)" },
2296 { "-deffnm", FALSE, etSTR, {&deffnm},
2297 "Set the default filenames (launch only)" },
2298 { "-resethway", FALSE, etBOOL, {&bResetCountersHalfWay},
2299 "HIDDENReset the cycle counters after half the number of steps or halfway [TT]-maxh[tt] (launch only)" }
2302 #define NFILE asize(fnm)
2304 seconds = gmx_gettime();
2306 if (!parse_common_args(&argc, argv, PCA_NOEXIT_ON_ARGS,
2307 NFILE, fnm, asize(pa), pa, asize(desc), desc,
2313 // procstring[0]Â is used inside two different conditionals further down
2314 GMX_RELEASE_ASSERT(procstring[0] != nullptr, "Options inconsistency; procstring[0]Â is NULL");
2316 /* Store the remaining unparsed command line entries in a string which
2317 * is then attached to the mdrun command line */
2319 ExtraArgs[0] = '\0';
2320 for (i = 1; i < argc; i++) /* argc will now be 1 if everything was understood */
2322 add_to_string(&ExtraArgs, argv[i]);
2323 add_to_string(&ExtraArgs, " ");
2326 if (opt2parg_bSet("-ntmpi", asize(pa), pa))
2329 if (opt2parg_bSet("-npstring", asize(pa), pa))
2331 fprintf(stderr, "WARNING: -npstring has no effect when using threads.\n");
2336 gmx_fatal(FARGS, "Can't run multi-threaded MPI simulation yet!");
2338 /* and now we just set this; a bit of an ugly hack*/
2341 /* Check for PME:PP ratio and whether tpr triggers additional output files */
2342 guessPMEratio = inspect_tpr(NFILE, fnm, &rcoulomb);
2344 /* Automatically set -beo options if -eo is set etc. */
2345 couple_files_options(NFILE, fnm);
2347 /* Construct the command line arguments for benchmark runs
2348 * as well as for the simulation run */
2351 sprintf(bbuf, " -ntmpi %d ", nthreads);
2355 /* This string will be used for MPI runs and will appear after the
2356 * mpirun command. */
2357 if (std::strcmp(procstring[0], "none") != 0)
2359 sprintf(bbuf, " -%s %d ", procstring[0], nnodes);
2369 create_command_line_snippets(bAppendFiles, bKeepAndNumCPT, bResetCountersHalfWay, presteps,
2370 NFILE, fnm, &cmd_args_bench, &cmd_args_launch, ExtraArgs, deffnm);
2372 /* Prepare to use checkpoint file if requested */
2374 if (opt2bSet("-cpi", NFILE, fnm))
2376 const char *filename = opt2fn("-cpi", NFILE, fnm);
2378 read_checkpoint_part_and_step(filename,
2379 &cpt_sim_part, &cpt_steps);
2380 if (cpt_sim_part == 0)
2382 gmx_fatal(FARGS, "Checkpoint file %s could not be read!", filename);
2384 /* tune_pme will run the next part of the simulation */
2385 sim_part = cpt_sim_part + 1;
2388 /* Open performance output file and write header info */
2389 fp = gmx_ffopen(opt2fn("-p", NFILE, fnm), "w");
2391 /* Make a quick consistency check of command line parameters */
2392 check_input(nnodes, repeats, &ntprs, &rmin, rcoulomb, &rmax,
2393 maxPMEfraction, minPMEfraction, npme_fixed,
2394 bench_nsteps, fnm, NFILE, sim_part, presteps,
2397 /* Determine the maximum and minimum number of PME nodes to test,
2398 * the actual list of settings is build in do_the_tests(). */
2399 if ((nnodes > 2) && (npme_fixed < -1))
2401 if (0 == std::strcmp(npmevalues_opt[0], "auto"))
2403 /* Determine the npme range automatically based on the PME:PP load guess */
2404 if (guessPMEratio > 1.0)
2406 /* More PME than PP work, probably we do not need separate PME nodes at all! */
2407 maxPMEnodes = nnodes/2;
2408 minPMEnodes = nnodes/2;
2412 /* PME : PP load is in the range 0..1, let's test around the guess */
2413 guessPMEnodes = static_cast<int>(nnodes/(1.0 + 1.0/guessPMEratio));
2414 minPMEnodes = static_cast<int>(std::floor(0.7*guessPMEnodes));
2415 maxPMEnodes = static_cast<int>(std::ceil(1.6*guessPMEnodes));
2416 maxPMEnodes = std::min(maxPMEnodes, nnodes/2);
2421 /* Determine the npme range based on user input */
2422 maxPMEnodes = static_cast<int>(std::floor(maxPMEfraction*nnodes));
2423 minPMEnodes = std::max(static_cast<int>(std::floor(minPMEfraction*nnodes)), 0);
2424 fprintf(stdout, "Will try runs with %d ", minPMEnodes);
2425 if (maxPMEnodes != minPMEnodes)
2427 fprintf(stdout, "- %d ", maxPMEnodes);
2429 fprintf(stdout, "PME-only ranks.\n Note that the automatic number of PME-only ranks and no separate PME ranks are always tested.\n");
2438 /* Get the commands we need to set up the runs from environment variables */
2439 get_program_paths(bThreads, &cmd_mpirun, &cmd_mdrun);
2440 if (bBenchmark && repeats > 0)
2442 check_mdrun_works(bThreads, cmd_mpirun, cmd_np, cmd_mdrun, nullptr != eligible_gpu_ids);
2445 /* Print some header info to file */
2447 fprintf(fp, "\n P E R F O R M A N C E R E S U L T S\n");
2449 fprintf(fp, "%s for GROMACS %s\n", output_env_get_program_display_name(oenv),
2453 fprintf(fp, "Number of ranks : %d\n", nnodes);
2454 fprintf(fp, "The mpirun command is : %s\n", cmd_mpirun);
2455 if (std::strcmp(procstring[0], "none") != 0)
2457 fprintf(fp, "Passing # of ranks via : -%s\n", procstring[0]);
2461 fprintf(fp, "Not setting number of ranks in system call\n");
2466 fprintf(fp, "Number of threads : %d\n", nnodes);
2469 fprintf(fp, "The mdrun command is : %s\n", cmd_mdrun);
2470 fprintf(fp, "mdrun args benchmarks : %s\n", cmd_args_bench);
2471 fprintf(fp, "Benchmark steps : ");
2472 fprintf(fp, "%" PRId64, bench_nsteps);
2474 fprintf(fp, "dlb equilibration steps : %d\n", presteps);
2477 fprintf(fp, "Checkpoint time step : ");
2478 fprintf(fp, "%" PRId64, cpt_steps);
2481 fprintf(fp, "mdrun args at launchtime: %s\n", cmd_args_launch);
2483 if (new_sim_nsteps >= 0)
2486 fprintf(stderr, "Note: Simulation input file %s will have ", opt2fn("-so", NFILE, fnm));
2487 fprintf(stderr, "%" PRId64, new_sim_nsteps+cpt_steps);
2488 fprintf(stderr, " steps.\n");
2489 fprintf(fp, "Simulation steps : ");
2490 fprintf(fp, "%" PRId64, new_sim_nsteps);
2495 fprintf(fp, "Repeats for each test : %d\n", repeats);
2498 if (npme_fixed >= -1)
2500 fprintf(fp, "Fixing -npme at : %d\n", npme_fixed);
2503 fprintf(fp, "Input file : %s\n", opt2fn("-s", NFILE, fnm));
2504 fprintf(fp, " PME/PP load estimate : %g\n", guessPMEratio);
2506 /* Allocate memory for the inputinfo struct: */
2508 info->nr_inputfiles = ntprs;
2509 for (i = 0; i < ntprs; i++)
2511 snew(info->rcoulomb, ntprs);
2512 snew(info->rvdw, ntprs);
2513 snew(info->rlist, ntprs);
2514 snew(info->nkx, ntprs);
2515 snew(info->nky, ntprs);
2516 snew(info->nkz, ntprs);
2517 snew(info->fsx, ntprs);
2518 snew(info->fsy, ntprs);
2519 snew(info->fsz, ntprs);
2521 /* Make alternative tpr files to test: */
2522 snew(tpr_names, ntprs);
2523 for (i = 0; i < ntprs; i++)
2525 snew(tpr_names[i], STRLEN);
2528 /* It can be that ntprs is reduced by make_benchmark_tprs if not enough
2529 * different grids could be found. */
2530 make_benchmark_tprs(opt2fn("-s", NFILE, fnm), tpr_names, bench_nsteps+presteps,
2531 cpt_steps, rmin, rmax, bScaleRvdw, &ntprs, info, fp);
2533 /********************************************************************************/
2534 /* Main loop over all scenarios we need to test: tpr files, PME nodes, repeats */
2535 /********************************************************************************/
2536 snew(perfdata, ntprs);
2539 GMX_RELEASE_ASSERT(npmevalues_opt[0] != nullptr, "Options inconsistency; npmevalues_opt[0] is NULL");
2540 do_the_tests(fp, tpr_names, maxPMEnodes, minPMEnodes, npme_fixed, npmevalues_opt[0], perfdata, &pmeentries,
2541 repeats, nnodes, ntprs, bThreads, cmd_mpirun, cmd_np, cmd_mdrun,
2542 cmd_args_bench, fnm, NFILE, presteps, cpt_steps, bCheck, eligible_gpu_ids);
2544 fprintf(fp, "\nTuning took%8.1f minutes.\n", (gmx_gettime()-seconds)/60.0);
2546 /* Analyse the results and give a suggestion for optimal settings: */
2547 bKeepTPR = analyze_data(fp, opt2fn("-p", NFILE, fnm), perfdata, nnodes, ntprs, pmeentries,
2548 repeats, info, &best_tpr, &best_npme);
2550 /* Take the best-performing tpr file and enlarge nsteps to original value */
2551 if (bKeepTPR && !bOverwrite)
2553 simulation_tpr = opt2fn("-s", NFILE, fnm);
2557 simulation_tpr = opt2fn("-so", NFILE, fnm);
2558 modify_PMEsettings(bOverwrite ? (new_sim_nsteps+cpt_steps) : info->orig_sim_steps,
2559 info->orig_init_step, tpr_names[best_tpr], simulation_tpr);
2562 /* Let's get rid of the temporary benchmark input files */
2563 for (i = 0; i < ntprs; i++)
2565 fprintf(stdout, "Deleting temporary benchmark input file %s\n", tpr_names[i]);
2566 remove(tpr_names[i]);
2569 /* Now start the real simulation if the user requested it ... */
2570 launch_simulation(bLaunch, fp, bThreads, cmd_mpirun, cmd_np, cmd_mdrun,
2571 cmd_args_launch, simulation_tpr, best_npme, eligible_gpu_ids);
2575 /* ... or simply print the performance results to screen: */
2578 finalize(opt2fn("-p", NFILE, fnm));