Parameters.h

Go to the documentation of this file.
00001 
00007 /*
00008    The Parameters class is used to read in and store all of the parameters
00009    from the parameter files.  These parameters are stored and used to assign
00010    constants to bonds and atoms as they are read in from the psf file
00011 */
00012 
00013 #ifndef PARAM_H
00014 
00015 #define PARAM_H
00016 
00017 #include "parm.h"
00018 
00019 #include "common.h"
00020 #include "structures.h"
00021 #include "strlib.h"
00022 #include "MStream.h"
00023 //****** BEGIN CHARMM/XPLOR type changes
00024 #include "SimParameters.h"
00025 //****** END CHARMM/XPLOR type changes
00026 
00027 #include "GromacsTopFile.h"
00028 
00029 #ifndef cbrt
00030   // cbrt() not in math.h on goneril
00031   #define cbrt(x)  pow((double)x,(double)(1.0/3.0))
00032 #endif
00033 
00034 class Communicate;
00035 class StringList;
00036 
00037 //  The class Parameters is used to store and query the parameters for
00038 //  bonds and atoms.  If the Parameter object resides on the master
00039 //  process, it is responsible for reading in all the parameters and
00040 //  then communicating them to the other processors.  To do this, first
00041 //  the routine read_paramter_file is called to read in each parameter
00042 //  file.  The information that is read in is stored in binary trees
00043 //  for vdw, bond, and angle information and in linked lists for 
00044 //  dihedrals and impropers.  Once all of the files have been read
00045 //  in, the routine done_reading_files is called.  At this point, all
00046 //  of the data that has been read in is copied to arrays.  This is
00047 //  so that each bond and atom can be assigned an index into these
00048 //  arrays to retreive the parameters in constant time.
00049 //
00050 //  Then the psf file is read in.  Each bond and atom is assigned an
00051 //  index into the parameter lists using the functions assign_*_index.
00052 //  Once the psf file has been read in, then the routine 
00053 //  done_reading_structre is called to tell the object that it no
00054 //  longer needs to store the binary trees and linked lists that were
00055 //  used to query the parameters based on atom type.  From this point
00056 //  on, only the indexes will be used.
00057 //
00058 //  The master node then uses send_parameters to send all of these
00059 //  parameters to the other nodes and the objects on all of the other
00060 //  nodes use receive_parameters to accept these parameters.
00061 //
00062 //  From this point on, all of the parameter data is static and the
00063 //  functions get_*_params are used to retrieve the parameter data
00064 //  that is desired.
00065 
00066 //  Define the number of Multiples that a Dihedral or Improper
00067 //  bond can have with the Charm22 parameter set
00068 #define MAX_MULTIPLICITY        6
00069 
00070 //  Number of characters maximum allowed for storing atom type names
00071 #define MAX_ATOMTYPE_CHARS      6
00072 
00073 //****** BEGIN CHARMM/XPLOR type changes
00074 //  Define the numbers associated with each possible parameter-file 
00075 //  type (format) NAMD can handle.
00076 #define paraXplor               0
00077 #define paraCharmm              1
00078 //****** END CHARMM/XPLOR type changes
00079 
00080 
00081 class BondValue {
00082 public:
00083         Real k;         //  Force constant for the bond
00084         Real x0;        //  Rest distance for the bond
00085 };
00086 
00087 class AngleValue {
00088 public:
00089         Real k;         //  Force constant for angle
00090         Real theta0;    //  Rest angle for angle
00091         Real k_ub;      //  Urey-Bradley force constant
00092         Real r_ub;      //  Urey-Bradley distance
00093   int normal; // Whether we use harmonic (0) or cos-based (1) angle terms
00094 };
00095 
00096 typedef struct four_body_consts
00097 {
00098         Real k;         //  Force constant
00099         int n;          //  Periodicity
00100         Real delta;     //  Phase shift
00101 } FourBodyConsts;
00102 
00103 class DihedralValue {
00104 public:
00105         int multiplicity;
00106         FourBodyConsts values[MAX_MULTIPLICITY];
00107 };
00108 
00109 class ImproperValue {
00110 public:
00111         int multiplicity;
00112         FourBodyConsts values[MAX_MULTIPLICITY];
00113 };
00114 
00115 struct CrosstermData { BigReal d00,d01,d10,d11; };
00116 
00117 class CrosstermValue {
00118 public:
00119         enum {dim=25};
00120         CrosstermData c[dim][dim];  // bicubic interpolation coefficients
00121 };
00122 
00123 // JLai
00124 class GromacsPairValue {
00125 public:
00126     BigReal pairC6;
00127     BigReal pairC12;
00128 /*    BigReal A;
00129     BigReal B;
00130     BigReal A14;
00131     BigReal B14;
00132     BigReal sqrt26;
00133     BigReal expo;*/
00134 };
00135 // End of JLai
00136 
00137 
00138 class NonbondedExclValue {
00139 public:
00140         // need to put parameters here...
00141         // for now, copy bond
00142         Real k;         //  Force constant for the bond
00143         Real x0;        //  Rest distance for the bond
00144 };
00145 
00146 typedef struct vdw_val
00147 {
00148         Real sigma;     //  Sigma value
00149         Real epsilon;   //  Epsilon value
00150         Real sigma14;   //  Sigma value for 1-4 interactions
00151         Real epsilon14; //  Epsilon value for 1-4 interactions
00152 } VdwValue;
00153 
00154 //  IndexedVdwPair is used to form a binary search tree that is
00155 //  indexed by vwd_type index.  This is the tree that will be
00156 //  used to search during the actual simulation
00157 
00158 typedef struct indexed_vdw_pair
00159 {
00160         Index ind1;             //  Index for first atom type
00161         Index ind2;             //  Index for second atom type
00162         Real A;                 //  Parameter A for this pair
00163         Real A14;               //  Parameter A for 1-4 interactions
00164         Real B;                 //  Parameter B for this pair
00165         Real B14;               //  Parameter B for 1-4 interactions
00166         struct indexed_vdw_pair *right;  //  Right child
00167         struct indexed_vdw_pair *left;   //  Left child
00168 } IndexedVdwPair;
00169 
00170 typedef struct indexed_nbthole_pair
00171 {
00172         Index ind1;             //  Index for first atom type
00173         Index ind2;             //  Index for second atom type
00174         Real alphai;           //  Parameter alpha for this pair
00175         Real alphaj;           //  Parameter alpha for this pair
00176         Real tholeij;          //  Parameter thole for this pair
00177         struct indexed_nbthole_pair *right;  //  Right child
00178         struct indexed_nbthole_pair *left;   //  Left child
00179 } IndexedNbtholePair;
00180 
00181 typedef struct nbthole_pair_value
00182 {
00183         Index ind1;             //  Index for first atom type
00184         Index ind2;             //  Index for second atom type
00185         Real alphai;           //  Parameter alpha for this pair
00186         Real alphaj;           //  Parameter alpha for this pair
00187         Real tholeij;          //  Parameter thole for this pair
00188 } NbtholePairValue;
00189 
00190 //  IndexedTablePair is used to form a binary search tree that is
00191 //  indexed by table_type index.  This is the tree that will be
00192 //  used to search during the actual simulation
00193 
00194 typedef struct indexed_table_pair
00195 {
00196         Index ind1;             //  Index for first atom type
00197         Index ind2;             //  Index for second atom type
00198         int K;                  //  number of the table type for this pair
00199         struct indexed_table_pair *right;        //  Right child
00200         struct indexed_table_pair *left;         //  Left child
00201 } IndexedTablePair;
00202 
00203 //  Structures that are defined in Parameters.C
00204 struct bond_params;
00205 struct angle_params;
00206 struct improper_params;
00207 struct dihedral_params;
00208 struct crossterm_params;
00209 struct vdw_params;
00210 struct vdw_pair_params;
00211 struct nbthole_pair_params;
00212 struct table_pair_params;
00213 
00214 class Parameters
00215 {
00216 private:
00217         void initialize();                      //  zeros out pointers
00218 
00219         char *atomTypeNames;                    //  Names of atom types
00220         Bool AllFilesRead;                      //  Flag TRUE imples that all
00221                                                 //  of the parameter files
00222                                                 //  have been read in and
00223                                                 //  the arrays have been
00224                                                 //  created.
00225 //****** BEGIN CHARMM/XPLOR type changes
00226         int paramType;                          //  Type (format) of parameter-file
00227 //****** END CHARMM/XPLOR type changes
00228   bool cosAngles; // True if some angles may be cos-based
00229         struct bond_params *bondp;              //  Binary tree of bond params
00230         struct angle_params *anglep;            //  Binary tree of angle params
00231         struct improper_params *improperp;      //  Linked list of improper par.
00232         struct dihedral_params *dihedralp;      //  Linked list of dihedral par.
00233         struct crossterm_params *crosstermp;    //  Linked list of cross-term par.
00234         struct vdw_params *vdwp;                //  Binary tree of vdw params
00235         struct vdw_pair_params *vdw_pairp;      //  Binary tree of vdw pairs
00236         struct nbthole_pair_params *nbthole_pairp;      //  Binary tree of nbthole pairs
00237         struct table_pair_params *table_pairp;  //  Binary tree of table pairs
00238 public:
00239         BondValue *bond_array;                  //  Array of bond params
00240         AngleValue *angle_array;                //  Array of angle params
00241         DihedralValue *dihedral_array;          //  Array of dihedral params
00242         ImproperValue *improper_array;          //  Array of improper params
00243         CrosstermValue *crossterm_array;        //  Array of crossterm params
00244         // JLai
00245         GromacsPairValue *gromacsPair_array;    //  Array of gromacsPair params
00246         // End of JLai
00247         VdwValue *vdw_array;                    //  Array of vdw params
00248         NbtholePairValue *nbthole_array;        //  Array of nbthole params
00249         int numenerentries;                     //  Number of entries for enertable
00250         int rowsize;
00251         int columnsize;
00252         BigReal* table_ener;                    //  Table for tabulated energies
00253         IndexedVdwPair *vdw_pair_tree;          //  Tree of vdw pair params
00254         IndexedNbtholePair *nbthole_pair_tree;      //  Tree of nbthole pair params
00255         IndexedTablePair *tab_pair_tree;                //  Tree of vdw pair params
00256   int tablenumtypes;
00257         int NumBondParams;                      //  Number of bond parameters
00258         int NumAngleParams;                     //  Number of angle parameters
00259   int NumCosAngles;       // Number of cosine-based angles
00260         int NumDihedralParams;                  //  Number of dihedral params
00261         int NumImproperParams;                  //  Number of improper params
00262         int NumCrosstermParams;                 //  Number of cross-term params
00263         // JLai
00264         int NumGromacsPairParams;               //  Number of gromacsPair params
00265         // End of JLai
00266         int NumVdwParams;                       //  Number of vdw parameters
00267         int NumTableParams;                     //  Number of table parameters
00268   int NumVdwParamsAssigned;               //  Number actually assigned
00269         int NumVdwPairParams;                   //  Number of vdw_pair params
00270         int NumNbtholePairParams;                   //  Number of nbthole_pair params
00271         int NumTablePairParams;                 //  Number of vdw_pair params
00272 private:
00273         ResizeArray<char *> error_msgs;         //  Avoids repeating warnings
00274 
00275         int *maxDihedralMults;                  //  Max multiplicity for
00276                                                 //  dihedral bonds
00277         int *maxImproperMults;                  //  Max multiplicity for
00278                                                 //  improper bonds
00279 
00280         void skip_stream_read(char *, FILE *);  // skip part of stream file
00281 
00282         void add_bond_param(char *);            //  Add a bond parameter
00283         struct bond_params *add_to_bond_tree(struct bond_params * , 
00284                                      struct bond_params *);
00285 
00286         void add_angle_param(char *);           //  Add an angle parameter
00287         struct angle_params *add_to_angle_tree(struct angle_params * , 
00288                                      struct angle_params *);
00289 
00290         void add_dihedral_param(char *, FILE *); //  Add a dihedral parameter
00291         void add_to_dihedral_list(struct dihedral_params *);
00292         void add_to_charmm_dihedral_list(struct dihedral_params *);
00293 
00294         void add_improper_param(char *, FILE *); //  Add an improper parameter
00295         void add_to_improper_list(struct improper_params *);
00296 
00297         void add_crossterm_param(char *, FILE *); //  Add an cross-term parameter
00298         void add_to_crossterm_list(struct crossterm_params *);
00299 
00300         void add_vdw_param(char *);             //  Add a vdw parameter
00301         struct vdw_params *add_to_vdw_tree(struct vdw_params *, 
00302                                      struct vdw_params *);
00303 
00304         void add_vdw_pair_param(char *);        //  Add a vdw pair parameter
00305         void add_nbthole_pair_param(char *);        //  Add a nbthole pair parameter        
00306         void add_table_pair_param(char *);      //  Add a table pair parameter
00307         void add_to_vdw_pair_list(struct vdw_pair_params *);
00308         void add_to_nbthole_pair_list(struct nbthole_pair_params *);
00309         void add_to_table_pair_list(struct table_pair_params *);
00310 
00311         void add_hb_pair_param(char *); //  Add a hydrogen bond pair parameter
00312 
00313         //  All of the traverse routines are used for debugging purposes
00314         //  to print out the appropriate list of parameters
00315         void traverse_vdw_pair_params(struct vdw_pair_params *);
00316         void traverse_nbthole_pair_params(struct nbthole_pair_params *);
00317         void traverse_vdw_params(struct vdw_params *);
00318         void traverse_dihedral_params(struct dihedral_params *);
00319         void traverse_improper_params(struct improper_params *);
00320         void traverse_crossterm_params(struct crossterm_params *);
00321         void traverse_angle_params(struct angle_params *);
00322         void traverse_bond_params(struct bond_params *);
00323 
00324         //  The index_* routines are used to index each of 
00325         //  the parameters and build the arrays that will be used
00326         //  for constant time access
00327         Index index_bonds(struct bond_params *, Index);
00328         Index index_angles(struct angle_params *, Index);
00329         Index index_vdw(struct vdw_params *, Index);
00330         void index_dihedrals();
00331         void index_impropers();
00332         void index_crossterms();
00333         
00334         void convert_vdw_pairs();
00335         void convert_nbthole_pairs();
00336         void convert_table_pairs();
00337         IndexedVdwPair *add_to_indexed_vdw_pairs(IndexedVdwPair *, IndexedVdwPair *);
00338         IndexedNbtholePair *add_to_indexed_nbthole_pairs(IndexedNbtholePair *, IndexedNbtholePair *);
00339         IndexedTablePair *add_to_indexed_table_pairs(IndexedTablePair *, IndexedTablePair *);
00340         
00341         int vdw_pair_to_arrays(int *, int *, Real *, Real *, Real *, Real *, 
00342                                int, IndexedVdwPair *);
00343 
00344         int nbthole_pair_to_arrays(int *, int *, Real *, Real *, Real *, int, IndexedNbtholePair *);
00345         
00346         int table_pair_to_arrays(int *, int *, int *, int, IndexedTablePair *);
00347 
00348         //  The free_* routines are used by the destructor to deallocate
00349         //  memory
00350         void free_bond_tree(struct bond_params *);
00351         void free_angle_tree(struct angle_params *);
00352         void free_dihedral_list(struct dihedral_params *);
00353         void free_improper_list(struct improper_params *);
00354         void free_crossterm_list(struct crossterm_params *);
00355         void free_vdw_tree(struct vdw_params *);
00356         void free_vdw_pair_tree(IndexedVdwPair *);
00357         void free_nbthole_pair_tree(IndexedNbtholePair *);
00358         void free_table_pair_tree(IndexedTablePair *);
00359         void free_vdw_pair_list();
00360         void free_nbthole_pair_list(); 
00361 
00362   BigReal interp_lin(BigReal, BigReal, BigReal, BigReal, BigReal); // Perform a linear interpolation for energy table 
00363 
00364         /* does the actual initialization, once the variables have all
00365            been given default values.  See Parameters() below */
00366         void read_parm(const GromacsTopFile *, Bool min);
00367 public:
00368         //****** BEGIN CHARMM/XPLOR type changes
00370         Parameters();
00371         Parameters(SimParameters *, StringList *f);
00372         //****** END CHARMM/XPLOR type changes
00373         
00374         Parameters(Ambertoppar *, BigReal);
00375         void read_parm(Ambertoppar *, BigReal);
00376 
00377         /* initializes this to hold the set of parameters in the
00378            GROMACS topology file <gf>.  If the parameter <min> is on,
00379            this assumes that we are going to do minimization and
00380            therefore can't have atoms with zero VDW - it will add a
00381            small repulsive term to these. */
00382         Parameters(const GromacsTopFile *gf, Bool min);
00383         
00384         ~Parameters();                          //  Destructor
00385 
00386         // return a string for the Nth atom type.  This can only be
00387         // called after all the param files have been read and the type
00388         // names have been indexed.  The Nth atom type refers to the same
00389         // index of the Nth vdw parameter (i.e. there are NumVdwParams names).
00390         char *atom_type_name(Index a) {
00391           return (atomTypeNames + (a * (MAX_ATOMTYPE_CHARS + 1)));
00392         }
00393 
00394         //  Read a parameter file
00395         void read_parameter_file(char *);
00396 
00397         //****** BEGIN CHARMM/XPLOR type changes
00398         void read_charmm_parameter_file(char *);
00399         //****** END CHARMM/XPLOR type changes
00400 
00401         //  Signal the parameter object that all of
00402         //  the parameter files have been read in
00403         void done_reading_files();
00404 
00405         //  Signal the parameter object that the
00406         //  structure file has been read in
00407         void done_reading_structure();
00408     
00409         //  The assign_*_index routines are used to assign
00410         //  an index to atoms or bonds.  If an specific atom
00411         //  or bond type can't be found, then the program 
00412         //  terminates
00413     #ifdef MEM_OPT_VERSION
00414     void assign_vdw_index(char *, AtomCstInfo *);
00415     #else
00416         void assign_vdw_index(char *, Atom *);  //  Assign a vdw index to
00417     #endif
00418                                                 //  an atom
00419         void assign_bond_index(char *, char *, Bond *); 
00420                                                 //  Assign a bond index
00421                                                 //  to a bond
00422         void assign_angle_index(char *, char *, char *, Angle *, int);
00423                                                 //  Assign an angle index
00424                                                 //  to an angle
00425         void assign_dihedral_index(char *, char*, char*, char *, Dihedral *, int, int);
00426                                                 //  Assign a dihedral index
00427                                                 //  to a dihedral
00428         void assign_improper_index(char *, char*, char*, char *, Improper *, int);
00429                                                 //  Assign an improper index
00430                                                 //  to an improper
00431         void assign_crossterm_index(char *, char*, char*, char *, char *, char*, char*, char *, Crossterm *);
00432 
00433         //  send_parameters is used by the master process to
00434         //  communicate the paramters to all the other processors
00435         void send_Parameters(MOStream *);
00436 
00437         //  receive_parameters is used by all the child processes
00438         //  to receive the parameters from the master process
00439         void receive_Parameters(MIStream *);
00440 
00441         //  The get_*_params routines are the routines that really
00442         //  do all the work for this object.  Given an index, they
00443         //  access the parameters and return the relevant information
00444         void get_bond_params(Real *k, Real *x0, Index index)
00445         {
00446                 *k = bond_array[index].k;
00447                 *x0 = bond_array[index].x0;
00448         }
00449 
00450         void get_angle_params(Real *k, Real *theta0, Real *k_ub, Real *r_ub,
00451                               Index index)
00452         {
00453                 *k = angle_array[index].k;
00454                 *theta0 = angle_array[index].theta0;
00455                 *k_ub = angle_array[index].k_ub;
00456                 *r_ub = angle_array[index].r_ub;
00457         }
00458 
00459         int get_improper_multiplicity(Index index)
00460         {
00461                 return(improper_array[index].multiplicity);
00462         }
00463 
00464         int get_dihedral_multiplicity(Index index)
00465         {
00466                 return(dihedral_array[index].multiplicity);
00467         }
00468 
00469         void get_improper_params(Real *k, int *n, Real *delta, 
00470                                  Index index, int mult)
00471         {
00472                 if ( (mult<0) || (mult>MAX_MULTIPLICITY) )
00473                 {
00474                         NAMD_die("Bad mult index in Parameters::get_improper_params");
00475                 }
00476 
00477                 *k = improper_array[index].values[mult].k;
00478                 *n = improper_array[index].values[mult].n;
00479                 *delta = improper_array[index].values[mult].delta;
00480         }
00481 
00482         void get_dihedral_params(Real *k, int *n, Real *delta, 
00483                                  Index index, int mult)
00484         {
00485                 if ( (mult<0) || (mult>MAX_MULTIPLICITY) )
00486                 {
00487                         NAMD_die("Bad mult index in Parameters::get_dihedral_params");
00488                 }
00489 
00490                 *k = dihedral_array[index].values[mult].k;
00491                 *n = dihedral_array[index].values[mult].n;
00492                 *delta = dihedral_array[index].values[mult].delta;
00493         }
00494 
00495         void get_vdw_params(Real *sigma, Real *epsilon, Real *sigma14, 
00496                             Real *epsilon14, Index index)
00497   {
00498     if ( vdw_array ) {
00499       *sigma = vdw_array[index].sigma;
00500       *epsilon = vdw_array[index].epsilon;
00501       *sigma14 = vdw_array[index].sigma14;
00502       *epsilon14 = vdw_array[index].epsilon14;
00503     } else {
00504       // sigma and epsilon from A and B for each vdw type's interaction with itself
00505       Real A; Real B; Real A14; Real B14;
00506       if (get_vdw_pair_params(index, index, &A, &B, &A14, &B14) ) {
00507         if (A && B) {
00508           *sigma = sqrt(cbrt(A)) / sqrt(cbrt(B));
00509           *epsilon = B*B / (4*A);
00510         }
00511         else {
00512           *sigma = 0; *epsilon = 0;
00513         }
00514         if (A14 && B14) {
00515           *sigma14 = sqrt(cbrt(A14)) / sqrt(cbrt(B14));
00516           *epsilon14 = B14*B14 / (4*A14);
00517         }
00518         else {
00519           *sigma14 = 0; *epsilon14 = 0;
00520         }
00521       }
00522       else {NAMD_die ("Function get_vdw_params failed to derive Lennard-Jones sigma and epsilon from A and B values\n");}
00523     }
00524   }
00525 
00526         int get_vdw_pair_params(Index ind1, Index ind2, Real *, Real *, Real *, Real *);
00527                                                 //  Find a vwd_pair parameter
00528 
00529         int get_num_vdw_params(void) { return NumVdwParamsAssigned; }
00530 
00531   int get_table_pair_params(Index, Index, int*);
00532 
00533         //  The print_*_params are provided for debugging purposes
00534         void print_bond_params();               //  Print bonds params
00535         void print_angle_params();              //  Print angle params
00536         void print_dihedral_params();           //  Print dihedral params
00537         void print_improper_params();           //  Print improper params
00538         void print_crossterm_params();          //  Print cross-term params
00539         void print_vdw_params();                //  Print vdw params
00540         void print_vdw_pair_params();           //  Print vdw_pair params
00541         void print_nbthole_pair_params();           //  Print nbthole_pair params
00542         void print_param_summary();             //  Print a summary of params
00543         void read_ener_table(SimParameters*); // Read an energy table file
00544   int get_int_table_type(char*); // Return the integer type for a named table interaction
00545 
00546   int read_energy_type(FILE*, const int, BigReal*, const float, const float);
00547   int read_energy_type_cubspline(FILE*, const int, BigReal*, const float, const float);
00548   int read_energy_type_bothcubspline(FILE*, const int, BigReal*, const float, const float);
00549 };
00550 
00551 #endif
00552 

Generated on Tue Nov 21 01:17:14 2017 for NAMD by  doxygen 1.4.7