Main Page | Namespace List | Class Hierarchy | Alphabetical List | Class List | File List | Class Members | File Members

Parameters Class Reference

#include <Parameters.h>

List of all members.

Public Member Functions

 Parameters ()
 Parameters (SimParameters *, StringList *f)
 Parameters (Ambertoppar *, BigReal)
void read_parm (Ambertoppar *, BigReal)
 Parameters (const GromacsTopFile *gf, Bool min)
 ~Parameters ()
char * atom_type_name (Index a)
void read_parameter_file (char *)
void read_charmm_parameter_file (char *)
void done_reading_files ()
void done_reading_structure ()
void assign_vdw_index (char *, Atom *)
void assign_bond_index (char *, char *, Bond *)
void assign_angle_index (char *, char *, char *, Angle *)
void assign_dihedral_index (char *, char *, char *, char *, Dihedral *, int)
void assign_improper_index (char *, char *, char *, char *, Improper *, int)
void assign_crossterm_index (char *, char *, char *, char *, char *, char *, char *, char *, Crossterm *)
void send_Parameters (MOStream *)
void receive_Parameters (MIStream *)
void get_bond_params (Real *k, Real *x0, Index index)
void get_angle_params (Real *k, Real *theta0, Real *k_ub, Real *r_ub, Index index)
int get_improper_multiplicity (Index index)
int get_dihedral_multiplicity (Index index)
void get_improper_params (Real *k, int *n, Real *delta, Index index, int mult)
void get_dihedral_params (Real *k, int *n, Real *delta, Index index, int mult)
void get_vdw_params (Real *sigma, Real *epsilon, Real *sigma14, Real *epsilon14, Index index)
int get_vdw_pair_params (Index ind1, Index ind2, Real *, Real *, Real *, Real *)
int get_num_vdw_params (void)
int get_table_pair_params (Index, Index, int *)
void print_bond_params ()
void print_angle_params ()
void print_dihedral_params ()
void print_improper_params ()
void print_crossterm_params ()
void print_vdw_params ()
void print_vdw_pair_params ()
void print_param_summary ()
void read_ener_table (SimParameters *)
int get_int_table_type (char *)
int read_energy_type (FILE *, const int, BigReal *, const float, const float)
int read_energy_type_cubspline (FILE *, const int, BigReal *, const float, const float)
int read_energy_type_bothcubspline (FILE *, const int, BigReal *, const float, const float)

Public Attributes

BondValuebond_array
AngleValueangle_array
DihedralValuedihedral_array
ImproperValueimproper_array
CrosstermValuecrossterm_array
VdwValuevdw_array
int numenerentries
int rowsize
int columnsize
BigRealtable_ener
IndexedVdwPairvdw_pair_tree
IndexedTablePairtab_pair_tree
int tablenumtypes
int NumBondParams
int NumAngleParams
int NumCosAngles
int NumDihedralParams
int NumImproperParams
int NumCrosstermParams
int NumVdwParams
int NumTableParams
int NumVdwParamsAssigned
int NumVdwPairParams
int NumTablePairParams


Constructor & Destructor Documentation

Parameters::Parameters  ) 
 

Definition at line 175 of file Parameters.C.

00175                        {
00176   initialize();
00177 }

Parameters::Parameters SimParameters ,
StringList f
 

Definition at line 227 of file Parameters.C.

References SimParameters::cosAngles, StringList::data, done_reading_files(), StringList::next, numenerentries, SimParameters::paraTypeCharmmOn, SimParameters::paraTypeXplorOn, read_charmm_parameter_file(), read_ener_table(), read_parameter_file(), simParams, table_ener, and SimParameters::tabulatedEnergies.

00228 {
00229   initialize();
00230 
00231   //****** BEGIN CHARMM/XPLOR type changes
00233   if (simParams->paraTypeXplorOn)
00234   {
00235     paramType = paraXplor;
00236   }
00237   else if (simParams->paraTypeCharmmOn)
00238   {
00239     paramType = paraCharmm;
00240   }
00241   //****** END CHARMM/XPLOR type changes
00242   //Test for cos-based angles
00243   if (simParams->cosAngles) {
00244     cosAngles = true;
00245   } else {
00246     cosAngles = false;
00247   }
00248 
00249   /* Set up AllFilesRead flag to FALSE.  Once all of the files    */
00250   /* have been read in, then this will be set to true and the     */
00251   /* arrays of parameters will be set up        */
00252   AllFilesRead = FALSE;
00253 
00254   numenerentries=0;
00255   table_ener = NULL;
00256   if (simParams->tabulatedEnergies) {
00257 
00258           fprintf(stdout,"Working on tables\n");
00259           read_ener_table(simParams);
00260   }
00261 
00262   if (NULL != f) 
00263   {
00264     do
00265     {
00266       //****** BEGIN CHARMM/XPLOR type changes
00267       if (paramType == paraXplor)
00268       {
00269         read_parameter_file(f->data);
00270       }
00271       else if (paramType == paraCharmm)
00272       {
00273         read_charmm_parameter_file(f->data);
00274       }
00275       //****** END CHARMM/XPLOR type changes
00276       f = f->next;
00277     } while ( f != NULL );
00278 
00279     done_reading_files();
00280   }
00281 
00282 }

Parameters::Parameters Ambertoppar ,
BigReal 
 

Definition at line 5749 of file Parameters.C.

References Ambertoppar, and read_parm().

05750 {
05751   initialize();
05752 
05753   // Read in parm parameters
05754   read_parm(amber_data,vdw14);
05755 }

Parameters::Parameters const GromacsTopFile gf,
Bool  min
 

Definition at line 5891 of file Parameters.C.

References read_parm().

05892 {
05893   initialize();
05894 
05895   // Read in parm parameters
05896   read_parm(gf,min);
05897 }

Parameters::~Parameters  ) 
 

Definition at line 294 of file Parameters.C.

References angle_array, bond_array, crossterm_array, dihedral_array, improper_array, ResizeArray< Elem >::resize(), ResizeArray< Elem >::size(), tab_pair_tree, vdw_array, and vdw_pair_tree.

00296 {
00297         if (atomTypeNames)
00298           delete [] atomTypeNames;
00299 
00300   if (bondp != NULL)
00301     free_bond_tree(bondp);
00302 
00303   if (anglep != NULL)
00304     free_angle_tree(anglep);
00305 
00306   if (dihedralp != NULL)
00307     free_dihedral_list(dihedralp);
00308 
00309   if (improperp != NULL)
00310     free_improper_list(improperp);
00311 
00312   if (crosstermp != NULL)
00313     free_crossterm_list(crosstermp);
00314 
00315   if (vdwp != NULL)
00316     free_vdw_tree(vdwp);
00317 
00318   if (vdw_pairp != NULL)
00319     free_vdw_pair_list();
00320 
00321   if (bond_array != NULL)
00322     delete [] bond_array;
00323 
00324   if (angle_array != NULL)
00325     delete [] angle_array;
00326 
00327   if (dihedral_array != NULL)
00328     delete [] dihedral_array;
00329 
00330   if (improper_array != NULL)
00331     delete [] improper_array;
00332 
00333   if (crossterm_array != NULL)
00334     delete [] crossterm_array;
00335 
00336   if (vdw_array != NULL)
00337     delete [] vdw_array;
00338   
00339   if (tab_pair_tree != NULL)
00340     free_table_pair_tree(tab_pair_tree);
00341 
00342   if (vdw_pair_tree != NULL)
00343     free_vdw_pair_tree(vdw_pair_tree);
00344 
00345   if (maxDihedralMults != NULL)
00346     delete [] maxDihedralMults;
00347 
00348   if (maxImproperMults != NULL)
00349     delete [] maxImproperMults;
00350 
00351   for( int i = 0; i < error_msgs.size(); ++i ) {
00352     delete [] error_msgs[i];
00353   }
00354   error_msgs.resize(0);
00355 }


Member Function Documentation

void Parameters::assign_angle_index char *  ,
char *  ,
char *  ,
Angle
 

Definition at line 3661 of file Parameters.C.

References Angle, angle::angle_type, angle_params::atom1name, angle_params::atom2name, angle_params::atom3name, angle_params::index, angle_params::left, NAMD_die(), and angle_params::right.

Referenced by getAngleData().

03664 {
03665   struct angle_params *ptr;  //  Current position in tree
03666   int comp_val;      //  value from strcasecmp
03667   int found=0;      //  flag 1->found a match
03668   char tmp_name[15];    //  Temporary atom name
03669 
03670   /*  Check to make sure the files have all been read    */
03671   if (!AllFilesRead)
03672   {
03673     NAMD_die("Tried to assign angle index before all parameter files were read");
03674   }
03675 
03676   /*  We need atom1 < atom3.  If that was not what we were   */
03677   /*  passed, switch them            */
03678   if (strcasecmp(atom1, atom3) > 0)
03679   {
03680     strcpy(tmp_name, atom1);
03681     strcpy(atom1, atom3);
03682     strcpy(atom3, tmp_name);
03683   }
03684 
03685   /*  Start at the top            */
03686   ptr=anglep;
03687 
03688   /*  While we don't have a match and we haven't reached the  */
03689   /*  bottom of the tree, compare values        */
03690   while (!found && (ptr != NULL))
03691   {
03692     comp_val = strcasecmp(atom1, ptr->atom1name);
03693 
03694     if (comp_val == 0)
03695     {
03696       /*  Atom 1 matches, so compare atom 2    */
03697       comp_val = strcasecmp(atom2, ptr->atom2name);
03698       
03699       if (comp_val == 0)
03700       {
03701         /*  Atoms 1&2 match, try atom 3    */
03702         comp_val = strcasecmp(atom3, ptr->atom3name);
03703       }
03704     }
03705 
03706     if (comp_val == 0)
03707     {
03708       /*  Found a match        */
03709       found = 1;
03710       angle_ptr->angle_type = ptr->index;
03711     }
03712     else if (comp_val < 0)
03713     {
03714       /*  Go left          */
03715       ptr=ptr->left;
03716     }
03717     else
03718     {
03719       /*  Go right          */
03720       ptr=ptr->right;
03721     }
03722   }
03723 
03724   /*  Make sure we found a match          */
03725   if (!found)
03726   {
03727     char err_msg[128];
03728 
03729     sprintf(err_msg, "UNABLE TO FIND ANGLE PARAMETERS FOR %s %s %s",
03730        atom1, atom2, atom3);
03731     NAMD_die(err_msg);
03732   }
03733 
03734   return;
03735 }

void Parameters::assign_bond_index char *  ,
char *  ,
Bond
 

Definition at line 3573 of file Parameters.C.

References bond_params::atom1name, bond_params::atom2name, Bond, bond::bond_type, bond_params::index, bond_params::left, NAMD_die(), and bond_params::right.

Referenced by getBondData().

03575 {
03576   struct bond_params *ptr;  //  Current location in tree
03577   int found=0;      //  Flag 1-> found a match
03578   int cmp_code;      //  return code from strcasecmp
03579   char tmp_name[15];    //  Temporary atom name
03580 
03581   /*  Check to make sure the files have all been read    */
03582   if (!AllFilesRead)
03583   {
03584     NAMD_die("Tried to assign bond index before all parameter files were read");
03585   }
03586 
03587   /*  We need atom1 < atom2, so if that's not the way they        */
03588   /*  were passed, flip them          */
03589   if (strcasecmp(atom1, atom2) > 0)
03590   {
03591     strcpy(tmp_name, atom1);
03592     strcpy(atom1, atom2);
03593     strcpy(atom2, tmp_name);
03594   }
03595 
03596   /*  Start at the top            */
03597   ptr=bondp;
03598 
03599   /*  While we haven't found a match and we're not at the end  */
03600   /*  of the tree, compare the bond passed in with the tree  */
03601   while (!found && (ptr!=NULL))
03602   {
03603     cmp_code=strcasecmp(atom1, ptr->atom1name);
03604 
03605     if (cmp_code == 0)
03606     {
03607       cmp_code=strcasecmp(atom2, ptr->atom2name);
03608     }
03609 
03610     if (cmp_code == 0)
03611     {
03612       /*  Found a match        */
03613       found=1;
03614       bond_ptr->bond_type = ptr->index;
03615     }
03616     else if (cmp_code < 0)
03617     {
03618       /*  Go left          */
03619       ptr=ptr->left;
03620     }
03621     else
03622     {
03623       /*  Go right          */
03624       ptr=ptr->right;
03625     }
03626   }
03627 
03628   /*  Check to see if we found anything        */
03629   if (!found)
03630   {
03631     char err_msg[128];
03632 
03633     sprintf(err_msg, "CAN'T FIND BOND PARAMETERS FOR BOND %s - %s IN PARAMETER FILES", atom1, atom2);
03634     NAMD_die(err_msg);
03635   }
03636 
03637   return;
03638 }

void Parameters::assign_crossterm_index char *  ,
char *  ,
char *  ,
char *  ,
char *  ,
char *  ,
char *  ,
char *  ,
Crossterm
 

Definition at line 3957 of file Parameters.C.

References crossterm_params::atom1name, crossterm_params::atom2name, crossterm_params::atom3name, crossterm_params::atom4name, crossterm_params::atom5name, crossterm_params::atom6name, crossterm_params::atom7name, crossterm_params::atom8name, Crossterm, crossterm::crossterm_type, crossterm_params::index, NAMD_die(), and crossterm_params::next.

Referenced by getCrosstermData().

03960 {
03961   struct crossterm_params *ptr;  //  Current position in list
03962   int found=0;      //  Flag 1->found a match
03963 
03964   /*  Start at the head of the list        */
03965   ptr=crosstermp;
03966 
03967   /*  While we haven't fuond a match and haven't reached the end  */
03968   /*  of the list, keep looking          */
03969   while (!found && (ptr!=NULL))
03970   {
03971     /*  Do a linear search through the linked list of   */
03972     /*  crossterm parameters.  Since the list is arranged    */
03973     /*  with wildcard paramters at the end of the list, we  */
03974     /*  can simply do a linear search and be guaranteed that*/
03975     /*  we will find exact matches before wildcard matches. */
03976     /*  Also, we must check for an exact match, and a match */
03977     /*  in reverse, since they are really the same          */
03978     /*  physically.            */
03979     if ( ( (strcasecmp(ptr->atom1name, atom1)==0) || 
03980            (strcasecmp(ptr->atom1name, "X")==0) ) &&
03981        ( (strcasecmp(ptr->atom2name, atom2)==0) || 
03982            (strcasecmp(ptr->atom2name, "X")==0) ) &&
03983        ( (strcasecmp(ptr->atom3name, atom3)==0) || 
03984            (strcasecmp(ptr->atom3name, "X")==0) ) &&
03985        ( (strcasecmp(ptr->atom4name, atom4)==0) || 
03986            (strcasecmp(ptr->atom4name, "X")==0) ) )
03987     {
03988       /*  Found an exact match      */
03989       found=1;
03990     }
03991     else if ( ( (strcasecmp(ptr->atom4name, atom1)==0) || 
03992            (strcasecmp(ptr->atom4name, "X")==0) ) &&
03993        ( (strcasecmp(ptr->atom3name, atom2)==0) || 
03994            (strcasecmp(ptr->atom3name, "X")==0) ) &&
03995        ( (strcasecmp(ptr->atom2name, atom3)==0) || 
03996            (strcasecmp(ptr->atom2name, "X")==0) ) &&
03997        ( (strcasecmp(ptr->atom1name, atom4)==0) || 
03998            (strcasecmp(ptr->atom1name, "X")==0) ) )
03999     {
04000       /*  Found a reverse match      */
04001       found=1;
04002     }
04003     if ( ! found ) {
04004       /*  Didn't find a match, go to the next node  */
04005       ptr=ptr->next;
04006       continue;
04007     }
04008     found = 0;
04009     if ( ( (strcasecmp(ptr->atom5name, atom5)==0) || 
04010            (strcasecmp(ptr->atom5name, "X")==0) ) &&
04011        ( (strcasecmp(ptr->atom6name, atom6)==0) || 
04012            (strcasecmp(ptr->atom6name, "X")==0) ) &&
04013        ( (strcasecmp(ptr->atom7name, atom7)==0) || 
04014            (strcasecmp(ptr->atom7name, "X")==0) ) &&
04015        ( (strcasecmp(ptr->atom8name, atom8)==0) || 
04016            (strcasecmp(ptr->atom8name, "X")==0) ) )
04017     {
04018       /*  Found an exact match      */
04019       found=1;
04020     }
04021     else if ( ( (strcasecmp(ptr->atom8name, atom5)==0) || 
04022            (strcasecmp(ptr->atom8name, "X")==0) ) &&
04023        ( (strcasecmp(ptr->atom7name, atom6)==0) || 
04024            (strcasecmp(ptr->atom7name, "X")==0) ) &&
04025        ( (strcasecmp(ptr->atom6name, atom7)==0) || 
04026            (strcasecmp(ptr->atom6name, "X")==0) ) &&
04027        ( (strcasecmp(ptr->atom5name, atom8)==0) || 
04028            (strcasecmp(ptr->atom5name, "X")==0) ) )
04029     {
04030       /*  Found a reverse match      */
04031       found=1;
04032     }
04033     if ( ! found ) {
04034       /*  Didn't find a match, go to the next node  */
04035       ptr=ptr->next;
04036     }
04037   }
04038 
04039   /*  Make sure we found a match          */
04040   if (!found)
04041   {
04042     char err_msg[128];
04043 
04044     sprintf(err_msg, "CAN'T FIND CROSSTERM PARAMETERS FOR %s  %s  %s  %s  %s  %s  %s  %s",
04045        atom1, atom2, atom3, atom4, atom5, atom6, atom7, atom8);
04046     
04047     NAMD_die(err_msg);
04048   }
04049 
04050   /*  Assign the constants          */
04051   crossterm_ptr->crossterm_type = ptr->index;
04052 
04053   return;
04054 }

