-/* -*- mode: c; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4; c-file-style: "stroustrup"; -*-
+/*
+ * 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) 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:
- * GROwing Monsters And Cloning Shrimps
+ * 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 "confio.h"
-#include "constr.h"
-#include "copyrite.h"
-#include "invblock.h"
-#include "main.h"
-#include "mdrun.h"
-#include "nrnb.h"
-#include "smalloc.h"
-#include "vec.h"
-#include "physics.h"
-#include "names.h"
-#include "txtdump.h"
-#include "domdec.h"
-#include "pdbio.h"
-#include "partdec.h"
-#include "splitter.h"
-#include "mtop_util.h"
-#include "gmxfio.h"
-#include "macros.h"
-#include "gmx_omp_nthreads.h"
+#include "gmxpre.h"
+
+#include "gromacs/legacyheaders/constr.h"
+
+#include <assert.h>
+#include <stdlib.h>
+
+#include "gromacs/essentialdynamics/edsam.h"
+#include "gromacs/fileio/confio.h"
+#include "gromacs/fileio/gmxfio.h"
+#include "gromacs/fileio/pdbio.h"
+#include "gromacs/legacyheaders/copyrite.h"
+#include "gromacs/legacyheaders/domdec.h"
+#include "gromacs/legacyheaders/gmx_omp_nthreads.h"
+#include "gromacs/legacyheaders/macros.h"
+#include "gromacs/legacyheaders/mdrun.h"
+#include "gromacs/legacyheaders/names.h"
+#include "gromacs/legacyheaders/nrnb.h"
+#include "gromacs/legacyheaders/splitter.h"
+#include "gromacs/legacyheaders/txtdump.h"
+#include "gromacs/legacyheaders/types/commrec.h"
+#include "gromacs/math/vec.h"
+#include "gromacs/pbcutil/pbc.h"
+#include "gromacs/pulling/pull.h"
+#include "gromacs/topology/block.h"
+#include "gromacs/topology/invblock.h"
+#include "gromacs/topology/mtop_util.h"
+#include "gromacs/utility/fatalerror.h"
+#include "gromacs/utility/smalloc.h"
typedef struct gmx_constr {
int ncon_tot; /* The total number of constraints */
int nblocks; /* The number of SHAKE blocks */
int *sblock; /* The SHAKE blocks */
int sblock_nalloc; /* The allocation size of sblock */
- real *lagr; /* Lagrange multipliers for SHAKE */
+ real *lagr; /* -2 times the Lagrange multipliers for SHAKE */
int lagr_nalloc; /* The allocation size of lagr */
int maxwarn; /* The maximum number of warnings */
int warncount_lincs;
atom_id blocknr;
} t_sortblock;
+/* delta_t should be used instead of ir->delta_t, to permit the time
+ step to be scaled by the calling code */
static void *init_vetavars(t_vetavars *vars,
gmx_bool constr_deriv,
- real veta, real vetanew, t_inputrec *ir, gmx_ekindata_t *ekind, gmx_bool bPscal)
+ real veta, real vetanew,
+ t_inputrec *ir, real delta_t,
+ gmx_ekindata_t *ekind, gmx_bool bPscal)
{
double g;
int i;
{
vars->alpha = 1.0;
}
- g = 0.5*veta*ir->delta_t;
+ g = 0.5*veta*delta_t;
vars->rscale = exp(g)*series_sinhx(g);
- g = -0.25*vars->alpha*veta*ir->delta_t;
+ g = -0.25*vars->alpha*veta*delta_t;
vars->vscale = exp(g)*series_sinhx(g);
vars->rvscale = vars->vscale*vars->rscale;
vars->veta = vetanew;
return nflexcon;
}
+static void clear_constraint_quantity_nonlocal(gmx_domdec_t *dd, rvec *q)
+{
+ int nonlocal_at_start, nonlocal_at_end, at;
+
+ dd_get_constraint_range(dd, &nonlocal_at_start, &nonlocal_at_end);
+
+ for (at = nonlocal_at_start; at < nonlocal_at_end; at++)
+ {
+ clear_rvec(q[at]);
+ }
+}
+
void too_many_constraint_warnings(int eConstrAlg, int warncount)
{
const char *abort = "- aborting to avoid logfile runaway.\n"
int start, int homenr, t_commrec *cr,
rvec x[], matrix box)
{
- char fname[STRLEN], format[STRLEN];
+ char fname[STRLEN];
FILE *out;
int dd_ac0 = 0, dd_ac1 = 0, i, ii, resnr;
gmx_domdec_t *dd;
char *anm, *resnm;
dd = NULL;
+ if (DOMAINDECOMP(cr))
+ {
+ dd = cr->dd;
+ dd_get_constraint_range(dd, &dd_ac0, &dd_ac1);
+ start = 0;
+ homenr = dd_ac1;
+ }
+
if (PAR(cr))
{
sprintf(fname, "%s_n%d.pdb", fn, cr->sim_nodeid);
- if (DOMAINDECOMP(cr))
- {
- dd = cr->dd;
- dd_get_constraint_range(dd, &dd_ac0, &dd_ac1);
- start = 0;
- homenr = dd_ac1;
- }
}
else
{
sprintf(fname, "%s.pdb", fn);
}
- sprintf(format, "%s\n", get_pdbformat());
out = gmx_fio_fopen(fname, "w");
ii = i;
}
gmx_mtop_atominfo_global(mtop, ii, &anm, &resnr, &resnm);
- fprintf(out, format, "ATOM", (ii+1)%100000,
- anm, resnm, ' ', resnr%10000, ' ',
- 10*x[i][XX], 10*x[i][YY], 10*x[i][ZZ]);
+ gmx_fprintf_pdb_atomline(out, epdbATOM, ii+1, anm, ' ', resnm, ' ', resnr, ' ',
+ 10*x[i][XX], 10*x[i][YY], 10*x[i][ZZ], 1.0, 0.0, "");
}
fprintf(out, "TER\n");
gmx_fio_fclose(out);
}
-static void dump_confs(FILE *fplog, gmx_large_int_t step, gmx_mtop_t *mtop,
+static void dump_confs(FILE *fplog, gmx_int64_t step, gmx_mtop_t *mtop,
int start, int homenr, t_commrec *cr,
rvec x[], rvec xprime[], matrix box)
{
struct gmx_constr *constr,
t_idef *idef, t_inputrec *ir, gmx_ekindata_t *ekind,
t_commrec *cr,
- gmx_large_int_t step, int delta_step,
+ gmx_int64_t step, int delta_step,
+ real step_scaling,
t_mdatoms *md,
rvec *x, rvec *xprime, rvec *min_proj,
gmx_bool bMolPBC, matrix box,
int ncons, settle_error;
tensor vir_r_m_dr;
rvec *vstor;
+ real scaled_delta_t;
real invdt, vir_fac, t;
t_ilist *settle;
int nsettle;
bOK = TRUE;
bDump = FALSE;
- start = md->start;
+ start = 0;
homenr = md->homenr;
nrend = start+homenr;
+ scaled_delta_t = step_scaling * ir->delta_t;
+
/* set constants for pressure control integration */
init_vetavars(&vetavar, econq != econqCoord,
- veta, vetanew, ir, ekind, bPscal);
+ veta, vetanew, ir, scaled_delta_t, ekind, bPscal);
+ /* Prepare time step for use in constraint implementations, and
+ avoid generating inf when ir->delta_t = 0. */
if (ir->delta_t == 0)
{
- invdt = 0;
+ invdt = 0.0;
}
else
{
- invdt = 1/ir->delta_t;
+ invdt = 1.0/scaled_delta_t;
}
if (ir->efep != efepNO && EI_DYNAMICS(ir->eI))
*/
if (cr->dd)
{
- dd_move_x_constraints(cr->dd, box, x, xprime);
- }
- else if (PARTDECOMP(cr))
- {
- pd_move_x_constraints(cr, x, xprime);
+ dd_move_x_constraints(cr->dd, box, x, xprime, econq == econqCoord);
+
+ if (v != NULL)
+ {
+ /* We need to initialize the non-local components of v.
+ * We never actually use these values, but we do increment them,
+ * so we should avoid uninitialized variables and overflows.
+ */
+ clear_constraint_quantity_nonlocal(cr->dd, v);
+ }
}
if (constr->lincsd != NULL)
}
else
{
- calcvir_atom_end = md->start + md->homenr;
+ calcvir_atom_end = md->homenr;
}
switch (econq)
{
char buf[256];
sprintf(buf,
- "\nstep " gmx_large_int_pfmt ": Water molecule starting at atom %d can not be "
+ "\nstep " "%"GMX_PRId64 ": Water molecule starting at atom %d can not be "
"settled.\nCheck for bad contacts and/or reduce the timestep if appropriate.\n",
step, ddglatnr(cr->dd, settle->iatoms[settle_error*(1+NRAL(F_SETTLE))+1]));
if (fplog)
if (vir != NULL)
{
+ /* The normal uses of constrain() pass step_scaling = 1.0.
+ * The call to constrain() for SD1 that passes step_scaling =
+ * 0.5 also passes vir = NULL, so cannot reach this
+ * assertion. This assertion should remain until someone knows
+ * that this path works for their intended purpose, and then
+ * they can use scaled_delta_t instead of ir->delta_t
+ * below. */
+ assert(gmx_within_tol(step_scaling, 1.0, GMX_REAL_EPS));
switch (econq)
{
case econqCoord:
}
}
-static void make_shake_sblock_pd(struct gmx_constr *constr,
- t_idef *idef, t_mdatoms *md)
+static void make_shake_sblock_serial(struct gmx_constr *constr,
+ t_idef *idef, t_mdatoms *md)
{
int i, j, m, ncons;
int bstart, bnr;
ncons = idef->il[F_CONSTR].nr/3;
init_blocka(&sblocks);
- gen_sblocks(NULL, md->start, md->start+md->homenr, idef, &sblocks, FALSE);
+ gen_sblocks(NULL, 0, md->homenr, idef, &sblocks, FALSE);
/*
bstart=(idef->nodeid > 0) ? blocks->multinr[idef->nodeid-1] : 0;
j, constr->nblocks, ncons);
for (i = 0; (i < ncons); i++)
{
- fprintf(stderr, "i: %5d sb[i].blocknr: %5u\n", i, sb[i].blocknr);
+ fprintf(stderr, "i: %5d sb[i].blocknr: %5d\n", i, sb[i].blocknr);
}
for (j = 0; (j <= constr->nblocks); j++)
{
}
else
{
- make_shake_sblock_pd(constr, idef, md);
+ make_shake_sblock_serial(constr, idef, md);
}
if (ncons > constr->lagr_nalloc)
{
biod.pnpi.spb.ru / alexxy / gromacs.git / blobdiff