NAMD
Molecule.h
Go to the documentation of this file.
1 
7 /*
8  This class is used to store all of the structural
9  information for a simulation. It reads in this information
10  from a .psf file, cross checks and obtains some information
11  from the Parameters object that is passed in, and then
12  stores all this information for later use.
13 */
14 
15 
16 #ifndef MOLECULE_H
17 
18 #define MOLECULE_H
19 
20 #include "parm.h"
21 
22 #include "common.h"
23 #include "NamdTypes.h"
24 #include "structures.h"
25 #include "ConfigList.h"
26 #include "Vector.h"
27 #include "UniqueSet.h"
28 #include "Hydrogen.h"
29 #include "GromacsTopFile.h"
30 /* BEGIN gf */
31 #include "GridForceGrid.h"
32 #include "Tensor.h"
33 /* END gf */
34 
35 // Go -- JLai
36 #define MAX_GO_CHAINS 10
37 #define MAX_RESTRICTIONS 10
38 // End of Go
39 
40 #include "molfile_plugin.h"
41 
42 #include <vector>
43 
44 class SimParameters;
45 class Parameters;
46 class ConfigList;
47 class PDB;
48 class MIStream;
49 class MOStream;
50 
51 class ExclElem;
52 class BondElem;
53 class AngleElem;
54 class DihedralElem;
55 class ImproperElem;
56 class TholeElem; // Drude model
57 class AnisoElem; // Drude model
58 class CrosstermElem;
59 // JLai
60 class GromacsPairElem;
61 // End of JLai
62 class ResidueLookupElem;
63 
64 struct OutputAtomRecord;
65 
66 template<class Type> class ObjectArena;
67 
69 public:
71  char *flags;
72 
74  min=0;
75  max=-1;
76  flags = NULL;
77  }
78 private:
79  ExclusionCheck(const ExclusionCheck& chk){
80  }
81  ExclusionCheck &operator=(const ExclusionCheck& chk){
82  return *this;
83  }
84 };
85 #define EXCHCK_FULL 1
86 #define EXCHCK_MOD 2
87 
89 {
90 public:
91  char mySegid[11];
92  ResidueLookupElem *next; // stored as a linked list
93  int firstResid; // valid resid is >= firstResid
94  int lastResid; // valid resid is <= lastResid
95  ResizeArray<int> atomIndex; // 0-based index for first atom in residue
96 
97  ResidueLookupElem(void) { next = 0; firstResid = -1; lastResid = -1; }
98  ~ResidueLookupElem(void) { delete next; }
99  int lookup(const char *segid, int resid, int *begin, int *end) const;
100  ResidueLookupElem* append(const char *segid, int resid, int aid);
101 };
102 
103 // Ported by JLai -- JE
104 typedef struct go_val
105 {
106  Real epsilon; // Epsilon
107  int exp_a; // First exponent for attractive L-J term
108  int exp_b; // Second exponent for attractive L-J term
109  int exp_rep; // Exponent for repulsive L-J term
110  Real sigmaRep; // Sigma for repulsive term
111  Real epsilonRep; // Epsilon for replusive term
112  Real cutoff; // Cutoff distance for Go calculation
113  int restrictions[MAX_RESTRICTIONS]; // List of residue ID differences to be excluded from Go calculation
114 } GoValue;
115 
116 typedef struct go_pair
117 {
118  int goIndxA; // indexA
119  int goIndxB; // indexB
120  double A; // double A for the LJ pair
121  double B; // double B for the LJ pair
122 } GoPair;
123 // End of port - JL
124 
125 #define QMLSSMODEDIST 1
126 #define QMLSSMODECOM 2
127 
128 #define QMFormatORCA 1
129 #define QMFormatMOPAC 2
130 #define QMFormatUSR 3
131 
132 #define QMCHRGNONE 0
133 #define QMCHRGMULLIKEN 1
134 #define QMCHRGCHELPG 2
135 
136 #define QMSCHEMECS 1
137 #define QMSCHEMERCD 2
138 #define QMSCHEMEZ1 3
139 #define QMSCHEMEZ2 4
140 #define QMSCHEMEZ3 5
141 
142 //only used for compressing the molecule information
143 typedef struct seg_resid
144 {
145  char segname[11];
146  int resid;
148 
149 // List maintaining the global atom indicies sorted by helix groups.
150 class Molecule
151 {
152  private:
153 typedef struct constraint_params
154 {
155  Real k; // Force constant
156  Vector refPos; // Reference position for restraint
157 } ConstraintParams;
158 
159 
160 
161 /* BEGIN gf */
162 typedef struct gridfrc_params
163 {
164  Real k; // force multiplier
165  Charge q; // charge
166 } GridforceParams;
167 /* END gf */
168 
169 
170 typedef struct stir_params
171 {
172  Real startTheta;
173  Vector refPos; // Reference position for stirring
174 } StirParams;
175 
176 typedef struct movdrag_params
177 {
178  Vector v; // Linear velocity (A/step)
179 } MovDragParams;
180 
181 
182 typedef struct rotdrag_params
183 {
184  Real v; // Angular velocity (deg/step)
185  Vector a; // Rotation axis
186  Vector p; // Rotation pivot point
187 } RotDragParams;
188 
189 typedef struct constorque_params
190 {
191  Real v; // "Torque" value (Kcal/(mol*A^2))
192  Vector a; // Rotation axis
193  Vector p; // Rotation pivot point
194 } ConsTorqueParams;
195 
196 #ifdef MEM_OPT_VERSION
197 //indicate a range of atoms from aid1 to aid2
198 struct AtomSet{
199  AtomID aid1;
200  AtomID aid2;
201 
202  int operator < (const AtomSet &a) const{
203  return aid1 < a.aid1;
204  }
205 };
206 typedef ResizeArray<AtomSet> AtomSetList;
207 
208  void load_atom_set(StringList *setfile, const char *setname,
209  int *numAtomsInSet, Molecule::AtomSetList **atomsSet) const;
210 
211 #endif
212 
213 friend class ExclElem;
214 friend class BondElem;
215 friend class AngleElem;
216 friend class DihedralElem;
217 friend class ImproperElem;
218 friend class TholeElem; // Drude model
219 friend class AnisoElem; // Drude model
220 friend class CrosstermElem;
221 // JLai
222 friend class GromacsPairElem;
223 // End of JLai
224 friend class WorkDistrib;
225 
226 private:
227 
228 #ifndef MEM_OPT_VERSION
229  Atom *atoms; // Array of atom structures
230  ObjectArena<char> *nameArena;
231  AtomNameInfo *atomNames;// Array of atom name info. Only maintained on node 0 for VMD interface
232  //replaced by atom signatures
233  Bond *bonds; // Array of bond structures
234  Angle *angles; // Array of angle structures
235  Dihedral *dihedrals; // Array of dihedral structures
236  Improper *impropers; // Array of improper structures
237  Crossterm *crossterms; // Array of cross-term structures
238  GromacsPair *gromacsPair; // Array of gromacs-pair structures
239 
240  //These will be replaced by exclusion signatures
241  Exclusion *exclusions; // Array of exclusion structures
242  UniqueSet<Exclusion> exclusionSet; // Used for building
243 
244  int32 *cluster; // first atom of connected cluster
245 
246  ObjectArena<int32> *tmpArena;
247  int32 **bondsWithAtom; // List of bonds involving each atom
248  ObjectArena<int32> *arena;
249 
250  //function is replaced by atom signatures
251  int32 **bondsByAtom; // List of bonds owned by each atom
252  int32 **anglesByAtom; // List of angles owned by each atom
253  int32 **dihedralsByAtom; // List of dihedrals owned by each atom
254  int32 **impropersByAtom; // List of impropers owned by each atom
255  int32 **crosstermsByAtom; // List of crossterms owned by each atom
256 
257  int32 **exclusionsByAtom; // List of exclusions owned by each atom
258  int32 **fullExclusionsByAtom; // List of atoms excluded for each atom
259  int32 **modExclusionsByAtom; // List of atoms modified for each atom
260 // JLai
261  int32 **gromacsPairByAtom; // gromacsPair exclusion list which by definition should not have any exclusions (still not sure if it should be empty or zeroed out)
262 // End of JLai
263 
264  ObjectArena<char> *exclArena;
265  ExclusionCheck *all_exclusions;
266 
267  // DRUDE
268  int32 **tholesByAtom; // list of Thole correction terms owned by each atom
269  int32 **anisosByAtom; // list of anisotropic terms owned by each atom
270  // DRUDE
271 
272 #else
273  //Indexing to constant pools to save space
275  Index *eachAtomMass; //after initialization, this could be freed (possibly)
276  Index *eachAtomCharge; //after initialization, this could be freed (possibly)
277  AtomNameIdx *atomNames;
278  ObjectArena<char> *nameArena; //the space for names to be allocated
279 
280  //A generally true assumption: all atoms are arranged in the order of clusters. In other words,
281  //all atoms for a cluster must appear before/after any atoms in other clusters
282  //The first atom in the cluster (which has the lowest atom id) stores the cluster size
283  //other atoms in the cluster stores -1
284  int32 *clusterSigs;
285  int32 numClusters;
286 
287  AtomSigID *eachAtomSig;
288  ExclSigID *eachAtomExclSig;
289 
290  AtomSetList *fixedAtomsSet;
291  AtomSetList *constrainedAtomsSet;
292 #endif
293 
294  ResidueLookupElem *resLookup; // find residues by name
295 
296  AtomSegResInfo *atomSegResids; //only used for generating compressed molecule info
297 
298  Bond *donors; // Array of hydrogen bond donor structures
299  Bond *acceptors; // Array of hydrogen bond acceptor
300 
301  // DRUDE
302  DrudeConst *drudeConsts; // supplement Atom data (length of Atom array)
303  Thole *tholes; // Thole (screened Coulomb) interactions
304  Aniso *anisos; // anisotropic terms
305  Lphost *lphosts; // lone pair hosts
306  int32 *lphostIndexes; // index for each LP into lphosts array
307  // DRUDE
308 
309  int32 *consIndexes; // Constraint indexes for each atom
310  ConstraintParams *consParams;
311 
312 /* BEGIN gf */
313  int32 **gridfrcIndexes;
314  GridforceParams **gridfrcParams;
315  GridforceGrid **gridfrcGrid;
316 /* END gf */
317 
318  // Parameters for each atom constrained
319  int32 *stirIndexes; //Stirring indexes for each atoms
320  StirParams *stirParams;
321  // Paramters for each atoms stirred
322  int32 *movDragIndexes; // Moving drag indexes for each atom
323  MovDragParams *movDragParams;
324  // Parameters for each atom moving-dragged
325  int32 *rotDragIndexes; // Rotating drag indexes for each atom
326  RotDragParams *rotDragParams;
327  // Parameters for each atom rotation-dragged
328 
329  Real *langevinParams; // b values for langevin dynamics
330  int32 *fixedAtomFlags; // 1 for fixed, -1 for fixed group, else 0
331  int32 *exPressureAtomFlags; // 1 for excluded, -1 for excluded group.
332 
333  //In the memory optimized version: it will be NULL if the general
334  //true assumption mentioned above holds. If not, its size is numClusters.
335  //In the ordinary version: its size is numAtoms, and indicates the size
336  //of connected cluster or 0.
337  int32 *clusterSize;
338 
339  Real *rigidBondLengths; // if H, length to parent or 0. or
340  // if not H, length between children or 0.
341 //fepb
342  unsigned char *fepAtomFlags;
343 //fepe
344 
345 //soluteScaling
346  unsigned char *ssAtomFlags;
347  unsigned int num_ss;
348 //soluteScaling
349 
350  //occupancy and bfactor data from plugin-based IO implementation of loading structures
351  float *occupancy;
352  float *bfactor;
353 
354 
355  // added during the trasition from 1x to 2
357  Parameters *params;
358 
359 private:
360  void initialize(SimParameters *, Parameters *param);
361  // Sets most data to zero
362 
363  //LCPO
364  int *lcpoParamType;
365 
366 #ifndef MEM_OPT_VERSION
367  void read_psf_file(char *, Parameters *);
368  // Read in a .psf file given
369  // the filename and the parameter
370  // object to use
371 
372  void read_atoms(FILE *, Parameters *);
373  // Read in atom info from .psf
374  void read_bonds(FILE *);
375  // Read in bond info from .psf
376  void process_bonds(Parameters *);
377  // Remove fake bonds after read_lphosts() to deal with TIP4 water
378  void read_angles(FILE *, Parameters *);
379  // Read in angle info from .psf
380  void read_dihedrals(FILE *, Parameters *);
381  // Read in dihedral info from .psf
382  void read_impropers(FILE *, Parameters *);
383  // Read in improper info from .psf
384  void read_crossterms(FILE *, Parameters *);
385  // Read in cross-term info from .psf
386  void read_donors(FILE *);
387  // Read in hydrogen bond donors from .psf
388  void read_acceptors(FILE *);
389  // Read in hydrogen bond acceptors from .psf
390  void read_exclusions(FILE *);
391  // Read in exclusion info from .psf
392  // JLai August 16th, 2012 Modifications
393  void read_exclusions(int*, int*, int);
394  /* Read in exclusion array and sort entries */
395  static bool goPairCompare (GoPair, GoPair);
396  // End of JLai August 16th, 2012 Modifications
397 
398 
399  // DRUDE: PSF reading
400  void read_lphosts(FILE *);
401  // Read in lone pair hosts from Drude PSF
402  void read_anisos(FILE *);
403  // Read in anisotropic terms from Drude PSF
404  // DRUDE
405 
406  //LCPO
407  //input type is Charmm/Amber/other
408  //0 - Charmm/Xplor
409  //1 - Amber TODO
410  //2 - Plugin TODO
411  //3 - Gromacs TODO
412  void assignLCPOTypes(int inputType);
413 
414  //pluginIO-based loading atoms' structure
415  void plgLoadAtomBasics(molfile_atom_t *atomarray);
416  void plgLoadBonds(int *from, int *to); //atom index is 1-based in the parameters
417  void plgLoadAngles(int *plgAngles);
418  void plgLoadDihedrals(int *plgDihedrals);
419  void plgLoadImpropers(int *plgImpropers);
420  void plgLoadCrossterms(int *plgCterms);
421 
422  // List of all exclusions, including
423  // explicit exclusions and those calculated
424  // from the bonded structure based on the
425  // exclusion policy. Also includes 1-4's.
426  void build_lists_by_atom();
427  // Build the list of structures by atom
428 
429  void build12excl(void);
430  void build13excl(void);
431  void build14excl(int);
432 
433  // DRUDE: extend exclusions for Drude and LP
434  void build_inherited_excl(int);
435  // DRUDE
436  void stripFepExcl(void);
437 
438  void build_exclusions();
439  // analyze the atoms, and determine which are oxygen, hb donors, etc.
440  // this is called after a molecule is sent our (or received in)
441  void build_atom_status(void);
442 
443 #else
444  //the method to load the signatures of atoms etc. (i.e. reading the file in
445  //text fomrat of the compressed psf file)
446  void read_mol_signatures(char *fname, Parameters *params, ConfigList *cfgList=0);
447  void load_one_inputatom(int aid, OutputAtomRecord *one, InputAtom *fAtom);
448  void build_excl_check_signatures();
449 #endif
450 
451  void read_parm(const GromacsTopFile *);
452 
453 public:
454 
455  // for solute scaling
458  int *ss_index;
459 
460  // DRUDE
461  int is_drude_psf; // flag for reading Drude PSF
462  int is_lonepairs_psf; // flag for reading lone pairs from PSF
463  // DRUDE
464 
465  // data for TIP4P
468 
469  // data for tail corrections
474  void compute_LJcorrection();
475  BigReal getEnergyTailCorr(const BigReal, const int);
477 
478  int const * getLcpoParamType() {
479  return lcpoParamType;
480  }
481 
482  BigReal GetAtomAlpha(int i) const {
483  return drudeConsts[i].alpha;
484  }
485 
486 #ifdef MEM_OPT_VERSION
487  AtomCstInfo *getAtoms() const { return atoms; }
488  AtomNameIdx *getAtomNames() const { return atomNames; }
489 #else
490  Atom *getAtoms () const { return atoms; }
491  AtomNameInfo *getAtomNames() const { return atomNames; }
492 #endif
493 
494  AtomSegResInfo *getAtomSegResInfo() const { return atomSegResids; }
495 
496  // return number of fixed atoms, guarded by SimParameters
497  int num_fixed_atoms() const {
498  // local variables prefixed by s_
499  int s_NumFixedAtoms = (simParams->fixedAtomsOn ? numFixedAtoms : 0);
500  return s_NumFixedAtoms; // value is "turned on" SimParameters
501  }
502 
503  int num_fixed_groups() const {
504  // local variables prefixed by s_
505  int s_NumFixedAtoms = num_fixed_atoms();
506  int s_NumFixedGroups = (s_NumFixedAtoms ? numFixedGroups : 0);
507  return s_NumFixedGroups; // value is "turned on" SimParameters
508  }
509 
510  int64_t num_group_deg_freedom() const {
511  // local variables prefixed by s_
512  int64_t s_NumGroupDegFreedom = 3 * (int64_t) numHydrogenGroups;
513  int64_t s_NumFixedAtoms = num_fixed_atoms();
514  int s_NumFixedGroups = num_fixed_groups();
515  if (s_NumFixedGroups) s_NumGroupDegFreedom -= 3 * s_NumFixedGroups;
516  if ( ! (s_NumFixedAtoms || numConstraints
517  || simParams->comMove || simParams->langevinOn) ) {
518  s_NumGroupDegFreedom -= 3;
519  }
520  return s_NumGroupDegFreedom;
521  }
522 
523  int64_t num_deg_freedom(int isInitialReport = 0) const {
524  // local variables prefixed by s_
525  int64_t s_NumDegFreedom = 3 * (int64_t) numAtoms;
526  int64_t s_NumFixedAtoms = num_fixed_atoms();
527  if (s_NumFixedAtoms) s_NumDegFreedom -= 3 * s_NumFixedAtoms;
528  if (numLonepairs) s_NumDegFreedom -= 3 * numLonepairs;
529  if ( ! (s_NumFixedAtoms || numConstraints
530  || simParams->comMove || simParams->langevinOn) ) {
531  s_NumDegFreedom -= 3;
532  }
533  if ( ! isInitialReport && simParams->pairInteractionOn) {
534  //
535  // DJH: a kludge? We want to initially report system degrees of freedom
536  //
537  // this doesn't attempt to deal with fixed atoms or constraints
538  s_NumDegFreedom = 3 * numFepInitial;
539  }
540  int s_NumFixedRigidBonds =
541  (simParams->fixedAtomsOn ? numFixedRigidBonds : 0);
542  if (simParams->watmodel == WAT_TIP4) {
543  // numLonepairs is subtracted here because all lonepairs have a rigid bond
544  // to oxygen, but all of the LP degrees of freedom are dealt with above
545  s_NumDegFreedom -= (numRigidBonds - s_NumFixedRigidBonds - numLonepairs);
546  }
547  else {
548  // Non-polarized systems don't have LPs.
549  // For Drude model, bonds that attach LPs are not counted as rigid;
550  // LPs have already been subtracted from degrees of freedom.
551  s_NumDegFreedom -= (numRigidBonds - s_NumFixedRigidBonds);
552  }
553  return s_NumDegFreedom;
554  }
555 
556  int numAtoms; // Number of atoms
557 
558  int numRealBonds; // Number of bonds for exclusion determination
559  int numBonds; // Number of bonds calculated, including extras
560  int numAngles; // Number of angles
561  int numDihedrals; // Number of dihedrals
562  int suspiciousAlchBonds; // angles dropped due to crossing FEP partitions
563  int alchDroppedAngles; // angles dropped due to crossing FEP partitions
564  int alchDroppedDihedrals; // dihedrals dropped due to crossing FEP partitions
565  int alchDroppedImpropers; // impropers dropped due to crossing FEP partitions
566  int numImpropers; // Number of impropers
567  int numCrossterms; // Number of cross-terms
568  int numDonors; // Number of hydrogen bond donors
569  int numAcceptors; // Number of hydrogen bond acceptors
570  int numExclusions; // Number of exclusions
571  int64 numTotalExclusions; // Real Total Number of Exclusions // hack
572  int num_alch_unpert_Bonds; // These are
573  int num_alch_unpert_Angles; // Shobana's
574  int num_alch_unpert_Dihedrals; // unperturbed bond
575  Bond *alch_unpert_bonds; // angle, dihedral
576  Angle *alch_unpert_angles; // terms to remove dummy
577  Dihedral *alch_unpert_dihedrals; //atom influence
578 
579  // DRUDE
582  int numTholes;
583  int numAnisos;
586  // DRUDE
587 
588  int numConstraints; // Number of atoms constrained
589 /* BEGIN gf */
590  int numGridforceGrids;// Number of gridforce grids
591  int *numGridforces; // Number of atoms in gridforce file (array, one per grid)
592 /* END gf */
593  int numMovDrag; // Number of atoms moving-dragged
594  int numRotDrag; // Number of atoms rotating-dragged
595  int numConsTorque; // Number of atoms "constant"-torqued
596  int numFixedAtoms; // Number of fixed atoms
597  int numStirredAtoms; // Number of stirred atoms
598  int numExPressureAtoms; // Number of atoms excluded from pressure
599  int numHydrogenGroups; // Number of hydrogen groups
600  int maxHydrogenGroupSize; // Max atoms per hydrogen group
601  int numMigrationGroups; // Number of migration groups
602  int maxMigrationGroupSize; // Max atoms per migration group
603  int numFixedGroups; // Number of totally fixed hydrogen groups
604  int numRigidBonds; // Number of rigid bonds
605  int numFixedRigidBonds; // Number of rigid bonds between fixed atoms
606 //fepb
607  int numFepInitial; // no. of fep atoms with initial flag
608  int numFepFinal; // no. of fep atoms with final flag
609 //fepe
610 
611  int numConsForce; // Number of atoms that have constant force applied
612  int32 *consForceIndexes;// Constant force indexes for each atom
613  Vector *consForce; // Constant force array
614 
615  int32 *consTorqueIndexes; // "Constant" torque indexes for each atom
616  ConsTorqueParams *consTorqueParams;
617  // Parameters for each atom "constant"-torqued
618 
619  // The following are needed for error checking because we
620  // eliminate bonds, etc. which involve only fixed atoms
621  int numCalcBonds; // Number of bonds requiring calculation
622  int numCalcAngles; // Number of angles requiring calculation
623  int numCalcDihedrals; // Number of dihedrals requiring calculation
624  int numCalcImpropers; // Number of impropers requiring calculation
625  int numCalcCrossterms; // Number of cross-terms requiring calculation
626  int64 numCalcExclusions; // Number of exclusions requiring calculation
627  int64 numCalcFullExclusions; // Number of full exclusions requiring calculation
628 
629  // DRUDE
630  int numCalcTholes; // Number of Thole correction terms requiring calculation
631  int numCalcAnisos; // Number of anisotropic terms requiring calculation
632  // DRUDE
633 
634  // Number of dihedrals with multiple periodicity
636  // Number of impropers with multiple periodicity
638  // indexes of "atoms" sorted by hydrogen groups
640 
641  // Ported by JLai -- JE - Go
642  int numGoAtoms; // Number of atoms subject to Go forces -- ported by JLai/ Original by JE
643  int32 *atomChainTypes; // Go chain type for each atom; from 1 to MAX_GO_CHAINS
644  int32 *goSigmaIndices; // Indices into goSigmas
645  int32 *goResidIndices; // Indices into goSigmas
646  Real *goSigmas; // Sigma values for Go forces L-J type formula
647  bool *goWithinCutoff; // Whether the reference atom-atom distance is within the Go cutoff
648  Real *goCoordinates; // Coordinates (x,y,z) for Go atoms in the native structure
649  int *goResids; // Residue ID from PDB
650  PDB *goPDB; // Pointer to PDB object to use
651  // NAMD-Go2 calculation code
652  int goNumLJPair; // Integer storing the total number of explicit pairs (LJ)
653  int *goIndxLJA; // Pointer to the array of atom indices for LJ atom A
654  int *goIndxLJB; // Pointer to the array of atom indices for LJ atom B
655  double *goSigmaPairA; // Pointer to the array of A LJ parameters
656  double *goSigmaPairB; // Pointer to the array of B LJ parameters
657  int *pointerToGoBeg; // Array of pointers to array
658  int *pointerToGoEnd; // Array of pointers to array
659  // Gromacs LJ Pair list calculation code
660  int numPair; // Integer storing the total number of explicit pairs (LJ + Gaussian)
661  int numLJPair; // Integer storing the number of explicit LJ pairs
662  int numCalcLJPair; // Number of explicit LJ pairs requiring calculation
663  int *pointerToLJBeg; // Array of pointers to array
664  int *pointerToLJEnd; // Array of pointers to array B
665  int *indxLJA; // Pointer to the array of atom indices for LJ atom A
666  int *indxLJB; // Pointer to the array of atom indices for LJ atom B
667  Real *pairC6; // Pointer to the array of C6 LJ parameters
668  Real *pairC12; // Pointer to the array of C12 LJ parameters
670  // Gromacs Gauss Pair list calculation code
671  int *pointerToGaussBeg; // Array of pointers to array B
672  int *pointerToGaussEnd; // Array of pointers to array B
673  int numGaussPair; // Integer storing the number of explicit Gaussian pairs
674  int *indxGaussA; // Pointer to the array of atom indices for Gaussian atom A
675  int *indxGaussB; // Pointer to the array of atom indices for Gaussian atom B
676  Real *gA; // Pointer to the array of force constants to the Gaussian potential
677  Real *gMu1; // Pointer to the array of mu (shifts Gaussian)
678  Real *giSigma1; // Pointer to the array of sigma (controls spread of Gaussian)
679  Real *gMu2; // Pointer to the array of mu (shifts Gaussian 2)
680  Real *giSigma2; // Pointer to the array of sigma (controls spread of Gaussian 2)
681  Real *gRepulsive; // Pointer to the a LJ-12 repulsive parameter that adds to the Gaussian
682 
683  // GO ENERGY CALCULATION CODE
684  BigReal energyNative; // Holds the energy value of the native structure
685  BigReal energyNonnative; // Holds the energy value of the nonnative structure
686  // GO ENERGY CALCULATION CODE
687  // End of port - JL
688 
689 
690 private:
692  // Begin QM
694 
695  int qmDroppedBonds, qmDroppedAngles, qmDroppedDihedrals;
696  int qmDroppedImpropers, qmDroppedCrossterms;
697  ResizeArray<Bond> qmExtraBonds;
698 
699  Bool qmReplaceAll; // Indicates if all forces should be replaced.
700  int qmNumQMAtoms; // Number of QM atoms, from all groups.
701 
702  // Array of abbreviated element names for all atoms.
703  char **qmElementArray;
704  // Array of elements for dummy atoms.
705  char **qmDummyElement;
706 
707  // Number of QM atoms in each QM group
708  int *qmGrpSizes;
709 
710  // QM group for each atom of the system. 0 means MM atom.
711  Real *qmAtomGroup;
712  // Number of different QM regions
713  int qmNumGrps ;
714 
715  // QM Atom charges.
716  Real *qmAtmChrg ;
717  // global indices of all QM atoms.
718  int *qmAtmIndx ;
719 
720  // Number of QM-MM bonds.
721  int qmNumBonds ;
722  // IDs of each group, that is, the value in the qmColumn, per group.
723  Real *qmGrpID ;
724  // The total charge of each QM group.
725  Real *qmGrpChrg;
726  // The multiplicity of QM groups
727  Real *qmGrpMult;
728  // Number of QM-MM bonds in each group.
729  int *qmGrpNumBonds ;
730  // Sequential list of bonds between a MM atom and a QM atom.
731  // (with the bonded atoms ins this order: MM first, QM second).
732  int **qmMMBond ;
733  // List of bond indexes for each QM group.
734  int **qmGrpBonds ;
735  // List of atom indexes for MM atoms in QM-MM bonds, per group.
736  // This is used to check if a point charge is actually an MM atom
737  // that need to be ignored, and will be substituted by a dummy atom (for example).
738  int **qmMMBondedIndx ;
739  // List of indices for MM atoms which will receive fractions of charges from
740  // the MM atom of a QM-MM bond. This is used in the Charge Shift scheme.
741  int **qmMMChargeTarget;
742  // Number of target MM atoms per QM-MM bond. Targets will receive the partial
743  // charge of the MM atom in a QM-MM bond. (in Charge shift scheme)
744  int *qmMMNumTargs;
745  // List of distances from thr QM atom where the dummy atom will be placed.
746  BigReal *qmDummyBondVal;
747  // Indicate if point charges will be used in QM calculations.
748  Bool qmNoPC;
749 
751  // These variables are used in case mechanichal embedding is requested (NoPC = TRUE)
752  // AND there are QM-MM bonds in the system.
753 
754  // Indicates the number of items in the arrays for Mechanical embedding with QM-MM bonds.
755  int qmMeNumBonds;
756  // Indicates the indices of MM atoms which participate in QM-MM bonds.
757  int *qmMeMMindx;
758  // Indicates the QM group to which the QM-MM bonds belongs.
759  Real *qmMeQMGrp;
761 
762  // Indicates frequency of selection of point charges.
763  int qmPCFreq;
764  // Custom PC array;
765  int *qmCustomPCIdxs;
766  // Number of Indexes associated with each QM Group.
767  int *qmCustPCSizes;
768  // Total number of custom PC charges.
769  int qmTotCustPCs;
770  // Number of solvent molecules that will be selected and updated by NAMD, per group.
771  int *qmLSSSize;
772  // Number of atoms per solvent molecule. We need all molecules to have the same size.
773  int qmLSSResidueSize;
774  // IDs of all atoms belonging to QM solvent molecules, from all groups.
775  int *qmLSSIdxs;
776  // MAss of each atom in LSS solvent molecules.
777  Mass *qmLSSMass;
778  // Atom IDs of QM atoms which will be used to select solvent molecules by distance.
779  int *qmLSSRefIDs;
780  // Mass of each QM atom used as reference for solvent selection.
781  Mass *qmLSSRefMass ;
782  // Number of atom IDs/Mass per group.
783  int *qmLSSRefSize;
784  int qmLSSFreq;
785  int qmLSSTotalNumAtms;
786  // This will map each classical solvent atom with a global index of solvent residues,
787  // so we don't depend on segment names and repeated residue indices.
788  std::map<int,int> qmClassicSolv ;
789 
791  // Conditional SMD
792 
793  void read_qm_csdm_file(std::map<Real,int> &qmGrpIDMap) ;
794 
795  Bool qmcSMD;
796  // Total number of cSMD instances.
797  int cSMDnumInst;
798  // Num instances per QM group
799  int* cSMDindxLen;
800  // Index of Conditional SMD guides per QM group
801  int** cSMDindex;
802  // Atom indices for Origin and Target of cSMD
803  int** cSMDpairs;
804  // Spring constants for cSMD
805  Real* cSMDKs;
806  // Speed of movement of guide particles for cSMD.
807  Real* cSMDVels;
808  // Distance cutoff for guide particles for cSMD.
809  Real* cSMDcoffs;
810 
812  // end QM
814 
815 public:
816  // QM
817  void set_qm_replaceAll(Bool newReplaceAll) { qmReplaceAll = newReplaceAll; } ;
818 
819  const Real * get_qmAtomGroup() const {return qmAtomGroup; } ;
820  Real get_qmAtomGroup(int indx) const {return qmAtomGroup[indx]; } ;
821 
822  Real *get_qmAtmChrg() {return qmAtmChrg; } ;
823  const int *get_qmAtmIndx() {return qmAtmIndx; } ;
824 
825  int get_numQMAtoms() {return qmNumQMAtoms; } ;
826  const char *const * get_qmElements() {return qmElementArray;} ;
827  int get_qmNumGrps() const {return qmNumGrps; } ;
828  const int * get_qmGrpSizes() {return qmGrpSizes; } ;
829  Real * get_qmGrpID() { return qmGrpID; } ;
830  Real *get_qmGrpChrg() {return qmGrpChrg; } ;
831  Real *get_qmGrpMult() {return qmGrpMult; } ;
832 
833  int get_qmNumBonds() { return qmNumBonds; } ;
834  const BigReal * get_qmDummyBondVal() { return qmDummyBondVal; } ;
835  const int * get_qmGrpNumBonds() { return qmGrpNumBonds; } ;
836  const int *const * get_qmMMBond() { return qmMMBond; } ;
837  const int *const * get_qmGrpBonds() { return qmGrpBonds; } ;
838  const int *const * get_qmMMBondedIndx() { return qmMMBondedIndx; } ;
839  const char *const * get_qmDummyElement() { return qmDummyElement; } ;
840 
841  const int *const * get_qmMMChargeTarget() { return qmMMChargeTarget; } ;
842  const int * get_qmMMNumTargs() { return qmMMNumTargs; } ;
843 
844  int* get_qmLSSSize() { return qmLSSSize;} ;
845  int* get_qmLSSIdxs() { return qmLSSIdxs;} ;
846  Mass *get_qmLSSMass() { return qmLSSMass; } ;
847  int get_qmLSSFreq() { return qmLSSFreq; } ;
848  int get_qmLSSResSize() { return qmLSSResidueSize; } ;
849  int *get_qmLSSRefIDs() { return qmLSSRefIDs ; } ;
850  int *get_qmLSSRefSize() { return qmLSSRefSize ; } ;
851  Mass *get_qmLSSRefMass() {return qmLSSRefMass; } ;
852  std::map<int,int> &get_qmMMSolv() { return qmClassicSolv;} ;
853 
854  Bool get_qmReplaceAll() {return qmReplaceAll; } ;
855 
856  Bool get_noPC() { return qmNoPC; } ;
857  int get_qmMeNumBonds() { return qmMeNumBonds; } ;
858  int *get_qmMeMMindx() { return qmMeMMindx; } ;
859  Real *get_qmMeQMGrp() { return qmMeQMGrp; } ;
860 
861  int get_qmPCFreq() { return qmPCFreq; } ;
862 
863  int get_qmTotCustPCs() { return qmTotCustPCs; } ;
864  int *get_qmCustPCSizes() { return qmCustPCSizes; } ;
865  int *get_qmCustomPCIdxs() { return qmCustomPCIdxs; } ;
866 
867  Bool get_qmcSMD() { return qmcSMD;} ;
868  int get_cSMDnumInst() { return cSMDnumInst;} ;
869  int const * get_cSMDindxLen() { return cSMDindxLen;} ;
870  int const * const * get_cSMDindex() {return cSMDindex;} ;
871  int const * const * get_cSMDpairs() {return cSMDpairs;} ;
872  const Real* get_cSMDKs() {return cSMDKs;} ;
873  const Real* get_cSMDVels() {return cSMDVels;} ;
874  const Real* get_cSMDcoffs() {return cSMDcoffs;} ;
875 
876  void prepare_qm(const char *pdbFileName, Parameters *params, ConfigList *cfgList) ;
877  void delete_qm_bonded() ;
878  // end QM
879 
880 
881 
882  Molecule(SimParameters *, Parameters *param);
883  Molecule(SimParameters *, Parameters *param, char *filename, ConfigList *cfgList=NULL);
884  Molecule(SimParameters *simParams, Parameters *param, molfile_plugin_t *pIOHdl, void *pIOFileHdl, int natoms);
885 
887  void read_parm(Ambertoppar *);
888 
890 
891  ~Molecule(); // Destructor
892 
893  void send_Molecule(MOStream *);
894  // send the molecular structure
895  // from the master to the clients
896 
897  void receive_Molecule(MIStream *);
898  // receive the molecular structure
899  // from the master on a client
900 
902  PDB *, char *);
903  // Build the set of harmonic constraint
904  // parameters
905 
906 /* BEGIN gf */
908  // Build the set of gridForce-style force pars
909 /* END gf */
910 
912  PDB *, char *);
913  // Build the set of moving drag pars
914 
917  PDB *, char *);
918  // Build the set of rotating drag pars
919 
922  PDB *, char *);
923  // Build the set of "constant" torque pars
924 
925 
926  void build_constant_forces(char *);
927  // Build the set of constant forces
928 
929  void build_langevin_params(BigReal coupling, BigReal drudeCoupling,
930  Bool doHydrogen);
931  void build_langevin_params(StringList *, StringList *, PDB *, char *);
932  // Build the set of langevin dynamics parameters
933 
934 #ifdef MEM_OPT_VERSION
935  void load_fixed_atoms(StringList *fixedFile);
936  void load_constrained_atoms(StringList *constrainedFile);
937 #endif
938 
939  void build_fixed_atoms(StringList *, StringList *, PDB *, char *);
940  // Determine which atoms are fixed (if any)
941 
942  void build_stirred_atoms(StringList *, StringList *, PDB *, char *);
943  // Determine which atoms are stirred (if any)
944 
945  void build_extra_bonds(Parameters *parameters, StringList *file);
946 
947 //fepb
948  void build_fep_flags(StringList *, StringList *, PDB *, char *, const char *);
949  // selection of the mutant atoms
950  void delete_alch_bonded(void);
951 
952  void build_alch_unpert_bond_lists(char *);
953  void read_alch_unpert_bonds(FILE *);
954  void read_alch_unpert_angles(FILE *);
955  void read_alch_unpert_dihedrals(FILE *);
956 //fepe
957 
967  void build_ss_flags(
968  const StringList *ssfile,
970  const StringList *sscol,
972  PDB *initial_pdb,
973  const char *cwd
974  );
975 
976  void build_exPressure_atoms(StringList *, StringList *, PDB *, char *);
977  // Determine which atoms are excluded from
978  // pressure (if any)
979 
980  // Ported by JLai -- Original JE - Go -- Change the unsigned int to ints
981  void print_go_sigmas(); // Print out Go sigma parameters
982  void build_go_sigmas(StringList *, char *);
983  // Determine which atoms have Go forces applied
984  // calculate sigmas from distances between Go atom pairs
985  void build_go_sigmas2(StringList *, char *);
986  // Determine which atoms have Go forces applied
987  // calculate sigmas from distances between Go atom pairs
988  void build_go_arrays(StringList *, char *);
989  // Determine which atoms have Go forces applied
990  BigReal get_gro_force(BigReal, BigReal, BigReal, int, int) const;
991  BigReal get_gro_force2(BigReal, BigReal, BigReal, int, int, BigReal *, BigReal *) const;
992  BigReal get_go_force(BigReal, int, int, BigReal *, BigReal *) const;
993  // Calculate the go force between a pair of atoms -- Modified to
994  // output Go energies
995  BigReal get_go_force_new(BigReal, int, int, BigReal *, BigReal *) const;
996  // Calculate the go force between a pair of atoms
997  BigReal get_go_force2(BigReal, BigReal, BigReal, int, int, BigReal *, BigReal *) const;
998  // Calculate the go force between a pair of atoms
999  Bool atoms_1to4(unsigned int, unsigned int);
1000 // End of port -- JL
1001 
1002  void reloadCharges(float charge[], int n);
1003 
1004  Bool is_lp(int); // return true if atom is a lone pair
1005  Bool is_drude(int); // return true if atom is a Drude particle
1006  Bool is_hydrogen(int); // return true if atom is hydrogen
1007  Bool is_oxygen(int); // return true if atom is oxygen
1008  Bool is_hydrogenGroupParent(int); // return true if atom is group parent
1009  Bool is_water(int); // return true if atom is part of water
1010  int get_groupSize(int); // return # atoms in (hydrogen) group
1011  int get_mother_atom(int); // return mother atom of a hydrogen
1012 
1013  #ifdef MEM_OPT_VERSION
1014  //the way to get the cluster size if the atom ids of the cluster are
1015  //contiguous. The input parameter is the atom id that leads the cluster.
1016  int get_cluster_size_con(int aid) const { return clusterSigs[aid]; }
1017  //the way to get the cluster size if the atoms ids of the cluster are
1018  //not contiguous. The input parameter is the cluster index.
1019  int get_cluster_size_uncon(int cIdx) const { return clusterSize[cIdx]; }
1020  int get_cluster_idx(int aid) const { return clusterSigs[aid]; }
1021  int get_num_clusters() const { return numClusters; }
1022  #else
1023  int get_cluster(int anum) const { return cluster[anum]; }
1024  int get_clusterSize(int anum) const { return clusterSize[anum]; }
1025  #endif
1026 
1027 #ifndef MEM_OPT_VERSION
1028  const float *getOccupancyData() { return (const float *)occupancy; }
1029  void setOccupancyData(molfile_atom_t *atomarray);
1030  void freeOccupancyData() { delete [] occupancy; occupancy=NULL; }
1031 
1032  const float *getBFactorData() { return (const float *)bfactor; }
1033  void setBFactorData(molfile_atom_t *atomarray);
1034  void freeBFactorData() { delete [] bfactor; bfactor=NULL; }
1035 #endif
1036 
1037  // Get the mass of an atom
1038  Real atommass(int anum) const
1039  {
1040  #ifdef MEM_OPT_VERSION
1041  return atomMassPool[eachAtomMass[anum]];
1042  #else
1043  return(atoms[anum].mass);
1044  #endif
1045  }
1046 
1047  // Get the charge of an atom
1048  Real atomcharge(int anum) const
1049  {
1050  #ifdef MEM_OPT_VERSION
1051  return atomChargePool[eachAtomCharge[anum]];
1052  #else
1053  return(atoms[anum].charge);
1054  #endif
1055  }
1056 
1057  // Get the vdw type of an atom
1058  Index atomvdwtype(int anum) const
1059  {
1060  return(atoms[anum].vdw_type);
1061  }
1062 
1063  #ifndef MEM_OPT_VERSION
1064  // Retrieve a bond structure
1065  Bond *get_bond(int bnum) const {return (&(bonds[bnum]));}
1066 
1067  // Retrieve an angle structure
1068  Angle *get_angle(int anum) const {return (&(angles[anum]));}
1069 
1070  // Retrieve an improper strutcure
1071  Improper *get_improper(int inum) const {return (&(impropers[inum]));}
1072 
1073  // Retrieve a dihedral structure
1074  Dihedral *get_dihedral(int dnum) const {return (&(dihedrals[dnum]));}
1075 
1076  // Retrieve a cross-term strutcure
1077  Crossterm *get_crossterm(int inum) const {return (&(crossterms[inum]));}
1078  #endif
1079 
1080  // DRUDE: retrieve lphost structure
1081  Lphost *get_lphost(int atomid) const {
1082  // don't call unless simParams->drudeOn == TRUE
1083  // otherwise lphostIndexes array doesn't exist!
1084  int index = lphostIndexes[atomid];
1085  return (index != -1 ? &(lphosts[index]) : NULL);
1086  }
1087  // DRUDE
1088 
1089  #ifndef MEM_OPT_VERSION
1090  Bond *getAllBonds() const {return bonds;}
1091  Angle *getAllAngles() const {return angles;}
1092  Improper *getAllImpropers() const {return impropers;}
1093  Dihedral *getAllDihedrals() const {return dihedrals;}
1094  Crossterm *getAllCrossterms() const {return crossterms;}
1095  #endif
1096 
1097  // DRUDE: retrieve entire lphosts array
1098  Lphost *getAllLphosts() const { return lphosts; }
1099  // DRUDE
1100 
1101  // Retrieve a hydrogen bond donor structure
1102  Bond *get_donor(int dnum) const {return (&(donors[dnum]));}
1103 
1104  // Retrieve a hydrogen bond acceptor structure
1105  Bond *get_acceptor(int dnum) const {return (&(acceptors[dnum]));}
1106 
1107  Bond *getAllDonors() const {return donors;}
1108  Bond *getAllAcceptors() const {return acceptors;}
1109 
1110  // Retrieve an exclusion structure
1111  #ifndef MEM_OPT_VERSION
1112  Exclusion *get_exclusion(int ex) const {return (&(exclusions[ex]));}
1113  #endif
1114 
1115  // Retrieve an atom type
1116  const char *get_atomtype(int anum) const
1117  {
1118  if (atomNames == NULL)
1119  {
1120  NAMD_die("Tried to find atom type on node other than node 0");
1121  }
1122 
1123  #ifdef MEM_OPT_VERSION
1124  return atomTypePool[atomNames[anum].atomtypeIdx];
1125  #else
1126  return(atomNames[anum].atomtype);
1127  #endif
1128  }
1129 
1130  // Lookup atom id from segment, residue, and name
1131  int get_atom_from_name(const char *segid, int resid, const char *aname) const;
1132 
1133  // Lookup number of atoms in residue from segment and residue
1134  int get_residue_size(const char *segid, int resid) const;
1135 
1136  // Lookup atom id from segment, residue, and index in residue
1137  int get_atom_from_index_in_residue(const char *segid, int resid, int index) const;
1138 
1139 
1140  // The following routines are used to get the list of bonds
1141  // for a given atom. This is used when creating the bond lists
1142  // for the force objects
1143 
1144  #ifndef MEM_OPT_VERSION
1146  { return bondsByAtom[anum]; }
1148  { return anglesByAtom[anum]; }
1150  { return dihedralsByAtom[anum]; }
1152  { return impropersByAtom[anum]; }
1154  { return crosstermsByAtom[anum]; }
1156  { return exclusionsByAtom[anum]; }
1157  const int32 *get_full_exclusions_for_atom(int anum) const
1158  { return fullExclusionsByAtom[anum]; }
1159  const int32 *get_mod_exclusions_for_atom(int anum) const
1160  { return modExclusionsByAtom[anum]; }
1161  #endif
1162 
1163  // Check for exclusions, either explicit or bonded.
1164  // Returns 1 for full, 2 for 1-4 exclusions.
1165  #ifdef MEM_OPT_VERSION
1166  int checkExclByIdx(int idx1, int atom1, int atom2) const;
1167  const ExclusionCheck *get_excl_check_for_idx(int idx) const{
1168  return &exclChkSigPool[idx];
1169  }
1170  #else
1171  int checkexcl(int atom1, int atom2) const;
1172 
1173  const ExclusionCheck *get_excl_check_for_atom(int anum) const{
1174  return &all_exclusions[anum];
1175  }
1176  #endif
1177 
1178 /* BEGIN gf */
1179  // Return true or false based on whether or not the atom
1180  // is subject to grid force
1181  Bool is_atom_gridforced(int atomnum, int gridnum) const
1182  {
1183  if (numGridforceGrids)
1184  {
1185  return(gridfrcIndexes[gridnum][atomnum] != -1);
1186  }
1187  else
1188  {
1189  return(FALSE);
1190  }
1191  }
1192 /* END gf */
1193 
1194 #ifndef MEM_OPT_VERSION
1195  // Return true or false based on whether the specified atom
1196  // is constrained or not.
1197  Bool is_atom_constrained(int atomnum) const
1198  {
1199  if (numConstraints)
1200  {
1201  // Check the index to see if it is constrained
1202  return(consIndexes[atomnum] != -1);
1203  }
1204  else
1205  {
1206  // No constraints at all, so just return FALSE
1207  return(FALSE);
1208  }
1209  }
1210 #endif
1211 
1212  // Return true or false based on whether the specified atom
1213  // is moving-dragged or not.
1214  Bool is_atom_movdragged(int atomnum) const
1215  {
1216  if (numMovDrag)
1217  {
1218  // Check the index to see if it is constrained
1219  return(movDragIndexes[atomnum] != -1);
1220  }
1221  else
1222  {
1223  // No constraints at all, so just return FALSE
1224  return(FALSE);
1225  }
1226  }
1227 
1228  // Return true or false based on whether the specified atom
1229  // is rotating-dragged or not.
1230  Bool is_atom_rotdragged(int atomnum) const
1231  {
1232  if (numRotDrag)
1233  {
1234  // Check the index to see if it is constrained
1235  return(rotDragIndexes[atomnum] != -1);
1236  }
1237  else
1238  {
1239  // No constraints at all, so just return FALSE
1240  return(FALSE);
1241  }
1242  }
1243 
1244  // Return true or false based on whether the specified atom
1245  // is "constant"-torqued or not.
1246  Bool is_atom_constorqued(int atomnum) const
1247  {
1248  if (numConsTorque)
1249  {
1250  // Check the index to see if it is constrained
1251  return(consTorqueIndexes[atomnum] != -1);
1252  }
1253  else
1254  {
1255  // No constraints at all, so just return FALSE
1256  return(FALSE);
1257  }
1258  }
1259 
1260 #ifndef MEM_OPT_VERSION
1261  // Get the harmonic constraints for a specific atom
1262  void get_cons_params(Real &k, Vector &refPos, int atomnum) const
1263  {
1264  k = consParams[consIndexes[atomnum]].k;
1265  refPos = consParams[consIndexes[atomnum]].refPos;
1266  }
1267 #endif
1268 
1269 /* BEGIN gf */
1270  void get_gridfrc_params(Real &k, Charge &q, int atomnum, int gridnum) const
1271  {
1272  k = gridfrcParams[gridnum][gridfrcIndexes[gridnum][atomnum]].k;
1273  q = gridfrcParams[gridnum][gridfrcIndexes[gridnum][atomnum]].q;
1274  }
1275 
1276  GridforceGrid* get_gridfrc_grid(int gridnum) const
1277  {
1278  GridforceGrid *result = NULL;
1279  if (gridnum >= 0 && gridnum < numGridforceGrids) {
1280  result = gridfrcGrid[gridnum];
1281  }
1282  return result;
1283  }
1284 
1285  int set_gridfrc_grid(int gridnum, GridforceGrid *grid)
1286  {
1287  if (grid && gridnum >= 0 && gridnum < numGridforceGrids) {
1288  gridfrcGrid[gridnum] = grid;
1289  return 0;
1290  } else {
1291  return -1;
1292  }
1293  }
1294 /* END gf */
1295 
1296  Real langevin_param(int atomnum) const
1297  {
1298  return(langevinParams ? langevinParams[atomnum] : 0.);
1299  }
1300 
1301  // Get the stirring constraints for a specific atom
1302  void get_stir_refPos(Vector &refPos, int atomnum) const
1303  {
1304  refPos = stirParams[stirIndexes[atomnum]].refPos;
1305  }
1306 
1307 
1308  void put_stir_startTheta(Real theta, int atomnum) const
1309  {
1310  stirParams[stirIndexes[atomnum]].startTheta = theta;
1311  }
1312 
1313 
1314  Real get_stir_startTheta(int atomnum) const
1315  {
1316  return stirParams[stirIndexes[atomnum]].startTheta;
1317  }
1318 
1319 
1320  // Get the moving drag factor for a specific atom
1321  void get_movdrag_params(Vector &v, int atomnum) const
1322  {
1323  v = movDragParams[movDragIndexes[atomnum]].v;
1324  }
1325 
1326  // Get the rotating drag pars for a specific atom
1328  int atomnum) const
1329  {
1330  v = rotDragParams[rotDragIndexes[atomnum]].v;
1331  a = rotDragParams[rotDragIndexes[atomnum]].a;
1332  p = rotDragParams[rotDragIndexes[atomnum]].p;
1333  }
1334 
1335  // Get the "constant" torque pars for a specific atom
1337  int atomnum) const
1338  {
1339  v = consTorqueParams[consTorqueIndexes[atomnum]].v;
1340  a = consTorqueParams[consTorqueIndexes[atomnum]].a;
1341  p = consTorqueParams[consTorqueIndexes[atomnum]].p;
1342  }
1343 
1344 //fepb
1345  unsigned char get_fep_type(int anum) const
1346  {
1347  return(fepAtomFlags[anum]);
1348  }
1349 
1350 //soluteScaling
1351  unsigned char get_ss_type(int anum) const
1352  {
1353  return(ssAtomFlags[anum]);
1354  }
1355 //soluteScaling
1356 
1357  /* BKR - Get the FEP type (i.e. 0, 1, or 2) of a bonded _interaction_ based
1358  on the atom indices of the atoms involved (internally converted to FEP
1359  types) and the order of the bonded interaction (i.e. 2, 3, or 4). This
1360  automatically accounts for whether or not purely alchemical interactions
1361  are being scaled (e.g. bonds between two atoms of fep_type 1).
1362 
1363  The logic here is admittedly a bit opaque. When the fep_type's are back
1364  mapped to -1,0,1, we can use the sum of all of the types to determine the
1365  type of interaction for all bonded term types:
1366 
1367  0 - only non-alchemical atoms involved
1368  >0 - _atleast_ one appearing atom
1369  <0 - _atleast_ one disappearing atom
1370 
1371  If the magnitude of the sum is equal to the interaction order (i.e. 2, 3,
1372  or 4), then it is a _purely_ alchemical interaction and it may be
1373  desireable to retain it for sampling purposes (note that this adds a
1374  constant to the free energy that will almost always differ at the
1375  endpoints). In order to avoid unexpected instability these interactions
1376  are retained by default, but can be scaled along with everything else by
1377  setting "alchBondDecouple on".
1378 
1379  NB: All pure alchemical interactions beyond order 2 are ALWAYS discarded
1380  by alchemify. This is good, because higher order interactions would break
1381  the above logic. For example, if we had an angle term between atoms of
1382  types (1,1,-1) the sum would be 1, but this term should receive no scaling
1383  because it involves groups -1 and 1 but not 0.
1384  */
1385  int get_fep_bonded_type(const int *atomID, unsigned int order) const
1386  {
1387  int typeSum = 0;
1388  for ( int i=0; i < order; ++i ) {
1389  typeSum += (fepAtomFlags[atomID[i]] == 2 ? -1 : fepAtomFlags[atomID[i]]);
1390  }
1391  // Increase the cutoff if scaling purely alchemical bonds.
1392  // This really only applies when order = 2.
1393  if ( simParams->alchBondDecouple ) order++;
1394 
1395  // The majority of interactions are type 0, so catch those first.
1396  if ( typeSum == 0 || abs(typeSum) == order ) return 0;
1397  else if ( 0 < typeSum && typeSum < order ) return 1;
1398  else if ( -order < typeSum && typeSum < 0 ) return 2;
1399 
1400  // Alchemify should always keep this from bombing, but just in case...
1401  NAMD_die("Unexpected alchemical bonded interaction!");
1402  return 0;
1403  }
1404 //fepe
1405 
1406 #ifndef MEM_OPT_VERSION
1407  Bool is_atom_fixed(int atomnum) const
1408  {
1409  return (numFixedAtoms && fixedAtomFlags[atomnum]);
1410  }
1411 #else
1412  //Since binary search is more expensive than direct array access,
1413  //and this function is usually called for consecutive atoms in this context,
1414  //the *listIdx returns the index to the range of atoms [aid1, aid2]
1415  //that are fixed. If the atom aid is fixed, then aid1=<aid<=aid2;
1416  //If the atom aid is not fixed, then aid1 indicates the smallest fixed atom
1417  //id that is larger than aid; so the listIdx could be equal the size of
1418  //fixedAtomsSet. --Chao Mei
1419  Bool is_atom_in_set(AtomSetList *localAtomsSet, int aid, int *listIdx) const;
1420  inline Bool is_atom_fixed(int aid, int *listIdx=NULL) const {
1421  return is_atom_in_set(fixedAtomsSet,aid,listIdx);
1422  }
1423  inline Bool is_atom_constrained(int aid, int *listIdx=NULL) const {
1424  return is_atom_in_set(constrainedAtomsSet,aid,listIdx);
1425  }
1426 #endif
1427 
1428  Bool is_atom_stirred(int atomnum) const
1429  {
1430  if (numStirredAtoms)
1431  {
1432  // Check the index to see if it is constrained
1433  return(stirIndexes[atomnum] != -1);
1434  }
1435  else
1436  {
1437  // No constraints at all, so just return FALSE
1438  return(FALSE);
1439  }
1440  }
1441 
1442 
1443  Bool is_group_fixed(int atomnum) const
1444  {
1445  return (numFixedAtoms && (fixedAtomFlags[atomnum] == -1));
1446  }
1447  Bool is_atom_exPressure(int atomnum) const
1448  {
1449  return (numExPressureAtoms && exPressureAtomFlags[atomnum]);
1450  }
1451  // 0 if not rigid or length to parent, for parent refers to H-H length
1452  // < 0 implies MOLLY but not SHAKE, > 1 implies both if MOLLY is on
1453  Real rigid_bond_length(int atomnum) const
1454  {
1455  return(rigidBondLengths[atomnum]);
1456  }
1457 
1458  void print_atoms(Parameters *);
1459  // Print out list of atoms
1460  void print_bonds(Parameters *);
1461  // Print out list of bonds
1462  void print_exclusions();// Print out list of exclusions
1463 
1464 public:
1466 
1467 #ifdef MEM_OPT_VERSION
1468  //read the per-atom file for the memory optimized version where the file
1469  //name already exists the range of atoms to be read are
1470  //[fromAtomID, toAtomID].
1471  void read_binary_atom_info(int fromAtomID, int toAtomID, InputAtomList &inAtoms);
1472 
1473  int64 getNumCalcExclusions(){return numCalcExclusions;}
1474  void setNumCalcExclusions(int64 x){numCalcExclusions= x;}
1475 
1476  Index getEachAtomMass(int i){return eachAtomMass[i];}
1477  Index getEachAtomCharge(int i){return eachAtomCharge[i];}
1478 
1479  ExclSigID getAtomExclSigId(int aid) const {
1480  return eachAtomExclSig[aid];
1481  }
1482 
1483  Real *getAtomMassPool(){return atomMassPool;}
1484  Real *getAtomChargePool(){return atomChargePool;}
1485  AtomCstInfo *getAtoms(){return atoms;}
1486 
1487  int atomSigPoolSize;
1489 
1490  /* All the following are temporary variables for reading the compressed psf file */
1491  //declarations for atoms' constant information
1492  int segNamePoolSize; //Its value is usually less than 5
1493  char **segNamePool; //This seems not to be important, but it only occupied very little space.
1494 
1495  int resNamePoolSize;
1496  char **resNamePool;
1497 
1498  int atomNamePoolSize;
1499  char **atomNamePool;
1500 
1501  int atomTypePoolSize;
1502  char **atomTypePool;
1503 
1504  int chargePoolSize;
1505  Real *atomChargePool;
1506 
1507  int massPoolSize;
1508  Real *atomMassPool;
1509 
1510  AtomSigID getAtomSigId(int aid) {
1511  return eachAtomSig[aid];
1512  }
1513 
1514  //Indicates the size of both exclSigPool and exclChkSigPool
1515  int exclSigPoolSize;
1516  //this will be deleted after build_lists_by_atom
1517  ExclusionSignature *exclSigPool;
1518  //This is the final data structure we want to store
1519  ExclusionCheck *exclChkSigPool;
1520 
1521  void addNewExclSigPool(const std::vector<ExclusionSignature>&);
1522 
1523  void delEachAtomSigs(){
1524  //for NAMD-smp version, only one Molecule object is held
1525  //on each node, therefore, only one deletion operation should
1526  //be taken on a node, otherwise, there possibly would be some
1527  //wierd memory problems. The same reason applies to other deletion
1528  //operations inside the Molecule object.
1529  if(CmiMyRank()) return;
1530 
1531  delete [] eachAtomSig;
1532  delete [] eachAtomExclSig;
1533  eachAtomSig = NULL;
1534  eachAtomExclSig = NULL;
1535  }
1536 
1537  void delChargeSpace(){
1538  if(CmiMyRank()) return;
1539 
1540  delete [] atomChargePool;
1541  delete [] eachAtomCharge;
1542  atomChargePool = NULL;
1543  eachAtomCharge = NULL;
1544  }
1545 
1546  void delMassSpace(){
1547  if(CmiMyRank()) return;
1548 
1549  delete [] atomMassPool;
1550  delete [] eachAtomMass;
1551  atomMassPool = NULL;
1552  eachAtomMass = NULL;
1553  }
1554 
1555  void delClusterSigs() {
1556  if(CmiMyRank()) return;
1557 
1558  delete [] clusterSigs;
1559  clusterSigs = NULL;
1560  }
1561 
1562  void delAtomNames(){
1563  if(CmiMyRank()) return;
1564  delete [] atomNamePool;
1565  delete [] atomNames;
1566  atomNamePool = NULL;
1567  atomNames = NULL;
1568  }
1569 
1570  void delFixedAtoms(){
1571  if(CmiMyRank()) return;
1572  delete fixedAtomsSet;
1573  fixedAtomsSet = NULL;
1574  }
1575 
1576 private:
1577  Index insert_new_mass(Real newMass);
1578 
1579 #endif
1580 
1581 // Go stuff
1582 public:
1583 
1584 GoValue go_array[MAX_GO_CHAINS*MAX_GO_CHAINS]; // Array of Go params -- JLai
1585 int go_indices[MAX_GO_CHAINS+1]; // Indices from chainIDS to go array -- JLai
1586 int NumGoChains; // Number of Go chain types -- JLai
1587 
1588 // Declares and initializes Go variables
1589 void goInit();
1590 
1591 // Builds explicit Gromacs pairs
1592 void build_gro_pair();
1593 
1594 // Builds the initial Go parameters
1595 void build_go_params(StringList *);
1596 
1597 // Read Go parameter file
1598 void read_go_file(char *);
1599 
1600 // Get Go cutoff for a given chain type pair
1601 Real get_go_cutoff(int chain1, int chain2) { return go_array[MAX_GO_CHAINS*chain1 + chain2].cutoff; };
1602 
1603 // Get Go epsilonRep for a given chain type pair
1604 Real get_go_epsilonRep(int chain1, int chain2) { return go_array[MAX_GO_CHAINS*chain1 + chain2].epsilonRep; };
1605 
1606 // Get Go sigmaRep for a given chain type pair
1607 Real get_go_sigmaRep(int chain1, int chain2) { return go_array[MAX_GO_CHAINS*chain1 + chain2].sigmaRep; };
1608 
1609 // Get Go epsilon for a given chain type pair
1610 Real get_go_epsilon(int chain1, int chain2) { return go_array[MAX_GO_CHAINS*chain1 + chain2].epsilon; };
1611 
1612 // Get Go exp_a for a given chain type pair
1613 int get_go_exp_a(int chain1, int chain2) { return go_array[MAX_GO_CHAINS*chain1 + chain2].exp_a; };
1614 
1615 // Get Go exp_b for a given chain type pair
1616 int get_go_exp_b(int chain1, int chain2) { return go_array[MAX_GO_CHAINS*chain1 + chain2].exp_b; };
1617 
1618 // Get Go exp_rep for a given chain type pair
1619 int get_go_exp_rep(int chain1, int chain2) { return go_array[MAX_GO_CHAINS*chain1 + chain2].exp_rep; };
1620 
1621 // Whether residue IDs with this difference are restricted
1622 Bool go_restricted(int, int, int);
1623 
1624 // Prints Go Params
1625 void print_go_params();
1626 
1627 void initialize();
1628 
1629 void send_GoMolecule(MOStream *);
1630 // send the molecular structure
1631 // from the master to the clients
1632 
1634 // receive the molecular structure
1635 // from the master on a client
1636 };
1637 
1638 #endif
1639 
int32 * get_exclusions_for_atom(int anum)
Definition: Molecule.h:1155
BigReal get_gro_force(BigReal, BigReal, BigReal, int, int) const
void print_go_params()
Definition: GoMolecule.C:548
int get_atom_from_name(const char *segid, int resid, const char *aname) const
Definition: Molecule.C:121
int numFixedGroups
Definition: Molecule.h:603
int exp_b
Definition: Molecule.h:108
int get_residue_size(const char *segid, int resid) const
Definition: Molecule.C:144
void build_gridforce_params(StringList *, StringList *, StringList *, StringList *, PDB *, char *)
Definition: Molecule.C:6373
BigReal tail_corr_virial_2
Definition: Molecule.h:473
int suspiciousAlchBonds
Definition: Molecule.h:562
Bool get_noPC()
Definition: Molecule.h:856
int * get_qmLSSRefSize()
Definition: Molecule.h:850
int num_alch_unpert_Dihedrals
Definition: Molecule.h:574
int get_go_exp_a(int chain1, int chain2)
Definition: Molecule.h:1613
Bond * getAllAcceptors() const
Definition: Molecule.h:1108
Real langevin_param(int atomnum) const
Definition: Molecule.h:1296
void get_stir_refPos(Vector &refPos, int atomnum) const
Definition: Molecule.h:1302
int * indxGaussA
Definition: Molecule.h:674
Real get_go_cutoff(int chain1, int chain2)
Definition: Molecule.h:1601
int32 * get_dihedrals_for_atom(int anum)
Definition: Molecule.h:1149
const int * get_qmMMNumTargs()
Definition: Molecule.h:842
Definition: PDB.h:35
int numCalcBonds
Definition: Molecule.h:621
Real * gMu2
Definition: Molecule.h:679
int alchDroppedDihedrals
Definition: Molecule.h:564
int32 * get_crossterms_for_atom(int anum)
Definition: Molecule.h:1153
int numCalcAnisos
Definition: Molecule.h:631
int numBonds
Definition: Molecule.h:559
double A
Definition: Molecule.h:120
Real * get_qmGrpMult()
Definition: Molecule.h:831
void print_bonds(Parameters *)
Definition: Molecule.C:5380
char segname[11]
Definition: Molecule.h:145
int get_qmNumGrps() const
Definition: Molecule.h:827
void send_GoMolecule(MOStream *)
Definition: GoMolecule.C:1635
int * pointerToGoBeg
Definition: Molecule.h:657
int numCalcLJPair
Definition: Molecule.h:662
Real * giSigma2
Definition: Molecule.h:680
struct go_pair GoPair
Bool is_atom_movdragged(int atomnum) const
Definition: Molecule.h:1214
PDB * goPDB
Definition: Molecule.h:650
short int32
Definition: dumpdcd.c:24
int AtomID
Definition: NamdTypes.h:29
const int * get_qmGrpNumBonds()
Definition: Molecule.h:835
BigReal get_go_force_new(BigReal, int, int, BigReal *, BigReal *) const
Definition: GoMolecule.C:1334
void build_extra_bonds(Parameters *parameters, StringList *file)
void setOccupancyData(molfile_atom_t *atomarray)
Definition: Molecule.C:3236
const int * get_qmAtmIndx()
Definition: Molecule.h:823
void print_exclusions()
Definition: Molecule.C:5417
Real epsilonRep
Definition: Molecule.h:111
BigReal get_gro_force2(BigReal, BigReal, BigReal, int, int, BigReal *, BigReal *) const
Definition: GoMolecule.C:1173
Bool is_hydrogen(int)
Real * goCoordinates
Definition: Molecule.h:648
BigReal tail_corr_virial_1
Definition: Molecule.h:473
int numHydrogenGroups
Definition: Molecule.h:599
void receive_GoMolecule(MIStream *)
Definition: GoMolecule.C:1744
void read_alch_unpert_dihedrals(FILE *)
Definition: Molecule.C:1931
int numAnisos
Number of anisotropic terms.
Definition: Molecule.h:583
Bool is_hydrogenGroupParent(int)
#define MAX_RESTRICTIONS
Definition: Molecule.h:37
Angle * alch_unpert_angles
Definition: Molecule.h:576
int get_cSMDnumInst()
Definition: Molecule.h:868
void send_Molecule(MOStream *)
Definition: Molecule.C:5448
Real get_go_epsilonRep(int chain1, int chain2)
Definition: Molecule.h:1604
Definition: Vector.h:64
BigReal tail_corr_dUdl_1
Definition: Molecule.h:472
Real * pairC6
Definition: Molecule.h:667
const int32 * get_mod_exclusions_for_atom(int anum) const
Definition: Molecule.h:1159
int get_numQMAtoms()
Definition: Molecule.h:825
Improper * get_improper(int inum) const
Definition: Molecule.h:1071
static __thread unsigned int * exclusions
static __thread atom * atoms
Mass * get_qmLSSMass()
Definition: Molecule.h:846
int num_alch_unpert_Bonds
Definition: Molecule.h:572
int goIndxB
Definition: Molecule.h:119
Atom * getAtoms() const
Definition: Molecule.h:490
int alchDroppedImpropers
Definition: Molecule.h:565
struct go_val GoValue
char * flags
Definition: Molecule.h:71
float Real
Definition: common.h:107
const int *const * get_qmMMBondedIndx()
Definition: Molecule.h:838
Real get_stir_startTheta(int atomnum) const
Definition: Molecule.h:1314
double * goSigmaPairA
Definition: Molecule.h:655
#define FALSE
Definition: common.h:116
static __thread int2 * exclusionsByAtom
Real rigid_bond_length(int atomnum) const
Definition: Molecule.h:1453
int get_atom_from_index_in_residue(const char *segid, int resid, int index) const
Definition: Molecule.C:158
int * pointerToLJEnd
Definition: Molecule.h:664
int numGridforceGrids
Definition: Molecule.h:590
void setBFactorData(molfile_atom_t *atomarray)
Definition: Molecule.C:3243
int numRealBonds
Definition: Molecule.h:558
Bond * getAllDonors() const
Definition: Molecule.h:1107
Real get_qmAtomGroup(int indx) const
Definition: Molecule.h:820
int num_alch_unpert_Angles
Definition: Molecule.h:573
Crossterm * getAllCrossterms() const
Definition: Molecule.h:1094
int ExclSigID
Definition: NamdTypes.h:83
bool * goWithinCutoff
Definition: Molecule.h:647
int get_cluster(int anum) const
Definition: Molecule.h:1023
void freeOccupancyData()
Definition: Molecule.h:1030
void reloadCharges(float charge[], int n)
int get_fep_bonded_type(const int *atomID, unsigned int order) const
Definition: Molecule.h:1385
int const * get_cSMDindxLen()
Definition: Molecule.h:869
int numMovDrag
Definition: Molecule.h:593
Bool get_qmcSMD()
Definition: Molecule.h:867
int num_fixed_atoms() const
Definition: Molecule.h:497
int * get_qmLSSRefIDs()
Definition: Molecule.h:849
BigReal tail_corr_ener
Definition: Molecule.h:470
void get_gridfrc_params(Real &k, Charge &q, int atomnum, int gridnum) const
Definition: Molecule.h:1270
struct seg_resid AtomSegResInfo
ResidueLookupElem * next
Definition: Molecule.h:92
int * get_qmLSSIdxs()
Definition: Molecule.h:845
std::map< int, int > & get_qmMMSolv()
Definition: Molecule.h:852
Dihedral * alch_unpert_dihedrals
Definition: Molecule.h:577
Real * gA
Definition: Molecule.h:676
#define MAX_GO_CHAINS
Definition: Molecule.h:36
const char *const * get_qmDummyElement()
Definition: Molecule.h:839
void read_go_file(char *)
Definition: GoMolecule.C:113
Bool is_oxygen(int)
Crossterm * get_crossterm(int inum) const
Definition: Molecule.h:1077
void build_go_arrays(StringList *, char *)
Definition: GoMolecule.C:950
int64_t num_deg_freedom(int isInitialReport=0) const
Definition: Molecule.h:523
int AtomSigID
Definition: NamdTypes.h:82
int exp_a
Definition: Molecule.h:107
int32 * get_angles_for_atom(int anum)
Definition: Molecule.h:1147
int get_qmTotCustPCs()
Definition: Molecule.h:863
int * indxLJB
Definition: Molecule.h:666
int numFepFinal
Definition: Molecule.h:608
Real * get_qmGrpChrg()
Definition: Molecule.h:830
Bond * get_donor(int dnum) const
Definition: Molecule.h:1102
const int *const * get_qmGrpBonds()
Definition: Molecule.h:837
Real get_go_epsilon(int chain1, int chain2)
Definition: Molecule.h:1610
int * pointerToGaussBeg
Definition: Molecule.h:671
int numGaussPair
Definition: Molecule.h:673
int32 * atomChainTypes
Definition: Molecule.h:643
int const *const * get_cSMDpairs()
Definition: Molecule.h:871
int numCalcCrossterms
Definition: Molecule.h:625
int isOccupancyValid
Definition: Molecule.h:1465
void get_constorque_params(BigReal &v, Vector &a, Vector &p, int atomnum) const
Definition: Molecule.h:1336
int numPair
Definition: Molecule.h:660
int32 * consTorqueIndexes
Definition: Molecule.h:615
HydrogenGroup hydrogenGroup
Definition: Molecule.h:639
void build_langevin_params(BigReal coupling, BigReal drudeCoupling, Bool doHydrogen)
int64 numTotalExclusions
Definition: Molecule.h:571
int Index
Definition: structures.h:26
Real * get_qmMeQMGrp()
Definition: Molecule.h:859
const Real * get_qmAtomGroup() const
Definition: Molecule.h:819
int get_qmMeNumBonds()
Definition: Molecule.h:857
ResidueLookupElem(void)
Definition: Molecule.h:97
const char * get_atomtype(int anum) const
Definition: Molecule.h:1116
Real r_ohc
Definition: Molecule.h:467
Real * gMu1
Definition: Molecule.h:677
Real * get_qmGrpID()
Definition: Molecule.h:829
#define order
Definition: PmeRealSpace.C:235
int numMultipleImpropers
Definition: Molecule.h:637
double B
Definition: Molecule.h:121
int * pointerToLJBeg
Definition: Molecule.h:663
Real * gRepulsive
Definition: Molecule.h:681
void goInit()
Definition: GoMolecule.C:54
Bool is_atom_rotdragged(int atomnum) const
Definition: Molecule.h:1230
int const * getLcpoParamType()
Definition: Molecule.h:478
Bool is_atom_stirred(int atomnum) const
Definition: Molecule.h:1428
int goIndxA
Definition: Molecule.h:118
int numMultipleDihedrals
Definition: Molecule.h:635
Bool is_atom_constorqued(int atomnum) const
Definition: Molecule.h:1246
BigReal getEnergyTailCorr(const BigReal, const int)
int numLonepairs
Number of lone pairs.
Definition: Molecule.h:580
int * pointerToGoEnd
Definition: Molecule.h:658
void build_gro_pair()
int * indxGaussB
Definition: Molecule.h:675
Lphost * get_lphost(int atomid) const
Definition: Molecule.h:1081
Bool is_drude(int)
int * get_qmMeMMindx()
Definition: Molecule.h:858
Real get_go_sigmaRep(int chain1, int chain2)
Definition: Molecule.h:1607
int64 numCalcFullExclusions
Definition: Molecule.h:627
BigReal get_go_force(BigReal, int, int, BigReal *, BigReal *) const
Definition: GoMolecule.C:1260
int * goIndxLJA
Definition: Molecule.h:653
Dihedral * get_dihedral(int dnum) const
Definition: Molecule.h:1074
~Molecule()
Definition: Molecule.C:558
int64_t num_group_deg_freedom() const
Definition: Molecule.h:510
int is_drude_psf
Definition: Molecule.h:461
const float * getBFactorData()
Definition: Molecule.h:1032
int goNumLJPair
Definition: Molecule.h:652
int checkexcl(int atom1, int atom2) const
int32 * get_impropers_for_atom(int anum)
Definition: Molecule.h:1151
int get_go_exp_rep(int chain1, int chain2)
Definition: Molecule.h:1619
int numFixedRigidBonds
Definition: Molecule.h:605
void build_constraint_params(StringList *, StringList *, StringList *, PDB *, char *)
int get_qmPCFreq()
Definition: Molecule.h:861
int32 * goResidIndices
Definition: Molecule.h:645
void build_constorque_params(StringList *, StringList *, StringList *, StringList *, StringList *, StringList *, PDB *, char *)
GoValue go_array[MAX_GO_CHAINS *MAX_GO_CHAINS]
Definition: Molecule.h:1584
void read_alch_unpert_bonds(FILE *)
Definition: Molecule.C:1660
int numFepInitial
Definition: Molecule.h:607
int Bool
Definition: common.h:131
int * pointerToGaussEnd
Definition: Molecule.h:672
unsigned char get_ss_type(int anum) const
Definition: Molecule.h:1351
int numAcceptors
Definition: Molecule.h:569
int numTholes
Number of Thole terms.
Definition: Molecule.h:582
Real * get_qmAtmChrg()
Definition: Molecule.h:822
int32 * get_bonds_for_atom(int anum)
Definition: Molecule.h:1145
int get_qmNumBonds()
Definition: Molecule.h:833
const int *const * get_qmMMBond()
Definition: Molecule.h:836
int numCalcDihedrals
Definition: Molecule.h:623
int * goResids
Definition: Molecule.h:649
Real * goSigmas
Definition: Molecule.h:646
int get_qmLSSFreq()
Definition: Molecule.h:847
int numExPressureAtoms
Definition: Molecule.h:598
void set_qm_replaceAll(Bool newReplaceAll)
Definition: Molecule.h:817
int lookup(const char *segid, int resid, int *begin, int *end) const
Definition: Molecule.C:76
Bond * alch_unpert_bonds
Definition: Molecule.h:575
BigReal energyNonnative
Definition: Molecule.h:685
void delete_qm_bonded()
Definition: MoleculeQM.C:1547
int numFixedAtoms
Definition: Molecule.h:596
void prepare_qm(const char *pdbFileName, Parameters *params, ConfigList *cfgList)
Definition: MoleculeQM.C:109
Bool go_restricted(int, int, int)
Definition: GoMolecule.C:525
Real cutoff
Definition: Molecule.h:112
Bond * get_bond(int bnum) const
Definition: Molecule.h:1065
int numCalcImpropers
Definition: Molecule.h:624
int numAngles
Definition: Molecule.h:560
void build_stirred_atoms(StringList *, StringList *, PDB *, char *)
const float * getOccupancyData()
Definition: Molecule.h:1028
Real * pairC12
Definition: Molecule.h:668
const Real * get_cSMDcoffs()
Definition: Molecule.h:874
Definition: parm.h:15
int numAtoms
Definition: Molecule.h:556
Real sigmaRep
Definition: Molecule.h:110
int * gromacsPair_type
Definition: Molecule.h:669
int numStirredAtoms
Definition: Molecule.h:597
GridforceGrid * get_gridfrc_grid(int gridnum) const
Definition: Molecule.h:1276
int numCrossterms
Definition: Molecule.h:567
void NAMD_die(const char *err_msg)
Definition: common.C:83
Bool is_lp(int)
int numConsForce
Definition: Molecule.h:611
HashPool< HashString > atomNamePool
Definition: CompressPsf.C:309
BigReal getVirialTailCorr(const BigReal)
void build_ss_flags(const StringList *ssfile, const StringList *sscol, PDB *initial_pdb, const char *cwd)
void build_rotdrag_params(StringList *, StringList *, StringList *, StringList *, StringList *, StringList *, PDB *, char *)
HashPool< HashString > resNamePool
Definition: CompressPsf.C:308
bool operator<(const intpair &lhs, const intpair &rhs)
Definition: ComputeGlobal.h:29
ConsTorqueParams * consTorqueParams
Definition: Molecule.h:616
AtomNameInfo * getAtomNames() const
Definition: Molecule.h:491
int get_groupSize(int)
int * ss_index
Definition: Molecule.h:458
int maxMigrationGroupSize
Definition: Molecule.h:602
char mySegid[11]
Definition: Molecule.h:91
void put_stir_startTheta(Real theta, int atomnum) const
Definition: Molecule.h:1308
int numDihedrals
Definition: Molecule.h:561
const int *const * get_qmMMChargeTarget()
Definition: Molecule.h:841
int * get_qmCustPCSizes()
Definition: Molecule.h:864
void build_fep_flags(StringList *, StringList *, PDB *, char *, const char *)
unsigned char get_fep_type(int anum) const
Definition: Molecule.h:1345
int numImpropers
Definition: Molecule.h:566
const char *const * get_qmElements()
Definition: Molecule.h:826
int numConsTorque
Definition: Molecule.h:595
ResidueLookupElem * append(const char *segid, int resid, int aid)
Definition: Molecule.C:89
double * goSigmaPairB
Definition: Molecule.h:656
int numConstraints
Definition: Molecule.h:588
Bool is_group_fixed(int atomnum) const
Definition: Molecule.h:1443
void build_alch_unpert_bond_lists(char *)
int64 numCalcExclusions
Definition: Molecule.h:626
Bool is_atom_constrained(int atomnum) const
Definition: Molecule.h:1197
void build_fixed_atoms(StringList *, StringList *, PDB *, char *)
void get_cons_params(Real &k, Vector &refPos, int atomnum) const
Definition: Molecule.h:1262
BigReal tail_corr_dUdl_2
Definition: Molecule.h:472
const Real * get_cSMDKs()
Definition: Molecule.h:872
void compute_LJcorrection()
long long int64
Definition: common.h:34
int num_fixed_groups() const
Definition: Molecule.h:503
int get_go_exp_b(int chain1, int chain2)
Definition: Molecule.h:1616
#define simParams
Definition: Output.C:127
Dihedral * getAllDihedrals() const
Definition: Molecule.h:1093
AtomSegResInfo * getAtomSegResInfo() const
Definition: Molecule.h:494
Exclusion * get_exclusion(int ex) const
Definition: Molecule.h:1112
Index atomvdwtype(int anum) const
Definition: Molecule.h:1058
int32 * consForceIndexes
Definition: Molecule.h:612
void print_atoms(Parameters *)
Definition: Molecule.C:5337
int maxHydrogenGroupSize
Definition: Molecule.h:600
BigReal get_go_force2(BigReal, BigReal, BigReal, int, int, BigReal *, BigReal *) const
Definition: GoMolecule.C:1456
void get_movdrag_params(Vector &v, int atomnum) const
Definition: Molecule.h:1321
int numMigrationGroups
Definition: Molecule.h:601
int go_indices[MAX_GO_CHAINS+1]
Definition: Molecule.h:1585
void build_go_sigmas(StringList *, char *)
Definition: GoMolecule.C:577
void get_rotdrag_params(BigReal &v, Vector &a, Vector &p, int atomnum) const
Definition: Molecule.h:1327
Real * giSigma1
Definition: Molecule.h:678
Real atomcharge(int anum) const
Definition: Molecule.h:1048
HashPool< HashString > segNamePool
Definition: CompressPsf.C:307
void delete_alch_bonded(void)
HashPool< HashString > atomTypePool
Definition: CompressPsf.C:310
BigReal tail_corr_virial
Definition: Molecule.h:471
void build_exPressure_atoms(StringList *, StringList *, PDB *, char *)
int numZeroMassAtoms
Number of atoms with zero mass.
Definition: Molecule.h:585
const ExclusionCheck * get_excl_check_for_atom(int anum) const
Definition: Molecule.h:1173
void build_constant_forces(char *)
float Mass
Definition: ComputeGBIS.inl:20
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 * get_qmCustomPCIdxs()
Definition: Molecule.h:865
Angle * getAllAngles() const
Definition: Molecule.h:1091
Vector * consForce
Definition: Molecule.h:613
#define WAT_TIP4
Definition: common.h:189
int * indxLJA
Definition: Molecule.h:665
~ResidueLookupElem(void)
Definition: Molecule.h:98
int set_gridfrc_grid(int gridnum, GridforceGrid *grid)
Definition: Molecule.h:1285
Real atommass(int anum) const
Definition: Molecule.h:1038
const Real * get_cSMDVels()
Definition: Molecule.h:873
BigReal energyNative
Definition: Molecule.h:684
void print_go_sigmas()
Definition: GoMolecule.C:1134
int * ss_vdw_type
Definition: Molecule.h:457
int * get_qmLSSSize()
Definition: Molecule.h:844
int * numGridforces
Definition: Molecule.h:591
int NumGoChains
Definition: Molecule.h:1586
int numDrudeAtoms
Number of Drude particles.
Definition: Molecule.h:581
Bool is_atom_exPressure(int atomnum) const
Definition: Molecule.h:1447
int restrictions[MAX_RESTRICTIONS]
Definition: Molecule.h:113
Real epsilon
Definition: Molecule.h:106
Bool atoms_1to4(unsigned int, unsigned int)
Definition: GoMolecule.C:1537
Bool is_atom_fixed(int atomnum) const
Definition: Molecule.h:1407
int is_lonepairs_psf
Definition: Molecule.h:462
int alchDroppedAngles
Definition: Molecule.h:563
int32 * goSigmaIndices
Definition: Molecule.h:644
int numCalcTholes
Definition: Molecule.h:630
BigReal GetAtomAlpha(int i) const
Definition: Molecule.h:482
gridSize x
void receive_Molecule(MIStream *)
Definition: Molecule.C:5806
void build_movdrag_params(StringList *, StringList *, StringList *, PDB *, char *)
int resid
Definition: Molecule.h:146
int numLphosts
Number of lone pair host records in PSF.
Definition: Molecule.h:584
Lphost * getAllLphosts() const
Definition: Molecule.h:1098
int * goIndxLJB
Definition: Molecule.h:654
Angle * get_angle(int anum) const
Definition: Molecule.h:1068
int numRigidBonds
Definition: Molecule.h:604
int const *const * get_cSMDindex()
Definition: Molecule.h:870
Mass * get_qmLSSRefMass()
Definition: Molecule.h:851
Bool is_water(int)
Bond * getAllBonds() const
Definition: Molecule.h:1090
void freeBFactorData()
Definition: Molecule.h:1034
const BigReal * get_qmDummyBondVal()
Definition: Molecule.h:834
int exp_rep
Definition: Molecule.h:109
int get_mother_atom(int)
void initialize()
int numDonors
Definition: Molecule.h:568
Bool get_qmReplaceAll()
Definition: Molecule.h:854
int numGoAtoms
Definition: Molecule.h:642
int isBFactorValid
Definition: Molecule.h:1465
int numCalcAngles
Definition: Molecule.h:622
float Charge
Definition: NamdTypes.h:32
const int32 * get_full_exclusions_for_atom(int anum) const
Definition: Molecule.h:1157
Improper * getAllImpropers() const
Definition: Molecule.h:1092
Real r_om
Definition: Molecule.h:466
int get_clusterSize(int anum) const
Definition: Molecule.h:1024
int ss_num_vdw_params
Definition: Molecule.h:456
int numRotDrag
Definition: Molecule.h:594
ResizeArray< int > atomIndex
Definition: Molecule.h:95
Molecule(SimParameters *, Parameters *param)
Definition: Molecule.C:429
double BigReal
Definition: common.h:112
int get_qmLSSResSize()
Definition: Molecule.h:848
Bool is_atom_gridforced(int atomnum, int gridnum) const
Definition: Molecule.h:1181
const int * get_qmGrpSizes()
Definition: Molecule.h:828
HashPool< AtomSigInfo > atomSigPool
Definition: CompressPsf.C:313
Bond * get_acceptor(int dnum) const
Definition: Molecule.h:1105
int numLJPair
Definition: Molecule.h:661
void build_go_params(StringList *)
Definition: GoMolecule.C:79
void build_go_sigmas2(StringList *, char *)
Definition: GoMolecule.C:747
int numExclusions
Definition: Molecule.h:570
void read_alch_unpert_angles(FILE *)
Definition: Molecule.C:1783