void Parameters::assign_dihedral_index char *  ,
char *  ,
char *  ,
char *  ,
Dihedral ,
int 
 

Definition at line 3759 of file Parameters.C.

References dihedral_params::atom1name, dihedral_params::atom1wild, dihedral_params::atom2name, dihedral_params::atom2wild, dihedral_params::atom3name, dihedral_params::atom3wild, dihedral_params::atom4name, dihedral_params::atom4wild, Dihedral, dihedral_array, dihedral::dihedral_type, dihedral_params::index, DihedralValue::multiplicity, NAMD_die(), and dihedral_params::next.

Referenced by getDihedralData().

03763 {
03764   struct dihedral_params *ptr;  //  Current position in list
03765   int found=0;      //  Flag 1->found a match
03766 
03767   /*  Start at the begining of the list        */
03768   ptr=dihedralp;
03769 
03770   /*  While we haven't found a match and we haven't reached       */
03771   /*  the end of the list, keep looking        */
03772   while (!found && (ptr!=NULL))
03773   {
03774     /*  Do a linear search through the linked list of   */
03775     /*  dihedral parameters.  Since the list is arranged    */
03776     /*  with wildcard paramters at the end of the list, we  */
03777     /*  can simply do a linear search and be guaranteed that*/
03778     /*  we will find exact matches before wildcard matches. */
03779     /*  Also, we must check for an exact match, and a match */
03780     /*  in reverse, since they are really the same          */
03781     /*  physically.            */
03782     if ( ( ptr->atom1wild || (strcasecmp(ptr->atom1name, atom1)==0) ) && 
03783          ( ptr->atom2wild || (strcasecmp(ptr->atom2name, atom2)==0) ) &&
03784          ( ptr->atom3wild || (strcasecmp(ptr->atom3name, atom3)==0) ) &&
03785          ( ptr->atom4wild || (strcasecmp(ptr->atom4name, atom4)==0) ) ) 
03786     {
03787       /*  Found an exact match      */
03788       found=1;
03789     }
03790     else if ( ( ptr->atom4wild || (strcasecmp(ptr->atom4name, atom1)==0) ) &&
03791               ( ptr->atom3wild || (strcasecmp(ptr->atom3name, atom2)==0) ) &&
03792               ( ptr->atom2wild || (strcasecmp(ptr->atom2name, atom3)==0) ) &&
03793               ( ptr->atom1wild || (strcasecmp(ptr->atom1name, atom4)==0) ) )
03794     {
03795       /*  Found a reverse match      */
03796       found=1;
03797     }
03798     else
03799     {
03800       /*  Didn't find a match, go to the next node  */
03801       ptr=ptr->next;
03802     }
03803   }
03804 
03805   /*  Make sure we found a match          */
03806   if (!found)
03807   {
03808     char err_msg[128];
03809 
03810     sprintf(err_msg, "CAN'T FIND DIHEDRAL PARAMETERS FOR %s  %s  %s  %s",
03811        atom1, atom2, atom3, atom4);
03812     
03813     NAMD_die(err_msg);
03814   }
03815 
03816   //  Check to make sure the number of multiples specified in the psf
03817   //  file doesn't exceed the number of parameters in the parameter
03818   //  files
03819   if (multiplicity > maxDihedralMults[ptr->index])
03820   {
03821     char err_msg[257];
03822 
03823     sprintf(err_msg, "Multiplicity of Paramters for diehedral bond %s %s %s %s of %d exceeded", atom1, atom2, atom3, atom4, maxDihedralMults[ptr->index]);
03824     NAMD_die(err_msg);
03825   }
03826 
03827   //  If the multiplicity from the current bond is larger than that
03828   //  seen in the past, increase the multiplicity for this bond
03829   if (multiplicity > dihedral_array[ptr->index].multiplicity)
03830   {
03831     dihedral_array[ptr->index].multiplicity = multiplicity;
03832   }
03833 
03834   dihedral_ptr->dihedral_type = ptr->index;
03835 
03836   return;
03837 }

void Parameters::assign_improper_index char *  ,
char *  ,
char *  ,
char *  ,
Improper ,
int 
 

Definition at line 3861 of file Parameters.C.

References improper_params::atom1name, improper_params::atom2name, improper_params::atom3name, improper_params::atom4name, Improper, improper_array, improper::improper_type, improper_params::index, ImproperValue::multiplicity, NAMD_die(), and improper_params::next.

Referenced by getImproperData().

03865 {
03866   struct improper_params *ptr;  //  Current position in list
03867   int found=0;      //  Flag 1->found a match
03868 
03869   /*  Start at the head of the list        */
03870   ptr=improperp;
03871 
03872   /*  While we haven't fuond a match and haven't reached the end  */
03873   /*  of the list, keep looking          */
03874   while (!found && (ptr!=NULL))
03875   {
03876     /*  Do a linear search through the linked list of   */
03877     /*  improper parameters.  Since the list is arranged    */
03878     /*  with wildcard paramters at the end of the list, we  */
03879     /*  can simply do a linear search and be guaranteed that*/
03880     /*  we will find exact matches before wildcard matches. */
03881     /*  Also, we must check for an exact match, and a match */
03882     /*  in reverse, since they are really the same          */
03883     /*  physically.            */
03884     if ( ( (strcasecmp(ptr->atom1name, atom1)==0) || 
03885            (strcasecmp(ptr->atom1name, "X")==0) ) &&
03886        ( (strcasecmp(ptr->atom2name, atom2)==0) || 
03887            (strcasecmp(ptr->atom2name, "X")==0) ) &&
03888        ( (strcasecmp(ptr->atom3name, atom3)==0) || 
03889            (strcasecmp(ptr->atom3name, "X")==0) ) &&
03890        ( (strcasecmp(ptr->atom4name, atom4)==0) || 
03891            (strcasecmp(ptr->atom4name, "X")==0) ) )
03892     {
03893       /*  Found an exact match      */
03894       found=1;
03895     }
03896     else if ( ( (strcasecmp(ptr->atom4name, atom1)==0) || 
03897            (strcasecmp(ptr->atom4name, "X")==0) ) &&
03898        ( (strcasecmp(ptr->atom3name, atom2)==0) || 
03899            (strcasecmp(ptr->atom3name, "X")==0) ) &&
03900        ( (strcasecmp(ptr->atom2name, atom3)==0) || 
03901            (strcasecmp(ptr->atom2name, "X")==0) ) &&
03902        ( (strcasecmp(ptr->atom1name, atom4)==0) || 
03903            (strcasecmp(ptr->atom1name, "X")==0) ) )
03904     {
03905       /*  Found a reverse match      */
03906       found=1;
03907     }
03908     else
03909     {
03910       /*  Didn't find a match, go to the next node  */
03911       ptr=ptr->next;
03912     }
03913   }
03914 
03915   /*  Make sure we found a match          */
03916   if (!found)
03917   {
03918     char err_msg[128];
03919 
03920     sprintf(err_msg, "CAN'T FIND IMPROPER PARAMETERS FOR %s  %s  %s  %s",
03921        atom1, atom2, atom3, atom4);
03922     
03923     NAMD_die(err_msg);
03924   }
03925 
03926   //  Check to make sure the number of multiples specified in the psf
03927   //  file doesn't exceed the number of parameters in the parameter
03928   //  files
03929   if (multiplicity > maxImproperMults[ptr->index])
03930   {
03931     char err_msg[257];
03932 
03933     sprintf(err_msg, "Multiplicity of Paramters for improper bond %s %s %s %s of %d exceeded", atom1, atom2, atom3, atom4, maxImproperMults[ptr->index]);
03934     NAMD_die(err_msg);
03935   }
03936 
03937   //  If the multiplicity from the current bond is larger than that
03938   //  seen in the past, increase the multiplicity for this bond
03939   if (multiplicity > improper_array[ptr->index].multiplicity)
03940   {
03941     improper_array[ptr->index].multiplicity = multiplicity;
03942   }
03943 
03944   /*  Assign the constants          */
03945   improper_ptr->improper_type = ptr->index;
03946 
03947   return;
03948 }

void Parameters::assign_vdw_index char *  ,
Atom
 

Definition at line 3285 of file Parameters.C.

References Atom, vdw_params::atomname, vdw_params::index, iout, iWARN(), vdw_params::left, NAMD_die(), and vdw_params::right.

03287 {
03288   struct vdw_params *ptr;    //  Current position in trees
03289   int found=0;      //  Flag 1->found match
03290   int comp_code;      //  return code from strcasecmp
03291 
03292   /*  Check to make sure the files have all been read    */
03293   if (!AllFilesRead)
03294   {
03295     NAMD_die("Tried to assign vdw index before all parameter files were read");
03296   }
03297 
03298   /*  Start at the top            */
03299   ptr=vdwp;
03300 
03301   /*  While we haven't found a match, and we haven't reached      */
03302   /*  the bottom of the tree, compare the atom passed in with     */
03303   /*  the current value and decide if we have a match, or if not, */
03304   /*  which way to go            */
03305   while (!found && (ptr!=NULL))
03306   {
03307     comp_code = strcasecmp(atomtype, ptr->atomname);
03308 
03309     if (comp_code == 0)
03310     {
03311       /*  Found a match!        */
03312       atom_ptr->vdw_type=ptr->index;
03313       found=1;
03314     }
03315     else if (comp_code < 0)
03316     {
03317       /*  Go to the left        */
03318       ptr=ptr->left;
03319     }
03320     else
03321     {
03322       /*  Go to the right        */
03323       ptr=ptr->right;
03324     }
03325   }
03326 
03327   //****** BEGIN CHARMM/XPLOR type changes
03328   if (!found)
03329   {
03330     // since CHARMM allows wildcards "*" in vdw typenames
03331     // we have to look again if necessary, this way, if
03332     // we already had an exact match, this is never executed
03333     size_t windx;                      //  wildcard index
03334 
03335     /*  Start again at the top                                */
03336     ptr=vdwp;
03337   
03338      while (!found && (ptr!=NULL))
03339      {
03340   
03341        // get index of wildcard wildcard, get index
03342        windx= strcspn(ptr->atomname,"*"); 
03343        if (windx == strlen(ptr->atomname))
03344        {
03345          // there is no wildcard here
03346          comp_code = strcasecmp(atomtype, ptr->atomname);   
03347        }
03348        else
03349        {
03350          comp_code = strncasecmp(atomtype, ptr->atomname, windx); 
03351        }  
03352 
03353        if (comp_code == 0)
03354        {
03355          /*  Found a match!                                */
03356          atom_ptr->vdw_type=ptr->index;
03357          found=1;
03358          char errbuf[100];
03359          sprintf(errbuf,"VDW TYPE NAME %s MATCHES PARAMETER TYPE NAME %s",
03360                         atomtype, ptr->atomname);
03361          int i;
03362          for(i=0; i<error_msgs.size(); i++) {
03363            if ( strcmp(errbuf,error_msgs[i]) == 0 ) break;
03364          }
03365          if ( i == error_msgs.size() ) {
03366            char *newbuf = new char[strlen(errbuf)+1];
03367            strcpy(newbuf,errbuf);
03368            error_msgs.add(newbuf);
03369            iout << iWARN << newbuf << "\n" << endi;
03370          }
03371        }
03372        else if (comp_code < 0)
03373        {
03374           /*  Go to the left                                */
03375                 ptr=ptr->left;
03376        }
03377        else
03378        {
03379          /*  Go to the right                                */
03380                 ptr=ptr->right;
03381        }
03382      
03383      }
03384                 
03385   }
03386   //****** END CHARMM/XPLOR type changes
03387 
03388   /*  Make sure we found it          */
03389   if (!found)
03390   {
03391     char err_msg[100];
03392 
03393     sprintf(err_msg, "DIDN'T FIND vdW PARAMETER FOR ATOM TYPE %s",
03394        atomtype);
03395     NAMD_die(err_msg);
03396   }
03397 
03398   return;
03399 }

char* Parameters::atom_type_name Index  a  )  [inline]
 

Definition at line 336 of file Parameters.h.

References MAX_ATOMTYPE_CHARS.

00336                                       {
00337           return (atomTypeNames + (a * (MAX_ATOMTYPE_CHARS + 1)));
00338         }

void Parameters::done_reading_files  ) 
 

Definition at line 2821 of file Parameters.C.

References angle_array, bond_array, crossterm_array, dihedral_array, improper_array, MAX_ATOMTYPE_CHARS, NAMD_die(), NumCrosstermParams, NumDihedralParams, NumImproperParams, NumVdwParams, NumVdwParamsAssigned, vdw_array, and VdwValue.

Referenced by Parameters().

02823 {
02824   AllFilesRead = TRUE;
02825 
02826   //  Allocate space for all of the arrays
02827   if (NumBondParams)
02828   {
02829     bond_array = new BondValue[NumBondParams];
02830 
02831     if (bond_array == NULL)
02832     {
02833       NAMD_die("memory allocation of bond_array failed!");
02834     }
02835   }
02836 
02837   if (NumAngleParams)
02838   {
02839     angle_array = new AngleValue[NumAngleParams];
02840 
02841     if (angle_array == NULL)
02842     {
02843       NAMD_die("memory allocation of angle_array failed!");
02844     }
02845   }
02846 
02847   if (NumDihedralParams)
02848   {
02849     dihedral_array = new DihedralValue[NumDihedralParams];
02850 
02851     if (dihedral_array == NULL)
02852     {
02853       NAMD_die("memory allocation of dihedral_array failed!");
02854     }
02855     memset(dihedral_array, 0, NumDihedralParams*sizeof(DihedralValue));
02856   }
02857 
02858   if (NumImproperParams)
02859   {
02860     improper_array = new ImproperValue[NumImproperParams];
02861 
02862     if (improper_array == NULL)
02863     {
02864       NAMD_die("memory allocation of improper_array failed!");
02865     }
02866     memset(improper_array, 0, NumImproperParams*sizeof(ImproperValue));
02867   }
02868 
02869   if (NumCrosstermParams)
02870   {
02871     crossterm_array = new CrosstermValue[NumCrosstermParams];
02872     memset(crossterm_array, 0, NumCrosstermParams*sizeof(CrosstermValue));
02873   }
02874 
02875   if (NumVdwParams)
02876   {
02877           atomTypeNames = new char[NumVdwParams*(MAX_ATOMTYPE_CHARS+1)];
02878     vdw_array = new VdwValue[NumVdwParams];
02879     
02880     if (vdw_array == NULL)
02881     {
02882       NAMD_die("memory allocation of vdw_array failed!");
02883     }
02884   }
02885 
02886   //  Assign indexes to each of the parameters and populate the
02887   //  arrays using the binary trees and linked lists that we have
02888   //  already read in
02889   index_bonds(bondp, 0);
02890   index_angles(anglep, 0);
02891   NumVdwParamsAssigned = index_vdw(vdwp, 0);
02892   index_dihedrals();
02893   index_impropers();
02894   index_crossterms();
02895   
02896   //  Convert the vdw pairs
02897   convert_vdw_pairs();
02898   convert_table_pairs();
02899 }

void Parameters::done_reading_structure  ) 
 

Definition at line 4686 of file Parameters.C.

04688 {
04689   if (bondp != NULL)
04690     free_bond_tree(bondp);
04691 
04692   if (anglep != NULL)
04693     free_angle_tree(anglep);
04694 
04695   if (dihedralp != NULL)
04696     free_dihedral_list(dihedralp);
04697 
04698   if (improperp != NULL)
04699     free_improper_list(improperp);
04700 
04701   if (crosstermp != NULL)
04702     free_crossterm_list(crosstermp);
04703 
04704   if (vdwp != NULL)
04705     free_vdw_tree(vdwp);
04706 
04707   //  Free the arrays used to track multiplicity for dihedrals
04708   //  and impropers
04709   if (maxDihedralMults != NULL)
04710     delete [] maxDihedralMults;
04711 
04712   if (maxImproperMults != NULL)
04713     delete [] maxImproperMults;
04714 
04715   bondp=NULL;
04716   anglep=NULL;
04717   dihedralp=NULL;
04718   improperp=NULL;
04719   crosstermp=NULL;
04720   vdwp=NULL;
04721   maxImproperMults=NULL;
04722   maxDihedralMults=NULL;
04723 }

