Update advice about grompp -t
[alexxy/gromacs.git] / src / gromacs / gmxpreprocess / grompp.c
index 04ee47de4053d37c9ffca49f81165448a4312224..0983650920dc88dc5cccafcd9ef21d8910183fa8 100644 (file)
@@ -3,7 +3,7 @@
  *
  * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
  * Copyright (c) 2001-2004, The GROMACS development team.
- * Copyright (c) 2013,2014, by the GROMACS development team, led by
+ * Copyright (c) 2013,2014,2015, 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.
@@ -1374,7 +1374,7 @@ static void set_verlet_buffer(const gmx_mtop_t *mtop,
                             &ls, &n_nonlin_vsite, &rlist_1x1);
 
     /* Set the pair-list buffer size in ir */
-    verletbuf_get_list_setup(FALSE, &ls);
+    verletbuf_get_list_setup(FALSE, FALSE, &ls);
     calc_verlet_buffer_size(mtop, det(box), ir, buffer_temp,
                             &ls, &n_nonlin_vsite, &ir->rlist);
 
@@ -1392,6 +1392,8 @@ static void set_verlet_buffer(const gmx_mtop_t *mtop,
            ls.cluster_size_i, ls.cluster_size_j,
            ir->rlist, ir->rlist-max(ir->rvdw, ir->rcoulomb));
 
+    printf("Note that mdrun will redetermine rlist based on the actual pair-list setup\n");
+
     if (sqr(ir->rlistlong) >= max_cutoff2(ir->ePBC, box))
     {
         gmx_fatal(FARGS, "The pair-list cut-off (%g nm) is longer than half the shortest box vector or longer than the smallest box diagonal element (%g nm). Increase the box size or decrease nstlist or increase verlet-buffer-tolerance.", ir->rlistlong, sqrt(max_cutoff2(ir->ePBC, box)));
@@ -1474,7 +1476,9 @@ int gmx_grompp(int argc, char *argv[])
         "with [TT]-cpi[tt]. If you wish to change the ensemble or things",
         "like output frequency, then supplying the checkpoint file to",
         "[THISMODULE] with [TT]-t[tt] along with a new [TT].mdp[tt] file",
-        "with [TT]-f[tt] is the recommended procedure.[PAR]",
+        "with [TT]-f[tt] is the recommended procedure. Actually preserving",
+        "the ensemble (if possible) still requires passing the checkpoint",
+        "file to [gmx-mdrun] [TT]-cpi[tt].[PAR]",
 
         "By default, all bonded interactions which have constant energy due to",
         "virtual site constructions will be removed. If this constant energy is",
@@ -1508,7 +1512,7 @@ int gmx_grompp(int argc, char *argv[])
     t_inputrec        *ir;
     int                natoms, nvsite, comb, mt;
     t_params          *plist;
-    t_state            state;
+    t_state           *state;
     matrix             box;
     real               max_spacing, fudgeQQ;
     double             reppow;
@@ -1627,8 +1631,9 @@ int gmx_grompp(int argc, char *argv[])
     {
         gmx_fatal(FARGS, "%s does not exist", fn);
     }
+    snew(state, 1);
     new_status(fn, opt2fn_null("-pp", NFILE, fnm), opt2fn("-c", NFILE, fnm),
-               opts, ir, bZero, bGenVel, bVerbose, &state,
+               opts, ir, bZero, bGenVel, bVerbose, state,
                atype, sys, &nmi, &mi, plist, &comb, &reppow, &fudgeQQ,
                opts->bMorse,
                wi);
@@ -1817,7 +1822,7 @@ int gmx_grompp(int argc, char *argv[])
     /* Check velocity for virtual sites and shells */
     if (bGenVel)
     {
-        check_vel(sys, state.v);
+        check_vel(sys, state->v);
     }
 
     /* check for shells and inpurecs */
@@ -1849,7 +1854,7 @@ int gmx_grompp(int argc, char *argv[])
     }
     do_index(mdparin, ftp2fn_null(efNDX, NFILE, fnm),
              sys, bVerbose, ir,
-             bGenVel ? state.v : NULL,
+             bGenVel ? state->v : NULL,
              wi);
 
     if (ir->cutoff_scheme == ecutsVERLET && ir->verletbuf_tol > 0 &&
@@ -1867,7 +1872,7 @@ int gmx_grompp(int argc, char *argv[])
                 }
                 else
                 {
-                    buffer_temp = calc_temp(sys, ir, state.v);
+                    buffer_temp = calc_temp(sys, ir, state->v);
                 }
                 if (buffer_temp > 0)
                 {
@@ -1918,13 +1923,13 @@ int gmx_grompp(int argc, char *argv[])
                     warning_note(wi, warn_buf);
                 }
 
