* 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 "gmxpre.h"
+
+#include "readir.h"
#include <ctype.h>
-#include <stdlib.h>
#include <limits.h>
-#include "sysstuff.h"
-#include "gromacs/utility/smalloc.h"
-#include "typedefs.h"
-#include "physics.h"
-#include "names.h"
-#include "gmx_fatal.h"
-#include "macros.h"
-#include "index.h"
-#include "symtab.h"
+#include <stdlib.h>
+
+#include "gromacs/gmxpreprocess/toputil.h"
+#include "gromacs/legacyheaders/chargegroup.h"
+#include "gromacs/legacyheaders/inputrec.h"
+#include "gromacs/legacyheaders/macros.h"
+#include "gromacs/legacyheaders/names.h"
+#include "gromacs/legacyheaders/network.h"
+#include "gromacs/legacyheaders/readinp.h"
+#include "gromacs/legacyheaders/typedefs.h"
+#include "gromacs/legacyheaders/warninp.h"
+#include "gromacs/math/units.h"
+#include "gromacs/math/vec.h"
+#include "gromacs/mdlib/calc_verletbuf.h"
+#include "gromacs/pbcutil/pbc.h"
+#include "gromacs/topology/block.h"
+#include "gromacs/topology/index.h"
+#include "gromacs/topology/mtop_util.h"
+#include "gromacs/topology/symtab.h"
#include "gromacs/utility/cstringutil.h"
-#include "readinp.h"
-#include "warninp.h"
-#include "readir.h"
-#include "toputil.h"
-#include "index.h"
-#include "network.h"
-#include "vec.h"
-#include "pbc.h"
-#include "mtop_util.h"
-#include "chargegroup.h"
-#include "inputrec.h"
-#include "calc_verletbuf.h"
+#include "gromacs/utility/fatalerror.h"
+#include "gromacs/utility/smalloc.h"
#define MAXPTR 254
#define NOGID 255
}
else if (simtemp->eSimTempScale == esimtempEXPONENTIAL)
{
- simtemp->temperatures[i] = simtemp->simtemp_low + (simtemp->simtemp_high-simtemp->simtemp_low)*((exp(temperature_lambdas[i])-1)/(exp(1.0)-1));
+ simtemp->temperatures[i] = simtemp->simtemp_low + (simtemp->simtemp_high-simtemp->simtemp_low)*(gmx_expm1(temperature_lambdas[i])/gmx_expm1(1.0));
}
else
{
{
warning_error(wi, "rlist should be >= 0");
}
+ sprintf(err_buf, "nstlist can not be smaller than 0. (If you were trying to use the heuristic neighbour-list update scheme for efficient buffering for improved energy conservation, please use the Verlet cut-off scheme instead.)");
+ CHECK(ir->nstlist < 0);
process_interaction_modifier(ir, &ir->coulomb_modifier);
process_interaction_modifier(ir, &ir->vdw_modifier);
{
warning_error(wi, "rlistlong can not be shorter than rlist");
}
- if (IR_TWINRANGE(*ir) && ir->nstlist <= 0)
+ if (IR_TWINRANGE(*ir) && ir->nstlist == 0)
{
- warning_error(wi, "Can not have nstlist<=0 with twin-range interactions");
+ warning_error(wi, "Can not have nstlist == 0 with twin-range interactions");
}
}
(ir->ePBC != epbcNONE) ||
(ir->rcoulomb != 0.0) || (ir->rvdw != 0.0));
- if (ir->nstlist < 0)
- {
- warning_error(wi, "Can not have heuristic neighborlist updates without cut-off");
- }
if (ir->nstlist > 0)
{
warning_note(wi, "Simulating without cut-offs can be (slightly) faster with nstlist=0, nstype=simple and only one MPI rank");
if (ir->cutoff_scheme == ecutsGROUP)
{
if (((ir->coulomb_modifier != eintmodNONE && ir->rcoulomb == ir->rlist) ||
- (ir->vdw_modifier != eintmodNONE && ir->rvdw == ir->rlist)) &&
- ir->nstlist != 1)
+ (ir->vdw_modifier != eintmodNONE && ir->rvdw == ir->rlist)))
{
warning_note(wi, "With exact cut-offs, rlist should be "
"larger than rcoulomb and rvdw, so that there "
warning_note(wi, "You have selected user tables with dispersion correction, the dispersion will be corrected to -C6/r^6 beyond rvdw_switch (the tabulated interaction between rvdw_switch and rvdw will not be double counted). Make sure that you really want dispersion correction to -C6/r^6.");
}
- if (ir->nstlist == -1)
- {
- sprintf(err_buf, "With nstlist=-1 rvdw and rcoulomb should be smaller than rlist to account for diffusion and possibly charge-group radii");
- CHECK(ir->rvdw >= ir->rlist || ir->rcoulomb >= ir->rlist);
- }
- sprintf(err_buf, "nstlist can not be smaller than -1");
- CHECK(ir->nstlist < -1);
-
if (ir->eI == eiLBFGS && (ir->coulombtype == eelCUT || ir->vdwtype == evdwCUT)
&& ir->rvdw != 0)
{
int np = 0;
char *copy0, *copy;
- copy0 = strdup(str);
+ copy0 = gmx_strdup(str);
copy = copy0;
ltrim(copy);
while (*copy != '\0')
}
for (i = 0; i < ir->nwall; i++)
{
- opts->wall_atomtype[i] = strdup(names[i]);
+ opts->wall_atomtype[i] = gmx_strdup(names[i]);
}
if (ir->wall_type == ewt93 || ir->wall_type == ewt104)
/* COM pulling */
CCTYPE("COM PULLING");
- CTYPE("Pull type: no, umbrella, constraint or constant-force");
- EETYPE("pull", ir->ePull, epull_names);
- if (ir->ePull != epullNO)
+ EETYPE("pull", ir->bPull, yesno_names);
+ if (ir->bPull)
{
snew(ir->pull, 1);
- is->pull_grp = read_pullparams(&ninp, &inp, ir->pull, &opts->pull_start, wi);
+ is->pull_grp = read_pullparams(&ninp, &inp, ir->pull, wi);
}
/* Enforced rotation */
{
if (ir->efep != efepNO)
{
- opts->couple_moltype = strdup(is->couple_moltype);
+ opts->couple_moltype = gmx_strdup(is->couple_moltype);
if (opts->couple_lam0 == opts->couple_lam1)
{
warning(wi, "The lambda=0 and lambda=1 states for coupling are identical");
} /* search_QMstring */
/* We would like gn to be const as well, but C doesn't allow this */
+/* TODO this is utility functionality (search for the index of a
+ string in a collection), so should be refactored and located more
+ centrally. */
int search_string(const char *s, int ng, char *gn[])
{
int i;
}
}
- if (ir->ePull == epullCONSTRAINT)
+ if (ir->bPull)
{
/* Correct nrdf for the COM constraints.
* We correct using the TC and VCM group of the first atom
for (i = 0; i < pull->ncoord; i++)
{
+ if (pull->coord[i].eType != epullCONSTRAINT)
+ {
+ continue;
+ }
+
imin = 1;
for (j = 0; j < 2; j++)
double a, phi;
int i;
- t = strdup(s);
+ t = gmx_strdup(s);
trim(t);
cosine->n = 0;
void make_IMD_group(t_IMD *IMDgroup, char *IMDgname, t_blocka *grps, char **gnames)
{
- int ig = -1, i;
+ int ig, i;
ig = search_string(IMDgname, grps->nr, gnames);
}
}
- if (ir->ePull != epullNO)
+ if (ir->bPull)
{
make_pull_groups(ir->pull, is->pull_grp, grps, gnames);
case efbposresSPHERE:
AbsRef[XX] = AbsRef[YY] = AbsRef[ZZ] = 1;
break;
+ case efbposresCYLINDERX:
+ AbsRef[YY] = AbsRef[ZZ] = 1;
+ break;
+ case efbposresCYLINDERY:
+ AbsRef[XX] = AbsRef[ZZ] = 1;
+ break;
case efbposresCYLINDER:
+ /* efbposres is a synonym for efbposresCYLINDERZ for backwards compatibility */
+ case efbposresCYLINDERZ:
AbsRef[XX] = AbsRef[YY] = 1;
break;
case efbposresX: /* d=XX */
gmx_fatal(FARGS, "Soft-core interactions are only supported with VdW repulsion power 12");
}
- if (ir->ePull != epullNO)
+ if (ir->bPull)
{
- gmx_bool bPullAbsoluteRef;
+ gmx_bool bWarned;
- bPullAbsoluteRef = FALSE;
- for (i = 0; i < ir->pull->ncoord; i++)
+ bWarned = FALSE;
+ for (i = 0; i < ir->pull->ncoord && !bWarned; i++)
{
- bPullAbsoluteRef = bPullAbsoluteRef ||
- ir->pull->coord[i].group[0] == 0 ||
- ir->pull->coord[i].group[1] == 0;
- }
- if (bPullAbsoluteRef)
- {
- absolute_reference(ir, sys, FALSE, AbsRef);
- for (m = 0; m < DIM; m++)
+ if (ir->pull->coord[i].group[0] == 0 ||
+ ir->pull->coord[i].group[1] == 0)
{
- if (ir->pull->dim[m] && !AbsRef[m])
+ absolute_reference(ir, sys, FALSE, AbsRef);
+ for (m = 0; m < DIM; m++)
{
- warning(wi, "You are using an absolute reference for pulling, but the rest of the system does not have an absolute reference. This will lead to artifacts.");
- break;
+ if (ir->pull->coord[i].dim[m] && !AbsRef[m])
+ {
+ warning(wi, "You are using an absolute reference for pulling, but the rest of the system does not have an absolute reference. This will lead to artifacts.");
+ bWarned = TRUE;
+ break;
+ }
}
}
}
- if (ir->pull->eGeom == epullgDIRPBC)
+ for (i = 0; i < 3; i++)
{
- for (i = 0; i < 3; i++)
+ for (m = 0; m <= i; m++)
{
- for (m = 0; m <= i; m++)
+ if ((ir->epc != epcNO && ir->compress[i][m] != 0) ||
+ ir->deform[i][m] != 0)
{
- if ((ir->epc != epcNO && ir->compress[i][m] != 0) ||
- ir->deform[i][m] != 0)
+ for (c = 0; c < ir->pull->ncoord; c++)
{
- for (c = 0; c < ir->pull->ncoord; c++)
+ if (ir->pull->coord[c].eGeom == epullgDIRPBC &&
+ ir->pull->coord[c].vec[m] != 0)
{
- if (ir->pull->coord[c].vec[m] != 0)
- {
- gmx_fatal(FARGS, "Can not have dynamic box while using pull geometry '%s' (dim %c)", EPULLGEOM(ir->pull->eGeom), 'x'+m);
- }
+ gmx_fatal(FARGS, "Can not have dynamic box while using pull geometry '%s' (dim %c)", EPULLGEOM(ir->pull->coord[c].eGeom), 'x'+m);
}
}
}
check_disre(sys);
}
-void double_check(t_inputrec *ir, matrix box, gmx_bool bConstr, warninp_t wi)
+void double_check(t_inputrec *ir, matrix box,
+ gmx_bool bHasNormalConstraints,
+ gmx_bool bHasAnyConstraints,
+ warninp_t wi)
{
real min_size;
gmx_bool bTWIN;
warning_error(wi, ptr);
}
- if (bConstr && ir->eConstrAlg == econtSHAKE)
+ if (bHasNormalConstraints && ir->eConstrAlg == econtSHAKE)
{
if (ir->shake_tol <= 0.0)
{
}
}
- if ( (ir->eConstrAlg == econtLINCS) && bConstr)
+ if ( (ir->eConstrAlg == econtLINCS) && bHasNormalConstraints)
{
/* If we have Lincs constraints: */
if (ir->eI == eiMD && ir->etc == etcNO &&
}
}
- if (bConstr && ir->epc == epcMTTK)
+ if (bHasAnyConstraints && ir->epc == epcMTTK)
{
- warning_note(wi, "MTTK with constraints is deprecated, and will be removed in GROMACS 5.1");
+ warning_error(wi, "Constraints are not implemented with MTTK pressure control.");
}
if (ir->LincsWarnAngle > 90.0)