void Parameters::get_angle_params Real k,
Real theta0,
Real k_ub,
Real r_ub,
Index  index
[inline]
 

Definition at line 396 of file Parameters.h.

References AngleValue::k, AngleValue::k_ub, AngleValue::r_ub, and AngleValue::theta0.

Referenced by getAngleData().

00398         {
00399                 *k = angle_array[index].k;
00400                 *theta0 = angle_array[index].theta0;
00401                 *k_ub = angle_array[index].k_ub;
00402                 *r_ub = angle_array[index].r_ub;
00403         }

void Parameters::get_bond_params Real k,
Real x0,
Index  index
[inline]
 

Definition at line 390 of file Parameters.h.

References BondValue::k, and BondValue::x0.

Referenced by getBondData(), and Molecule::print_bonds().

00391         {
00392                 *k = bond_array[index].k;
00393                 *x0 = bond_array[index].x0;
00394         }

int Parameters::get_dihedral_multiplicity Index  index  )  [inline]
 

Definition at line 410 of file Parameters.h.

References DihedralValue::multiplicity.

00411         {
00412                 return(dihedral_array[index].multiplicity);
00413         }

void Parameters::get_dihedral_params Real k,
int *  n,
Real delta,
Index  index,
int  mult
[inline]
 

Definition at line 428 of file Parameters.h.

References four_body_consts::delta, four_body_consts::k, four_body_consts::n, NAMD_die(), and DihedralValue::values.

00430         {
00431                 if ( (mult<0) || (mult>MAX_MULTIPLICITY) )
00432                 {
00433                         NAMD_die("Bad mult index in Parameters::get_dihedral_params");
00434                 }
00435 
00436                 *k = dihedral_array[index].values[mult].k;
00437                 *n = dihedral_array[index].values[mult].n;
00438                 *delta = dihedral_array[index].values[mult].delta;
00439         }

int Parameters::get_improper_multiplicity Index  index  )  [inline]
 

Definition at line 405 of file Parameters.h.

References ImproperValue::multiplicity.

00406         {
00407                 return(improper_array[index].multiplicity);
00408         }

void Parameters::get_improper_params Real k,
int *  n,
Real delta,
Index  index,
int  mult
[inline]
 

Definition at line 415 of file Parameters.h.

References four_body_consts::delta, four_body_consts::k, four_body_consts::n, NAMD_die(), and ImproperValue::values.

00417         {
00418                 if ( (mult<0) || (mult>MAX_MULTIPLICITY) )
00419                 {
00420                         NAMD_die("Bad mult index in Parameters::get_improper_params");
00421                 }
00422 
00423                 *k = improper_array[index].values[mult].k;
00424                 *n = improper_array[index].values[mult].n;
00425                 *delta = improper_array[index].values[mult].delta;
00426         }

int Parameters::get_int_table_type char *   ) 
 

Definition at line 6939 of file Parameters.C.

References table_types.

06939                                                   {
06940   for (int i=0; i<tablenumtypes; i++) {
06941     if (!strncmp(tabletype, table_types[i], 5)) {
06942       return i;
06943     }
06944   }
06945 
06946   return -1;
06947 }

int Parameters::get_num_vdw_params void   )  [inline]
 

Definition at line 475 of file Parameters.h.

Referenced by LJTable::LJTable().

00475 { return NumVdwParamsAssigned; }

int Parameters::get_table_pair_params Index  ,
Index  ,
int * 
 

Definition at line 3420 of file Parameters.C.

References indexed_table_pair::ind1, indexed_table_pair::ind2, Index, IndexedTablePair, indexed_table_pair::K, indexed_table_pair::left, and indexed_table_pair::right.

03420                                                                     {
03421   IndexedTablePair *ptr;
03422   Index temp;
03423   int found=FALSE;
03424 
03425   ptr=tab_pair_tree;
03426   //
03427   //  We need the smaller type in ind1, so if it isn't already that 
03428   //  way, switch them        */
03429   if (ind1 > ind2)
03430   {
03431     temp = ind1;
03432     ind1 = ind2;
03433     ind2 = temp;
03434   }
03435 
03436   /*  While we haven't found a match and we're not at the end  */
03437   /*  of the tree, compare the bond passed in with the tree  */
03438   while (!found && (ptr!=NULL))
03439   {
03440 //    printf("Comparing %i with %i and %i with %i\n", ind1, ptr->ind1, ind2, ptr->ind2);
03441     if ( (ind1 == ptr->ind1) && (ind2 == ptr->ind2) )
03442     {
03443        found = TRUE;
03444     }
03445     else if ( (ind1 < ptr->ind1) || 
03446         ( (ind1==ptr->ind1) && (ind2 < ptr->ind2) ) )
03447     {
03448       /*  Go left          */
03449       ptr=ptr->left;
03450     }
03451     else
03452     {
03453       /*  Go right          */
03454       ptr=ptr->right;
03455     }
03456   }
03457 
03458   /*  If we found a match, assign the values      */
03459   if (found)
03460   {
03461     *K = ptr->K;
03462     return(TRUE);
03463   }
03464   else
03465   {
03466     return(FALSE);
03467   }
03468 }

int Parameters::get_vdw_pair_params Index  ind1,
Index  ind2,
Real ,
Real ,
Real ,
Real
 

Definition at line 3495 of file Parameters.C.

References indexed_vdw_pair::A, indexed_vdw_pair::A14, indexed_vdw_pair::B, indexed_vdw_pair::B14, indexed_vdw_pair::ind1, indexed_vdw_pair::ind2, Index, IndexedVdwPair, indexed_vdw_pair::left, and indexed_vdw_pair::right.

03498 {
03499   IndexedVdwPair *ptr;    //  Current location in tree
03500   Index temp;      //  Temporary value for swithcing
03501           // values
03502   int found=FALSE;    //  Flag 1-> found a match
03503 
03504   ptr=vdw_pair_tree;
03505 
03506   //  We need the smaller type in ind1, so if it isn't already that 
03507   //  way, switch them        */
03508   if (ind1 > ind2)
03509   {
03510     temp = ind1;
03511     ind1 = ind2;
03512     ind2 = temp;
03513   }
03514 
03515   /*  While we haven't found a match and we're not at the end  */
03516   /*  of the tree, compare the bond passed in with the tree  */
03517   while (!found && (ptr!=NULL))
03518   {
03519     if ( (ind1 == ptr->ind1) && (ind2 == ptr->ind2) )
03520     {
03521        found = TRUE;
03522     }
03523     else if ( (ind1 < ptr->ind1) || 
03524         ( (ind1==ptr->ind1) && (ind2 < ptr->ind2) ) )
03525     {
03526       /*  Go left          */
03527       ptr=ptr->left;
03528     }
03529     else
03530     {
03531       /*  Go right          */
03532       ptr=ptr->right;
03533     }
03534   }
03535 
03536   /*  If we found a match, assign the values      */
03537   if (found)
03538   {
03539     *A = ptr->A;
03540     *B = ptr->B;
03541     *A14 = ptr->A14;
03542     *B14 = ptr->B14;
03543 
03544     return(TRUE);
03545   }
03546   else
03547   {
03548     return(FALSE);
03549   }
03550 }

void Parameters::get_vdw_params Real sigma,
Real epsilon,
Real sigma14,
Real epsilon14,
Index  index
[inline]
 

Definition at line 441 of file Parameters.h.

References cbrt, vdw_val::epsilon, vdw_val::epsilon14, NAMD_die(), Real, vdw_val::sigma, and vdw_val::sigma14.

Referenced by ComputeNonbondedCUDA::doWork(), and Molecule::print_atoms().

00443   {
00444     if ( vdw_array ) {
00445       *sigma = vdw_array[index].sigma;
00446       *epsilon = vdw_array[index].epsilon;
00447       *sigma14 = vdw_array[index].sigma14;
00448       *epsilon14 = vdw_array[index].epsilon14;
00449     } else {
00450       // sigma and epsilon from A and B for each vdw type's interaction with itself
00451       Real A; Real B; Real A14; Real B14;
00452       if (get_vdw_pair_params(index, index, &A, &B, &A14, &B14) ) {
00453         if (A && B) {
00454           *sigma = sqrt(cbrt(A)) / sqrt(cbrt(B));
00455           *epsilon = B*B / (4*A);
00456         }
00457         else {
00458           *sigma = 0; *epsilon = 0;
00459         }
00460         if (A14 && B14) {
00461           *sigma14 = sqrt(cbrt(A14)) / sqrt(cbrt(B14));
00462           *epsilon14 = B14*B14 / (4*A14);
00463         }
00464         else {
00465           *sigma14 = 0; *epsilon14 = 0;
00466         }
00467       }
00468       else {NAMD_die ("Function get_vdw_params failed to derive Lennard-Jones sigma and epsilon from A and B values\n");}
00469     }
00470   }

void Parameters::print_angle_params  ) 
 

Definition at line 4570 of file Parameters.C.

References DebugM, and NumAngleParams.

04571 {
04572   DebugM(3,NumAngleParams << " ANGLE PARAMETERS\n"
04573       << "*****************************************" );
04574   traverse_angle_params(anglep);
04575 }

void Parameters::print_bond_params  ) 
 

Definition at line 4552 of file Parameters.C.

References DebugM, and NumBondParams.

04553 {
04554   DebugM(3,NumBondParams << " BOND PARAMETERS\n" \
04555       << "*****************************************"  \
04556       );
04557 
04558   traverse_bond_params(bondp);
04559 }

void Parameters::print_crossterm_params  ) 
 

void Parameters::print_dihedral_params  ) 
 

Definition at line 4586 of file Parameters.C.

References DebugM, and NumDihedralParams.

04587 {
04588   DebugM(3,NumDihedralParams << " DIHEDRAL PARAMETERS\n" \
04589       << "*****************************************" );
04590 
04591   traverse_dihedral_params(dihedralp);
04592 }

void Parameters::print_improper_params  ) 
 

Definition at line 4603 of file Parameters.C.

References DebugM, and NumImproperParams.

04604 {
04605   DebugM(3,NumImproperParams << " IMPROPER PARAMETERS\n" \
04606       << "*****************************************" );
04607 
04608   traverse_improper_params(improperp);
04609 }

void Parameters::print_param_summary  ) 
 

Definition at line 4654 of file Parameters.C.

References iINFO(), iout, NumAngleParams, NumBondParams, NumCosAngles, NumCrosstermParams, NumDihedralParams, NumImproperParams, NumVdwPairParams, and NumVdwParams.

Referenced by NamdState::configListInit().

04655 {
04656   iout << iINFO << "SUMMARY OF PARAMETERS:\n" 
04657        << iINFO << NumBondParams << " BONDS\n" 
04658        << iINFO << NumAngleParams << " ANGLES\n" << endi;
04659        if (cosAngles) {
04660          iout << iINFO << "  " << (NumAngleParams - NumCosAngles) << " HARMONIC\n"
04661               << iINFO << "  " << NumCosAngles << " COSINE-BASED\n" << endi;
04662        }
04663   iout << iINFO << NumDihedralParams << " DIHEDRAL\n"
04664        << iINFO << NumImproperParams << " IMPROPER\n"
04665        << iINFO << NumCrosstermParams << " CROSSTERM\n"
04666        << iINFO << NumVdwParams << " VDW\n"
04667        << iINFO << NumVdwPairParams << " VDW_PAIRS\n" << endi;
04668 }

void Parameters::print_vdw_pair_params  ) 
 

Definition at line 4637 of file Parameters.C.

References DebugM, and NumVdwPairParams.

04638 {
04639   DebugM(3,NumVdwPairParams << " vdW PAIR PARAMETERS\n" \
04640       << "*****************************************" );
04641 
04642   traverse_vdw_pair_params(vdw_pairp);
04643 }

void Parameters::print_vdw_params  ) 
 

Definition at line 4620 of file Parameters.C.

References DebugM, and NumVdwParams.

04621 {
04622   DebugM(3,NumVdwParams << " vdW PARAMETERS\n" \
04623       << "*****************************************" );
04624 
04625   traverse_vdw_params(vdwp);
04626 }

void Parameters::read_charmm_parameter_file char *   ) 
 

Definition at line 509 of file Parameters.C.

References iout, iWARN(), NAMD_blank_string(), NAMD_die(), NAMD_find_first_word(), NAMD_read_line(), NumAngleParams, NumBondParams, NumCrosstermParams, NumDihedralParams, NumImproperParams, NumTablePairParams, NumVdwPairParams, NumVdwParams, and table_ener.

Referenced by Parameters().

