NAMD
GromacsTopFile.h
Go to the documentation of this file.
1 #ifndef GROMACSTOPFILE_H
2 #define GROMACSTOPFILE_H
3 
4 typedef float Real;
5 
6 #define NAMESIZE 5
7 #define LONGNAMESIZE 20
8 
9 /* A GenericBond represents a bond between two atoms in a template
10  that can be used one or more times within a simulation. This is an
11  immutable type.
12  Here, atomi and atomj start at 0. */
13 class GenericBond {
14  private:
15  int atomi; /* the first atom of the bond */
16  int atomj; /* the second atom of the bond */
17  int type; /* an index into the bond's parameters, stored elsewhere */
18 
19  public:
20  /* initializes this to be a bond between atoms <i> and <j> of type
21  <type> */
22  GenericBond(int i, int j, int theType);
23 
24  int getAtomi() const {return atomi;}
25  int getAtomj() const {return atomj;}
26  int getType() const {return type;}
27 
28 };
29 
30 /* A GenericAngle represents a angle between two atoms in a template
31  that can be used one or more times within a simulation. This is an
32  immutable type.
33  Here, atomi and atomj start at 0. */
34 class GenericAngle {
35  private:
36  int atomi; /* the first atom of the angle */
37  int atomj; /* the second atom of the angle */
38  int atomk; /* the third atom of the angle */
39  int type; /* an index into the angle's parameters, stored elsewhere */
40 
41  public:
42  /* initializes this to be a angle between atoms <i>, <j>, and <k> of
43  type <type> */
44  GenericAngle(int i, int j, int k, int theType);
45 
46  int getAtomi() const {return atomi;}
47  int getAtomj() const {return atomj;}
48  int getAtomk() const {return atomk;}
49  int getType() const {return type;}
50 
51 };
52 
53 /* A GenericDihedral represents a dihedral angle between four atoms in
54  a template that can be used one or more times within a simulation.
55  This is an immutable type.
56  Here, atomi and atomj start at 0. */
58  private:
59  int atomi; /* the first atom of the angle */
60  int atomj; /* the second atom of the angle */
61  int atomk; /* the third atom of the angle */
62  int atoml; /* the fourth atom of the angle */
63  int type; /* an index into the parameters, stored elsewhere */
64 
65  public:
66  /* initializes this to be a angle between atoms <i>, <j>, <k>, and
67  <l> of type <type> */
68  GenericDihedral(int i, int j, int k, int l, int theType);
69 
70  int getAtomi() const {return atomi;}
71  int getAtomj() const {return atomj;}
72  int getAtomk() const {return atomk;}
73  int getAtoml() const {return atoml;}
74  int getType() const {return type;}
75 
76 };
77 
78 /* A GenericAtom represents an atom in a template that can be used one
79  or more times within a simulation. This is an immutable type. */
80 class GenericAtom {
81  private:
82  char type[NAMESIZE+1];
83  int typeNum;
84  int resNum;
85  char resType[NAMESIZE+1];
86  char atomName[NAMESIZE+1];
87  Real charge;
88  Real mass;
89 
90  public:
91 
92  /* Initializes this to be the atom specified by the given
93  parameters. All the strings will be copied, and must be at most
94  four characters long. */
95  GenericAtom(const char *theType, int theTypeNum, int theResNum,
96  const char *theResType,
97  const char *theAtomName, Real theCharge, Real theMass);
98 
99  /* The following just return the parameters of the atom. */
100  const char *getType() const { return type; }
101  int getTypeNum() const { return typeNum; }
102  int getResNum() const {return resNum;}
103  const char *getResType() const {return resType;}
104  const char *getAtomName() const {return atomName;}
105  Real getCharge() const {return charge;}
106  Real getMass() const {return mass;}
107 
108 };
109 
110 /* A GenericMol represents a complete molecule template. This will be
111  a mutable type, for convenience's sake. The atoms here are
112  numbered starting at ZERO */
113 class GenericMol {
114  private:
115  /* the name of this Molecule, so it can be identified later */
116  const char *name;
117 
118  /* the rest of the list */
119  const GenericMol *next;
120 
121  /* the list of GenericAtoms */
123 
124  /* the list of GenericBonds */
126 
127  /* the list of GenericAngles */
129 
130  /* the list of GenericDihedrals */
132 
133  public:
134  /* initializes this to be the molecule with name <name> */
135  GenericMol(const char *theName);
136 
137  /* getNumAtoms returns the number of atoms stored here. */
138  int getNumAtoms() const { return atomList.size(); }
139  int getNumBonds() const { return bondList.size(); }
140  int getNumAngles() const { return angleList.size(); }
141  int getNumDihedrals() const { return dihedralList.size(); }
142  int getNumRes() const { return atomList[atomList.size()-1]->getResNum(); }
143 
144  /* adds an atom to the end of list. The residue numbers must start
145  at zero and count up in unit steps */
146  void addAtom(const char *theType, int theTypeNum, int theResNum,
147  const char *theResType,
148  const char *theAtomName, Real theCharge, Real theMass);
149 
150  /* adds a bond/angle/dihedral to the molecule */
151  void addBond(int atomi, int atomj, int type);
152  void addAngle(int atomi, int atomj, int atomk, int type);
153  void addDihedral(int atomi, int atomj, int atomk, int atoml, int type);
154 
155  /* gets the name */
156  const char *getName() const { return name; }
157 
158  /* gets atom <n>, where atom 0 is the first atom in the list. */
159  const GenericAtom *getAtom(int n) const;
160 
161  /* gets bond <n>, where bond 0 is the first bond in the list. */
162  const GenericBond *getBond(int n) const;
163 
164  /* gets angle <n>, where angle 0 is the first angle in the list. */
165  const GenericAngle *getAngle(int n) const;
166 
167  /* gets dihedral <n>, where angle 0 is the first one in the list. */
168  const GenericDihedral *getDihedral(int n) const;
169 
170 };
171 
172 /* A MolInst represents a number of instances of a molecule. These
173  objects are meant to be assembled (in order) into a list, to
174  represent all the instances of molecules in a system. */
175 class MolInst {
176  int num;
177  const GenericMol *mol;
178  public:
179  /* initializes this to represent <theNum> copies of <theMol> */
180  MolInst(const GenericMol *theMol, int theNum);
181 
182  /* get the number */
183  int getNum() const { return num; }
184 
185  /* get the total number of X here */
186  int getNumAtoms() const;
187  int getNumRes() const;
188  int getNumBonds() const;
189  int getNumAngles() const;
190  int getNumDihedrals() const;
191 
192  /* get the molecule */
193  const GenericMol *getMol() const { return mol; }
194 };
195 
196 /* An AtomTable stores all the known types of atoms.
197  mass: the mass in amu
198  charge: the charge in e
199  c6: the LJ parameter c6 in kcal/mol A^6
200  c12: the LJ parameter c12 in kcal/mol A^12 */
201 class AtomTable {
202  private:
203  /* implementation of the extensible array */
204  ResizeArray<Real> mArray;
205  ResizeArray<Real> qArray;
206  ResizeArray<Real> c6Array;
207  ResizeArray<Real> c12Array;
208  ResizeArray<const char *> typeArray; /* should never be NULL */
209 
210  public:
211  /* this adds a entry for an atom type to the table */
212  void addType(const char *type, Real m, Real q, Real c6, Real c12);
213 
214  /* This looks up the atom by type - if it cannot be found, this
215  returns -1, otherwise this returns the index of the atomtype */
216  int getParams(const char *type, Real *m, Real *q,
217  Real *c6, Real *c12) const;
218 
219  /* This looks up the atom type by number, returning it in the string
220  <type>, which must be at least 11 characters long. */
221  void getType(int num, char *type) const;
222 
223 
224  /* returns the number of types in the table */
225  int size() const { return mArray.size(); }
226 };
227 
228 /* A BondTable stores all the different types of bonds, so that we
229  only have to have one copy of each.
230  b0: the natural bond length in A.
231  kB: the spring constant in kcal/mol/A^2, where E=kx^2 to 1st order.
232  funct: 1 for harmonic potential, 2 for fourth-order GROMACS96
233  potential. */
234 class BondTable {
235  private:
236  /* implementation of the extensible array */
237  ResizeArray<Real> b0Array;
238  ResizeArray<Real> kBArray;
239  ResizeArray<int> functArray;
240  ResizeArray<const char *> typeaArray; /* NULL if no type is set */
241  ResizeArray<const char *> typebArray; /* NULL if no type is set */
242 
243  public:
244  /* returns the number of types in the table */
245  int size() const { return b0Array.size(); }
246 
247  /* this gets the index of a bond in the table (adding an entry if
248  none exists). */
249  int getIndex(Real b0, Real kB, int funct);
250 
251  /* this adds a entry for a bond type to the table, including the two
252  atoms involved in the bond */
253  void addType(const char *typea, const char *typeb, Real b0,
254  Real kB, int funct);
255 
256  /* getParams puts the parameters for bond-type <num> into the
257  spaces pointed to by the other arguments. */
258  void getParams(int num, Real *b0, Real *kB, int *funct) const;
259 
260  /* This version looks up the bond by atom type - the order of the
261  types doesn't matter! If the specified bond/function combination
262  cannot be found, it returns -1, otherwise it returns the index of
263  the bondtype. */
264  int getParams(const char *typea, const char *typeb, int funct,
265  Real *b0, Real *kB) const;
266 };
267 
268 /* A AngleTable stores all the different types of angles, so that we
269  only have to have one copy of each. Units:
270  th0 - degrees
271  kth - kcal/mol/rad^2 (funct 1) OR kcal/mol (funct 2)
272  funct 1: harmonic angle potential
273  funct 2: harmonic G96 cosine potential
274  The atoms are in standard geometry order - A--B--C for angle /ABC
275 */
276 class AngleTable {
277  private:
278  /* implementation of the extensible array */
279  ResizeArray<Real> th0Array;
280  ResizeArray<Real> kthArray;
281  ResizeArray<int> functArray;
282  ResizeArray<const char *> typeaArray; /* NULL if no type is set */
283  ResizeArray<const char *> typebArray; /* NULL if no type is set */
284  ResizeArray<const char *> typecArray; /* NULL if no type is set */
285 
286  public:
287  /* returns the number of types in the table */
288  int size() const { return th0Array.size(); }
289 
290  /* this gets the index of a angle in the table (adding an entry if
291  none exists).
292  */
293  int getIndex(Real th0, Real kth, int funct);
294 
295  /* this adds a entry for a angle type to the table, including the
296  three atoms involved in the angle */
297  void addType(const char *typea, const char *typeb,
298  const char *typec, Real th0,
299  Real kth, int funct);
300 
301  /* getAngleParams puts the parameters for angle-type <num> into the
302  spaces pointed to by the other arguments.
303  */
304  void getParams(int num, Real *th0, Real *kth, int *funct) const;
305 
306  /* This version looks up the angle by atom type - the direction of
307  the types doesn't matter (it can be A--B--C or C--B--A)! If the
308  specified angle/function combination cannot be found, it returns
309  -1, otherwise it returns the index of the angletype. */
310  int getParams(const char *typea, const char *typeb,
311  const char *typec, int funct,
312  Real *th0, Real *kth) const;
313 };
314 
315 /* A DihedralTable stores all the different types of angles, so that we
316  only have to have one copy of each. Units:
317  funct 1: proper dihedral
318  c0 - thmax (deg)
319  c1 - fc (kcal/mol)
320  a,b are the two inner atoms
321  funct 2: improper dihedral
322  c0 - th0 (deg)
323  c1 - fc (kcal/mol)
324  a,b are the two outer atoms
325  funct 3: RB dihedral
326  c0-c5 - C0-C5 (kcal/mol)
327  a,b are the two inner atoms
328 
329  The atoms are in standard geometry order - A--B--C--D for the
330  dihedral between the planes containing ABC and BCD.
331 */
333  private:
334  /* implementation of the extensible array */
335  ResizeArray<int> functArray;
336  ResizeArray<const char *> typeaArray; /* NULL if no type is set */
337  ResizeArray<const char *> typebArray; /* NULL if no type is set */
338  ResizeArray<int> multArray; /* 0 if funct is not 1 */
339 
340  /* cArray is where the actual paramters live - each Real * is a
341  pointer to a list of Reals that is at least as long as the
342  number of parameters required for the given type. */
344 
345  public:
346  /* returns the number of types in the table */
347  int size() const { return functArray.size(); }
348 
349  /* this function gets the index of a dihedral angle in the table
350  (adding an entry if none exists). The required number of
351  parameters (see notes on class DihedralTable) must be stored in
352  the array <c> */
353  int getIndex(const Real *c, int mult, int funct);
354 
355  /* this adds a entry for a angle type to the table, including the
356  two atoms involved in the angle. The required number of
357  parameters (see notes on class DihedralTable) must be stored in
358  the array <c>. Note that these two angles really refer to either
359  atoms A and D or B and C depending on the dihedral type. */
360  void addType(const char *typea, const char *typeb,
361  const Real *c, int mult, int funct);
362 
363  /* getParams puts the parameters for angle-type <num> into the
364  spaces pointed to by the other arguments. The required number of
365  parameters (see notes on class DihedralTable) will be stored in
366  the array <c>, so make sure that c has size >= 6 */
367  void getParams(int num, Real *c, int *mult, int *funct) const;
368 
369  /* This version looks up the angle by atom type - the direction of
370  the types doesn't matter (it can be ABCD or DCBA)! If the
371  specified angle/function combination cannot be found, it returns
372  -1, otherwise it returns the index of the dihedraltye. This is a
373  little strange because the different dihedral types use different
374  atoms - but since we always will know all four AND the funct when
375  we need to look up a dihedral, this will work. The required
376  number of parameters (see notes on class DihedralTable) will be
377  stored in the array <c>, so make sure that c is big enough for
378  the specified function. */
379  int getParams(const char *typea, const char *typeb,
380  const char *typec, const char *typed, int funct,
381  Real *c, int *mult) const;
382 };
383 
384 /* this stores the VDW interaction parameters as a function of the
385  atom types. */
386 class VDWTable {
387  private:
388  /* atom types a and b */
389  ResizeArray<const char *> typeAArray; /* NULL if no type is set */
390  ResizeArray<const char *> typeBArray; /* NULL if no type is set */
391 
392  /* the VDW parameters */
393  ResizeArray<Real> c6Array;
394  ResizeArray<Real> c12Array;
395  ResizeArray<Real> c6PairArray; // these two are
396  ResizeArray<Real> c12PairArray; // for 1-4 pairs
397 
398  /* finds the index from the two interacting types, or returns -1 */
399  int getIndex(const char *typea, const char *typeb) const;
400 
401  public:
402  /* these add VDW parameters to the list */
403  void addType(const char *typea, const char *typeb, Real c6, Real c12);
404  void add14Type(const char *typea, const char *typeb,
405  Real c6pair, Real c12pair);
406 
407  /* this returns the VDW interaction parameters that were added to
408  the list earlier, and returns the index or the parameters (a
409  useless number) or -1 if it can't find the specified types. */
410  int getParams(const char *typea, const char *typeb,
411  Real *c6, Real *c12, Real *c6pair, Real *c12pair) const;
412 };
413 
414 /* this stores the pair interaction parameters as a function of the
415  atom types.
416  JLai August 16th, 2012
417 */
418 struct GroLJPair {
419  int indxLJA;
420  int indxLJB;
423 };
424 
425 struct GroGaussPair {
434 };
435 
436 class PairTable {
437  private:
438  // Data structure for LJ Structure-based potential
439  ResizeArray<int> indxLJA;
440  ResizeArray<int> indxLJB;
441  ResizeArray<Real> c6pair;
442  ResizeArray<Real> c12pair;
443  ResizeArray<GroLJPair> pairlistLJ;
444 
445  // Data structure for Gaussian Structure-based potential
446  ResizeArray<int> indxGaussA;
447  ResizeArray<int> indxGaussB;
449  ResizeArray<Real> gMu1;
450  ResizeArray<Real> giSigma1;
451  ResizeArray<Real> gMu2;
452  ResizeArray<Real> giSigma2;
453  ResizeArray<Real> gRepulsive;
454  ResizeArray<GroGaussPair> pairlistGauss;
455 
456  public:
457  /* these add pair parameters to the list */
458  int addPairLJType2(int typea, int typeb, Real c6, Real c12);
459  int addPairGaussType2(int typea, int typeb, Real gA,
460  Real gMu1, Real gSigma1);
461  int addPairGaussType2(int typea, int typeb, Real gA,
462  Real gMu1, Real gSigma1, Real gRepulsive);
463  int addPairGaussType2(int typea, int typeb, Real gA,
464  Real gMu1, Real gSigma1, Real gMu2, Real gSigma2, Real gRepulsive);
465  void getPairLJArrays2(int *indexA, int *indexB, Real *pairC6, Real *pairC12);
466  void getPairGaussArrays2(int *indexA, int *indexB, Real *gaussA, Real *gaussMu1,
467  Real *gaussSigma1, Real *gaussMu2, Real *gaussSigma2,
468  Real *gaussRepulsive);
469 
470  static bool GroLJCompare (GroLJPair, GroLJPair);
472 
473 };
474 
475 /* A GromacsTopFile represents the information stored in a GROMACS
476  topology file. This is an immutable type. */
478  private:
479  /* The system name */
480  char *systemName;
481 
482  /* 1-4 parameter scaling values */
483  Real fudgeLJ, fudgeQQ;
484 
485  /* whether or not to generate the 1-4 interactions automatically */
486  int genPairs;
487 
488  /* a list of molecule templates */
489  ResizeArray <GenericMol *> genericMols;
490 
491  /* the list of molecule instances */
493 
494  /* the table of bonds/angles/atoms */
495  AtomTable atomTable;
496  BondTable bondTable;
497  AngleTable angleTable;
498  DihedralTable dihedralTable;
499  VDWTable vdwTable;
500  PairTable pairTable;
501 
502  public:
503 
504  /* initializes this to represent the data stored in the topology
505  file <filename>, or exits on error. */
506  GromacsTopFile(char *filename);
507 
508  /* getSystemName returns the name of the system */
509  char *getSystemName() const { return systemName; }
510 
511  /* returns the number of atoms/bonds/angles/dihedrals in the file */
512  int getNumAtoms() const;
513  int getNumBonds() const;
514  int getNumAngles() const;
515  int getNumDihedrals() const;
516 
517  /* returns the number of atoms/bonds/angles/dihedrals params in the file */
518  int getNumAtomParams() const { return atomTable.size(); }
519  int getNumBondParams() const { return bondTable.size(); }
520  int getNumAngleParams() const { return angleTable.size(); }
521  int getNumDihedralParams() const { return dihedralTable.size(); }
522 
523  /* getAtom puts the information about atom number <num> into the
524  spaces pointed to by the other arguments. The string buffers
525  must be at least 11 characters long. Atom number 0 is the first
526  atom in the list, so that it corresponds with our normal idea of
527  atom numbers.
528  charge - e
529  mass - amu
530  */
531  void getAtom(int num, int *residue_number, char *residue_name,
532  char *atom_name, char *atom_type, int *atom_typenum,
533  Real *charge, Real *mass)
534  const;
535 
536  /* getAtomParams puts the parameters for atom-type <num> into the
537  spaces pointed to by the other arguements. Currently, all anyone
538  cares about is the atom type string, so that's all we're
539  returning! This may change later, but for now we get everything
540  else with the getAtom() method. The string buffer must be at
541  least 11 characters long. */
542  void getAtomParams(int num, char *type) const {
543  atomTable.getType(num,type);
544  }
545 
546  // JLai
547  int getNumPair() const;
548  int getNumLJPair() const;
549  int getNumGaussPair() const;
550  int getNumExclusions() const;
551  void getExclusions(int *, int *) const;
552  // End of JLai
553 
554  /* getBond puts the information about bond number <num> into the
555  spaces pointed to by the other arguments. Bond number 0 is the
556  first bond in the list.
557  */
558  void getBond(int num, int *atomi, int *atomj, int *bondtype) const;
559 
560  /* getBondParams puts the parameters for bond-type <num> into the
561  spaces pointed to by the other arguments.
562  b0 - natural length in A
563  kB - spring constant in kcal/mol/A^2 - E=kx^2
564  funct - 1 for harmonic, 2 for special fourth-order GROMOS96 thing */
565  void getBondParams(int num, Real *b0, Real *kB, int *funct) const;
566 
567  /* getAngle puts the information about angle number <num> into the
568  spaces pointed to by the other arguments. Angle number 0 is the
569  first angle in the list.
570  */
571  void getAngle(int num, int *atomi, int *atomj, int *atomk,
572  int *angletype) const;
573 
574  /* getAngleParams puts the parameters for angle-type <num> into the
575  spaces pointed to by the other arguments.
576  th0 - natural angle in degrees
577  kth - spring constant in kcal/rad - E=kx^2
578  funct - 1 for harmonic, 2 for special GROMOS96 thing */
579  void getAngleParams(int num, Real *th0, Real *kth, int *funct) const;
580 
581  /* getDihedral puts the information about dihedral number <num> into
582  the spaces pointed to by the other arguments. Dihedral number 0
583  is the first angle in the list. */
584  void getDihedral(int num, int *atomi, int *atomj, int *atomk,
585  int *atoml, int *type) const;
586 
587  /* getDihedralParams puts the parameters for dihedral-type <num>
588  into the spaces pointed to by the other arguments. The array of
589  Reals <c> must be of size >= 6, since up to 6 parameters will be
590  set there. These parameters are as follows:
591  funct 1: proper dihedral
592  c0 - thmax (deg)
593  c1 - fc (kcal/mol)
594  mult is the multiplicity
595  funct 2: improper dihedral
596  c0 - th0 (deg)
597  c1 - fc (kcal/mol/rad2)
598  funct 3: RB dihedral
599  c0-c5 - C0-C5 (kcal/mol)
600  */
601  void getDihedralParams(int num, Real *c, int *mult, int *funct) const;
602 
603  // JLai
604  void getPairLJArrays2(int *indexA, int *indexB, Real *pairC6, Real *pairC12);
605  void getPairGaussArrays2(int *indexA, int *indexB, Real *gaussA, Real *gaussMu1,
606  Real *gaussSigma1, Real *gaussMu2, Real *gaussSigma2,
607  Real *gaussRepulsive);
608 
609  /* getVDWParams returs the Lennard-Jones bonding parameters for the
610  specified two atom types, and the modified bonding parameters
611  for 1-4 L-J interactons (c6pair, c12pair).
612  Note that unlike the other types of parameters, you don't refer
613  to this one by number - any two atom types define a VDWParam */
614  void getVDWParams(int typea, int typeb,
615  Real *c6, Real *c12, Real *c6pair, Real *c7) const;
616 
617 };
618 
619 #endif // GROMACSTOPFILE_H
int getType() const
Real getCharge() const
void getBond(int num, int *atomi, int *atomj, int *bondtype) const
GenericMol(const char *theName)
int getNumBonds() const
int addPairGaussType2(int typea, int typeb, Real gA, Real gMu1, Real gSigma1)
GenericBond(int i, int j, int theType)
int getParams(const char *typea, const char *typeb, Real *c6, Real *c12, Real *c6pair, Real *c12pair) const
void getPairLJArrays2(int *indexA, int *indexB, Real *pairC6, Real *pairC12)
static bool GroGaussCompare(GroGaussPair, GroGaussPair)
int getNumDihedrals() const
int getNumDihedrals() const
const GenericBond * getBond(int n) const
const char * getAtomName() const
void add14Type(const char *typea, const char *typeb, Real c6pair, Real c12pair)
int getAtoml() const
int getResNum() const
void getVDWParams(int typea, int typeb, Real *c6, Real *c12, Real *c6pair, Real *c7) const
void getType(int num, char *type) const
void addType(const char *type, Real m, Real q, Real c6, Real c12)
int getParams(const char *type, Real *m, Real *q, Real *c6, Real *c12) const
float Real
Definition: common.h:109
int getNumRes() const
int getNumAngles() const
const char * getType() const
int getNumAtoms() const
int size() const
int getNumExclusions() const
int size() const
const GenericAngle * getAngle(int n) const
void addType(const char *typea, const char *typeb, const char *typec, Real th0, Real kth, int funct)
int getNum() const
void getBondParams(int num, Real *b0, Real *kB, int *funct) const
int getIndex(Real b0, Real kB, int funct)
int getNumLJPair() const
void getDihedral(int num, int *atomi, int *atomj, int *atomk, int *atoml, int *type) const
char * getSystemName() const
int getIndex(Real th0, Real kth, int funct)
int getNumAngles() const
void addType(const char *typea, const char *typeb, const Real *c, int mult, int funct)
int getType() const
int getNumGaussPair() const
int getAtomk() const
const char * getResType() const
int getType() const
void getParams(int num, Real *b0, Real *kB, int *funct) const
void getParams(int num, Real *c, int *mult, int *funct) const
GromacsTopFile(char *filename)
void getExclusions(int *, int *) const
int addPairLJType2(int typea, int typeb, Real c6, Real c12)
void getParams(int num, Real *th0, Real *kth, int *funct) const
void addAngle(int atomi, int atomj, int atomk, int type)
int getNumBonds() const
const GenericAtom * getAtom(int n) const
int getNumDihedralParams() const
#define NAMESIZE
Definition: GromacsTopFile.h:6
const GenericMol * getMol() const
void addType(const char *typea, const char *typeb, Real c6, Real c12)
void addType(const char *typea, const char *typeb, Real b0, Real kB, int funct)
int getNumAtomParams() const
int getNumAngles() const
int getTypeNum() const
void getPairLJArrays2(int *indexA, int *indexB, Real *pairC6, Real *pairC12)
const GenericDihedral * getDihedral(int n) const
int getIndex(const Real *c, int mult, int funct)
void getAngle(int num, int *atomi, int *atomj, int *atomk, int *angletype) const
int getNumAtoms() const
int getAtomj() const
void getDihedralParams(int num, Real *c, int *mult, int *funct) const
static bool GroLJCompare(GroLJPair, GroLJPair)
GenericAngle(int i, int j, int k, int theType)
int getAtomj() const
int getAtomi() const
int getNumAngleParams() const
void getAtomParams(int num, char *type) const
int size() const
k< npairi;++k){TABENERGY(const int numtypes=simParams->tableNumTypes;const float table_spacing=simParams->tableSpacing;const int npertype=(int)(namdnearbyint(simParams->tableMaxDist/simParams->tableSpacing)+1);) int table_i=(r2iilist[2 *k] >> 14)+r2_delta_expc;const int j=pairlisti[k];#define p_j BigReal diffa=r2list[k]-r2_table[table_i];#define table_four_i TABENERGY(register const int tabtype=-1-(lj_pars->A< 0?lj_pars->A:0);) BigReal kqq=kq_i *p_j-> charge
int getNumBondParams() const
int getNumRes() const
Real getMass() const
void getPairGaussArrays2(int *indexA, int *indexB, Real *gaussA, Real *gaussMu1, Real *gaussSigma1, Real *gaussMu2, Real *gaussSigma2, Real *gaussRepulsive)
int getNumAtoms() const
int getNumDihedrals() const
void getPairGaussArrays2(int *indexA, int *indexB, Real *gaussA, Real *gaussMu1, Real *gaussSigma1, Real *gaussMu2, Real *gaussSigma2, Real *gaussRepulsive)
int size(void) const
Definition: ResizeArray.h:127
const char * getName() const
int getAtomi() const
int getAtomk() const
int getAtomi() const
int size() const
void addBond(int atomi, int atomj, int type)
int getAtomj() const
void addDihedral(int atomi, int atomj, int atomk, int atoml, int type)
int getNumPair() const
GenericAtom(const char *theType, int theTypeNum, int theResNum, const char *theResType, const char *theAtomName, Real theCharge, Real theMass)
GenericDihedral(int i, int j, int k, int l, int theType)
void addAtom(const char *theType, int theTypeNum, int theResNum, const char *theResType, const char *theAtomName, Real theCharge, Real theMass)
void getAtom(int num, int *residue_number, char *residue_name, char *atom_name, char *atom_type, int *atom_typenum, Real *charge, Real *mass) const
void getAngleParams(int num, Real *th0, Real *kth, int *funct) const
MolInst(const GenericMol *theMol, int theNum)
int getNumBonds() const