-                set_verlet_buffer(sys, ir, buffer_temp, state.box, wi);
+                set_verlet_buffer(sys, ir, buffer_temp, state->box, wi);
             }
         }
     }
 
     /* Init the temperature coupling state */
-    init_gtc_state(&state, ir->opts.ngtc, 0, ir->opts.nhchainlength); /* need to add nnhpres here? */
+    init_gtc_state(state, ir->opts.ngtc, 0, ir->opts.nhchainlength); /* need to add nnhpres here? */
 
     if (bVerbose)
     {
@@ -1964,24 +1969,24 @@ int gmx_grompp(int argc, char *argv[])
             fprintf(stderr, "getting data from old trajectory ...\n");
         }
         cont_status(ftp2fn(efTRN, NFILE, fnm), ftp2fn_null(efEDR, NFILE, fnm),
-                    bNeedVel, bGenVel, fr_time, ir, &state, sys, oenv);
+                    bNeedVel, bGenVel, fr_time, ir, state, sys, oenv);
     }
 
     if (ir->ePBC == epbcXY && ir->nwall != 2)
     {
-        clear_rvec(state.box[ZZ]);
+        clear_rvec(state->box[ZZ]);
     }
 
     if (ir->cutoff_scheme != ecutsVERLET && ir->rlist > 0)
     {
         set_warning_line(wi, mdparin, -1);
-        check_chargegroup_radii(sys, ir, state.x, wi);
+        check_chargegroup_radii(sys, ir, state->x, wi);
     }
 
     if (EEL_FULL(ir->coulombtype) || EVDW_PME(ir->vdwtype))
     {
         /* Calculate the optimal grid dimensions */
-        copy_mat(state.box, box);
+        copy_mat(state->box, box);
         if (ir->ePBC == epbcXY && ir->nwall == 2)
         {
             svmul(ir->wall_ewald_zfac, box[ZZ], box[ZZ]);
@@ -2005,13 +2010,13 @@ int gmx_grompp(int argc, char *argv[])
        potentially conflict if not handled correctly. */
     if (ir->efep != efepNO)
     {
-        state.fep_state = ir->fepvals->init_fep_state;
+        state->fep_state = ir->fepvals->init_fep_state;
         for (i = 0; i < efptNR; i++)
         {
             /* init_lambda trumps state definitions*/
             if (ir->fepvals->init_lambda >= 0)
             {
-                state.lambda[i] = ir->fepvals->init_lambda;
+                state->lambda[i] = ir->fepvals->init_lambda;
             }
             else
             {
@@ -2021,7 +2026,7 @@ int gmx_grompp(int argc, char *argv[])
                 }
                 else
                 {
-                    state.lambda[i] = ir->fepvals->all_lambda[i][state.fep_state];
+                    state->lambda[i] = ir->fepvals->all_lambda[i][state->fep_state];
                 }
             }
         }
@@ -2029,12 +2034,12 @@ int gmx_grompp(int argc, char *argv[])
 
     if (ir->ePull != epullNO)
     {
-        set_pull_init(ir, sys, state.x, state.box, state.lambda[efptMASS], oenv, opts->pull_start);
+        set_pull_init(ir, sys, state->x, state->box, state->lambda[efptMASS], oenv, opts->pull_start);
     }
 
     if (ir->bRot)
     {
-        set_reference_positions(ir->rot, state.x, state.box,
+        set_reference_positions(ir->rot, state->x, state->box,
                                 opt2fn("-ref", NFILE, fnm), opt2bSet("-ref", NFILE, fnm),
                                 wi);
     }
@@ -2043,7 +2048,7 @@ int gmx_grompp(int argc, char *argv[])
 
     if (EEL_PME(ir->coulombtype))
     {
-        float ratio = pme_load_estimate(sys, ir, state.box);
+        float ratio = pme_load_estimate(sys, ir, state->box);
         fprintf(stderr, "Estimate for the relative computational load of the PME mesh part: %.2f\n", ratio);
         /* With free energy we might need to do PME both for the A and B state
          * charges. This will double the cost, but the optimal performance will
@@ -2085,11 +2090,13 @@ int gmx_grompp(int argc, char *argv[])
     }
 
     done_warning(wi, FARGS);
-    write_tpx_state(ftp2fn(efTPX, NFILE, fnm), ir, &state, sys);
+    write_tpx_state(ftp2fn(efTPX, NFILE, fnm), ir, state, sys);
 
     /* Output IMD group, if bIMD is TRUE */
-    write_IMDgroup_to_file(ir->bIMD, ir, &state, sys, NFILE, fnm);
+    write_IMDgroup_to_file(ir->bIMD, ir, state, sys, NFILE, fnm);
 
+    done_state(state);
+    sfree(state);
     done_atomtype(atype);
     done_mtop(sys, TRUE);
     done_inputrec_strings();