00511 {
00512   int  par_type=0;         //  What type of parameter are we currently
00513                            //  dealing with? (vide infra)
00514   int  skipline;           //  skip this line?
00515   int  skipall = 0;        //  skip rest of file;
00516   char buffer[512];           //  Buffer to store each line of the file
00517   char first_word[512];           //  First word of the current line
00518   FILE *pfile;                   //  File descriptor for the parameter file
00519 
00520   /*  Check to make sure that we haven't previously been told     */
00521   /*  that all the files were read                                */
00522   if (AllFilesRead)
00523   {
00524     NAMD_die("Tried to read another parameter file after being told that all files were read!");
00525   }
00526 
00527   /*  Try and open the file                                        */
00528   if ( (pfile = fopen(fname, "r")) == NULL)
00529   {
00530     char err_msg[256];
00531 
00532     sprintf(err_msg, "UNABLE TO OPEN CHARMM PARAMETER FILE %s\n", fname);
00533     NAMD_die(err_msg);
00534   }
00535 
00536   /*  Keep reading in lines until we hit the EOF                        */
00537   while (NAMD_read_line(pfile, buffer) != -1)
00538   {
00539     /*  Get the first word of the line                        */
00540     NAMD_find_first_word(buffer, first_word);
00541     skipline=0;
00542 
00543     /*  First, screen out things that we ignore.                   */   
00544     /*  blank lines, lines that start with '!' or '*', lines that  */
00545     /*  start with "END".                                          */
00546     if (!NAMD_blank_string(buffer) &&
00547         (strncmp(first_word, "!", 1) != 0) &&
00548          (strncmp(first_word, "*", 1) != 0) &&
00549         (strncasecmp(first_word, "END", 3) != 0))
00550     {
00551       if ( skipall ) {
00552         iout << iWARN << "SKIPPING PART OF PARAMETER FILE AFTER RETURN STATEMENT\n" << endi;
00553         break;
00554       }
00555       /*  Now, determine the apropriate parameter type.   */
00556       if (strncasecmp(first_word, "bond", 4)==0)
00557       {
00558         par_type=1; skipline=1;
00559       }
00560       else if (strncasecmp(first_word, "angl", 4)==0)
00561       {
00562         par_type=2; skipline=1;
00563       }
00564       else if (strncasecmp(first_word, "thet", 4)==0)
00565       {
00566         par_type=2; skipline=1;
00567       }
00568       else if (strncasecmp(first_word, "dihe", 4)==0)
00569       {
00570         par_type=3; skipline=1;
00571       }
00572       else if (strncasecmp(first_word, "phi", 3)==0)
00573       {
00574         par_type=3; skipline=1;
00575       }
00576       else if (strncasecmp(first_word, "impr", 4)==0)
00577       {
00578         par_type=4; skipline=1;
00579       }
00580       else if (strncasecmp(first_word, "imph", 4)==0)
00581       {
00582         par_type=4; skipline=1;
00583       }
00584       else if (strncasecmp(first_word, "nbon", 4)==0)
00585       {
00586         par_type=5; skipline=1;
00587       }
00588       else if (strncasecmp(first_word, "nonb", 4)==0)
00589       {
00590         par_type=5; skipline=1;
00591       }
00592       else if (strncasecmp(first_word, "nbfi", 4)==0)
00593       {
00594         par_type=6; skipline=1;
00595       }
00596       else if (strncasecmp(first_word, "hbon", 4)==0)
00597       {
00598         par_type=7; skipline=1;
00599       }
00600       else if (strncasecmp(first_word, "cmap", 4)==0)
00601       {
00602         par_type=8; skipline=1;
00603       }
00604       else if (strncasecmp(first_word, "nbta", 4)==0)
00605       {
00606         par_type=9; skipline=1;
00607       }
00608       else if (strncasecmp(first_word, "read", 4)==0)
00609       {
00610         skip_stream_read(buffer, pfile);  skipline=1;
00611       }
00612       else if (strncasecmp(first_word, "return", 4)==0)
00613       {
00614         skipall=1;  skipline=1;
00615       }
00616       else if ((strncasecmp(first_word, "nbxm", 4) == 0) ||
00617                (strncasecmp(first_word, "grou", 4) == 0) ||
00618                (strncasecmp(first_word, "cdie", 4) == 0) ||
00619                (strncasecmp(first_word, "shif", 4) == 0) ||
00620                (strncasecmp(first_word, "vgro", 4) == 0) ||
00621                (strncasecmp(first_word, "vdis", 4) == 0) ||
00622                (strncasecmp(first_word, "vswi", 4) == 0) ||
00623                (strncasecmp(first_word, "cutn", 4) == 0) ||
00624                (strncasecmp(first_word, "ctof", 4) == 0) ||
00625                (strncasecmp(first_word, "cton", 4) == 0) ||
00626                (strncasecmp(first_word, "eps", 3) == 0) ||
00627                (strncasecmp(first_word, "e14f", 4) == 0) ||
00628                (strncasecmp(first_word, "wmin", 4) == 0) ||
00629                (strncasecmp(first_word, "aexp", 4) == 0) ||
00630                (strncasecmp(first_word, "rexp", 4) == 0) ||
00631                (strncasecmp(first_word, "haex", 4) == 0) ||
00632                (strncasecmp(first_word, "aaex", 4) == 0) ||
00633                (strncasecmp(first_word, "noac", 4) == 0) ||
00634                (strncasecmp(first_word, "hbno", 4) == 0) ||
00635                (strncasecmp(first_word, "cuth", 4) == 0) ||
00636                (strncasecmp(first_word, "ctof", 4) == 0) ||
00637                (strncasecmp(first_word, "cton", 4) == 0) ||
00638                (strncasecmp(first_word, "cuth", 4) == 0) ||
00639                (strncasecmp(first_word, "ctof", 4) == 0) ||
00640                (strncasecmp(first_word, "cton", 4) == 0) ) 
00641       {
00642         if ((par_type != 5) && (par_type != 6) && (par_type != 7) && (par_type != 9))
00643         {
00644           char err_msg[512];
00645 
00646           sprintf(err_msg, "ERROR IN CHARMM PARAMETER FILE %s\nLINE=*%s*",fname, buffer);
00647           NAMD_die(err_msg);
00648         }
00649         else 
00650         {
00651           skipline = 1;
00652         }
00653       }        
00654       else if (par_type == 0)
00655       {
00656         /*  This is an unknown paramter.        */
00657         /*  This is BAD                                */
00658         char err_msg[512];
00659 
00660         sprintf(err_msg, "UNKNOWN PARAMETER IN CHARMM PARAMETER FILE %s\nLINE=*%s*",fname, buffer);
00661         NAMD_die(err_msg);
00662       }
00663     }
00664     else
00665     {
00666       skipline=1;
00667     }
00668 
00669     if ( (par_type != 0) && (!skipline) )
00670     {
00671       /*  Now, call the appropriate function based    */
00672       /*  on the type of parameter we have                */
00673       /*  I know, this should really be a switch ...  */
00674       if (par_type == 1)
00675       {
00676         add_bond_param(buffer);
00677         NumBondParams++;
00678       }
00679       else if (par_type == 2)
00680       {
00681         add_angle_param(buffer);
00682         NumAngleParams++;
00683       }
00684       else if (par_type == 3)
00685       {
00686         add_dihedral_param(buffer, pfile);
00687         NumDihedralParams++;
00688       }
00689       else if (par_type == 4)
00690       {
00691         add_improper_param(buffer, pfile);
00692         NumImproperParams++;
00693       }
00694       else if (par_type == 5)
00695       {
00696         add_vdw_param(buffer);
00697         NumVdwParams++;
00698       }
00699       else if (par_type == 6)
00700       {
00701         add_vdw_pair_param(buffer);
00702         NumVdwPairParams++; 
00703       }
00704       else if (par_type == 7)
00705       {
00706         add_hb_pair_param(buffer);                  
00707       }
00708       else if (par_type == 8)
00709       {
00710         add_crossterm_param(buffer, pfile);                  
00711         NumCrosstermParams++;
00712       }
00713       else if (par_type == 9)
00714       {
00715 
00716         if (table_ener == NULL) {
00717           continue;
00718         }
00719 
00720         add_table_pair_param(buffer);                  
00721         NumTablePairParams++;
00722       }
00723       else
00724       {
00725         /*  This really should not occour!      */
00726         /*  This is an internal error.          */
00727         /*  This is VERY BAD                        */
00728         char err_msg[512];
00729 
00730         sprintf(err_msg, "INTERNAL ERROR IN CHARMM PARAMETER FILE %s\nLINE=*%s*",fname, buffer);
00731         NAMD_die(err_msg);
00732       }
00733     }
00734   }
00735 
00736   /*  Close the file                                                */
00737   fclose(pfile);
00738 
00739   return;
00740 }

void Parameters::read_ener_table SimParameters  ) 
 

Definition at line 6103 of file Parameters.C.

References BigReal, SimParameters::cutoff, iout, mynearbyint, NAMD_die(), numenerentries, read_energy_type(), read_energy_type_bothcubspline(), Real, simParams, table_ener, table_spacing, table_types, SimParameters::tableInterpType, SimParameters::tableMaxDist, tablenumtypes, SimParameters::tableNumTypes, SimParameters::tableSpacing, and SimParameters::tabulatedEnergiesFile.

Referenced by Parameters().

06103                                                          {
06104         char* table_file = simParams->tabulatedEnergiesFile;
06105   char* interp_type = simParams->tableInterpType;
06106         FILE* enertable;
06107         enertable = fopen(table_file, "r");
06108 
06109         if (enertable == NULL) {
06110                 NAMD_die("ERROR: Could not open energy table file!\n");
06111         }
06112 
06113         char tableline[121];
06114   char* newtypename;
06115   int numtypes;
06116         int atom1;
06117         int atom2;
06118         int distbin;
06119   int readcount;
06120         Real dist;
06121         Real energy;
06122         Real force;
06123         Real table_spacing;
06124         Real maxdist;
06125 
06126 /* First find the header line and set the variables we need */
06127         iout << "Beginning table read\n" << endi;
06128         while(fgets(tableline,120,enertable)) {
06129                 if (strncmp(tableline,"#",1)==0) {
06130                         continue;
06131                 }
06132     readcount = sscanf(tableline, "%i %f %f", &numtypes, &table_spacing, &maxdist);
06133     if (readcount != 3) {
06134       NAMD_die("ERROR: Couldn't parse table header line\n");
06135     }
06136     break;
06137   }
06138 
06139   if (maxdist < simParams->cutoff) {
06140     NAMD_die("Tabulated energies must at least extend to the cutoff distance\n");
06141   }
06142 
06143         if (maxdist > simParams->cutoff) {
06144                 iout << "Info: Reducing max table size to match nonbond cutoff\n";
06145                 maxdist = ceil(simParams->cutoff);
06146         }
06147 
06148 /* Now allocate memory for the table; we know what we should be getting */
06149         numenerentries = 2 * numtypes * int (mynearbyint(maxdist/table_spacing) + 1);
06150         //Set up a full energy lookup table from a file
06151         //Allocate the table; layout is atom1 x atom2 x distance energy force
06152         fprintf(stdout, "Table has %i entries\n",numenerentries);
06153         //iout << "Allocating original energy table\n" << endi;
06154         if ( table_ener ) delete [] table_ener;
06155         table_ener = new BigReal[numenerentries];
06156   if ( table_types ) delete [] table_types;
06157   table_types = new char*[numtypes];
06158 
06159 /* Finally, start reading the table */
06160   int numtypesread = 0; //Number of types read so far
06161 
06162         while(fgets(tableline,120,enertable)) {
06163                 if (strncmp(tableline,"#",1)==0) {
06164                         continue;
06165                 }
06166     if (strncmp(tableline,"TYPE",4)==0) {
06167       // Read a new type
06168       newtypename = new char[5];
06169       int readcount = sscanf(tableline, "%*s %s", newtypename);
06170       if (readcount != 1) {
06171         NAMD_die("Couldn't read interaction type from TYPE line\n");
06172       }
06173 //      printf("Setting table_types[%i] to %s\n", numtypesread, newtypename);
06174       table_types[numtypesread] = newtypename;
06175 
06176       if (numtypesread == numtypes) {
06177         NAMD_die("Error: Number of types in table doesn't match table header\n");
06178       }
06179 
06180       // Read the current energy type with the proper interpolation
06181       if (!strncasecmp(interp_type, "linear", 6)) {
06182         if (read_energy_type(enertable, numtypesread, table_ener, table_spacing, maxdist) != 0) {
06183           char err_msg[512];
06184           sprintf(err_msg, "Failed to read table energy (linear) type %s\n", newtypename);
06185           NAMD_die(err_msg);
06186         }
06187       } else if (!strncasecmp(interp_type, "cubic", 5)) {
06188         if (read_energy_type_bothcubspline(enertable, numtypesread, table_ener, table_spacing, maxdist) != 0) {
06189           char err_msg[512];
06190           sprintf(err_msg, "Failed to read table energy (cubic) type %s\n", newtypename);
06191           NAMD_die(err_msg);
06192         }
06193       } else {
06194         NAMD_die("Table interpolation type must be linear or cubic\n");
06195       }
06196 
06197       numtypesread++;
06198       continue;
06199     }
06200     //char err_msg[512];
06201     //sprintf(err_msg, "Unrecognized line in energy table file: %s\n", tableline);
06202     //NAMD_die(err_msg);
06203   }
06204 
06205   // See if we got what we expected
06206   if (numtypesread != numtypes) {
06207     char err_msg[512];
06208     sprintf(err_msg, "ERROR: Expected %i tabulated energy types but got %i\n", numtypes, numtypesread);
06209     NAMD_die(err_msg);
06210   }
06211 
06212   // Move the necessary information into simParams
06213   simParams->tableNumTypes = numtypes;
06214   simParams->tableSpacing = table_spacing;
06215   simParams->tableMaxDist = maxdist;
06216   tablenumtypes = numtypes;
06217 
06218   /*
06219 xxxxxx
06220         int numtypes = simParams->tableNumTypes;
06221 
06222         //This parameter controls the table spacing
06223         BigReal table_spacing = 0.01;
06224         BigReal maxdist = 20.0;
06225         if (maxdist > simParams->cutoff) {
06226                 iout << "Info: Reducing max table size to match nonbond cutoff\n";
06227                 maxdist = ceil(simParams->cutoff);
06228         }
06229 
06230         numenerentries = (numtypes + 1) * numtypes * int (ceil(maxdist/table_spacing));
06231 //      fprintf(stdout, "Table arithmetic: maxdist %f, table_spacing %f, numtypes %i, numentries %i\n", maxdist, table_spacing, numtypes, numenerentries);
06232         columnsize = 2 * mynearbyint(maxdist/table_spacing);
06233         rowsize = numtypes * columnsize;
06234         //Set up a full energy lookup table from a file
06235         //Allocate the table; layout is atom1 x atom2 x distance energy force
06236         fprintf(stdout, "Table has %i entries\n",numenerentries);
06237         //iout << "Allocating original energy table\n" << endi;
06238         if ( table_ener ) delete [] table_ener;
06239         table_ener = new Real[numenerentries];
06240         //
06241         //Set sentinel values for uninitialized table energies
06242         for (int i=0 ; i<numenerentries ; i++) {
06243                 table_ener[i] = 1337.0;
06244         }
06245         Real compval = 1337.0;
06246 
06247         //    iout << "Entering table reading\n" << endi;
06248         //iout << "Finished allocating table\n" << endi;
06249 
06250         //Counter to make sure we read the right number of energy entries
06251         int numentries = 0;
06252 
06253         //Now, start reading from the file
06254         char* table_file = simParams->tabulatedEnergiesFile;
06255         FILE* enertable;
06256         enertable = fopen(table_file, "r");
06257 
06258         if (enertable == NULL) {
06259                 NAMD_die("ERROR: Could not open energy table file!\n");
06260         }
06261 
06262         char tableline[121];
06263         int atom1;
06264         int atom2;
06265         int distbin;
06266         Real dist;
06267         Real energy;
06268         Real force;
06269 
06270         iout << "Beginning table read\n" << endi;
06271         //Read the table entries
06272         while(fgets(tableline,120,enertable)) {
06273 //              IOut << "Processing line " << tableline << "\n" << endi;
06274                 if (strncmp(tableline,"#",1)==0) {
06275                         continue;
06276                 }
06277 
06278 
06279                 sscanf(tableline, "%i %i %f %f %f\n", &atom1, &atom2, &dist, &energy, &force);
06280                 distbin = int(mynearbyint(dist/table_spacing));
06281 //              fprintf(stdout, "Working on atoms %i and %i at distance %f\n",atom1,atom2,dist);
06282                 if ((atom2 > atom1) || (distbin > int(mynearbyint(maxdist/table_spacing)))) {
06283 //                      fprintf(stdout,"Rejected\n");
06284 //                      fprintf(stdout, "Error: Throwing out energy line beyond bounds\n");
06285         //              fprintf(stdout, "Bad line: Atom 1: %i Atom 2: %i Distance Bin: %i Max Distance Bin: %i \n", atom1, atom2, distbin, int(mynearbyint(maxdist/table_spacing)));
06286                 } else {
06287                         //The magic formula for the number of columns from previous rows
06288                         //in the triangular matrix is (2ni+i-i**2)/2
06289                         //Here n is the number of types, and i is atom2
06290 //                      fprintf(stdout, "Input: atom1 %f atom2 %f columnsize %f \n", float(atom1), float(atom2), float(columnsize));
06291 //                      fprintf(stdout, "Testing arithmetic: Part1: %i Part2: %i Part3: %i Total: %i\n", columnsize*((2*numtypes*atom2 + atom2 - atom2*atom2)/2), columnsize*(atom1-atom2), 2*distbin, columnsize*((2*numtypes*atom2 + atom2 - atom2*atom2)/2) + columnsize*(atom1-atom2) + 2*distbin - 2);
06292                         int eneraddress = columnsize*((2*numtypes*atom2 + atom2 - atom2*atom2)/2) + columnsize*(atom1-atom2) + 2*distbin - 2;
06293                         int forceaddress = eneraddress + 1;
06294 //                              fprintf(stdout, "Tableloc: %i Atom 1: %i Atom 2: %i Distance Bin: %i Energy: %f Force: %f\n", eneraddress, atom1, atom2, distbin, table_ener[eneraddress], table_ener[forceaddress]);
06295                 fflush(stdout);
06296 //                      fprintf(stdout, "Checking for dupes: Looking at: %f %f \n", table_ener[eneraddress], table_ener[forceaddress]);
06297                         if ((table_ener[eneraddress] == compval && table_ener[forceaddress] == compval)) {
06298                                 numentries++;
06299                                 table_ener[eneraddress] = energy;
06300                                 table_ener[forceaddress] = force;
06301 //                              fprintf(stdout, "Tableloc: %i Atom 1: %i Atom 2: %i Distance Bin: %i Energy: %f Force: %f\n", eneraddress, atom1, atom2, distbin, table_ener[eneraddress], table_ener[forceaddress]);
06302                                 //table_ener[rowsize*atom2 + columnsize*atom1 + 2*distbin] = energy;
06303                                 //table_ener[rowsize*atom2 + columnsize*atom1 + 2*distbin + 1] = force;
06304 //                              fflush(stdout);
06305                         } else {
06306 //                              fprintf(stdout,"Rejecting duplicate entry\n");
06307                         }
06308                 }
06309                 //      iout << "Done with line\n"<< endi;
06310         }
06311 
06312         //    int should = numtypes * numtypes * (maxdist/table_spacing);
06313         //    cout << "Entries: " << numentries << " should be " << should << "\n" << endi;
06314 //      int exptypes = ceil((numtypes+1) * numtypes * (maxdist/table_spacing));
06315 //fprintf(stdout,"numtypes: %i maxdist: %f table_spacing: %f exptypes: %i\n",numtypes,maxdist,table_spacing);
06316         if (numentries != int(numenerentries/2)) {
06317                 fprintf(stdout,"ERROR: Expected %i entries but got %i\n",numenerentries/2, numentries);
06318                 NAMD_die("Number of energy table entries did not match expected value\n");
06319         }
06320         //      iout << "Done with table\n"<< endi;
06321         fprintf(stdout, "Numenerentries: %i\n",numenerentries/2);
06322   */
06323 } /* END of function read_ener_table */ 

int Parameters::read_energy_type FILE *  ,
const   int,
BigReal ,
const   float,
const   float
 

Definition at line 6821 of file Parameters.C.

References BigReal, mynearbyint, NAMD_die(), and table_spacing.

Referenced by read_ener_table().

