-/* -*- 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:
- * 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 "gmxpre.h"
+#include "toppush.h"
+
+#include <assert.h>
#include <ctype.h>
#include <math.h>
-
-#include "sysstuff.h"
-#include "smalloc.h"
-#include "macros.h"
-#include "string2.h"
-#include "names.h"
-#include "toputil.h"
-#include "toppush.h"
-#include "topdirs.h"
-#include "readir.h"
-#include "symtab.h"
-#include "gmx_fatal.h"
-#include "warninp.h"
-#include "gpp_atomtype.h"
-#include "gpp_bond_atomtype.h"
+#include <stdlib.h>
+
+#include "gromacs/gmxpreprocess/gpp_atomtype.h"
+#include "gromacs/gmxpreprocess/gpp_bond_atomtype.h"
+#include "gromacs/gmxpreprocess/readir.h"
+#include "gromacs/gmxpreprocess/topdirs.h"
+#include "gromacs/gmxpreprocess/toputil.h"
+#include "gromacs/legacyheaders/macros.h"
+#include "gromacs/legacyheaders/names.h"
+#include "gromacs/legacyheaders/warninp.h"
+#include "gromacs/topology/symtab.h"
+#include "gromacs/utility/cstringutil.h"
+#include "gromacs/utility/fatalerror.h"
+#include "gromacs/utility/smalloc.h"
void generate_nbparams(int comb, int ftype, t_params *plist, gpp_atomtype_t atype,
warninp_t wi)
break;
case eCOMB_ARITHMETIC:
- /* c0 and c1 are epsilon and sigma */
+ /* c0 and c1 are sigma and epsilon */
for (i = k = 0; (i < nr); i++)
{
for (j = 0; (j < nr); j++, k++)
cj0 = get_atomtype_nbparam(j, 0, atype);
ci1 = get_atomtype_nbparam(i, 1, atype);
cj1 = get_atomtype_nbparam(j, 1, atype);
- plist->param[k].c[0] = (ci0+cj0)*0.5;
+ plist->param[k].c[0] = (fabs(ci0) + fabs(cj0))*0.5;
+ /* Negative sigma signals that c6 should be set to zero later,
+ * so we need to propagate that through the combination rules.
+ */
+ if (ci0 < 0 || cj0 < 0)
+ {
+ plist->param[k].c[0] *= -1;
+ }
plist->param[k].c[1] = sqrt(ci1*cj1);
}
}
break;
case eCOMB_GEOM_SIG_EPS:
- /* c0 and c1 are epsilon and sigma */
+ /* c0 and c1 are sigma and epsilon */
for (i = k = 0; (i < nr); i++)
{
for (j = 0; (j < nr); j++, k++)
cj0 = get_atomtype_nbparam(j, 0, atype);
ci1 = get_atomtype_nbparam(i, 1, atype);
cj1 = get_atomtype_nbparam(j, 1, atype);
- plist->param[k].c[0] = sqrt(ci0*cj0);
+ plist->param[k].c[0] = sqrt(fabs(ci0*cj0));
+ /* Negative sigma signals that c6 should be set to zero later,
+ * so we need to propagate that through the combination rules.
+ */
+ if (ci0 < 0 || cj0 < 0)
+ {
+ plist->param[k].c[0] *= -1;
+ }
plist->param[k].c[1] = sqrt(ci1*cj1);
}
}
{
sprintf(errbuf, "Overriding %s parameters.%s",
interaction_function[ftype].longname,
- (ftype == F_PDIHS) ? "\nUse dihedraltype 4 to allow several multiplicity terms." : "");
+ (ftype == F_PDIHS) ?
+ "\nUse dihedraltype 9 to allow several multiplicity terms. Only consecutive lines are combined. Non-consective lines overwrite each other."
+ : "");
warning(wi, errbuf);
- fprintf(stderr, " old:");
+ fprintf(stderr, " old: ");
for (j = 0; (j < nrfp); j++)
{
fprintf(stderr, " %g", bt->param[i].c[j]);
/* When the B topology parameters are not set,
* copy them from topology A
*/
+ assert(nrfp <= 4);
for (i = n; i < nrfp; i++)
{
c[i] = c[i-2];
nrfp = nrfpA+nrfpB;
/* Allocate memory for the CMAP grid */
- bt->ncmap += nrfp;
- srenew(bt->cmap, bt->ncmap);
+ bt[F_CMAP].ncmap += nrfp;
+ srenew(bt[F_CMAP].cmap, bt[F_CMAP].ncmap);
/* Read in CMAP parameters */
sl = 0;
}
nn = sscanf(line+start+sl, " %s ", s);
sl += strlen(s);
- bt->cmap[i+(bt->ncmap)-nrfp] = strtod(s, NULL);
+ bt[F_CMAP].cmap[i+(bt[F_CMAP].ncmap)-nrfp] = strtod(s, NULL);
if (nn == 1)
{
{
for (i = 0; i < ncmap; i++)
{
- bt->cmap[i+ncmap] = bt->cmap[i];
+ bt[F_CMAP].cmap[i+ncmap] = bt[F_CMAP].cmap[i];
}
}
else
/* Set grid spacing and the number of grids (we assume these numbers to be the same for all grids
* so we can safely assign them each time
*/
- bt->grid_spacing = nxcmap; /* Or nycmap, they need to be equal */
- bt->nc = bt->nc + 1; /* Since we are incrementing here, we need to subtract later, see (*****) */
- nct = (nral+1) * bt->nc;
+ bt[F_CMAP].grid_spacing = nxcmap; /* Or nycmap, they need to be equal */
+ bt[F_CMAP].nc = bt[F_CMAP].nc + 1; /* Since we are incrementing here, we need to subtract later, see (*****) */
+ nct = (nral+1) * bt[F_CMAP].nc;
/* Allocate memory for the cmap_types information */
- srenew(bt->cmap_types, nct);
+ srenew(bt[F_CMAP].cmap_types, nct);
for (i = 0; (i < nral); i++)
{
}
/* Assign a grid number to each cmap_type */
- bt->cmap_types[bt->nct++] = get_bond_atomtype_type(alc[i], bat);
+ bt[F_CMAP].cmap_types[bt[F_CMAP].nct++] = get_bond_atomtype_type(alc[i], bat);
}
/* Assign a type number to this cmap */
- bt->cmap_types[bt->nct++] = bt->nc-1; /* Since we inremented earlier, we need to subtrac here, to get the types right (****) */
+ bt[F_CMAP].cmap_types[bt[F_CMAP].nct++] = bt[F_CMAP].nc-1; /* Since we inremented earlier, we need to subtrac here, to get the types right (****) */
/* Check for the correct number of atoms (again) */
- if (bt->nct != nct)
+ if (bt[F_CMAP].nct != nct)
{
- gmx_fatal(FARGS, "Incorrect number of atom types (%d) in cmap type %d\n", nct, bt->nc);
+ gmx_fatal(FARGS, "Incorrect number of atom types (%d) in cmap type %d\n", nct, bt[F_CMAP].nc);
}
/* Is this correct?? */
static void push_atom_now(t_symtab *symtab, t_atoms *at, int atomnr,
int atomicnumber,
int type, char *ctype, int ptype,
- char *resnumberic, int cgnumber,
+ char *resnumberic,
char *resname, char *name, real m0, real q0,
int typeB, char *ctypeB, real mB, real qB)
{
push_cg(cgs, lastcg, cgnumber, nr);
push_atom_now(symtab, at, atomnr, get_atomtype_atomnumber(type, atype),
- type, ctype, ptype, resnumberic, cgnumber,
+ type, ctype, ptype, resnumberic,
resname, name, m0, q0, typeB,
typeB == type ? ctype : ctypeB, mB, qB);
}
return bFound;
}
-static gmx_bool default_cmap_params(int ftype, t_params bondtype[],
+static gmx_bool default_cmap_params(t_params bondtype[],
t_atoms *at, gpp_atomtype_t atype,
t_param *p, gmx_bool bB,
int *cmap_type, int *nparam_def)
ct = 0;
/* Match the current cmap angle against the list of cmap_types */
- for (i = 0; i < bondtype->nct && !bFound; i += 6)
+ for (i = 0; i < bondtype[F_CMAP].nct && !bFound; i += 6)
{
if (bB)
{
else
{
if (
- (get_atomtype_batype(at->atom[p->a[0]].type, atype) == bondtype->cmap_types[i]) &&
- (get_atomtype_batype(at->atom[p->a[1]].type, atype) == bondtype->cmap_types[i+1]) &&
- (get_atomtype_batype(at->atom[p->a[2]].type, atype) == bondtype->cmap_types[i+2]) &&
- (get_atomtype_batype(at->atom[p->a[3]].type, atype) == bondtype->cmap_types[i+3]) &&
- (get_atomtype_batype(at->atom[p->a[4]].type, atype) == bondtype->cmap_types[i+4]))
+ (get_atomtype_batype(at->atom[p->a[0]].type, atype) == bondtype[F_CMAP].cmap_types[i]) &&
+ (get_atomtype_batype(at->atom[p->a[1]].type, atype) == bondtype[F_CMAP].cmap_types[i+1]) &&
+ (get_atomtype_batype(at->atom[p->a[2]].type, atype) == bondtype[F_CMAP].cmap_types[i+2]) &&
+ (get_atomtype_batype(at->atom[p->a[3]].type, atype) == bondtype[F_CMAP].cmap_types[i+3]) &&
+ (get_atomtype_batype(at->atom[p->a[4]].type, atype) == bondtype[F_CMAP].cmap_types[i+4]))
{
/* Found cmap torsion */
bFound = TRUE;
- ct = bondtype->cmap_types[i+5];
+ ct = bondtype[F_CMAP].cmap_types[i+5];
nparam_found = 1;
}
}
bFoundA = default_params(ftype, bondtype, at, atype, ¶m, FALSE, ¶m_defA, &nparam_defA);
if (bFoundA)
{
- /* Copy the A-state and B-state default parameters */
+ /* Copy the A-state and B-state default parameters. */
+ assert(NRFPA(ftype)+NRFPB(ftype) <= MAXFORCEPARAM);
for (j = 0; (j < NRFPA(ftype)+NRFPB(ftype)); j++)
{
param.c[j] = param_defA->c[j];
void push_cmap(directive d, t_params bondtype[], t_params bond[],
t_atoms *at, gpp_atomtype_t atype, char *line,
- gmx_bool *bWarn_copy_A_B,
warninp_t wi)
{
const char *aaformat[MAXATOMLIST+1] =
}
/* Get the cmap type for this cmap angle */
- bFound = default_cmap_params(ftype, bondtype, at, atype, ¶m, FALSE, &cmap_type, &ncmap_params);
+ bFound = default_cmap_params(bondtype, at, atype, ¶m, FALSE, &cmap_type, &ncmap_params);
/* We want exactly one parameter (the cmap type in state A (currently no state B) back */
if (bFound && ncmap_params == 1)
-void push_vsitesn(directive d, t_params bondtype[], t_params bond[],
- t_atoms *at, gpp_atomtype_t atype, char *line,
+void push_vsitesn(directive d, t_params bond[],
+ t_atoms *at, char *line,
warninp_t wi)
{
char *ptr;
static void convert_pairs_to_pairsQ(t_params *plist,
real fudgeQQ, t_atoms *atoms)
{
- t_param *param;
- int i;
- real v, w;
+ t_param *paramp1, *paramp2, *paramnew;
+ int i, j, p1nr, p2nr, p2newnr;
+
+ /* Add the pair list to the pairQ list */
+ p1nr = plist[F_LJ14].nr;
+ p2nr = plist[F_LJC14_Q].nr;
+ p2newnr = p1nr + p2nr;
+ snew(paramnew, p2newnr);
+
+ paramp1 = plist[F_LJ14].param;
+ paramp2 = plist[F_LJC14_Q].param;
+
+ /* Fill in the new F_LJC14_Q array with the old one. NOTE:
+ it may be possible to just ADD the converted F_LJ14 array
+ to the old F_LJC14_Q array, but since we have to create
+ a new sized memory structure, better just to deep copy it all.
+ */
- /* Copy the pair list to the pairQ list */
- plist[F_LJC14_Q] = plist[F_LJ14];
- /* Empty the pair list */
- plist[F_LJ14].nr = 0;
- plist[F_LJ14].param = NULL;
- param = plist[F_LJC14_Q].param;
- for (i = 0; i < plist[F_LJC14_Q].nr; i++)
+ for (i = 0; i < p2nr; i++)
{
- v = param[i].c[0];
- w = param[i].c[1];
- param[i].c[0] = fudgeQQ;
- param[i].c[1] = atoms->atom[param[i].a[0]].q;
- param[i].c[2] = atoms->atom[param[i].a[1]].q;
- param[i].c[3] = v;
- param[i].c[4] = w;
+ /* Copy over parameters */
+ for (j = 0; j < 5; j++) /* entries are 0-4 for F_LJC14_Q */
+ {
+ paramnew[i].c[j] = paramp2[i].c[j];
+ }
+
+ /* copy over atoms */
+ for (j = 0; j < 2; j++)
+ {
+ paramnew[i].a[j] = paramp2[i].a[j];
+ }
}
+
+ for (i = p2nr; i < p2newnr; i++)
+ {
+ j = i-p2nr;
+
+ /* Copy over parameters */
+ paramnew[i].c[0] = fudgeQQ;
+ paramnew[i].c[1] = atoms->atom[paramp1[j].a[0]].q;
+ paramnew[i].c[2] = atoms->atom[paramp1[j].a[1]].q;
+ paramnew[i].c[3] = paramp1[j].c[0];
+ paramnew[i].c[4] = paramp1[j].c[1];
+
+ /* copy over atoms */
+ paramnew[i].a[0] = paramp1[j].a[0];
+ paramnew[i].a[1] = paramp1[j].a[1];
+ }
+
+ /* free the old pairlists */
+ sfree(plist[F_LJC14_Q].param);
+ sfree(plist[F_LJ14].param);
+
+ /* now assign the new data to the F_LJC14_Q structure */
+ plist[F_LJC14_Q].param = paramnew;
+ plist[F_LJC14_Q].nr = p2newnr;
+
+ /* Empty the LJ14 pairlist */
+ plist[F_LJ14].nr = 0;
+ plist[F_LJ14].param = NULL;
}
static void generate_LJCpairsNB(t_molinfo *mol, int nb_funct, t_params *nbp)
biod.pnpi.spb.ru / alexxy / gromacs.git / blobdiff