/*
+ * This file is part of the GROMACS molecular simulation package.
*
- * 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
+ * Copyright (c) 2001-2004, The GROMACS development team.
+ * Copyright (c) 2011,2012,2013,2014, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and 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.
*
- * 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.
+ * 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.
*
- * 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.
+ * 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.
*
- * For more info, check our website at http://www.gromacs.org
+ * 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.
*
- * And Hey:
- * Gallium Rubidium Oxygen Manganese Argon Carbon Silicon
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
*/
#ifdef HAVE_CONFIG_H
#include <config.h>
#endif
#include "typedefs.h"
-#include "smalloc.h"
+#include "gromacs/utility/smalloc.h"
#include "sysstuff.h"
#include "vec.h"
-#include "statutil.h"
#include "vcm.h"
#include "mdebin.h"
#include "nrnb.h"
#include "vsite.h"
#include "update.h"
#include "ns.h"
-#include "trnio.h"
-#include "xtcio.h"
#include "mdrun.h"
#include "md_support.h"
#include "md_logging.h"
-#include "confio.h"
#include "network.h"
-#include "pull.h"
#include "xvgr.h"
#include "physics.h"
#include "names.h"
#include "pme.h"
#include "mdatoms.h"
#include "repl_ex.h"
+#include "deform.h"
#include "qmmm.h"
#include "domdec.h"
#include "domdec_network.h"
-#include "partdec.h"
-#include "topsort.h"
+#include "gromacs/gmxlib/topsort.h"
#include "coulomb.h"
#include "constr.h"
#include "shellfc.h"
-#include "compute_io.h"
-#include "mvdata.h"
+#include "gromacs/gmxpreprocess/compute_io.h"
#include "checkpoint.h"
#include "mtop_util.h"
#include "sighandler.h"
#include "txtdump.h"
-#include "string2.h"
+#include "gromacs/utility/cstringutil.h"
#include "pme_loadbal.h"
#include "bondf.h"
#include "membed.h"
#include "nbnxn_cuda_data_mgmt.h"
#include "gromacs/utility/gmxmpi.h"
+#include "gromacs/fileio/confio.h"
+#include "gromacs/fileio/trajectory_writing.h"
+#include "gromacs/fileio/trnio.h"
+#include "gromacs/fileio/trxio.h"
+#include "gromacs/fileio/xtcio.h"
+#include "gromacs/timing/wallcycle.h"
+#include "gromacs/timing/walltime_accounting.h"
+#include "gromacs/pulling/pull.h"
+#include "gromacs/swap/swapcoords.h"
+#include "gromacs/imd/imd.h"
#ifdef GMX_FAHCORE
#include "corewrap.h"
#endif
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_int64_t step,
+ gmx_int64_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];
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", gmx_gettime());
}
double do_md(FILE *fplog, t_commrec *cr, int nfile, const t_filenm fnm[],
int repl_ex_nst, int repl_ex_nex, int repl_ex_seed, gmx_membed_t membed,
real cpt_period, real max_hours,
const char gmx_unused *deviceOptions,
+ int imdport,
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 run_time;
+ gmx_mdoutf_t outf = NULL;
+ gmx_int64_t step, step_rel;
+ double elapsed_time;
double t, t0, lam0[efptNR];
gmx_bool bGStatEveryStep, bGStat, bCalcVir, bCalcEner;
gmx_bool bNS, bNStList, bSimAnn, bStopCM, bRerunMD, bNotLastFrame = FALSE,
gmx_repl_ex_t repl_ex = NULL;
int nchkpt = 1;
gmx_localtop_t *top;
- t_mdebin *mdebin = NULL;
+ t_mdebin *mdebin = NULL;
t_state *state = NULL;
rvec *f_global = NULL;
gmx_enerdata_t *enerd;
gmx_update_t upd = NULL;
t_graph *graph = NULL;
globsig_t gs;
- gmx_rng_t mcrng = NULL;
gmx_groups_t *groups;
gmx_ekindata_t *ekind, *ekind_save;
gmx_shellfc_t shellfc;
int count, nconverged = 0;
- real timestep = 0;
- double tcount = 0;
- gmx_bool bConverged = TRUE, bOK, bSumEkinhOld, bExchanged;
+ real timestep = 0;
+ double tcount = 0;
+ gmx_bool bConverged = TRUE, bOK, bSumEkinhOld, bDoReplEx, bExchanged, bNeedRepartition;
gmx_bool bAppend;
gmx_bool bResetCountersHalfMaxH = FALSE;
gmx_bool bVV, bIterativeCase, bFirstIterate, bTemp, bPres, bTrotter;
gmx_bool bUpdateDoLR;
real dvdl_constr;
- int a0, a1;
rvec *cbuf = NULL;
matrix lastbox;
real veta_save, scalevir, tracevir;
char sbuf[STEPSTRSIZE], sbuf2[STEPSTRSIZE];
int handled_stop_condition = gmx_stop_cond_none; /* compare to get_stop_condition*/
gmx_iterate_t iterate;
- gmx_large_int_t multisim_nsteps = -1; /* number of steps to do before first multisim
- simulation stops. If equal to zero, don't
- communicate any more between multisims.*/
+ gmx_int64_t multisim_nsteps = -1; /* number of steps to do before first multisim
+ simulation stops. If equal to zero, don't
+ communicate any more between multisims.*/
/* PME load balancing data for GPU kernels */
pme_load_balancing_t pme_loadbal = NULL;
double cycles_pmes;
gmx_bool bPMETuneTry = FALSE, bPMETuneRunning = FALSE;
+ /* Interactive MD */
+ gmx_bool bIMDstep = FALSE;
+
#ifdef GMX_FAHCORE
/* Temporary addition for FAHCORE checkpointing */
int chkpt_ret;
if (bRerunMD)
{
- ir->nstxtcout = 0;
+ ir->nstxout_compressed = 0;
}
groups = &top_global->groups;
&(state_global->fep_state), lam0,
nrnb, top_global, &upd,
nfile, fnm, &outf, &mdebin,
- force_vir, shake_vir, mu_tot, &bSimAnn, &vcm, Flags);
+ force_vir, shake_vir, mu_tot, &bSimAnn, &vcm, Flags, wcycle);
clear_mat(total_vir);
clear_mat(pres);
(ir->bContinuation ||
(DOMAINDECOMP(cr) && !MASTER(cr))) ?
NULL : state_global->x);
+ if (shellfc && ir->nstcalcenergy != 1)
+ {
+ gmx_fatal(FARGS, "You have nstcalcenergy set to a value (%d) that is different from 1.\nThis is not supported in combinations with shell particles.\nPlease make a new tpr file.", ir->nstcalcenergy);
+ }
+ if (shellfc && DOMAINDECOMP(cr))
+ {
+ gmx_fatal(FARGS, "Shell particles are not implemented with domain decomposition, use a single rank");
+ }
+ if (shellfc && ir->eI == eiNM)
+ {
+ /* Currently shells don't work with Normal Modes */
+ gmx_fatal(FARGS, "Normal Mode analysis is not supported with shells.\nIf you'd like to help with adding support, we have an open discussion at http://redmine.gromacs.org/issues/879\n");
+ }
+
+ if (vsite && ir->eI == eiNM)
+ {
+ /* Currently virtual sites don't work with Normal Modes */
+ gmx_fatal(FARGS, "Normal Mode analysis is not supported with virtual sites.\nIf you'd like to help with adding support, we have an open discussion at http://redmine.gromacs.org/issues/879\n");
+ }
if (DEFORM(*ir))
{
-#ifdef GMX_THREAD_MPI
tMPI_Thread_mutex_lock(&deform_init_box_mutex);
-#endif
set_deform_reference_box(upd,
deform_init_init_step_tpx,
deform_init_box_tpx);
-#ifdef GMX_THREAD_MPI
tMPI_Thread_mutex_unlock(&deform_init_box_mutex);
-#endif
}
{
}
else
{
- if (PAR(cr))
- {
- /* Initialize the particle decomposition and split the topology */
- top = split_system(fplog, top_global, ir, cr);
-
- pd_cg_range(cr, &fr->cg0, &fr->hcg);
- pd_at_range(cr, &a0, &a1);
- }
- else
- {
- top = gmx_mtop_generate_local_top(top_global, ir);
+ top = gmx_mtop_generate_local_top(top_global, ir);
- a0 = 0;
- a1 = top_global->natoms;
- }
+ forcerec_set_excl_load(fr, top);
- forcerec_set_excl_load(fr, top, cr);
-
- state = partdec_init_local_state(cr, state_global);
+ state = serial_init_local_state(state_global);
f_global = f;
- atoms2md(top_global, ir, 0, NULL, a0, a1-a0, mdatoms);
+ atoms2md(top_global, ir, 0, NULL, top_global->natoms, mdatoms);
if (vsite)
{
}
setup_bonded_threading(fr, &top->idef);
-
- if (ir->pull && PAR(cr))
- {
- dd_make_local_pull_groups(NULL, ir->pull, mdatoms);
- }
}
+ /* Set up interactive MD (IMD) */
+ init_IMD(ir, cr, top_global, fplog, ir->nstcalcenergy, state_global->x,
+ nfile, fnm, oenv, imdport, Flags);
+
if (DOMAINDECOMP(cr))
{
/* Distribute the charge groups over the nodes from the master node */
if (ir->bExpanded)
{
- init_expanded_ensemble(bStateFromCP,ir,&mcrng,&state->dfhist);
+ init_expanded_ensemble(bStateFromCP, ir, &state->dfhist);
}
if (MASTER(cr))
update_energyhistory(&state_global->enerhist, mdebin);
}
- if ((state->flags & (1<<estLD_RNG)) && (Flags & MD_READ_RNG))
- {
- /* Set the random state if we read a checkpoint file */
- set_stochd_state(upd, state);
- }
-
- if (state->flags & (1<<estMC_RNG))
- {
- set_mc_state(mcrng, state);
- }
-
/* Initialize constraints */
- if (constr)
+ if (constr && !DOMAINDECOMP(cr))
{
- if (!DOMAINDECOMP(cr))
- {
- set_constraints(constr, top, ir, mdatoms, cr);
- }
+ set_constraints(constr, top, ir, mdatoms, cr);
}
- if (repl_ex_nst > 0)
- {
- /* We need to be sure replica exchange can only occur
- * when the energies are current */
- check_nst_param(fplog, cr, "nstcalcenergy", ir->nstcalcenergy,
- "repl_ex_nst", &repl_ex_nst);
- /* This check needs to happen before inter-simulation
- * signals are initialized, too */
- }
if (repl_ex_nst > 0 && MASTER(cr))
{
repl_ex = init_replica_exchange(fplog, cr->ms, state_global, ir,
}
/* PME tuning is only supported with GPUs or PME nodes and not with rerun.
- * With perturbed charges with soft-core we should not change the cut-off.
+ * PME tuning is not supported with PME only for LJ and not for Coulomb.
*/
if ((Flags & MD_TUNEPME) &&
EEL_PME(fr->eeltype) &&
( (fr->cutoff_scheme == ecutsVERLET && fr->nbv->bUseGPU) || !(cr->duty & DUTY_PME)) &&
- !(ir->efep != efepNO && mdatoms->nChargePerturbed > 0 && ir->fepvals->bScCoul) &&
!bRerunMD)
{
pme_loadbal_init(&pme_loadbal, ir, state->box, fr->ic, fr->pmedata);
if (mdatoms->cFREEZE && (state->flags & (1<<estV)))
{
/* Set the velocities of frozen particles to zero */
- for (i = mdatoms->start; i < mdatoms->start+mdatoms->homenr; i++)
+ for (i = 0; i < mdatoms->homenr; i++)
{
for (m = 0; m < DIM; m++)
{
/* Construct the virtual sites for the initial configuration */
construct_vsites(vsite, state->x, ir->delta_t, NULL,
top->idef.iparams, top->idef.il,
- fr->ePBC, fr->bMolPBC, graph, cr, state->box);
+ fr->ePBC, fr->bMolPBC, cr, state->box);
}
}
fprintf(fplog, "\n");
}
- /* Set and write start time */
- runtime_start(runtime);
- print_date_and_time(fplog, cr->nodeid, "Started mdrun", runtime);
+ walltime_accounting_start(walltime_accounting);
wallcycle_start(wcycle, ewcRUN);
- if (fplog)
- {
- fprintf(fplog, "\n");
- }
+ print_start(fplog, cr, walltime_accounting, "mdrun");
/* safest point to do file checkpointing is here. More general point would be immediately before integrator call */
#ifdef GMX_FAHCORE
bStartingFromCpt = (Flags & MD_STARTFROMCPT) && bInitStep;
bLastStep = FALSE;
bSumEkinhOld = FALSE;
+ bDoReplEx = FALSE;
bExchanged = FALSE;
+ bNeedRepartition = FALSE;
init_global_signals(&gs, cr, ir, repl_ex_nst);
&& (ir->bExpanded) && (step > 0) && (!bStartingFromCpt));
}
+ bDoReplEx = ((repl_ex_nst > 0) && (step > 0) && !bLastStep &&
+ do_per_step(step, repl_ex_nst));
+
if (bSimAnn)
{
update_annealing_target_temp(&(ir->opts), t);
if (bRerunMD)
{
- if (!(DOMAINDECOMP(cr) && !MASTER(cr)))
+ if (!DOMAINDECOMP(cr) || MASTER(cr))
{
for (i = 0; i < state_global->natoms; i++)
{
{
if (DOMAINDECOMP(cr))
{
- gmx_fatal(FARGS, "Vsite recalculation with -rerun is not implemented for domain decomposition, use particle decomposition");
+ gmx_fatal(FARGS, "Vsite recalculation with -rerun is not implemented with domain decomposition, use a single rank");
}
if (graph)
{
}
construct_vsites(vsite, state->x, ir->delta_t, state->v,
top->idef.iparams, top->idef.il,
- fr->ePBC, fr->bMolPBC, graph, cr, state->box);
+ fr->ePBC, fr->bMolPBC, cr, state->box);
if (graph)
{
unshift_self(graph, state->box, state->x);
/* Determine whether or not to do Neighbour Searching and LR */
bNStList = (ir->nstlist > 0 && step % ir->nstlist == 0);
- bNS = (bFirstStep || bExchanged || bNStList || bDoFEP ||
+ bNS = (bFirstStep || bExchanged || bNeedRepartition || bNStList || bDoFEP ||
(ir->nstlist == -1 && nlh.nabnsb > 0));
if (bNS && ir->nstlist == -1)
{
- set_nlistheuristics(&nlh, bFirstStep || bExchanged || bDoFEP, step);
+ set_nlistheuristics(&nlh, bFirstStep || bExchanged || bNeedRepartition || bDoFEP, step);
}
}
do_ene = (do_per_step(step, ir->nstenergy) || bLastStep);
- if (do_ene || do_log)
+ if (do_ene || do_log || bDoReplEx)
{
bCalcVir = TRUE;
bCalcEner = TRUE;
nrnb, wcycle, graph, groups,
shellfc, fr, bBornRadii, t, mu_tot,
&bConverged, vsite,
- outf->fp_field);
+ mdoutf_get_fp_field(outf));
tcount += count;
if (bConverged)
state->box, state->x, &state->hist,
f, force_vir, mdatoms, enerd, fcd,
state->lambda, graph,
- fr, vsite, mu_tot, t, outf->fp_field, ed, bBornRadii,
+ fr, vsite, mu_tot, t, mdoutf_get_fp_field(outf), ed, bBornRadii,
(bNS ? GMX_FORCE_NS : 0) | force_flags);
}
if (bVV && !bStartingFromCpt && !bRerunMD)
/* ############### START FIRST UPDATE HALF-STEP FOR VV METHODS############### */
{
+ wallcycle_start(wcycle, ewcUPDATE);
if (ir->eI == eiVV && bInitStep)
{
/* if using velocity verlet with full time step Ekin,
* step to combine the long-range forces on these steps.
* For nstcalclr=1 this is not done, since the forces would have been added
* directly to the short-range forces already.
+ *
+ * TODO Remove various aspects of VV+twin-range in master
+ * branch, because VV integrators did not ever support
+ * twin-range multiple time stepping with constraints.
*/
bUpdateDoLR = (fr->bTwinRange && do_per_step(step, ir->nstcalclr));
update_coords(fplog, step, ir, mdatoms, state, fr->bMolPBC,
- f, bUpdateDoLR, fr->f_twin, fcd,
+ f, bUpdateDoLR, fr->f_twin, bCalcVir ? &fr->vir_twin_constr : NULL, fcd,
ekind, M, upd, bInitStep, etrtVELOCITY1,
cr, nrnb, constr, &top->idef);
bOK = TRUE;
if (!bRerunMD || rerun_fr.bV || bForceUpdate) /* Why is rerun_fr.bV here? Unclear. */
{
+ wallcycle_stop(wcycle, ewcUPDATE);
update_constraints(fplog, step, NULL, ir, ekind, mdatoms,
state, fr->bMolPBC, graph, f,
&top->idef, shake_vir,
cr, nrnb, wcycle, upd, constr,
TRUE, bCalcVir, vetanew);
+ wallcycle_start(wcycle, ewcUPDATE);
+
+ if (bCalcVir && bUpdateDoLR && ir->nstcalclr > 1)
+ {
+ /* Correct the virial for multiple time stepping */
+ m_sub(shake_vir, fr->vir_twin_constr, shake_vir);
+ }
if (!bOK)
{
So we need information from the last step in the first half of the integration */
if (bGStat || do_per_step(step-1, nstglobalcomm))
{
+ wallcycle_stop(wcycle, ewcUPDATE);
compute_globals(fplog, gstat, cr, ir, fr, ekind, state, state_global, mdatoms, nrnb, vcm,
wcycle, enerd, force_vir, shake_vir, total_vir, pres, mu_tot,
constr, NULL, FALSE, state->box,
time step kinetic energy for the pressure (always true now, since we want accurate statistics).
b) If we are using EkinAveEkin for the kinetic energy for the temperature control, we still feed in
EkinAveVel because it's needed for the pressure */
+ wallcycle_start(wcycle, ewcUPDATE);
}
/* temperature scaling and pressure scaling to produce the extended variables at t+dt */
if (!bInitStep)
{
if (bExchanged)
{
-
+ wallcycle_stop(wcycle, ewcUPDATE);
/* We need the kinetic energy at minus the half step for determining
* the full step kinetic energy and possibly for T-coupling.*/
/* This may not be quite working correctly yet . . . . */
constr, NULL, FALSE, state->box,
top_global, &bSumEkinhOld,
CGLO_RERUNMD | CGLO_GSTAT | CGLO_TEMPERATURE);
+ wallcycle_start(wcycle, ewcUPDATE);
}
}
}
{
copy_rvecn(cbuf, state->v, 0, state->natoms);
}
+ wallcycle_stop(wcycle, ewcUPDATE);
}
/* MRS -- now done iterating -- compute the conserved quantity */
statistics, but if performing simulated tempering, we
do update the velocities and the tau_t. */
- lamnew = ExpandedEnsembleDynamics(fplog, ir, enerd, state, &MassQ, state->fep_state, &state->dfhist, step, mcrng, state->v, mdatoms);
+ lamnew = ExpandedEnsembleDynamics(fplog, ir, enerd, state, &MassQ, state->fep_state, &state->dfhist, step, state->v, mdatoms);
/* history is maintained in state->dfhist, but state_global is what is sent to trajectory and log output */
- copy_df_history(&state_global->dfhist,&state->dfhist);
+ copy_df_history(&state_global->dfhist, &state->dfhist);
}
/* Now we have the energies and forces corresponding to the
* coordinates at time t. We must output all of this before
* the update.
*/
- do_trajectory_writing(fplog, cr, nfile, fnm, step, step_rel, t,
- ir, state, state_global, top_global, fr, upd,
- outf, mdebin, ekind, f, f_global,
- wcycle, mcrng, &nchkpt,
- bCPT, bRerunMD, bLastStep, (Flags & MD_CONFOUT),
- bSumEkinhOld);
+ do_md_trajectory_writing(fplog, cr, nfile, fnm, step, step_rel, t,
+ ir, state, state_global, top_global, fr,
+ outf, mdebin, ekind, f, f_global,
+ &nchkpt,
+ bCPT, bRerunMD, bLastStep, (Flags & MD_CONFOUT),
+ bSumEkinhOld);
+ /* Check if IMD step and do IMD communication, if bIMD is TRUE. */
+ bIMDstep = do_IMD(ir->bIMD, step, cr, bNS, state->box, state->x, ir, t, wcycle);
/* kludge -- virial is lost with restart for NPT control. Must restart */
if (bStartingFromCpt && bVV)
copy_mat(state->fvir_prev, force_vir);
}
- /* Determine the wallclock run time up till now */
- run_time = gmx_gettime() - (double)runtime->real;
+ 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)
handled_stop_condition = (int)gmx_get_stop_condition();
}
else if (MASTER(cr) && (bNS || ir->nstlist <= 0) &&
- (max_hours > 0 && 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 */
}
if (bResetCountersHalfMaxH && MASTER(cr) &&
- run_time > max_hours*60.0*60.0*0.495)
+ elapsed_time > max_hours*60.0*60.0*0.495)
{
gs.sig[eglsRESETCOUNTERS] = 1;
}
if (MASTER(cr) && ((bGStat || !PAR(cr)) &&
cpt_period >= 0 &&
(cpt_period == 0 ||
- run_time >= nchkpt*cpt_period*60.0)) &&
+ elapsed_time >= nchkpt*cpt_period*60.0)) &&
gs.set[eglsCHKPT] == 0)
{
gs.sig[eglsCHKPT] = 1;
{
if (!bInitStep)
{
- update_tcouple(step, ir, state, ekind, upd, &MassQ, mdatoms);
+ update_tcouple(step, ir, state, ekind, &MassQ, mdatoms);
}
if (ETC_ANDERSEN(ir->etc)) /* keep this outside of update_tcouple because of the extra info required to pass */
{
gmx_bool bIfRandomize;
- bIfRandomize = update_randomize_velocities(ir, step, mdatoms, state, upd, &top->idef, constr);
+ bIfRandomize = update_randomize_velocities(ir, step, cr, mdatoms, state, upd, constr);
/* if we have constraints, we have to remove the kinetic energy parallel to the bonds */
if (constr && bIfRandomize)
{
}
else
{
- update_tcouple(step, ir, state, ekind, upd, &MassQ, mdatoms);
+ update_tcouple(step, ir, state, ekind, &MassQ, mdatoms);
update_pcouple(fplog, step, ir, state, pcoupl_mu, M, bInitStep);
}
/* velocity half-step update */
update_coords(fplog, step, ir, mdatoms, state, fr->bMolPBC, f,
- bUpdateDoLR, fr->f_twin, fcd,
+ bUpdateDoLR, fr->f_twin, bCalcVir ? &fr->vir_twin_constr : NULL, fcd,
ekind, M, upd, FALSE, etrtVELOCITY2,
cr, nrnb, constr, &top->idef);
}
bUpdateDoLR = (fr->bTwinRange && do_per_step(step, ir->nstcalclr));
update_coords(fplog, step, ir, mdatoms, state, fr->bMolPBC, f,
- bUpdateDoLR, fr->f_twin, fcd,
+ bUpdateDoLR, fr->f_twin, bCalcVir ? &fr->vir_twin_constr : NULL, fcd,
ekind, M, upd, bInitStep, etrtPOSITION, cr, nrnb, constr, &top->idef);
wallcycle_stop(wcycle, ewcUPDATE);
cr, nrnb, wcycle, upd, constr,
FALSE, bCalcVir, state->veta);
+ if (bCalcVir && bUpdateDoLR && ir->nstcalclr > 1)
+ {
+ /* Correct the virial for multiple time stepping */
+ m_sub(shake_vir, fr->vir_twin_constr, shake_vir);
+ }
+
if (ir->eI == eiVVAK)
{
/* erase F_EKIN and F_TEMP here? */
bUpdateDoLR = (fr->bTwinRange && do_per_step(step, ir->nstcalclr));
update_coords(fplog, step, ir, mdatoms, state, fr->bMolPBC, f,
- bUpdateDoLR, fr->f_twin, fcd,
+ bUpdateDoLR, fr->f_twin, bCalcVir ? &fr->vir_twin_constr : NULL, fcd,
ekind, M, upd, bInitStep, etrtPOSITION, cr, nrnb, constr, &top->idef);
wallcycle_stop(wcycle, ewcUPDATE);
}
construct_vsites(vsite, state->x, ir->delta_t, state->v,
top->idef.iparams, top->idef.il,
- fr->ePBC, fr->bMolPBC, graph, cr, state->box);
+ fr->ePBC, fr->bMolPBC, cr, state->box);
if (graph != NULL)
{
}
/* ######### END PREPARING EDR OUTPUT ########### */
- /* Time for performance */
- if (((step % stepout) == 0) || bLastStep)
- {
- runtime_upd_proc(runtime);
- }
-
/* Output stuff */
if (MASTER(cr))
{
do_dr = do_per_step(step, ir->nstdisreout);
do_or = do_per_step(step, ir->nstorireout);
- print_ebin(outf->fp_ene, do_ene, do_dr, do_or, do_log ? fplog : NULL,
+ print_ebin(mdoutf_get_fp_ene(outf), do_ene, do_dr, do_or, do_log ? fplog : NULL,
step, t,
eprNORMAL, bCompact, mdebin, fcd, groups, &(ir->opts));
}
/* 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);
+ }
+
+ /* Ion/water position swapping.
+ * Not done in last step since trajectory writing happens before this call
+ * in the MD loop and exchanges would be lost anyway. */
+ bNeedRepartition = FALSE;
+ if ((ir->eSwapCoords != eswapNO) && (step > 0) && !bLastStep &&
+ do_per_step(step, ir->swap->nstswap))
+ {
+ bNeedRepartition = do_swapcoords(cr, step, t, ir, wcycle,
+ bRerunMD ? rerun_fr.x : state->x,
+ bRerunMD ? rerun_fr.box : state->box,
+ top_global, MASTER(cr) && bVerbose, bRerunMD);
+
+ if (bNeedRepartition && DOMAINDECOMP(cr))
+ {
+ dd_collect_state(cr->dd, state, state_global);
+ }
}
/* Replica exchange */
bExchanged = FALSE;
- if ((repl_ex_nst > 0) && (step > 0) && !bLastStep &&
- do_per_step(step, repl_ex_nst))
+ if (bDoReplEx)
{
bExchanged = replica_exchange(fplog, cr, repl_ex,
state_global, enerd,
state, step, t);
+ }
- if (bExchanged && DOMAINDECOMP(cr))
- {
- dd_partition_system(fplog, step, cr, TRUE, 1,
- state_global, top_global, ir,
- state, &f, mdatoms, top, fr,
- vsite, shellfc, constr,
- nrnb, wcycle, FALSE);
- }
+ if ( (bExchanged || bNeedRepartition) && DOMAINDECOMP(cr) )
+ {
+ dd_partition_system(fplog, step, cr, TRUE, 1,
+ state_global, top_global, ir,
+ state, &f, mdatoms, top, fr,
+ vsite, shellfc, constr,
+ nrnb, wcycle, FALSE);
}
bFirstStep = FALSE;
}
dd_bcast(cr->dd, sizeof(gmx_bool), &bPMETuneRunning);
+ if (bPMETuneRunning &&
+ fr->nbv->bUseGPU && DOMAINDECOMP(cr) &&
+ !(cr->duty & DUTY_PME))
+ {
+ /* Lock DLB=auto to off (does nothing when DLB=yes/no).
+ * With GPUs + separate PME ranks, we don't want DLB.
+ * This could happen when we scan coarse grids and
+ * it would then never be turned off again.
+ * This would hurt performance at the final, optimal
+ * grid spacing, where DLB almost never helps.
+ * Also, DLB can limit the cut-off for PME tuning.
+ */
+ dd_dlb_set_lock(cr->dd, TRUE);
+ }
+
if (bPMETuneRunning || step_rel > ir->nstlist*50)
{
bPMETuneTry = FALSE;
step);
/* Update constants in forcerec/inputrec to keep them in sync with fr->ic */
- fr->ewaldcoeff = fr->ic->ewaldcoeff;
- fr->rlist = fr->ic->rlist;
- fr->rlistlong = fr->ic->rlistlong;
- fr->rcoulomb = fr->ic->rcoulomb;
- fr->rvdw = fr->ic->rvdw;
+ fr->ewaldcoeff_q = fr->ic->ewaldcoeff_q;
+ fr->ewaldcoeff_lj = fr->ic->ewaldcoeff_lj;
+ fr->rlist = fr->ic->rlist;
+ fr->rlistlong = fr->ic->rlistlong;
+ fr->rcoulomb = fr->ic->rcoulomb;
+ fr->rvdw = fr->ic->rvdw;
+
+ if (ir->eDispCorr != edispcNO)
+ {
+ calc_enervirdiff(NULL, ir->eDispCorr, fr);
+ }
+
+ if (!bPMETuneRunning &&
+ DOMAINDECOMP(cr) &&
+ dd_dlb_is_locked(cr->dd))
+ {
+ /* Unlock the DLB=auto, DLB is allowed to activate
+ * (but we don't expect it to activate in most cases).
+ */
+ dd_dlb_set_lock(cr->dd, FALSE);
+ }
}
cycles_pmes = 0;
}
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))
gmx_pme_send_resetcounters(cr, step);
}
/* Correct max_hours for the elapsed time */
- max_hours -= run_time/(60.0*60.0);
+ max_hours -= elapsed_time/(60.0*60.0);
bResetCountersHalfMaxH = FALSE;
gs.set[eglsRESETCOUNTERS] = 0;
}
+ /* If bIMD is TRUE, the master updates the IMD energy record and sends positions to VMD client */
+ IMD_prep_energies_send_positions(ir->bIMD && MASTER(cr), bIMDstep, ir->imd, enerd, step, bCalcEner, wcycle);
+
}
/* End of main MD loop */
debug_gmx();
- /* Stop the time */
- runtime_end(runtime);
+ /* Closing TNG files can include compressing data. Therefore it is good to do that
+ * before stopping the time measurements. */
+ mdoutf_tng_close(outf);
+
+ /* Stop measuring walltime */
+ walltime_accounting_end(walltime_accounting);
if (bRerunMD && MASTER(cr))
{
{
if (ir->nstcalcenergy > 0 && !bRerunMD)
{
- print_ebin(outf->fp_ene, FALSE, FALSE, FALSE, fplog, step, t,
+ print_ebin(mdoutf_get_fp_ene(outf), FALSE, FALSE, FALSE, fplog, step, t,
eprAVER, FALSE, mdebin, fcd, groups, &(ir->opts));
}
}
done_mdoutf(outf);
-
debug_gmx();
if (ir->nstlist == -1 && nlh.nns > 0 && fplog)
print_replica_exchange_statistics(fplog, repl_ex);
}
- runtime->nsteps_done = step_rel;
+ /* IMD cleanup, if bIMD is TRUE. */
+ IMD_finalize(ir->bIMD, ir->imd);
+
+ walltime_accounting_set_nsteps_done(walltime_accounting, step_rel);
return 0;
}
biod.pnpi.spb.ru / alexxy / gromacs.git / blobdiff