06821                                                                                                                                           {
06822 
06823   char tableline[120];
06824 
06825   /* Last values read from table */
06826   BigReal readdist;
06827   BigReal readener;
06828   BigReal readforce;
06829   BigReal lastdist;
06830   BigReal lastener;
06831   BigReal lastforce;
06832   readdist = -1.0;
06833   readener = 0.0;
06834   readforce = 0.0;
06835 
06836   /* Keep track of where in the table we are */
06837   float currdist;
06838   int distbin;
06839   currdist = -1.0;
06840   distbin = -1;
06841 
06842         while(fgets(tableline,120,enertable) && distbin <= (int) (mynearbyint(maxdist / table_spacing) + 1)) {
06843     printf("At distance %f + %f vs. %f\n", currdist, table_spacing, maxdist);
06844                 if (strncmp(tableline,"#",1)==0) {
06845                         continue;
06846                 }
06847     if (strncmp(tableline,"TYPE",4)==0) {
06848       break;
06849     }
06850 
06851     // Read an energy line from the table
06852     lastdist = readdist;
06853     lastener = readener;
06854     lastforce = readforce;
06855     int readcount = sscanf(tableline, "%lf %lf %lf", &readdist, &readener, &readforce);
06856     if (distbin == -1) {
06857       if (readdist != 0.0) {
06858         NAMD_die("ERROR: Energy/force tables must start at d=0.0\n");
06859       } else {
06860         distbin = 0;
06861         continue;
06862       }
06863     }
06864    // printf("Read an energy line: %f %f %f\n", readdist, readener, readforce);
06865     if (readcount != 3) {
06866       char err_msg[512];
06867       sprintf(err_msg, "ERROR: Failed to parse table line %s!\n", tableline);
06868       NAMD_die(err_msg);
06869     }
06870 
06871     //Sanity check the current entry
06872     if (readdist < lastdist) {
06873       NAMD_die("ERROR: Encountered badly ordered entries in energy table!\n");
06874     }
06875 
06876     currdist = lastdist;
06877 
06878     while (currdist <= readdist && distbin <= (int) (mynearbyint(maxdist / table_spacing))) {
06879       distbin = (int) (mynearbyint(currdist / table_spacing));
06880       int table_loc = 2 * (distbin + (typeindex * (1 + (int) mynearbyint(maxdist / table_spacing))));
06881       printf("Doing interpolation for energy between %f %f and %f %f: Dist %f\n", readener, readdist, lastener, lastdist, currdist);
06882       table_ener[table_loc] = interp_lin(readener, lastener, readdist, lastdist, currdist);
06883       table_ener[table_loc + 1] = interp_lin(readforce, lastforce, readdist, lastdist, currdist);
06884       printf("Adding energy/force entry: %f %f in distbin %i (distance %f) to address %i/%i\n", table_ener[table_loc], table_ener[table_loc + 1], distbin, currdist, table_loc, table_loc + 1);
06885       currdist += table_spacing;
06886       distbin++;
06887     }
06888   }
06889 
06890   // Go back one line, since we should now be into the next TYPE block
06891   fseek(enertable, -1 * strlen(tableline), SEEK_CUR); 
06892 
06893   // Clean up and make sure everything worked ok
06894   distbin--;
06895   printf("Testing: %i vs %i (from %f / %f)\n", distbin, (int) (mynearbyint(maxdist / table_spacing)), maxdist, table_spacing);
06896   if (distbin != (int) (mynearbyint(maxdist / table_spacing))) return 1;
06897   return 0;
06898 }

int Parameters::read_energy_type_bothcubspline FILE *  ,
const   int,
BigReal ,
const   float,
const   float
 

Definition at line 6346 of file Parameters.C.

References BigReal, INDEX, j, mynearbyint, NAMD_die(), and table_spacing.

Referenced by read_ener_table().

06346                                                                                                                                                         {
06347 
06348   char tableline[120];
06349   int i,j;
06350 
06351   /* Last values read from table */
06352   BigReal readdist;
06353   BigReal readener;
06354   BigReal readforce;
06355   BigReal spacing;
06356 //  BigReal readforce;
06357   BigReal lastdist;
06358 //  BigReal lastener;
06359 //  BigReal lastforce;
06360 //  readdist = -1.0;
06361 //  readener = 0.0;
06362 //  readforce = 0.0;
06363 
06364   /* Create arrays for holding the input values */
06365   std::vector<BigReal>  dists;
06366   std::vector<BigReal> enervalues;
06367   std::vector<BigReal> forcevalues;
06368   int numentries = 0;
06369 
06370 
06371   /* Keep track of where in the table we are */
06372   BigReal currdist;
06373   int distbin;
06374   currdist = 0.0;
06375   lastdist = -1.0;
06376   distbin = 0;
06377 
06378   // Read all of the values first -- we'll interpolate later
06379         while(fgets(tableline,120,enertable) && distbin <= (int) (mynearbyint(maxdist / table_spacing) + 1)) {
06380                 if (strncmp(tableline,"#",1)==0) {
06381                         continue;
06382                 }
06383     if (strncmp(tableline,"TYPE",4)==0) {
06384       fseek(enertable, -1 * strlen(tableline), SEEK_CUR); 
06385       break;
06386     }
06387 
06388     // Read an energy line from the table
06389     int readcount = sscanf(tableline, "%lf %lf %lf", &readdist, &readener, &readforce);
06390 
06391     //printf("Read an energy line: %g %g %g\n", readdist, readener, readforce);
06392     if (readcount != 3) {
06393       char err_msg[512];
06394       sprintf(err_msg, "ERROR: Failed to parse table line %s!\n", tableline);
06395       NAMD_die(err_msg);
06396     }
06397 
06398     //Sanity check the current entry
06399     if (readdist < lastdist) {
06400       NAMD_die("ERROR: Encountered badly ordered entries in energy table!\n");
06401     }
06402 
06403     lastdist = readdist;
06404     dists.push_back(readdist);
06405     enervalues.push_back(readener);
06406     forcevalues.push_back(readforce);
06407     numentries++;
06408   }
06409 
06410   // Check the spacing and make sure it is uniform
06411   if (dists[0] != 0.0) {
06412     NAMD_die("ERROR: First data point for energy table must be at r=0\n");
06413   }
06414   spacing = dists[1] - dists[0];
06415   for (i=2; i<(numentries - 1); i++) {
06416     BigReal myspacing;
06417     myspacing = dists[i] - dists[i-1];
06418     if (fabs(myspacing - spacing) > 0.00001) {
06419       printf("Bad spacing in table: %f should be %f (between distances %f and %f)\n", myspacing, spacing, dists[i-1], dists[i]);
06420       NAMD_die("ERROR: Nonuniform table encountered on cubic interpolation. Use a uniformly spaced table or switch to linear interpolation.\n");
06421     }
06422   }
06423 
06424 // Perform cubic spline interpolation to get the energies and forces
06425 
06426   /* allocate spline coefficient matrix */
06427   // xe and xf / be and bf for energies and forces, respectively
06428   double* m = new double[numentries*numentries];
06429   double* xe = new double[numentries];
06430   double* xf = new double[numentries];
06431   double* be = new double[numentries];
06432   double* bf = new double[numentries];
06433 
06434   be[0] = 3 * (enervalues[1] - enervalues[0]);
06435   for (i=1; i < (numentries - 1); i++) {
06436 //    printf("Control point %i at %f\n", i, enervalues[i]);
06437     be[i] = 3 * (enervalues[i+1] - enervalues[i-1]);
06438 //    printf("be is %f\n", be[i]);
06439   }
06440   be[numentries - 1] = 3 * (enervalues[numentries - 1] - enervalues[numentries - 2]);
06441 
06442   bf[0] = 3 * (forcevalues[1] - forcevalues[0]);
06443   for (i=1; i < (numentries - 1); i++) {
06444 //    printf("Control point %i at %f\n", i, forcevalues[i]);
06445     bf[i] = 3 * (forcevalues[i+1] - forcevalues[i-1]);
06446 //    printf("bf is %f\n", bf[i]);
06447   }
06448   bf[numentries - 1] = 3 * (forcevalues[numentries - 1] - forcevalues[numentries - 2]);
06449 
06450   memset(m,0,numentries*numentries*sizeof(double));
06451 
06452   /* initialize spline coefficient matrix */
06453   m[0] = 2;
06454   for (i = 1;  i < numentries;  i++) {
06455     m[INDEX(numentries,i-1,i)] = 1;
06456     m[INDEX(numentries,i,i-1)] = 1;
06457     m[INDEX(numentries,i,i)] = 4;
06458   }
06459   m[INDEX(numentries,numentries-1,numentries-1)] = 2;
06460 
06461   /* Now we need to solve the equation M X = b for X */
06462   // Do this simultaneously for energy and force -- ONLY because we use the same matrix
06463 
06464   //Put m in upper triangular form and apply corresponding operations to b
06465   for (i=0; i<numentries; i++) {
06466     // zero the ith column in all rows below i
06467     const BigReal baseval = m[INDEX(numentries,i,i)];
06468     for (j=i+1; j<numentries; j++) {
06469       const BigReal myval = m[INDEX(numentries,j,i)];
06470       const BigReal factor = myval / baseval;
06471 
06472       for (int k=i; k<numentries; k++) {
06473         const BigReal subval = m[INDEX(numentries,i,k)];
06474         m[INDEX(numentries,j,k)] -= (factor * subval);
06475       }
06476 
06477       be[j] -= (factor * be[i]);
06478       bf[j] -= (factor * bf[i]);
06479 
06480     }
06481   }
06482 
06483   //Now work our way diagonally up from the bottom right to assign values to X
06484   for (i=numentries-1; i>=0; i--) {
06485 
06486     //Subtract the effects of previous columns
06487     for (j=i+1; j<numentries; j++) {
06488       be[i] -= ( m[INDEX(numentries,i,j)] * xe[j] );
06489       bf[i] -= ( m[INDEX(numentries,i,j)] * xf[j] );
06490     }
06491 
06492     xe[i] = be[i] / m[INDEX(numentries,i,i)];
06493     xf[i] = bf[i] / m[INDEX(numentries,i,i)];
06494 
06495   }
06496 
06497   // We now have the coefficient information we need to make the table
06498   // Now interpolate on each interval we want
06499 
06500   distbin = 0;
06501   int entriesperseg = (int) ceil(spacing / table_spacing);
06502   int table_prefix = 2 * typeindex * (int) (mynearbyint(maxdist / table_spacing) + 1);
06503 
06504   for (i=0; i<numentries-1; i++) {
06505     BigReal Ae,Be,Ce,De;
06506     BigReal Af,Bf,Cf,Df;
06507     currdist = dists[i];
06508 
06509 //    printf("Interpolating on interval %i\n", i);
06510 
06511     // Set up the coefficients for this section
06512     Ae = enervalues[i];
06513     Be = xe[i];
06514     Ce = 3 * (enervalues[i+1] - enervalues[i]) - (2 * xe[i]) - xe[i+1];
06515     De = 2 * (enervalues[i] - enervalues[i+1]) + xe[i] + xe[i+1];
06516 
06517     Af = forcevalues[i];
06518     Bf = xf[i];
06519     Cf = 3 * (forcevalues[i+1] - forcevalues[i]) - (2 * xf[i]) - xf[i+1];
06520     Df = 2 * (forcevalues[i] - forcevalues[i+1]) + xf[i] + xf[i+1];
06521 
06522     // Go over the region of interest and fill in the table
06523     for (j=0; j<entriesperseg; j++) {
06524       const BigReal mydist = currdist + (j * table_spacing);
06525       const BigReal mydistfrac = (float) j / (entriesperseg - 1);
06526       distbin = (int) mynearbyint(mydist / table_spacing);
06527       if (distbin >= (int) mynearbyint(maxdist / table_spacing)) break;
06528       BigReal energy;
06529       BigReal force;
06530 
06531       energy = Ae + (Be * mydistfrac) + (Ce * mydistfrac * mydistfrac) + (De * mydistfrac * mydistfrac * mydistfrac);
06532       force = Af + (Bf * mydistfrac) + (Cf * mydistfrac * mydistfrac) + (Df * mydistfrac * mydistfrac * mydistfrac);
06533 
06534 //      printf("Adding energy/force entry %f / %f in bins %i / %i for distance %f (%f)\n", energy, force, (table_prefix + 2 * distbin), (table_prefix + 2 * distbin + 1), mydist, mydistfrac);
06535       table_ener[table_prefix + 2 * distbin] = energy;
06536       table_ener[table_prefix + 2 * distbin + 1] = force;
06537       distbin++;
06538     }
06539     if (currdist >= maxdist) break;
06540   }
06541 
06542   //The procedure above leaves out the last entry -- add it explicitly
06543   distbin = (int) mynearbyint(maxdist / table_spacing);
06544 //  printf("Adding energy/force entry %f / %f in bins %i / %i\n", enervalues[numentries - 1], 0.0, (table_prefix + 2 * distbin), (table_prefix + 2 * distbin + 1));
06545   table_ener[table_prefix + 2 * distbin] = enervalues[numentries - 1];
06546   table_ener[table_prefix + 2 * distbin + 1] = 0.0;
06547   distbin++;
06548 
06549 
06550   // Clean up and make sure everything worked ok
06551   delete m;
06552   delete xe;
06553   delete xf;
06554   delete be;
06555   delete bf;
06556   distbin--;
06557   printf("Testing: %i vs %i (from %f / %f)\n", distbin, (int) (mynearbyint(maxdist / table_spacing)), maxdist, table_spacing);
06558   if (distbin != (int) (mynearbyint(maxdist / table_spacing))) return 1;
06559   return 0;
06560 } /* end read_energy_type_bothcubspline */

int Parameters::read_energy_type_cubspline FILE *  ,
const   int,
BigReal ,
const   float,
const   float
 

Definition at line 6582 of file Parameters.C.

References BigReal, INDEX, j, mynearbyint, NAMD_die(), and table_spacing.

06582                                                                                                                                                     {
06583 
06584   char tableline[120];
06585   int i,j;
06586 
06587   /* Last values read from table */
06588   BigReal readdist;
06589   BigReal readener;
06590   BigReal spacing;
06591 //  BigReal readforce;
06592   BigReal lastdist;
06593 //  BigReal lastener;
06594 //  BigReal lastforce;
06595 //  readdist = -1.0;
06596 //  readener = 0.0;
06597 //  readforce = 0.0;
06598 
06599   /* Create arrays for holding the input values */
06600   std::vector<BigReal>  dists;
06601   std::vector<BigReal> enervalues;
06602   int numentries = 0;
06603 
06604 
06605   /* Keep track of where in the table we are */
06606   BigReal currdist;
06607   int distbin;
06608   currdist = 0.0;
06609   lastdist = -1.0;
06610   distbin = 0;
06611 
06612   // Read all of the values first -- we'll interpolate later
06613         while(fgets(tableline,120,enertable) && distbin <= (int) (mynearbyint(maxdist / table_spacing) + 1)) {
06614                 if (strncmp(tableline,"#",1)==0) {
06615                         continue;
06616                 }
06617     if (strncmp(tableline,"TYPE",4)==0) {
06618       fseek(enertable, -1 * strlen(tableline), SEEK_CUR); 
06619       break;
06620     }
06621 
06622     // Read an energy line from the table
06623     int readcount = sscanf(tableline, "%lf %lf", &readdist, &readener);
06624 
06625    // printf("Read an energy line: %f %f %f\n", readdist, readener, readforce);
06626     if (readcount != 2) {
06627       char err_msg[512];
06628       sprintf(err_msg, "ERROR: Failed to parse table line %s!\n", tableline);
06629       NAMD_die(err_msg);
06630     }
06631 
06632     //Sanity check the current entry
06633     if (readdist < lastdist) {
06634       NAMD_die("ERROR: Encountered badly ordered entries in energy table!\n");
06635     }
06636 
06637     lastdist = readdist;
06638     dists.push_back(readdist);
06639     enervalues.push_back(readener);
06640     numentries++;
06641   }
06642 
06643   // Check the spacing and make sure it is uniform
06644   if (dists[0] != 0.0) {
06645     NAMD_die("ERROR: First data point for energy table must be at r=0\n");
06646   }
06647   spacing = dists[1] - dists[0];
06648   for (i=2; i<(numentries - 1); i++) {
06649     BigReal myspacing;
06650     myspacing = dists[i] - dists[i-1];
06651     if (fabs(myspacing - spacing) > 0.00001) {
06652       printf("Bad spacing in table: %f should be %f (between distances %f and %f)\n", myspacing, spacing, dists[i-1], dists[i]);
06653       NAMD_die("ERROR: Nonuniform table encountered on cubic interpolation. Use a uniformly spaced table or switch to linear interpolation.\n");
06654     }
06655   }
06656 
06657 // Perform cubic spline interpolation to get the energies and forces
06658 
06659   /* allocate spline coefficient matrix */
06660   double* m = new double[numentries*numentries];
06661   double* x = new double[numentries];
06662   double* b = new double[numentries];
06663 
06664   b[0] = 3 * (enervalues[1] - enervalues[0]);
06665   for (i=1; i < (numentries - 1); i++) {
06666     printf("Control point %i at %f\n", i, enervalues[i]);
06667     b[i] = 3 * (enervalues[i+1] - enervalues[i-1]);
06668     printf("b is %f\n", b[i]);
06669   }
06670   b[numentries - 1] = 3 * (enervalues[numentries - 1] - enervalues[numentries - 2]);
06671 
06672   memset(m,0,numentries*numentries*sizeof(double));
06673 
06674   /* initialize spline coefficient matrix */
06675   m[0] = 2;
06676   for (i = 1;  i < numentries;  i++) {
06677     m[INDEX(numentries,i-1,i)] = 1;
06678     m[INDEX(numentries,i,i-1)] = 1;
06679     m[INDEX(numentries,i,i)] = 4;
06680   }
06681   m[INDEX(numentries,numentries-1,numentries-1)] = 2;
06682 
06683   /* Now we need to solve the equation M X = b for X */
06684 
06685   printf("Solving the matrix equation: \n");
06686 
06687   for (i=0; i<numentries; i++) {
06688     printf(" ( ");
06689     for (j=0; j<numentries; j++) {
06690       printf(" %6.3f,", m[INDEX(numentries, i, j)]);
06691     }
06692     printf(" )  ( D%-3i )  =  ( %6.3f )\n", i, b[i]);
06693   }
06694 
06695   //Put m in upper triangular form and apply corresponding operations to b
06696   for (i=0; i<numentries; i++) {
06697     // zero the ith column in all rows below i
06698     const BigReal baseval = m[INDEX(numentries,i,i)];
06699     for (j=i+1; j<numentries; j++) {
06700       const BigReal myval = m[INDEX(numentries,j,i)];
06701       const BigReal factor = myval / baseval;
06702 
06703       for (int k=i; k<numentries; k++) {
06704         const BigReal subval = m[INDEX(numentries,i,k)];
06705         m[INDEX(numentries,j,k)] -= (factor * subval);
06706       }
06707 
06708       b[j] -= (factor * b[i]);
06709 
06710     }
06711   }
06712 
06713   printf(" In upper diagonal form, equation is:\n");
06714   for (i=0; i<numentries; i++) {
06715     printf(" ( ");
06716     for (j=0; j<numentries; j++) {
06717       printf(" %6.3f,", m[INDEX(numentries, i, j)]);
06718     }
06719     printf(" )  ( D%-3i )  =  ( %6.3f )\n", i, b[i]);
06720   }
06721 
06722   //Now work our way diagonally up from the bottom right to assign values to X
06723   for (i=numentries-1; i>=0; i--) {
06724 
06725     //Subtract the effects of previous columns
06726     for (j=i+1; j<numentries; j++) {
06727       b[i] -= ( m[INDEX(numentries,i,j)] * x[j] );
06728     }
06729 
06730     x[i] = b[i] / m[INDEX(numentries,i,i)];
06731 
06732   }
06733 
06734   printf(" Solution vector is:\n\t(");
06735   for (i=0; i<numentries; i++) {
06736     printf(" %6.3f ", x[i]);
06737   }
06738   printf(" ) \n");
06739 
06740   // We now have the coefficient information we need to make the table
06741   // Now interpolate on each interval we want
06742 
06743   distbin = 0;
06744   int entriesperseg = (int) ceil(spacing / table_spacing);
06745   int table_prefix = 2 * typeindex * (int) (mynearbyint(maxdist / table_spacing) + 1);
06746 
06747   for (i=0; i<numentries-1; i++) {
06748     BigReal A,B,C,D;
06749     currdist = dists[i];
06750 
06751     printf("Interpolating on interval %i\n", i);
06752 
06753     // Set up the coefficients for this section
06754     A = enervalues[i];
06755     B = x[i];
06756     C = 3 * (enervalues[i+1] - enervalues[i]) - (2 * x[i]) - x[i+1];
06757     D = 2 * (enervalues[i] - enervalues[i+1]) + x[i] + x[i+1];
06758 
06759     printf("Coefficients for this interval: %f %f %f %f\n", A, B, C, D);
06760 
06761     // Go over the region of interest and fill in the table
06762     for (j=0; j<entriesperseg; j++) {
06763       const BigReal mydist = currdist + (j * table_spacing);
06764       const BigReal mydistfrac = (float) j / (entriesperseg - 1);
06765       distbin = (int) mynearbyint(mydist / table_spacing);
06766       if (distbin >= (int) mynearbyint(maxdist / table_spacing)) break;
06767       BigReal energy;
06768       BigReal force;
06769 
06770       energy = A + (B * mydistfrac) + (C * mydistfrac * mydistfrac) + (D * mydistfrac * mydistfrac * mydistfrac);
06771       force = B + (2 * C * mydistfrac) + (3 * D * mydistfrac * mydistfrac);
06772       // Multiply force by 1 / (interval length)
06773       force *= (1.0 / spacing);
06774 
06775       printf("Adding energy/force entry %f / %f in bins %i / %i for distance %f (%f)\n", energy, force, (table_prefix + 2 * distbin), (table_prefix + 2 * distbin + 1), mydist, mydistfrac);
06776       table_ener[table_prefix + 2 * distbin] = energy;
06777       table_ener[table_prefix + 2 * distbin + 1] = force;
06778       distbin++;
06779     }
06780     if (currdist >= maxdist) break;
06781   }
06782 
06783   //The procedure above leaves out the last entry -- add it explicitly
06784   distbin = (int) mynearbyint(maxdist / table_spacing);
06785   printf("Adding energy/force entry %f / %f in bins %i / %i\n", enervalues[numentries - 1], 0.0, (table_prefix + 2 * distbin), (table_prefix + 2 * distbin + 1));
06786   table_ener[table_prefix + 2 * distbin] = enervalues[numentries - 1];
06787   table_ener[table_prefix + 2 * distbin + 1] = 0.0;
06788   distbin++;
06789 
06790 
06791   // Clean up and make sure everything worked ok
06792   delete m;
06793   delete x;
06794   delete b;
06795   distbin--;
06796   printf("Testing: %i vs %i (from %f / %f)\n", distbin, (int) (mynearbyint(maxdist / table_spacing)), maxdist, table_spacing);
06797   if (distbin != (int) (mynearbyint(maxdist / table_spacing))) return 1;
06798   return 0;
06799 } /* end read_energy_type_cubspline */

void Parameters::read_parameter_file char *   ) 
 

Definition at line 374 of file Parameters.C.

References Fclose(), Fopen(), NAMD_blank_string(), NAMD_die(), NAMD_find_first_word(), NAMD_read_line(), NumAngleParams, NumBondParams, NumDihedralParams, NumImproperParams, NumTablePairParams, NumVdwPairParams, NumVdwParams, and table_ener.

Referenced by Parameters().

00376 {
00377   char buffer[512];  //  Buffer to store each line of the file
00378   char first_word[512];  //  First word of the current line
00379   FILE *pfile;    //  File descriptor for the parameter file
00380 
00381   /*  Check to make sure that we haven't previously been told     */
00382   /*  that all the files were read        */
00383   if (AllFilesRead)
00384   {
00385     NAMD_die("Tried to read another parameter file after being told that all files were read!");
00386   }
00387 
00388   /*  Try and open the file          */
00389   if ( (pfile = Fopen(fname, "r")) == NULL)
00390   {
00391     char err_msg[256];
00392 
00393     sprintf(err_msg, "UNABLE TO OPEN XPLOR PARAMETER FILE %s\n", fname);
00394     NAMD_die(err_msg);
00395   }
00396 
00397   /*  Keep reading in lines until we hit the EOF      */
00398   while (NAMD_read_line(pfile, buffer) != -1)
00399   {
00400     /*  Get the first word of the line      */
00401     NAMD_find_first_word(buffer, first_word);
00402 
00403     /*  First, screen out things that we ignore, such as    */
00404     /*  blank lines, lines that start with '!', lines that  */
00405     /*  start with "REMARK", lines that start with set",    */
00406     /*  and most of the HBOND parameters which include      */
00407     /*  AEXP, REXP, HAEX, AAEX, but not the HBOND statement */
00408     /*  which is parsed.                                    */
00409     if ((buffer[0] != '!') && 
00410         !NAMD_blank_string(buffer) &&
00411         (strncasecmp(first_word, "REMARK", 6) != 0) &&
00412         (strcasecmp(first_word, "set")!=0) &&
00413         (strncasecmp(first_word, "AEXP", 4) != 0) &&
00414         (strncasecmp(first_word, "REXP", 4) != 0) &&
00415         (strncasecmp(first_word, "HAEX", 4) != 0) &&
00416         (strncasecmp(first_word, "AAEX", 4) != 0) &&
00417         (strncasecmp(first_word, "NBOND", 5) != 0) &&
00418         (strncasecmp(first_word, "CUTNB", 5) != 0) &&
00419         (strncasecmp(first_word, "END", 3) != 0) &&
00420         (strncasecmp(first_word, "CTONN", 5) != 0) &&
00421         (strncasecmp(first_word, "EPS", 3) != 0) &&
00422         (strncasecmp(first_word, "VSWI", 4) != 0) &&
00423         (strncasecmp(first_word, "NBXM", 4) != 0) &&
00424         (strncasecmp(first_word, "INHI", 4) != 0) )
00425     {
00426       /*  Now, call the appropriate function based    */
00427       /*  on the type of parameter we have    */
00428       if (strncasecmp(first_word, "bond", 4)==0)
00429       {
00430         add_bond_param(buffer);
00431         NumBondParams++;
00432       }
00433       else if (strncasecmp(first_word, "angl", 4)==0)
00434       {
00435         add_angle_param(buffer);
00436         NumAngleParams++;
00437       }
00438       else if (strncasecmp(first_word, "dihe", 4)==0)
00439       {
00440         add_dihedral_param(buffer, pfile);
00441         NumDihedralParams++;
00442       }
00443       else if (strncasecmp(first_word, "impr", 4)==0)
00444       {
00445         add_improper_param(buffer, pfile);
00446         NumImproperParams++;
00447       }
00448       else if (strncasecmp(first_word, "nonb", 4)==0)
00449       {
00450         add_vdw_param(buffer);
00451         NumVdwParams++; 
00452       }
00453       else if (strncasecmp(first_word, "nbfi", 4)==0)
00454       {
00455         add_vdw_pair_param(buffer);
00456         NumVdwPairParams++; 
00457       }
00458       else if (strncasecmp(first_word, "nbta", 4)==0)
00459       {
00460 
00461         if (table_ener == NULL) {
00462           continue;
00463         }
00464 
00465         add_table_pair_param(buffer);
00466         NumTablePairParams++; 
00467       }
00468       else if (strncasecmp(first_word, "hbon", 4)==0)
00469       {
00470         add_hb_pair_param(buffer);
00471       }
00472       else
00473       {
00474         /*  This is an unknown paramter.        */
00475         /*  This is BAD        */
00476         char err_msg[512];
00477 
00478         sprintf(err_msg, "UNKNOWN PARAMETER IN XPLOR PARAMETER FILE %s\nLINE=*%s*",
00479            fname, buffer);
00480         NAMD_die(err_msg);
00481       }
00482     }
00483   }
00484 
00485   /*  Close the file            */
00486   Fclose(pfile);
00487 
00488   return;
00489 }

void Parameters::read_parm Ambertoppar ,
BigReal 
 

Definition at line 5771 of file Parameters.C.

References Ambertoppar, angle_array, bond_array, parm::data_read, dihedral_array, IndexedVdwPair, four_body_consts::k, AngleValue::k, BondValue::k, AngleValue::k_ub, MAX_MULTIPLICITY, DihedralValue::multiplicity, four_body_consts::n, NAMD_die(), AngleValue::normal, parm::Nptra, parm::Numang, NumAngleParams, parm::Numbnd, NumBondParams, NumDihedralParams, parm::Phase, parm::Pk, parm::Pn, AngleValue::r_ub, parm::Req, parm::Rk, parm::Teq, AngleValue::theta0, parm::Tk, DihedralValue::values, and BondValue::x0.

Referenced by Parameters().

05772 { 
05773   int i,j,ico,numtype,mul;
05774   IndexedVdwPair *new_node;
05775 
05776   if (!amber_data->data_read)
05777     NAMD_die("No data read from parm file yet!");
05778 
05779   // Copy bond parameters
05780   NumBondParams = amber_data->Numbnd;
05781   if (NumBondParams)
05782   { bond_array = new BondValue[NumBondParams];
05783     if (bond_array == NULL)
05784       NAMD_die("memory allocation of bond_array failed!");
05785   }
05786   for (i=0; i<NumBondParams; ++i)
05787   { bond_array[i].k = amber_data->Rk[i];
05788     bond_array[i].x0 = amber_data->Req[i];
05789   }
05790 
05791   // Copy angle parameters
05792   NumAngleParams = amber_data->Numang;
05793   if (NumAngleParams)
05794   { angle_array = new AngleValue[NumAngleParams];
05795     if (angle_array == NULL)
05796       NAMD_die("memory allocation of angle_array failed!");
05797   }
05798   for (i=0; i<NumAngleParams; ++i)
05799   { angle_array[i].k = amber_data->Tk[i];
05800     angle_array[i].theta0 = amber_data->Teq[i];
05801     // Amber has no k_ub and r_ub for angle parameters, so they're set to 0
05802     angle_array[i].k_ub = angle_array[i].r_ub = 0;
05803     // All angles are harmonic
05804     angle_array[i].normal = 1;
05805   }
05806 
05807   // Copy dihedral parameters
05808   // Note: If the periodicity is negative, it means the following
05809   //  entry is another term in a multitermed dihedral, and the
05810   //  absolute value is the true periodicity; in this case the
05811   //  following entry in "dihedral_array" should be left empty,
05812   //  NOT be skipped, in order to keep the next dihedral's index
05813   //  unchanged.
05814   NumDihedralParams = amber_data->Nptra;
05815   if (NumDihedralParams)
05816   { dihedral_array = new DihedralValue[amber_data->Nptra];
05817     if (dihedral_array == NULL)
05818       NAMD_die("memory allocation of dihedral_array failed!");
05819   }
05820   mul = 0;
05821   for (i=0; i<NumDihedralParams; ++i)
05822   { dihedral_array[i-mul].values[mul].k = amber_data->Pk[i];
05823     dihedral_array[i-mul].values[mul].n = int(fabs(amber_data->Pn[i])+0.5);
05824     if (dihedral_array[i-mul].values[mul].n == 0)
05825     { char err_msg[128];
05826       sprintf(err_msg, "The periodicity of dihedral # %d is zero!", i+1);
05827       NAMD_die(err_msg);
05828     }
05829     dihedral_array[i-mul].values[mul].delta = amber_data->Phase[i];
05830     // If the periodicity is positive, it means the following
05831     // entry is a new dihedral term.
05832     if (amber_data->Pn[i] > 0)
05833     { dihedral_array[i-mul].multiplicity = mul+1;
05834       mul = 0;
05835     }
05836     else if (++mul >= MAX_MULTIPLICITY)
05837     { char err_msg[181];
05838       sprintf(err_msg, "Multiple dihedral with multiplicity of %d greater than max of %d",
05839          mul+1, MAX_MULTIPLICITY);
05840       NAMD_die(err_msg);
05841     }
05842   }
05843   if (mul > 0)
05844     NAMD_die("Negative periodicity in the last dihedral entry!");
05845 
05846   // Copy VDW parameters: AMBER explicitly gives VDW parameters between every
05847   // 2 atom types, so the data are copied to vdw_pair_tree
05848   // In AMBER, all 1-4 VDW interactions are divided by factor vdw14
05849   NumVdwParamsAssigned = numtype = amber_data->Ntypes; // Number of atom types
05850   NumVdwPairParams = numtype * (numtype+1) / 2;
05851   for (i=0; i<numtype; ++i)
05852     for (j=i; j<numtype; ++j)
05853     { new_node = new IndexedVdwPair;
05854       if (new_node == NULL)
05855         NAMD_die("memory allocation of vdw_pair_tree failed!");
05856       new_node->ind1 = i;
05857       new_node->ind2 = j;
05858       new_node->left = new_node->right = NULL;
05859       // ico is the index of interaction between atom type i and j into
05860       // the parameter arrays. If it's positive, the interaction is
05861       // 6-12 VDW; otherwise it's 10-12 H-bond interaction. NAMD doesn't
05862       // have 10-12 term, so if ico is negative, then the 10-12
05863       // coefficients must be 0, otherwise die.
05864       ico = amber_data->Cno[numtype*i+j];
05865       if (ico>0)
05866       { new_node->A = amber_data->Cn1[ico-1];
05867         new_node->A14 = new_node->A / vdw14;
05868         new_node->B = amber_data->Cn2[ico-1];
05869         new_node->B14 = new_node->B / vdw14;
05870       }
05871       else if (amber_data->HB12[abs(ico)-1]==0.0 && amber_data->HB6[abs(ico)-1]==0.0)
05872       { new_node->A = new_node->A14 = new_node->B = new_node->B14 = 0.0;
05873         iout << iWARN << "Encounter 10-12 H-bond term\n";
05874       }
05875       else
05876         NAMD_die("Encounter non-zero 10-12 H-bond term!");
05877       // Add the node to the binary tree
05878       vdw_pair_tree = add_to_indexed_vdw_pairs(new_node, vdw_pair_tree);
05879     }
05880 }

void Parameters::receive_Parameters MIStream  ) 
 

Definition at line 5053 of file Parameters.C.

References indexed_vdw_pair::A, indexed_vdw_pair::A14, angle_array, indexed_vdw_pair::B, indexed_vdw_pair::B14, BigReal, bond_array, crossterm_array, four_body_consts::delta, dihedral_array, vdw_val::epsilon, vdw_val::epsilon14, MIStream::get(), improper_array, indexed_table_pair::ind1, indexed_vdw_pair::ind1, indexed_table_pair::ind2, indexed_vdw_pair::ind2, IndexedTablePair, IndexedVdwPair, j, indexed_table_pair::K, four_body_consts::k, AngleValue::k, BondValue::k, AngleValue::k_ub, indexed_table_pair::left, indexed_vdw_pair::left, MAX_ATOMTYPE_CHARS, ImproperValue::multiplicity, DihedralValue::multiplicity, four_body_consts::n, NAMD_die(), AngleValue::normal, NumAngleParams, NumBondParams, NumCrosstermParams, NumDihedralParams, numenerentries, NumImproperParams, NumTablePairParams, NumVdwPairParams, NumVdwParams, NumVdwParamsAssigned, AngleValue::r_ub, Real, indexed_table_pair::right, indexed_vdw_pair::right, vdw_val::sigma, vdw_val::sigma14, tab_pair_tree, table_ener, AngleValue::theta0, ImproperValue::values, DihedralValue::values, vdw_array, vdw_pair_tree, VdwValue, and BondValue::x0.

05055 {
05056   int i, j;      //  Loop counters
05057   Real *a1, *a2, *a3, *a4;  //  Temporary arrays to get data from message in
05058   int *i1, *i2, *i3;      //  Temporary int array to get data from message in
05059   IndexedVdwPair *new_node;  //  New node for vdw pair param tree
05060   IndexedTablePair *new_tab_node;
05061   Real **kvals;      //  Force constant values for dihedrals and impropers
05062   int **nvals;      //  Periodicity values for dihedrals and impropers
05063   Real **deltavals;    //  Phase shift values for dihedrals and impropers
05064 
05065   //  Get the bonded parameters
05066   msg->get(NumBondParams);
05067 
05068   if (NumBondParams)
05069   {
05070     bond_array = new BondValue[NumBondParams];
05071     a1 = new Real[NumBondParams];
05072     a2 = new Real[NumBondParams];
05073 
05074     if ( (bond_array == NULL) || (a1 == NULL) || (a2 == NULL) )
05075     {
05076       NAMD_die("memory allocation failed in Parameters::receive_Parameters");
05077     }
05078 
05079     msg->get(NumBondParams, a1);
05080     msg->get(NumBondParams, a2);
05081 
05082     for (i=0; i<NumBondParams; i++)
05083     {
05084       bond_array[i].k = a1[i];
05085       bond_array[i].x0 = a2[i];
05086     }
05087 
05088     delete [] a1;
05089     delete [] a2;
05090   }
05091 
05092   //  Get the angle parameters
05093   msg->get(NumAngleParams);
05094 
05095   if (NumAngleParams)
05096   {
05097     angle_array = new AngleValue[NumAngleParams];
05098     a1 = new Real[NumAngleParams];
05099     a2 = new Real[NumAngleParams];
05100     a3 = new Real[NumAngleParams];
05101     a4 = new Real[NumAngleParams];
05102     i1 = new int[NumAngleParams];
05103 
05104     if ( (angle_array == NULL) || (a1 == NULL) || (a2 == NULL) ||
05105          (a3 == NULL) || (a4 == NULL) || (i1 == NULL))
05106     {
05107       NAMD_die("memory allocation failed in Parameters::receive_Parameters");
05108     }
05109 
05110     msg->get(NumAngleParams, a1);
05111     msg->get(NumAngleParams, a2);
05112     msg->get(NumAngleParams, a3);
05113     msg->get(NumAngleParams, a4);
05114     msg->get(NumAngleParams, i1);
05115 
05116     for (i=0; i<NumAngleParams; i++)
05117     {
05118       angle_array[i].k = a1[i];
05119       angle_array[i].theta0 = a2[i];
05120       angle_array[i].k_ub = a3[i];
05121       angle_array[i].r_ub = a4[i];
05122       angle_array[i].normal = i1[i];
05123     }
05124 
05125     delete [] a1;
05126     delete [] a2;
05127     delete [] a3;
05128     delete [] a4;
05129     delete [] i1;
05130   }
05131 
05132   //  Get the dihedral parameters
05133   msg->get(NumDihedralParams);
05134 
05135   if (NumDihedralParams)
05136   {
05137     dihedral_array = new DihedralValue[NumDihedralParams];
05138 
05139     i1 = new int[NumDihedralParams];
05140     kvals = new Real *[MAX_MULTIPLICITY];
05141     nvals = new int *[MAX_MULTIPLICITY];
05142     deltavals = new Real *[MAX_MULTIPLICITY];
05143 
05144     if ( (i1 == NULL) || (kvals == NULL) || (nvals == NULL) || 
05145          (deltavals == NULL) || (dihedral_array == NULL) )
05146     {
05147       NAMD_die("memory allocation failed in Parameters::send_Parameters");
05148     }
05149 
05150     for (i=0; i<MAX_MULTIPLICITY; i++)
05151     {
05152       kvals[i] = new Real[NumDihedralParams];
05153       nvals[i] = new int[NumDihedralParams];
05154       deltavals[i] = new Real[NumDihedralParams];
05155 
05156       if ( (kvals[i] == NULL) || (nvals[i] == NULL) || (deltavals[i] == NULL) )
05157       {
05158         NAMD_die("memory allocation failed in Parameters::send_Parameters");
05159       }
05160     }
05161 
05162     msg->get(NumDihedralParams, i1);
05163 
05164     for (i=0; i<MAX_MULTIPLICITY; i++)
05165     {
05166       msg->get(NumDihedralParams, kvals[i]);
05167       msg->get(NumDihedralParams, nvals[i]);
05168       msg->get(NumDihedralParams, deltavals[i]);
05169     }
05170 
05171     for (i=0; i<NumDihedralParams; i++)
05172     {
05173       dihedral_array[i].multiplicity = i1[i];
05174 
05175       for (j=0; j<MAX_MULTIPLICITY; j++)
05176       {
05177         dihedral_array[i].values[j].k = kvals[j][i];
05178         dihedral_array[i].values[j].n = nvals[j][i];
05179         dihedral_array[i].values[j].delta = deltavals[j][i];
05180       }
05181     }
05182 
05183     for (i=0; i<MAX_MULTIPLICITY; i++)
05184     {
05185       delete [] kvals[i];
05186       delete [] nvals[i];
05187       delete [] deltavals[i];
05188     }
05189 
05190     delete [] i1;
05191     delete [] kvals;
05192     delete [] nvals;
05193     delete [] deltavals;
05194   }
05195 
05196   //  Get the improper parameters
05197   msg->get(NumImproperParams);
05198 
05199   if (NumImproperParams)
05200   {
05201     improper_array = new ImproperValue[NumImproperParams];
05202     i1 = new int[NumImproperParams];
05203     kvals = new Real *[MAX_MULTIPLICITY];
05204     nvals = new int *[MAX_MULTIPLICITY];
05205     deltavals = new Real *[MAX_MULTIPLICITY];
05206 
05207     if ( (i1 == NULL) || (kvals == NULL) || (nvals == NULL) || 
05208          (deltavals == NULL) || (improper_array==NULL) )
05209     {
05210       NAMD_die("memory allocation failed in Parameters::send_Parameters");
05211     }
05212 
05213     for (i=0; i<MAX_MULTIPLICITY; i++)
05214     {
05215       kvals[i] = new Real[NumImproperParams];
05216       nvals[i] = new int[NumImproperParams];
05217       deltavals[i] = new Real[NumImproperParams];
05218 
05219       if ( (kvals[i] == NULL) || (nvals[i] == NULL) || (deltavals[i] == NULL) )
05220       {
05221         NAMD_die("memory allocation failed in Parameters::send_Parameters");
05222       }
05223     }
05224 
05225     msg->get(NumImproperParams,i1);
05226 
05227     for (i=0; i<MAX_MULTIPLICITY; i++)
05228     {
05229       msg->get(NumImproperParams,kvals[i]);
05230       msg->get(NumImproperParams,nvals[i]);
05231       msg->get(NumImproperParams,deltavals[i]);
05232     }
05233 
05234     for (i=0; i<NumImproperParams; i++)
05235     {
05236       improper_array[i].multiplicity = i1[i];
05237 
05238       for (j=0; j<MAX_MULTIPLICITY; j++)
05239       {
05240         improper_array[i].values[j].k = kvals[j][i];
05241         improper_array[i].values[j].n = nvals[j][i];
05242         improper_array[i].values[j].delta = deltavals[j][i];
05243       }
05244     }
05245 
05246     for (i=0; i<MAX_MULTIPLICITY; i++)
05247     {
05248       delete [] kvals[i];
05249       delete [] nvals[i];
05250       delete [] deltavals[i];
05251     }
05252 
05253     delete [] i1;
05254     delete [] kvals;
05255     delete [] nvals;
05256     delete [] deltavals;
05257   }
05258 
05259   //  Get the crossterm parameters
05260   msg->get(NumCrosstermParams);
05261 
05262   if (NumCrosstermParams)
05263   {
05264     crossterm_array = new CrosstermValue[NumCrosstermParams];
05265 
05266     for (i=0; i<NumCrosstermParams; ++i) {
05267       int nvals = CrosstermValue::dim * CrosstermValue::dim * 2 * 2;
05268       msg->get(nvals,&crossterm_array[i].c[0][0].d00);
05269     }
05270   }
05271   
05272   //Get the energy table
05273   msg->get(numenerentries);
05274   if (numenerentries > 0) {
05275     //fprintf(stdout, "Getting tables\n");
05276     //fprintf(infofile, "Trying to allocate table\n");
05277     table_ener = new BigReal[numenerentries];
05278     //fprintf(infofile, "Finished table allocation\n");
05279     if (table_ener==NULL)
05280     {
05281       NAMD_die("memory allocation failed in Parameters::receive_Parameters");
05282     }
05283 
05284     msg->get(numenerentries, table_ener);
05285   }
05286 
05287   //  Get the vdw parameters
05288   msg->get(NumVdwParams);
05289   msg->get(NumVdwParamsAssigned);
05290 
05291   if (NumVdwParams)
05292   {
05293           atomTypeNames = new char[NumVdwParams*(MAX_ATOMTYPE_CHARS+1)];
05294     vdw_array = new VdwValue[NumVdwParams];
05295     a1 = new Real[NumVdwParams];
05296     a2 = new Real[NumVdwParams];
05297     a3 = new Real[NumVdwParams];
05298     a4 = new Real[NumVdwParams];
05299 
05300     if ( (vdw_array==NULL) || (a1==NULL) || (a2==NULL) || (a3==NULL)
05301              || (a4==NULL) || (atomTypeNames==NULL))
05302     {
05303       NAMD_die("memory allocation failed in Parameters::receive_Parameters");
05304     }
05305 
05306     msg->get(NumVdwParams * (MAX_ATOMTYPE_CHARS+1), atomTypeNames);
05307     msg->get(NumVdwParams, a1);
05308     msg->get(NumVdwParams, a2);
05309     msg->get(NumVdwParams, a3);
05310     msg->get(NumVdwParams, a4);
05311 
05312     for (i=0; i<NumVdwParams; i++)
05313     {
05314       vdw_array[i].sigma = a1[i];
05315       vdw_array[i].epsilon = a2[i];
05316       vdw_array[i].sigma14 = a3[i];
05317       vdw_array[i].epsilon14 = a4[i];
05318     }
05319 
05320     delete [] a1;
05321     delete [] a2;
05322     delete [] a3;
05323     delete [] a4;
05324   }
05325   
05326   //  Get the vdw_pair_parameters
05327   msg->get(NumVdwPairParams);
05328   
05329   if (NumVdwPairParams)
05330   {
05331     a1 = new Real[NumVdwPairParams];
05332     a2 = new Real[NumVdwPairParams];
05333     a3 = new Real[NumVdwPairParams];
05334     a4 = new Real[NumVdwPairParams];
05335     i1 = new int[NumVdwPairParams];
05336     i2 = new int[NumVdwPairParams];
05337 
05338     if ( (a1 == NULL) || (a2 == NULL) || (a3 == NULL) || (a4 == NULL) || 
05339          (i1 == NULL) || (i2 == NULL) )
05340     {
05341       NAMD_die("memory allocation failed in Parameters::send_Parameters");
05342     }
05343     
05344     msg->get(NumVdwPairParams, i1);
05345     msg->get(NumVdwPairParams, i2);
05346     msg->get(NumVdwPairParams, a1);
05347     msg->get(NumVdwPairParams, a2);
05348     msg->get(NumVdwPairParams, a3);
05349     msg->get(NumVdwPairParams, a4);
05350     
05351     for (i=0; i<NumVdwPairParams; i++)
05352     {
05353       new_node = (IndexedVdwPair *) malloc(sizeof(IndexedVdwPair));
05354       
05355       if (new_node == NULL)
05356       {
05357          NAMD_die("memory allocation failed in Parameters::receive_Parameters");
05358       }
05359       
05360       new_node->ind1 = i1[i];
05361       new_node->ind2 = i2[i];
05362       new_node->A = a1[i];
05363       new_node->A14 = a2[i];
05364       new_node->B = a3[i];
05365       new_node->B14 = a4[i];
05366       new_node->left = NULL;
05367       new_node->right = NULL;
05368       
05369       vdw_pair_tree = add_to_indexed_vdw_pairs(new_node, vdw_pair_tree);
05370     }
05371     
05372     delete [] i1;
05373     delete [] i2;
05374     delete [] a1;
05375     delete [] a2;
05376     delete [] a3;
05377     delete [] a4;
05378   }
05379  
05380   //  Get the table_pair_parameters
05381   msg->get(NumTablePairParams);
05382   
05383   if (NumTablePairParams)
05384   {
05385     i1 = new int[NumTablePairParams];
05386     i2 = new int[NumTablePairParams];
05387     i3 = new int[NumTablePairParams];
05388 
05389     if ( (i3 == NULL) || (i1 == NULL) || (i2 == NULL) )
05390     {
05391       NAMD_die("memory allocation failed in Parameters::send_Parameters");
05392     }
05393     
05394     msg->get(NumTablePairParams, i1);
05395     msg->get(NumTablePairParams, i2);
05396     msg->get(NumTablePairParams, i3);
05397     
05398     for (i=0; i<NumTablePairParams; i++)
05399     {
05400       new_tab_node = (IndexedTablePair *) malloc(sizeof(IndexedTablePair));
05401       
05402       if (new_tab_node == NULL)
05403       {
05404          NAMD_die("memory allocation failed in Parameters::receive_Parameters");
05405       }
05406       
05407 //      printf("Adding new table pair with ind1 %i ind2 %i k %i\n", i1[i], i2[i],i3[i]);
05408       new_tab_node->ind1 = i1[i];
05409       new_tab_node->ind2 = i2[i];
05410       new_tab_node->K = i3[i];
05411       new_tab_node->left = NULL;
05412       new_tab_node->right = NULL;
05413       
05414       tab_pair_tree = add_to_indexed_table_pairs(new_tab_node, tab_pair_tree);
05415     }
05416     
05417     delete [] i1;
05418     delete [] i2;
05419     delete [] i3;
05420   }
05421   
05422   // receive the hydrogen bond parameters
05423   // hbondParams.receive_message(msg);
05424 
05425   AllFilesRead = TRUE;
05426 
05427   delete msg;
05428 }

void Parameters::send_Parameters MOStream  ) 
 

Definition at line 4735 of file Parameters.C.

References angle_array, bond_array, crossterm_array, four_body_consts::delta, dihedral_array, MOStream::end(), vdw_val::epsilon, vdw_val::epsilon14, improper_array, j, four_body_consts::k, AngleValue::k, BondValue::k, AngleValue::k_ub, MAX_ATOMTYPE_CHARS, ImproperValue::multiplicity, DihedralValue::multiplicity, four_body_consts::n, NAMD_die(), AngleValue::normal, NumAngleParams, NumBondParams, NumCrosstermParams, NumDihedralParams, numenerentries, NumImproperParams, NumTablePairParams, NumVdwPairParams, NumVdwParams, NumVdwParamsAssigned, MOStream::put(), AngleValue::r_ub, Real, vdw_val::sigma, vdw_val::sigma14, tab_pair_tree, table_ener, AngleValue::theta0, ImproperValue::values, DihedralValue::values, vdw_array, vdw_pair_tree, and BondValue::x0.

04736 {
04737   Real *a1, *a2, *a3, *a4;        //  Temporary arrays for sending messages
04738   int *i1, *i2, *i3;      //  Temporary int array
04739   int i, j;      //  Loop counters
04740   Real **kvals;      //  Force constant values for dihedrals and impropers
04741   int **nvals;      //  Periodicity values for  dihedrals and impropers
04742   Real **deltavals;    //  Phase shift values for  dihedrals and impropers
04743   /*MOStream *msg=comm_obj->newOutputStream(ALLBUTME, STATICPARAMSTAG, BUFSIZE);
04744   if ( msg == NULL )
04745   {
04746     NAMD_die("memory allocation failed in Parameters::send_Parameters");
04747   }*/
04748 
04749   //  Send the bond parameters
04750   msg->put(NumBondParams);
04751 
04752   if (NumBondParams)
04753   {
04754     a1 = new Real[NumBondParams];
04755     a2 = new Real[NumBondParams];
04756 
04757     if ( (a1 == NULL) || (a2 == NULL) )
04758     {
04759       NAMD_die("memory allocation failed in Parameters::send_Parameters");
04760     }
04761 
04762     for (i=0; i<NumBondParams; i++)
04763     {
04764       a1[i] = bond_array[i].k;
04765       a2[i] = bond_array[i].x0;
04766     }
04767 
04768     msg->put(NumBondParams, a1)->put(NumBondParams, a2);
04769 
04770     delete [] a1;
04771     delete [] a2;
04772   }
04773 
04774   //  Send the angle parameters
04775   msg->put(NumAngleParams);
04776 
04777   if (NumAngleParams)
04778   {
04779     a1 = new Real[NumAngleParams];
04780     a2 = new Real[NumAngleParams];
04781     a3 = new Real[NumAngleParams];
04782     a4 = new Real[NumAngleParams];
04783     i1 = new int[NumAngleParams];
04784 
04785     if ( (a1 == NULL) || (a2 == NULL) || (a3 == NULL) ||
04786          (a4 == NULL) || (i1 == NULL))
04787     {
04788       NAMD_die("memory allocation failed in Parameters::send_Parameters");
04789     }
04790 
04791     for (i=0; i<NumAngleParams; i++)
04792     {
04793       a1[i] = angle_array[i].k;
04794       a2[i] = angle_array[i].theta0;
04795       a3[i] = angle_array[i].k_ub;
04796       a4[i] = angle_array[i].r_ub;
04797       i1[i] = angle_array[i].normal;
04798     }
04799 
04800     msg->put(NumAngleParams, a1)->put(NumAngleParams, a2);
04801     msg->put(NumAngleParams, a3)->put(NumAngleParams, a4);
04802     msg->put(NumAngleParams, i1);
04803 
04804     delete [] a1;
04805     delete [] a2;
04806     delete [] a3;
04807     delete [] a4;
04808     delete [] i1;
04809   }
04810 
04811   //  Send the dihedral parameters
04812   msg->put(NumDihedralParams);
04813 
04814   if (NumDihedralParams)
04815   {
04816     i1 = new int[NumDihedralParams];
04817     kvals = new Real *[MAX_MULTIPLICITY];
04818     nvals = new int *[MAX_MULTIPLICITY];
04819     deltavals = new Real *[MAX_MULTIPLICITY];
04820 
04821     if ( (i1 == NULL) || (kvals == NULL) || (nvals == NULL) || 
04822          (deltavals == NULL) )
04823     {
04824       NAMD_die("memory allocation failed in Parameters::send_Parameters");
04825     }
04826 
04827     for (i=0; i<MAX_MULTIPLICITY; i++)
04828     {
04829       kvals[i] = new Real[NumDihedralParams];
04830       nvals[i] = new int[NumDihedralParams];
04831       deltavals[i] = new Real[NumDihedralParams];
04832 
04833       if ( (kvals[i] == NULL) || (nvals[i] == NULL) || (deltavals[i] == NULL) )
04834       {
04835         NAMD_die("memory allocation failed in Parameters::send_Parameters");
04836       }
04837     }
04838 
04839     for (i=0; i<NumDihedralParams; i++)
04840     {
04841       i1[i] = dihedral_array[i].multiplicity;
04842 
04843       for (j=0; j<MAX_MULTIPLICITY; j++)
04844       {
04845         kvals[j][i] = dihedral_array[i].values[j].k;
04846         nvals[j][i] = dihedral_array[i].values[j].n;
04847         deltavals[j][i] = dihedral_array[i].values[j].delta;
04848       }
04849     }
04850 
04851     msg->put(NumDihedralParams, i1);
04852 
04853     for (i=0; i<MAX_MULTIPLICITY; i++)
04854     {
04855       msg->put(NumDihedralParams, kvals[i]);
04856       msg->put(NumDihedralParams, nvals[i]);
04857       msg->put(NumDihedralParams, deltavals[i]);
04858 
04859       delete [] kvals[i];
04860       delete [] nvals[i];
04861       delete [] deltavals[i];
04862     }
04863 
04864     delete [] i1;
04865     delete [] kvals;
04866     delete [] nvals;
04867     delete [] deltavals;
04868   }
04869 
04870   //  Send the improper parameters
04871   msg->put(NumImproperParams);
04872 
04873   if (NumImproperParams)
04874   {
04875     i1 = new int[NumImproperParams];
04876     kvals = new Real *[MAX_MULTIPLICITY];
04877     nvals = new int *[MAX_MULTIPLICITY];
04878     deltavals = new Real *[MAX_MULTIPLICITY];
04879 
04880     if ( (i1 == NULL) || (kvals == NULL) || (nvals == NULL) || 
04881          (deltavals == NULL) )
04882     {
04883       NAMD_die("memory allocation failed in Parameters::send_Parameters");
04884     }
04885 
04886     for (i=0; i<MAX_MULTIPLICITY; i++)
04887     {
04888       kvals[i] = new Real[NumImproperParams];
04889       nvals[i] = new int[NumImproperParams];
04890       deltavals[i] = new Real[NumImproperParams];
04891 
04892       if ( (kvals[i] == NULL) || (nvals[i] == NULL) || (deltavals[i] == NULL) )
04893       {
04894         NAMD_die("memory allocation failed in Parameters::send_Parameters");
04895       }
04896     }
04897 
04898     for (i=0; i<NumImproperParams; i++)
04899     {
04900       i1[i] = improper_array[i].multiplicity;
04901 
04902       for (j=0; j<MAX_MULTIPLICITY; j++)
04903       {
04904         kvals[j][i] = improper_array[i].values[j].k;
04905         nvals[j][i] = improper_array[i].values[j].n;
04906         deltavals[j][i] = improper_array[i].values[j].delta;
04907       }
04908     }
04909 
04910     msg->put(NumImproperParams, i1);
04911 
04912     for (i=0; i<MAX_MULTIPLICITY; i++)
04913     {
04914       msg->put(NumImproperParams, kvals[i]);
04915       msg->put(NumImproperParams, nvals[i]);
04916       msg->put(NumImproperParams, deltavals[i]);
04917 
04918       delete [] kvals[i];
04919       delete [] nvals[i];
04920       delete [] deltavals[i];
04921     }
04922 
04923     delete [] i1;
04924     delete [] kvals;
04925     delete [] nvals;
04926     delete [] deltavals;
04927   }
04928 
04929   //  Send the crossterm parameters
04930   msg->put(NumCrosstermParams);
04931 
04932   if (NumCrosstermParams)
04933   {
04934     for (i=0; i<NumCrosstermParams; ++i) {
04935       int nvals = CrosstermValue::dim * CrosstermValue::dim * 2 * 2;
04936       msg->put(nvals,&crossterm_array[i].c[0][0].d00);
04937     }
04938   }
04939   //
04940   //Send the energy table parameters
04941   msg->put(numenerentries);
04942 
04943   if (numenerentries) {
04944           /*
04945     b1 = new Real[numenerentries];
04946     if (b1 == NULL) 
04947     {
04948       NAMD_die("memory allocation failed in Parameters::send_Parameters");
04949     }
04950     */
04951     
04952     msg->put(numenerentries, table_ener);
04953   }
04954 
04955   //  Send the vdw parameters
04956   msg->put(NumVdwParams);
04957   msg->put(NumVdwParamsAssigned);
04958 
04959   if (NumVdwParams)
04960   {
04961     a1 = new Real[NumVdwParams];
04962     a2 = new Real[NumVdwParams];
04963     a3 = new Real[NumVdwParams];
04964     a4 = new Real[NumVdwParams];
04965 
04966     if ( (a1 == NULL) || (a2 == NULL) || (a3 == NULL) )
04967     {
04968       NAMD_die("memory allocation failed in Parameters::send_Parameters");
04969     }
04970 
04971     for (i=0; i<NumVdwParams; i++)
04972     {
04973       a1[i] = vdw_array[i].sigma;
04974       a2[i] = vdw_array[i].epsilon;
04975       a3[i] = vdw_array[i].sigma14;
04976       a4[i] = vdw_array[i].epsilon14;
04977     }
04978 
04979     msg->put(NumVdwParams * (MAX_ATOMTYPE_CHARS+1), atomTypeNames);
04980     msg->put(NumVdwParams, a1);
04981     msg->put(NumVdwParams, a2);
04982     msg->put(NumVdwParams, a3);
04983     msg->put(NumVdwParams, a4);
04984 
04985     delete [] a1;
04986     delete [] a2;
04987     delete [] a3;
04988     delete [] a4;
04989   }
04990   
04991   //  Send the vdw pair parameters
04992   msg->put(NumVdwPairParams);
04993   
04994   if (NumVdwPairParams)
04995   {
04996     a1 = new Real[NumVdwPairParams];
04997     a2 = new Real[NumVdwPairParams];
04998     a3 = new Real[NumVdwPairParams];
04999     a4 = new Real[NumVdwPairParams];
05000     i1 = new int[NumVdwPairParams];
05001     i2 = new int[NumVdwPairParams];    
05002 
05003     if ( (a1 == NULL) || (a2 == NULL) || (a3 == NULL) || (a4 == NULL) || 
05004          (i1 == NULL) || (i2 == NULL) )
05005     {
05006       NAMD_die("memory allocation failed in Parameters::send_Parameters");
05007     }
05008     
05009     vdw_pair_to_arrays(i1, i2, a1, a2, a3, a4, 0, vdw_pair_tree);
05010     
05011     msg->put(NumVdwPairParams, i1)->put(NumVdwPairParams, i2);
05012     msg->put(NumVdwPairParams, a1);
05013     msg->put(NumVdwPairParams, a2)->put(NumVdwPairParams, a3);
05014     msg->put(NumVdwPairParams, a4);
05015   }
05016   
05017   //  Send the table pair parameters
05018   //printf("Pairs: %i\n", NumTablePairParams);
05019   msg->put(NumTablePairParams);
05020   
05021   if (NumTablePairParams)
05022   {
05023     i1 = new int[NumTablePairParams];
05024     i2 = new int[NumTablePairParams];    
05025     i3 = new int[NumTablePairParams];
05026 
05027     if ( (i3 == NULL) || (i1 == NULL) || (i2 == NULL) )
05028     {
05029       NAMD_die("memory allocation failed in Parameters::send_Parameters");
05030     }
05031     
05032     table_pair_to_arrays(i1, i2, i3, 0, tab_pair_tree);
05033     
05034     msg->put(NumTablePairParams, i1)->put(NumTablePairParams, i2);
05035     msg->put(NumTablePairParams, i3);
05036   }
05037 
05038   // send the hydrogen bond parameters
05039   // hbondParams.create_message(msg);
05040   msg->end();
05041   delete msg;
05042 }


Member Data Documentation

AngleValue* Parameters::angle_array
 

Definition at line 204 of file Parameters.h.

Referenced by done_reading_files(), AngleElem::getParameterPointers(), read_parm(), receive_Parameters(), send_Parameters(), and ~Parameters().

BondValue* Parameters::bond_array
 

Definition at line 203 of file Parameters.h.

Referenced by done_reading_files(), BondElem::getParameterPointers(), read_parm(), receive_Parameters(), send_Parameters(), and ~Parameters().

int Parameters::columnsize
 

Definition at line 211 of file Parameters.h.

Referenced by ComputeNonbondedUtil::select().

CrosstermValue* Parameters::crossterm_array
 

Definition at line 207 of file Parameters.h.

Referenced by done_reading_files(), CrosstermElem::getParameterPointers(), receive_Parameters(), send_Parameters(), and ~Parameters().

DihedralValue* Parameters::dihedral_array
 

Definition at line 205 of file Parameters.h.

Referenced by assign_dihedral_index(), done_reading_files(), DihedralElem::getParameterPointers(), outputCompressedFile(), outputPsfFile(), read_parm(), receive_Parameters(), send_Parameters(), and ~Parameters().

ImproperValue* Parameters::improper_array
 

Definition at line 206 of file Parameters.h.

Referenced by assign_improper_index(), done_reading_files(), ImproperElem::getParameterPointers(), outputCompressedFile(), outputPsfFile(), receive_Parameters(), send_Parameters(), and ~Parameters().

int Parameters::NumAngleParams
 

Definition at line 217 of file Parameters.h.

Referenced by getExtraBonds(), print_angle_params(), print_param_summary(), read_charmm_parameter_file(), read_parameter_file(), read_parm(), receive_Parameters(), and send_Parameters().

int Parameters::NumBondParams
 

Definition at line 216 of file Parameters.h.

Referenced by getExtraBonds(), print_bond_params(), print_param_summary(), read_charmm_parameter_file(), read_parameter_file(), read_parm(), receive_Parameters(), and send_Parameters().

int Parameters::NumCosAngles
 

Definition at line 218 of file Parameters.h.

Referenced by print_param_summary().

int Parameters::NumCrosstermParams
 

Definition at line 221 of file Parameters.h.

Referenced by done_reading_files(), print_param_summary(), read_charmm_parameter_file(), receive_Parameters(), and send_Parameters().

int Parameters::NumDihedralParams
 

Definition at line 219 of file Parameters.h.

Referenced by done_reading_files(), getExtraBonds(), outputCompressedFile(), outputPsfFile(), print_dihedral_params(), print_param_summary(), read_charmm_parameter_file(), read_parameter_file(), read_parm(), receive_Parameters(), and send_Parameters().

int Parameters::numenerentries
 

Definition at line 209 of file Parameters.h.

Referenced by Parameters(), read_ener_table(), receive_Parameters(), and send_Parameters().

int Parameters::NumImproperParams
 

Definition at line 220 of file Parameters.h.

Referenced by done_reading_files(), getExtraBonds(), outputCompressedFile(), outputPsfFile(), print_improper_params(), print_param_summary(), read_charmm_parameter_file(), read_parameter_file(), receive_Parameters(), and send_Parameters().

int Parameters::NumTablePairParams
 

Definition at line 226 of file Parameters.h.

Referenced by read_charmm_parameter_file(), read_parameter_file(), receive_Parameters(), and send_Parameters().

int Parameters::NumTableParams
 

Definition at line 223 of file Parameters.h.

int Parameters::NumVdwPairParams
 

Definition at line 225 of file Parameters.h.

Referenced by print_param_summary(), print_vdw_pair_params(), read_charmm_parameter_file(), read_parameter_file(), receive_Parameters(), and send_Parameters().

int Parameters::NumVdwParams
 

Definition at line 222 of file Parameters.h.

Referenced by done_reading_files(), print_param_summary(), print_vdw_params(), read_charmm_parameter_file(), read_parameter_file(), receive_Parameters(), and send_Parameters().

int Parameters::NumVdwParamsAssigned
 

Definition at line 224 of file Parameters.h.

Referenced by done_reading_files(), receive_Parameters(), and send_Parameters().

int Parameters::rowsize
 

Definition at line 210 of file Parameters.h.

Referenced by ComputeNonbondedUtil::select().

IndexedTablePair* Parameters::tab_pair_tree
 

Definition at line 214 of file Parameters.h.

Referenced by receive_Parameters(), send_Parameters(), and ~Parameters().

BigReal* Parameters::table_ener
 

Definition at line 212 of file Parameters.h.

Referenced by Parameters(), read_charmm_parameter_file(), read_ener_table(), read_parameter_file(), receive_Parameters(), ComputeNonbondedUtil::select(), and send_Parameters().

int Parameters::tablenumtypes
 

Definition at line 215 of file Parameters.h.

Referenced by read_ener_table().

VdwValue* Parameters::vdw_array
 

Definition at line 208 of file Parameters.h.

Referenced by done_reading_files(), receive_Parameters(), send_Parameters(), and ~Parameters().

IndexedVdwPair* Parameters::vdw_pair_tree
 

Definition at line 213 of file Parameters.h.

Referenced by receive_Parameters(), send_Parameters(), and ~Parameters().


The documentation for this class was generated from the following files:
Generated on Mon Nov 23 05:00:04 2009 for NAMD by  doxygen 1.3.9.1