NAMD
Public Member Functions | Public Attributes | List of all members
Parameters Class Reference

#include <Parameters.h>

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 (Bool)
 
void done_reading_structure ()
 
void assign_vdw_index (const char *, Atom *)
 
void assign_bond_index (const char *, const char *, Bond *)
 
void assign_angle_index (const char *, const char *, const char *, Angle *, int)
 
void assign_dihedral_index (const char *, const char *, const char *, const char *, Dihedral *, int, int)
 
void assign_improper_index (const char *, const char *, const char *, const char *, Improper *, int)
 
void assign_crossterm_index (const char *, const char *, const char *, const char *, const char *, const char *, const char *, const 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_nbthole_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
 
GromacsPairValuegromacsPair_array
 
VdwValuevdw_array
 
NbtholePairValuenbthole_array
 
int numenerentries
 
int rowsize
 
int columnsize
 
BigRealtable_ener
 
IndexedVdwPairvdw_pair_tree
 
IndexedNbtholePairnbthole_pair_tree
 
IndexedTablePairtab_pair_tree
 
int tablenumtypes
 
int NumBondParams
 
int NumAngleParams
 
int NumCosAngles
 
int NumDihedralParams
 
int NumImproperParams
 
int NumCrosstermParams
 
int NumGromacsPairParams
 
int NumVdwParams
 
int NumTableParams
 
int NumVdwParamsAssigned
 
int NumVdwPairParams
 
int NumNbtholePairParams
 
int NumTablePairParams
 

Detailed Description

Definition at line 216 of file Parameters.h.

Constructor & Destructor Documentation

Parameters::Parameters ( )

Definition at line 186 of file Parameters.C.

186  {
187  initialize();
188 }
Parameters::Parameters ( SimParameters simParams,
StringList f 
)

Definition at line 249 of file Parameters.C.

References SimParameters::cosAngles, StringList::data, done_reading_files(), SimParameters::drudeOn, FALSE, StringList::next, paraCharmm, SimParameters::paraTypeCharmmOn, SimParameters::paraTypeXplorOn, paraXplor, read_charmm_parameter_file(), read_ener_table(), read_parameter_file(), and SimParameters::tabulatedEnergies.

250 {
251  initialize();
252 
254  if (simParams->paraTypeXplorOn)
255  {
256  paramType = paraXplor;
257  }
258  else if (simParams->paraTypeCharmmOn)
259  {
260  paramType = paraCharmm;
261  }
262  //****** END CHARMM/XPLOR type changes
263  //Test for cos-based angles
264  if (simParams->cosAngles) {
265  cosAngles = true;
266  } else {
267  cosAngles = false;
268  }
269 
270  if (simParams->tabulatedEnergies) {
271  CkPrintf("Working on tables\n");
272  read_ener_table(simParams);
273  }
274 
275  //****** BEGIN CHARMM/XPLOR type changes
276  /* Set up AllFilesRead flag to FALSE. Once all of the files */
277  /* have been read in, then this will be set to true and the */
278  /* arrays of parameters will be set up */
279  AllFilesRead = FALSE;
280 
281  if (NULL != f)
282  {
283  do
284  {
285  //****** BEGIN CHARMM/XPLOR type changes
286  if (paramType == paraXplor)
287  {
289  }
290  else if (paramType == paraCharmm)
291  {
293  }
294  //****** END CHARMM/XPLOR type changes
295  f = f->next;
296  } while ( f != NULL );
297 
298  done_reading_files(simParams->drudeOn && paramType == paraCharmm);
299  }
300 
301 }
#define paraCharmm
Definition: Parameters.h:77
void read_charmm_parameter_file(char *)
Definition: Parameters.C:539
#define FALSE
Definition: common.h:116
void read_parameter_file(char *)
Definition: Parameters.C:404
void done_reading_files(Bool)
Definition: Parameters.C:2959
StringList * next
Definition: ConfigList.h:49
char * data
Definition: ConfigList.h:48
#define paraXplor
Definition: Parameters.h:76
void read_ener_table(SimParameters *)
Definition: Parameters.C:6728
Bool tabulatedEnergies
Parameters::Parameters ( Ambertoppar amber_data,
BigReal  vdw14 
)

Definition at line 6344 of file Parameters.C.

6345 {
6346  initialize();
6347 
6348  // Read in parm parameters
6349  read_parm(amber_data,vdw14);
6350 }
Parameters::Parameters ( const GromacsTopFile gf,
Bool  min 
)

Definition at line 6486 of file Parameters.C.

6487 {
6488  initialize();
6489 
6490  // Read in parm parameters
6491  read_parm(gf,min);
6492 }
Parameters::~Parameters ( )

Definition at line 313 of file Parameters.C.

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

315 {
316  if (atomTypeNames)
317  delete [] atomTypeNames;
318 
319  if (bondp != NULL)
320  free_bond_tree(bondp);
321 
322  if (anglep != NULL)
323  free_angle_tree(anglep);
324 
325  if (dihedralp != NULL)
326  free_dihedral_list(dihedralp);
327 
328  if (improperp != NULL)
329  free_improper_list(improperp);
330 
331  if (crosstermp != NULL)
332  free_crossterm_list(crosstermp);
333 
334  if (vdwp != NULL)
335  free_vdw_tree(vdwp);
336 
337  if (vdw_pairp != NULL)
338  free_vdw_pair_list();
339 
340  if (nbthole_pairp != NULL)
341  free_nbthole_pair_list();
342 
343  if (bond_array != NULL)
344  delete [] bond_array;
345 
346  if (angle_array != NULL)
347  delete [] angle_array;
348 
349  if (dihedral_array != NULL)
350  delete [] dihedral_array;
351 
352  if (improper_array != NULL)
353  delete [] improper_array;
354 
355  if (crossterm_array != NULL)
356  delete [] crossterm_array;
357 
358  // JLai
359  if (gromacsPair_array != NULL)
360  delete [] gromacsPair_array;
361  // End of JLai
362 
363  if (vdw_array != NULL)
364  delete [] vdw_array;
365 
366  if (tab_pair_tree != NULL)
367  free_table_pair_tree(tab_pair_tree);
368 
369  if (vdw_pair_tree != NULL)
370  free_vdw_pair_tree(vdw_pair_tree);
371 
372  if (nbthole_pair_tree != NULL)
373  free_nbthole_pair_tree(nbthole_pair_tree);
374 
375  if (maxDihedralMults != NULL)
376  delete [] maxDihedralMults;
377 
378  if (maxImproperMults != NULL)
379  delete [] maxImproperMults;
380 
381  for( int i = 0; i < error_msgs.size(); ++i ) {
382  delete [] error_msgs[i];
383  }
384  error_msgs.resize(0);
385 }
GromacsPairValue * gromacsPair_array
Definition: Parameters.h:247
DihedralValue * dihedral_array
Definition: Parameters.h:243
CrosstermValue * crossterm_array
Definition: Parameters.h:245
IndexedVdwPair * vdw_pair_tree
Definition: Parameters.h:255
IndexedNbtholePair * nbthole_pair_tree
Definition: Parameters.h:256
AngleValue * angle_array
Definition: Parameters.h:242
ImproperValue * improper_array
Definition: Parameters.h:244
VdwValue * vdw_array
Definition: Parameters.h:249
void resize(int i)
Definition: ResizeArray.h:84
BondValue * bond_array
Definition: Parameters.h:241
int size(void) const
Definition: ResizeArray.h:127
IndexedTablePair * tab_pair_tree
Definition: Parameters.h:257

Member Function Documentation

void Parameters::assign_angle_index ( const char *  atom1,
const char *  atom2,
const char *  atom3,
Angle angle_ptr,
int  notFoundIndex 
)

Definition at line 3854 of file Parameters.C.

References angle::angle_type, angle::atom1, angle_params::atom1name, angle::atom2, angle_params::atom2name, angle::atom3, angle_params::atom3name, endi(), angle_params::index, iout, iWARN(), angle_params::left, NAMD_die(), and angle_params::right.

3857 {
3858  struct angle_params *ptr; // Current position in tree
3859  int comp_val; // value from strcasecmp
3860  int found=0; // flag 1->found a match
3861 
3862  /* Check to make sure the files have all been read */
3863  if (!AllFilesRead)
3864  {
3865  NAMD_die("Tried to assign angle index before all parameter files were read");
3866  }
3867 
3868  /* We need atom1 < atom3. If that was not what we were */
3869  /* passed, switch them */
3870  if (strcasecmp(atom1, atom3) > 0)
3871  {
3872  const char *tmp_name = atom1;
3873  atom1 = atom3;
3874  atom3 = tmp_name;
3875  }
3876 
3877  /* Start at the top */
3878  ptr=anglep;
3879 
3880  /* While we don't have a match and we haven't reached the */
3881  /* bottom of the tree, compare values */
3882  while (!found && (ptr != NULL))
3883  {
3884  comp_val = strcasecmp(atom1, ptr->atom1name);
3885 
3886  if (comp_val == 0)
3887  {
3888  /* Atom 1 matches, so compare atom 2 */
3889  comp_val = strcasecmp(atom2, ptr->atom2name);
3890 
3891  if (comp_val == 0)
3892  {
3893  /* Atoms 1&2 match, try atom 3 */
3894  comp_val = strcasecmp(atom3, ptr->atom3name);
3895  }
3896  }
3897 
3898  if (comp_val == 0)
3899  {
3900  /* Found a match */
3901  found = 1;
3902  angle_ptr->angle_type = ptr->index;
3903  }
3904  else if (comp_val < 0)
3905  {
3906  /* Go left */
3907  ptr=ptr->left;
3908  }
3909  else
3910  {
3911  /* Go right */
3912  ptr=ptr->right;
3913  }
3914  }
3915 
3916  /* Make sure we found a match */
3917  if (!found)
3918  {
3919  char err_msg[128];
3920 
3921  sprintf(err_msg, "UNABLE TO FIND ANGLE PARAMETERS FOR %s %s %s (ATOMS %i %i %i)",
3922  atom1, atom2, atom3, angle_ptr->atom1+1, angle_ptr->atom2+1, angle_ptr->atom3+1);
3923 
3924  if ( notFoundIndex ) {
3925  angle_ptr->angle_type = notFoundIndex;
3926  iout << iWARN << err_msg << "\n" << endi;
3927  return;
3928  } else NAMD_die(err_msg);
3929  }
3930 
3931  return;
3932 }
std::ostream & iWARN(std::ostream &s)
Definition: InfoStream.C:108
#define iout
Definition: InfoStream.h:87
char atom3name[11]
Definition: Parameters.C:64
int32 atom2
Definition: structures.h:56
int32 atom1
Definition: structures.h:55
int32 atom3
Definition: structures.h:57
char atom1name[11]
Definition: Parameters.C:62
void NAMD_die(const char *err_msg)
Definition: common.C:83
Index angle_type
Definition: structures.h:58
struct angle_params * left
Definition: Parameters.C:71
Index index
Definition: Parameters.C:70
char atom2name[11]
Definition: Parameters.C:63
infostream & endi(infostream &s)
Definition: InfoStream.C:38
struct angle_params * right
Definition: Parameters.C:72
void Parameters::assign_bond_index ( const char *  atom1,
const char *  atom2,
Bond bond_ptr 
)

Definition at line 3758 of file Parameters.C.

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

3760 {
3761  struct bond_params *ptr; // Current location in tree
3762  int found=0; // Flag 1-> found a match
3763  int cmp_code; // return code from strcasecmp
3764 
3765  /* Check to make sure the files have all been read */
3766  if (!AllFilesRead)
3767  {
3768  NAMD_die("Tried to assign bond index before all parameter files were read");
3769  }
3770 
3771  /* We need atom1 < atom2, so if that's not the way they */
3772  /* were passed, flip them */
3773  if (strcasecmp(atom1, atom2) > 0)
3774  {
3775  const char *tmp_name = atom1;
3776  atom1 = atom2;
3777  atom2 = tmp_name;
3778  }
3779 
3780  /* Start at the top */
3781  ptr=bondp;
3782 
3783  /* While we haven't found a match and we're not at the end */
3784  /* of the tree, compare the bond passed in with the tree */
3785  while (!found && (ptr!=NULL))
3786  {
3787  cmp_code=strcasecmp(atom1, ptr->atom1name);
3788 
3789  if (cmp_code == 0)
3790  {
3791  cmp_code=strcasecmp(atom2, ptr->atom2name);
3792  }
3793 
3794  if (cmp_code == 0)
3795  {
3796  /* Found a match */
3797  found=1;
3798  bond_ptr->bond_type = ptr->index;
3799  }
3800  else if (cmp_code < 0)
3801  {
3802  /* Go left */
3803  ptr=ptr->left;
3804  }
3805  else
3806  {
3807  /* Go right */
3808  ptr=ptr->right;
3809  }
3810  }
3811 
3812  /* Check to see if we found anything */
3813  if (!found)
3814  {
3815  if ((strcmp(atom1, "DRUD")==0 || strcmp(atom2, "DRUD")==0)
3816  && (strcmp(atom1, "X")!=0 && strcmp(atom2, "X")!=0)) {
3817  /* try a wildcard DRUD X match for this Drude bond */
3818  char a1[8] = "DRUD", a2[8] = "X";
3819  return assign_bond_index(a1, a2, bond_ptr); /* recursive call */
3820  }
3821  else {
3822  char err_msg[128];
3823 
3824  sprintf(err_msg, "UNABLE TO FIND BOND PARAMETERS FOR %s %s (ATOMS %i %i)",
3825  atom1, atom2, bond_ptr->atom1+1, bond_ptr->atom2+1);
3826  NAMD_die(err_msg);
3827  }
3828  }
3829 
3830  return;
3831 }
Index index
Definition: Parameters.C:51
void assign_bond_index(const char *, const char *, Bond *)
Definition: Parameters.C:3758
int32 atom1
Definition: structures.h:48
struct bond_params * right
Definition: Parameters.C:53
char atom1name[11]
Definition: Parameters.C:47
void NAMD_die(const char *err_msg)
Definition: common.C:83
struct bond_params * left
Definition: Parameters.C:52
char atom2name[11]
Definition: Parameters.C:48
Index bond_type
Definition: structures.h:50
int32 atom2
Definition: structures.h:49
void Parameters::assign_crossterm_index ( const char *  atom1,
const char *  atom2,
const char *  atom3,
const char *  atom4,
const char *  atom5,
const char *  atom6,
const char *  atom7,
const char *  atom8,
Crossterm crossterm_ptr 
)

Definition at line 4164 of file Parameters.C.

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

4167 {
4168  struct crossterm_params *ptr; // Current position in list
4169  int found=0; // Flag 1->found a match
4170 
4171  /* Start at the head of the list */
4172  ptr=crosstermp;
4173 
4174  /* While we haven't fuond a match and haven't reached the end */
4175  /* of the list, keep looking */
4176  while (!found && (ptr!=NULL))
4177  {
4178  /* Do a linear search through the linked list of */
4179  /* crossterm parameters. Since the list is arranged */
4180  /* with wildcard paramters at the end of the list, we */
4181  /* can simply do a linear search and be guaranteed that*/
4182  /* we will find exact matches before wildcard matches. */
4183  /* Also, we must check for an exact match, and a match */
4184  /* in reverse, since they are really the same */
4185  /* physically. */
4186  if ( ( (strcasecmp(ptr->atom1name, atom1)==0) ||
4187  (strcasecmp(ptr->atom1name, "X")==0) ) &&
4188  ( (strcasecmp(ptr->atom2name, atom2)==0) ||
4189  (strcasecmp(ptr->atom2name, "X")==0) ) &&
4190  ( (strcasecmp(ptr->atom3name, atom3)==0) ||
4191  (strcasecmp(ptr->atom3name, "X")==0) ) &&
4192  ( (strcasecmp(ptr->atom4name, atom4)==0) ||
4193  (strcasecmp(ptr->atom4name, "X")==0) ) )
4194  {
4195  /* Found an exact match */
4196  found=1;
4197  }
4198  else if ( ( (strcasecmp(ptr->atom4name, atom1)==0) ||
4199  (strcasecmp(ptr->atom4name, "X")==0) ) &&
4200  ( (strcasecmp(ptr->atom3name, atom2)==0) ||
4201  (strcasecmp(ptr->atom3name, "X")==0) ) &&
4202  ( (strcasecmp(ptr->atom2name, atom3)==0) ||
4203  (strcasecmp(ptr->atom2name, "X")==0) ) &&
4204  ( (strcasecmp(ptr->atom1name, atom4)==0) ||
4205  (strcasecmp(ptr->atom1name, "X")==0) ) )
4206  {
4207  /* Found a reverse match */
4208  found=1;
4209  }
4210  if ( ! found ) {
4211  /* Didn't find a match, go to the next node */
4212  ptr=ptr->next;
4213  continue;
4214  }
4215  found = 0;
4216  if ( ( (strcasecmp(ptr->atom5name, atom5)==0) ||
4217  (strcasecmp(ptr->atom5name, "X")==0) ) &&
4218  ( (strcasecmp(ptr->atom6name, atom6)==0) ||
4219  (strcasecmp(ptr->atom6name, "X")==0) ) &&
4220  ( (strcasecmp(ptr->atom7name, atom7)==0) ||
4221  (strcasecmp(ptr->atom7name, "X")==0) ) &&
4222  ( (strcasecmp(ptr->atom8name, atom8)==0) ||
4223  (strcasecmp(ptr->atom8name, "X")==0) ) )
4224  {
4225  /* Found an exact match */
4226  found=1;
4227  }
4228  else if ( ( (strcasecmp(ptr->atom8name, atom5)==0) ||
4229  (strcasecmp(ptr->atom8name, "X")==0) ) &&
4230  ( (strcasecmp(ptr->atom7name, atom6)==0) ||
4231  (strcasecmp(ptr->atom7name, "X")==0) ) &&
4232  ( (strcasecmp(ptr->atom6name, atom7)==0) ||
4233  (strcasecmp(ptr->atom6name, "X")==0) ) &&
4234  ( (strcasecmp(ptr->atom5name, atom8)==0) ||
4235  (strcasecmp(ptr->atom5name, "X")==0) ) )
4236  {
4237  /* Found a reverse match */
4238  found=1;
4239  }
4240  if ( ! found ) {
4241  /* Didn't find a match, go to the next node */
4242  ptr=ptr->next;
4243  }
4244  }
4245 
4246  /* Make sure we found a match */
4247  if (!found)
4248  {
4249  char err_msg[128];
4250 
4251  sprintf(err_msg, "UNABLE TO FIND CROSSTERM PARAMETERS FOR %s %s %s %s %s %s %s %s\n"
4252  "(ATOMS %i %i %i %i %i %i %i %i)",
4253  atom1, atom2, atom3, atom4, atom5, atom6, atom7, atom8,
4254  crossterm_ptr->atom1+1, crossterm_ptr->atom1+2, crossterm_ptr->atom1+3, crossterm_ptr->atom4+1,
4255  crossterm_ptr->atom1+5, crossterm_ptr->atom1+6, crossterm_ptr->atom1+7, crossterm_ptr->atom8+1);
4256 
4257  NAMD_die(err_msg);
4258  }
4259 
4260  /* Assign the constants */
4261  crossterm_ptr->crossterm_type = ptr->index;
4262 
4263  return;
4264 }
char atom8name[11]
Definition: Parameters.C:131
Index crossterm_type
Definition: structures.h:89
char atom5name[11]
Definition: Parameters.C:128
char atom1name[11]
Definition: Parameters.C:124
int32 atom8
Definition: structures.h:88
struct crossterm_params * next
Definition: Parameters.C:135
int32 atom1
Definition: structures.h:81
int32 atom4
Definition: structures.h:84
char atom7name[11]
Definition: Parameters.C:130
char atom3name[11]
Definition: Parameters.C:126
void NAMD_die(const char *err_msg)
Definition: common.C:83
char atom4name[11]
Definition: Parameters.C:127
char atom2name[11]
Definition: Parameters.C:125
char atom6name[11]
Definition: Parameters.C:129
void Parameters::assign_dihedral_index ( const char *  atom1,
const char *  atom2,
const char *  atom3,
const char *  atom4,
Dihedral dihedral_ptr,
int  multiplicity,
int  notFoundIndex 
)

Definition at line 3956 of file Parameters.C.

References dihedral::atom1, dihedral_params::atom1name, dihedral_params::atom1wild, dihedral::atom2, dihedral_params::atom2name, dihedral_params::atom2wild, dihedral::atom3, dihedral_params::atom3name, dihedral_params::atom3wild, dihedral::atom4, dihedral_params::atom4name, dihedral_params::atom4wild, dihedral_array, dihedral::dihedral_type, endi(), dihedral_params::index, iout, iWARN(), dihedral_params::multiplicity, DihedralValue::multiplicity, NAMD_die(), dihedral_params::next, and paraXplor.

3960 {
3961  struct dihedral_params *ptr; // Current position in list
3962  int found=0; // Flag 1->found a match
3963 
3964  /* Start at the begining of the list */
3965  ptr=dihedralp;
3966 
3967  /* While we haven't found a match and we haven't reached */
3968  /* the end of the list, keep looking */
3969  while (!found && (ptr!=NULL))
3970  {
3971  /* Do a linear search through the linked list of */
3972  /* dihedral parameters. Since the list is arranged */
3973  /* with wildcard paramters at the end of the list, we */
3974  /* can simply do a linear search and be guaranteed that*/
3975  /* we will find exact matches before wildcard matches. */
3976  /* Also, we must check for an exact match, and a match */
3977  /* in reverse, since they are really the same */
3978  /* physically. */
3979  if ( ( ptr->atom1wild || (strcasecmp(ptr->atom1name, atom1)==0) ) &&
3980  ( ptr->atom2wild || (strcasecmp(ptr->atom2name, atom2)==0) ) &&
3981  ( ptr->atom3wild || (strcasecmp(ptr->atom3name, atom3)==0) ) &&
3982  ( ptr->atom4wild || (strcasecmp(ptr->atom4name, atom4)==0) ) )
3983  {
3984  /* Found an exact match */
3985  found=1;
3986  }
3987  else if ( ( ptr->atom4wild || (strcasecmp(ptr->atom4name, atom1)==0) ) &&
3988  ( ptr->atom3wild || (strcasecmp(ptr->atom3name, atom2)==0) ) &&
3989  ( ptr->atom2wild || (strcasecmp(ptr->atom2name, atom3)==0) ) &&
3990  ( ptr->atom1wild || (strcasecmp(ptr->atom1name, atom4)==0) ) )
3991  {
3992  /* Found a reverse match */
3993  found=1;
3994  }
3995  else
3996  {
3997  /* Didn't find a match, go to the next node */
3998  ptr=ptr->next;
3999  }
4000  }
4001 
4002  /* Make sure we found a match */
4003  if (!found)
4004  {
4005  char err_msg[128];
4006 
4007  sprintf(err_msg, "UNABLE TO FIND DIHEDRAL PARAMETERS FOR %s %s %s %s (ATOMS %i %i %i %i)",
4008  atom1, atom2, atom3, atom4, dihedral_ptr->atom1+1, dihedral_ptr->atom2+1,
4009  dihedral_ptr->atom3+1, dihedral_ptr->atom4+1);
4010 
4011  if ( notFoundIndex ) {
4012  dihedral_ptr->dihedral_type = notFoundIndex;
4013  iout << iWARN << err_msg << "\n" << endi;
4014  return;
4015  } else NAMD_die(err_msg);
4016  }
4017 
4018  if (paramType == paraXplor) {
4019  // Check to make sure the number of multiples specified in the psf
4020  // file doesn't exceed the number of parameters in the parameter
4021  // files
4022  if (multiplicity > maxDihedralMults[ptr->index])
4023  {
4024  char err_msg[257];
4025 
4026  sprintf(err_msg, "Multiplicity of Paramters for diehedral bond %s %s %s %s of %d exceeded", atom1, atom2, atom3, atom4, maxDihedralMults[ptr->index]);
4027  NAMD_die(err_msg);
4028  }
4029 
4030  // If the multiplicity from the current bond is larger than that
4031  // seen in the past, increase the multiplicity for this bond
4033  {
4035  }
4036  }
4037 
4038  dihedral_ptr->dihedral_type = ptr->index;
4039 
4040  return;
4041 }
char atom4name[11]
Definition: Parameters.C:86
int32 atom3
Definition: structures.h:65
std::ostream & iWARN(std::ostream &s)
Definition: InfoStream.C:108
DihedralValue * dihedral_array
Definition: Parameters.h:243
char atom3name[11]
Definition: Parameters.C:85
#define iout
Definition: InfoStream.h:87
int32 atom4
Definition: structures.h:66
Index dihedral_type
Definition: structures.h:67
void NAMD_die(const char *err_msg)
Definition: common.C:83
char atom1name[11]
Definition: Parameters.C:83
struct dihedral_params * next
Definition: Parameters.C:94
int32 atom2
Definition: structures.h:64
#define paraXplor
Definition: Parameters.h:76
char atom2name[11]
Definition: Parameters.C:84
infostream & endi(infostream &s)
Definition: InfoStream.C:38
int32 atom1
Definition: structures.h:63
void Parameters::assign_improper_index ( const char *  atom1,
const char *  atom2,
const char *  atom3,
const char *  atom4,
Improper improper_ptr,
int  multiplicity 
)

Definition at line 4065 of file Parameters.C.

References improper::atom1, improper_params::atom1name, improper::atom2, improper_params::atom2name, improper::atom3, improper_params::atom3name, improper::atom4, improper_params::atom4name, improper_array, improper::improper_type, improper_params::index, improper_params::multiplicity, ImproperValue::multiplicity, NAMD_die(), improper_params::next, and paraXplor.

4069 {
4070  struct improper_params *ptr; // Current position in list
4071  int found=0; // Flag 1->found a match
4072 
4073  /* Start at the head of the list */
4074  ptr=improperp;
4075 
4076  /* While we haven't fuond a match and haven't reached the end */
4077  /* of the list, keep looking */
4078  while (!found && (ptr!=NULL))
4079  {
4080  /* Do a linear search through the linked list of */
4081  /* improper parameters. Since the list is arranged */
4082  /* with wildcard paramters at the end of the list, we */
4083  /* can simply do a linear search and be guaranteed that*/
4084  /* we will find exact matches before wildcard matches. */
4085  /* Also, we must check for an exact match, and a match */
4086  /* in reverse, since they are really the same */
4087  /* physically. */
4088  if ( ( (strcasecmp(ptr->atom1name, atom1)==0) ||
4089  (strcasecmp(ptr->atom1name, "X")==0) ) &&
4090  ( (strcasecmp(ptr->atom2name, atom2)==0) ||
4091  (strcasecmp(ptr->atom2name, "X")==0) ) &&
4092  ( (strcasecmp(ptr->atom3name, atom3)==0) ||
4093  (strcasecmp(ptr->atom3name, "X")==0) ) &&
4094  ( (strcasecmp(ptr->atom4name, atom4)==0) ||
4095  (strcasecmp(ptr->atom4name, "X")==0) ) )
4096  {
4097  /* Found an exact match */
4098  found=1;
4099  }
4100  else if ( ( (strcasecmp(ptr->atom4name, atom1)==0) ||
4101  (strcasecmp(ptr->atom4name, "X")==0) ) &&
4102  ( (strcasecmp(ptr->atom3name, atom2)==0) ||
4103  (strcasecmp(ptr->atom3name, "X")==0) ) &&
4104  ( (strcasecmp(ptr->atom2name, atom3)==0) ||
4105  (strcasecmp(ptr->atom2name, "X")==0) ) &&
4106  ( (strcasecmp(ptr->atom1name, atom4)==0) ||
4107  (strcasecmp(ptr->atom1name, "X")==0) ) )
4108  {
4109  /* Found a reverse match */
4110  found=1;
4111  }
4112  else
4113  {
4114  /* Didn't find a match, go to the next node */
4115  ptr=ptr->next;
4116  }
4117  }
4118 
4119  /* Make sure we found a match */
4120  if (!found)
4121  {
4122  char err_msg[128];
4123 
4124  sprintf(err_msg, "UNABLE TO FIND IMPROPER PARAMETERS FOR %s %s %s %s (ATOMS %i %i %i %i)",
4125  atom1, atom2, atom3, atom4, improper_ptr->atom1+1, improper_ptr->atom2+1,
4126  improper_ptr->atom3+1, improper_ptr->atom4+1);
4127 
4128  NAMD_die(err_msg);
4129  }
4130 
4131  if (paramType == paraXplor) {
4132  // Check to make sure the number of multiples specified in the psf
4133  // file doesn't exceed the number of parameters in the parameter
4134  // files
4135  if (multiplicity > maxImproperMults[ptr->index])
4136  {
4137  char err_msg[257];
4138 
4139  sprintf(err_msg, "Multiplicity of Paramters for improper bond %s %s %s %s of %d exceeded", atom1, atom2, atom3, atom4, maxImproperMults[ptr->index]);
4140  NAMD_die(err_msg);
4141  }
4142 
4143  // If the multiplicity from the current bond is larger than that
4144  // seen in the past, increase the multiplicity for this bond
4146  {
4148  }
4149  }
4150 
4151  /* Assign the constants */
4152  improper_ptr->improper_type = ptr->index;
4153 
4154  return;
4155 }
struct improper_params * next
Definition: Parameters.C:113
int32 atom4
Definition: structures.h:75
Index improper_type
Definition: structures.h:76
char atom2name[11]
Definition: Parameters.C:107
char atom3name[11]
Definition: Parameters.C:108
ImproperValue * improper_array
Definition: Parameters.h:244
void NAMD_die(const char *err_msg)
Definition: common.C:83
char atom4name[11]
Definition: Parameters.C:109
int32 atom1
Definition: structures.h:72
#define paraXplor
Definition: Parameters.h:76
int32 atom2
Definition: structures.h:73
char atom1name[11]
Definition: Parameters.C:106
int32 atom3
Definition: structures.h:74
void Parameters::assign_vdw_index ( const char *  atomtype,
Atom atom_ptr 
)

Definition at line 3470 of file Parameters.C.

References vdw_params::atomname, endi(), vdw_params::index, iout, iWARN(), vdw_params::left, NAMD_die(), vdw_params::right, and atom_constants::vdw_type.

Referenced by Molecule::prepare_qm().

3472 {
3473  struct vdw_params *ptr; // Current position in trees
3474  int found=0; // Flag 1->found match
3475  int comp_code; // return code from strcasecmp
3476 
3477  /* Check to make sure the files have all been read */
3478  if (!AllFilesRead)
3479  {
3480  NAMD_die("Tried to assign vdw index before all parameter files were read");
3481  }
3482 
3483  /* Start at the top */
3484  ptr=vdwp;
3485 
3486  /* While we haven't found a match, and we haven't reached */
3487  /* the bottom of the tree, compare the atom passed in with */
3488  /* the current value and decide if we have a match, or if not, */
3489  /* which way to go */
3490  while (!found && (ptr!=NULL))
3491  {
3492  comp_code = strcasecmp(atomtype, ptr->atomname);
3493 
3494  if (comp_code == 0)
3495  {
3496  /* Found a match! */
3497  atom_ptr->vdw_type=ptr->index;
3498  found=1;
3499  }
3500  else if (comp_code < 0)
3501  {
3502  /* Go to the left */
3503  ptr=ptr->left;
3504  }
3505  else
3506  {
3507  /* Go to the right */
3508  ptr=ptr->right;
3509  }
3510  }
3511 
3512  //****** BEGIN CHARMM/XPLOR type changes
3513  if (!found)
3514  {
3515  // since CHARMM allows wildcards "*" in vdw typenames
3516  // we have to look again if necessary, this way, if
3517  // we already had an exact match, this is never executed
3518  size_t windx; // wildcard index
3519 
3520  /* Start again at the top */
3521  ptr=vdwp;
3522 
3523  while (!found && (ptr!=NULL))
3524  {
3525 
3526  // get index of wildcard wildcard, get index
3527  windx= strcspn(ptr->atomname,"*");
3528  if (windx == strlen(ptr->atomname))
3529  {
3530  // there is no wildcard here
3531  comp_code = strcasecmp(atomtype, ptr->atomname);
3532  }
3533  else
3534  {
3535  comp_code = strncasecmp(atomtype, ptr->atomname, windx);
3536  }
3537 
3538  if (comp_code == 0)
3539  {
3540  /* Found a match! */
3541  atom_ptr->vdw_type=ptr->index;
3542  found=1;
3543  char errbuf[100];
3544  sprintf(errbuf,"VDW TYPE NAME %s MATCHES PARAMETER TYPE NAME %s",
3545  atomtype, ptr->atomname);
3546  int i;
3547  for(i=0; i<error_msgs.size(); i++) {
3548  if ( strcmp(errbuf,error_msgs[i]) == 0 ) break;
3549  }
3550  if ( i == error_msgs.size() ) {
3551  char *newbuf = new char[strlen(errbuf)+1];
3552  strcpy(newbuf,errbuf);
3553  error_msgs.add(newbuf);
3554  iout << iWARN << newbuf << "\n" << endi;
3555  }
3556  }
3557  else if (comp_code < 0)
3558  {
3559  /* Go to the left */
3560  ptr=ptr->left;
3561  }
3562  else
3563  {
3564  /* Go to the right */
3565  ptr=ptr->right;
3566  }
3567 
3568  }
3569 
3570  }
3571  //****** END CHARMM/XPLOR type changes
3572 
3573  /* Make sure we found it */
3574  if (!found)
3575  {
3576  char err_msg[100];
3577 
3578  sprintf(err_msg, "DIDN'T FIND vdW PARAMETER FOR ATOM TYPE %s",
3579  atomtype);
3580  NAMD_die(err_msg);
3581  }
3582 
3583  return;
3584 }
std::ostream & iWARN(std::ostream &s)
Definition: InfoStream.C:108
#define iout
Definition: InfoStream.h:87
char atomname[11]
Definition: Parameters.C:143
struct vdw_params * left
Definition: Parameters.C:149
void NAMD_die(const char *err_msg)
Definition: common.C:83
int add(const Elem &elem)
Definition: ResizeArray.h:97
Index vdw_type
Definition: structures.h:40
Index index
Definition: Parameters.C:148
int size(void) const
Definition: ResizeArray.h:127
infostream & endi(infostream &s)
Definition: InfoStream.C:38
struct vdw_params * right
Definition: Parameters.C:150
char* Parameters::atom_type_name ( Index  a)
inline

Definition at line 392 of file Parameters.h.

References MAX_ATOMTYPE_CHARS.

392  {
393  return (atomTypeNames + (a * (MAX_ATOMTYPE_CHARS + 1)));
394  }
#define MAX_ATOMTYPE_CHARS
Definition: Parameters.h:71
void Parameters::done_reading_files ( Bool  addDrudeBond)

Definition at line 2959 of file Parameters.C.

References angle_array, bond_array, crossterm_array, dihedral_array, FALSE, gromacsPair_array, improper_array, MAX_ATOMTYPE_CHARS, NAMD_die(), nbthole_array, AngleValue::normal, NumAngleParams, NumBondParams, NumCrosstermParams, NumDihedralParams, NumGromacsPairParams, NumImproperParams, NumNbtholePairParams, NumVdwParams, NumVdwParamsAssigned, TRUE, and vdw_array.

Referenced by Parameters().

2961 {
2962  AllFilesRead = TRUE;
2963 
2964  if (addDrudeBond) {
2965  // default definition for Drude bonds if none given
2966  NumBondParams++;
2967  add_bond_param("X DRUD 500.0 0.0\n", FALSE);
2968  }
2969  // Allocate space for all of the arrays
2970  if (NumBondParams)
2971  {
2973 
2974  if (bond_array == NULL)
2975  {
2976  NAMD_die("memory allocation of bond_array failed!");
2977  }
2978  memset(bond_array, 0, NumBondParams*sizeof(BondValue));
2979  }
2980 
2981  if (NumAngleParams)
2982  {
2984 
2985  if (angle_array == NULL)
2986  {
2987  NAMD_die("memory allocation of angle_array failed!");
2988  }
2989  memset(angle_array, 0, NumAngleParams*sizeof(AngleValue));
2990  for ( Index i=0; i<NumAngleParams; ++i ) {
2991  angle_array[i].normal = 1;
2992  }
2993  }
2994 
2995  if (NumDihedralParams)
2996  {
2998 
2999  if (dihedral_array == NULL)
3000  {
3001  NAMD_die("memory allocation of dihedral_array failed!");
3002  }
3003  memset(dihedral_array, 0, NumDihedralParams*sizeof(DihedralValue));
3004  }
3005 
3006  if (NumImproperParams)
3007  {
3009 
3010  if (improper_array == NULL)
3011  {
3012  NAMD_die("memory allocation of improper_array failed!");
3013  }
3014  memset(improper_array, 0, NumImproperParams*sizeof(ImproperValue));
3015  }
3016 
3017  if (NumCrosstermParams)
3018  {
3021  }
3022 
3023  // JLai
3025  {
3028  }
3029  // End of JLai
3030 
3031  if (NumVdwParams)
3032  {
3033  atomTypeNames = new char[NumVdwParams*(MAX_ATOMTYPE_CHARS+1)];
3035 
3036  if (vdw_array == NULL)
3037  {
3038  NAMD_die("memory allocation of vdw_array failed!");
3039  }
3040  }
3042  {
3044 
3045  if(nbthole_array == NULL)
3046  {
3047  NAMD_die("memory allocation of nbthole_array failed!");
3048  }
3049  }
3050  // Assign indexes to each of the parameters and populate the
3051  // arrays using the binary trees and linked lists that we have
3052  // already read in
3053 
3054  // Note that if parameters have been overwritten (matching
3055  // atom patterns but different parameter values) the tree
3056  // contains fewer elements than Num...Params would suggest.
3057  // The arrays are initialized above because the end values
3058  // may not be occupied. Modifying the Num...Params values
3059  // would break backwards compatibility of memopt extraBonds.
3060 
3061  index_bonds(bondp, 0);
3062  index_angles(anglep, 0);
3063  NumVdwParamsAssigned = index_vdw(vdwp, 0);
3064  index_dihedrals();
3065  index_impropers();
3066  index_crossterms();
3067  convert_nbthole_pairs();
3068  // Convert the vdw pairs
3069  convert_vdw_pairs();
3070  convert_table_pairs();
3071 }
GromacsPairValue * gromacsPair_array
Definition: Parameters.h:247
int NumBondParams
Definition: Parameters.h:259
int NumVdwParams
Definition: Parameters.h:268
#define FALSE
Definition: common.h:116
DihedralValue * dihedral_array
Definition: Parameters.h:243
int NumAngleParams
Definition: Parameters.h:260
CrosstermValue * crossterm_array
Definition: Parameters.h:245
int NumDihedralParams
Definition: Parameters.h:262
AngleValue * angle_array
Definition: Parameters.h:242
int Index
Definition: structures.h:26
int NumNbtholePairParams
Definition: Parameters.h:272
int normal
Definition: Parameters.h:95
int NumImproperParams
Definition: Parameters.h:263
int NumCrosstermParams
Definition: Parameters.h:264
ImproperValue * improper_array
Definition: Parameters.h:244
void NAMD_die(const char *err_msg)
Definition: common.C:83
int NumVdwParamsAssigned
Definition: Parameters.h:270
NbtholePairValue * nbthole_array
Definition: Parameters.h:250
VdwValue * vdw_array
Definition: Parameters.h:249
BondValue * bond_array
Definition: Parameters.h:241
int NumGromacsPairParams
Definition: Parameters.h:266
#define TRUE
Definition: common.h:117
#define MAX_ATOMTYPE_CHARS
Definition: Parameters.h:71
void Parameters::done_reading_structure ( )

Definition at line 4999 of file Parameters.C.

5001 {
5002  if (bondp != NULL)
5003  free_bond_tree(bondp);
5004 
5005  if (anglep != NULL)
5006  free_angle_tree(anglep);
5007 
5008  if (dihedralp != NULL)
5009  free_dihedral_list(dihedralp);
5010 
5011  if (improperp != NULL)
5012  free_improper_list(improperp);
5013 
5014  if (crosstermp != NULL)
5015  free_crossterm_list(crosstermp);
5016 
5017  if (vdwp != NULL)
5018  free_vdw_tree(vdwp);
5019 
5020  // Free the arrays used to track multiplicity for dihedrals
5021  // and impropers
5022  if (maxDihedralMults != NULL)
5023  delete [] maxDihedralMults;
5024 
5025  if (maxImproperMults != NULL)
5026  delete [] maxImproperMults;
5027 
5028  bondp=NULL;
5029  anglep=NULL;
5030  dihedralp=NULL;
5031  improperp=NULL;
5032  crosstermp=NULL;
5033  vdwp=NULL;
5034  maxImproperMults=NULL;
5035  maxDihedralMults=NULL;
5036 }
void Parameters::get_angle_params ( Real k,
Real theta0,
Real k_ub,
Real r_ub,
Index  index 
)
inline

Definition at line 452 of file Parameters.h.

References angle_array, vdw_params::index, AngleValue::k, AngleValue::k_ub, AngleValue::r_ub, and AngleValue::theta0.

454  {
455  *k = angle_array[index].k;
456  *theta0 = angle_array[index].theta0;
457  *k_ub = angle_array[index].k_ub;
458  *r_ub = angle_array[index].r_ub;
459  }
Real r_ub
Definition: Parameters.h:94
AngleValue * angle_array
Definition: Parameters.h:242
Index index
Definition: Parameters.C:148
Real theta0
Definition: Parameters.h:92
Real k_ub
Definition: Parameters.h:93
void Parameters::get_bond_params ( Real k,
Real x0,
Index  index 
)
inline

Definition at line 446 of file Parameters.h.

References bond_array, vdw_params::index, BondValue::k, and BondValue::x0.

Referenced by Molecule::print_bonds().

447  {
448  *k = bond_array[index].k;
449  *x0 = bond_array[index].x0;
450  }
Index index
Definition: Parameters.C:148
Real x0
Definition: Parameters.h:85
BondValue * bond_array
Definition: Parameters.h:241
Real k
Definition: Parameters.h:84
int Parameters::get_dihedral_multiplicity ( Index  index)
inline

Definition at line 466 of file Parameters.h.

References dihedral_array.

467  {
468  return(dihedral_array[index].multiplicity);
469  }
DihedralValue * dihedral_array
Definition: Parameters.h:243
Index index
Definition: Parameters.C:148
void Parameters::get_dihedral_params ( Real k,
int *  n,
Real delta,
Index  index,
int  mult 
)
inline

Definition at line 484 of file Parameters.h.

References four_body_consts::delta, dihedral_array, vdw_params::index, four_body_consts::k, MAX_MULTIPLICITY, four_body_consts::n, NAMD_die(), and DihedralValue::values.

486  {
487  if ( (mult<0) || (mult>MAX_MULTIPLICITY) )
488  {
489  NAMD_die("Bad mult index in Parameters::get_dihedral_params");
490  }
491 
492  *k = dihedral_array[index].values[mult].k;
493  *n = dihedral_array[index].values[mult].n;
494  *delta = dihedral_array[index].values[mult].delta;
495  }
DihedralValue * dihedral_array
Definition: Parameters.h:243
FourBodyConsts values[MAX_MULTIPLICITY]
Definition: Parameters.h:108
#define MAX_MULTIPLICITY
Definition: Parameters.h:68
void NAMD_die(const char *err_msg)
Definition: common.C:83
Index index
Definition: Parameters.C:148
int Parameters::get_improper_multiplicity ( Index  index)
inline

Definition at line 461 of file Parameters.h.

References improper_array.

462  {
463  return(improper_array[index].multiplicity);
464  }
ImproperValue * improper_array
Definition: Parameters.h:244
Index index
Definition: Parameters.C:148
void Parameters::get_improper_params ( Real k,
int *  n,
Real delta,
Index  index,
int  mult 
)
inline

Definition at line 471 of file Parameters.h.

References four_body_consts::delta, improper_array, vdw_params::index, four_body_consts::k, MAX_MULTIPLICITY, four_body_consts::n, NAMD_die(), and ImproperValue::values.

473  {
474  if ( (mult<0) || (mult>MAX_MULTIPLICITY) )
475  {
476  NAMD_die("Bad mult index in Parameters::get_improper_params");
477  }
478 
479  *k = improper_array[index].values[mult].k;
480  *n = improper_array[index].values[mult].n;
481  *delta = improper_array[index].values[mult].delta;
482  }
#define MAX_MULTIPLICITY
Definition: Parameters.h:68
ImproperValue * improper_array
Definition: Parameters.h:244
void NAMD_die(const char *err_msg)
Definition: common.C:83
FourBodyConsts values[MAX_MULTIPLICITY]
Definition: Parameters.h:114
Index index
Definition: Parameters.C:148
int Parameters::get_int_table_type ( char *  tabletype)

Definition at line 7564 of file Parameters.C.

References table_types, and tablenumtypes.

7564  {
7565  for (int i=0; i<tablenumtypes; i++) {
7566  if (!strncmp(tabletype, table_types[i], 5)) {
7567  return i;
7568  }
7569  }
7570 
7571  return -1;
7572 }
static char ** table_types
Definition: Parameters.C:39
int tablenumtypes
Definition: Parameters.h:258
int Parameters::get_num_vdw_params ( void  )
inline
int Parameters::get_table_pair_params ( Index  ind1,
Index  ind2,
int *  K 
)

Definition at line 3605 of file Parameters.C.

References FALSE, indexed_table_pair::ind1, indexed_table_pair::ind2, indexed_table_pair::K, indexed_table_pair::left, indexed_table_pair::right, tab_pair_tree, and TRUE.

3605  {
3606  IndexedTablePair *ptr;
3607  Index temp;
3608  int found=FALSE;
3609 
3610  ptr=tab_pair_tree;
3611  //
3612  // We need the smaller type in ind1, so if it isn't already that
3613  // way, switch them */
3614  if (ind1 > ind2)
3615  {
3616  temp = ind1;
3617  ind1 = ind2;
3618  ind2 = temp;
3619  }
3620 
3621  /* While we haven't found a match and we're not at the end */
3622  /* of the tree, compare the bond passed in with the tree */
3623  while (!found && (ptr!=NULL))
3624  {
3625 // printf("Comparing %i with %i and %i with %i\n", ind1, ptr->ind1, ind2, ptr->ind2);
3626  if ( (ind1 == ptr->ind1) && (ind2 == ptr->ind2) )
3627  {
3628  found = TRUE;
3629  }
3630  else if ( (ind1 < ptr->ind1) ||
3631  ( (ind1==ptr->ind1) && (ind2 < ptr->ind2) ) )
3632  {
3633  /* Go left */
3634  ptr=ptr->left;
3635  }
3636  else
3637  {
3638  /* Go right */
3639  ptr=ptr->right;
3640  }
3641  }
3642 
3643  /* If we found a match, assign the values */
3644  if (found)
3645  {
3646  *K = ptr->K;
3647  return(TRUE);
3648  }
3649  else
3650  {
3651  return(FALSE);
3652  }
3653 }
struct indexed_table_pair * right
Definition: Parameters.h:201
#define FALSE
Definition: common.h:116
int Index
Definition: structures.h:26
struct indexed_table_pair * left
Definition: Parameters.h:202
IndexedTablePair * tab_pair_tree
Definition: Parameters.h:257
#define TRUE
Definition: common.h:117
int Parameters::get_vdw_pair_params ( Index  ind1,
Index  ind2,
Real A,
Real B,
Real A14,
Real B14 
)

Definition at line 3680 of file Parameters.C.

References indexed_vdw_pair::A, indexed_vdw_pair::A14, indexed_vdw_pair::B, indexed_vdw_pair::B14, FALSE, indexed_vdw_pair::ind1, indexed_vdw_pair::ind2, indexed_vdw_pair::left, indexed_vdw_pair::right, TRUE, and vdw_pair_tree.

Referenced by get_vdw_params(), and ComputeNonbondedUtil::select().

3683 {
3684  IndexedVdwPair *ptr; // Current location in tree
3685  Index temp; // Temporary value for swithcing
3686  // values
3687  int found=FALSE; // Flag 1-> found a match
3688 
3689  ptr=vdw_pair_tree;
3690 
3691  // We need the smaller type in ind1, so if it isn't already that
3692  // way, switch them */
3693  if (ind1 > ind2)
3694  {
3695  temp = ind1;
3696  ind1 = ind2;
3697  ind2 = temp;
3698  }
3699 
3700  /* While we haven't found a match and we're not at the end */
3701  /* of the tree, compare the bond passed in with the tree */
3702  while (!found && (ptr!=NULL))
3703  {
3704  if ( (ind1 == ptr->ind1) && (ind2 == ptr->ind2) )
3705  {
3706  found = TRUE;
3707  }
3708  else if ( (ind1 < ptr->ind1) ||
3709  ( (ind1==ptr->ind1) && (ind2 < ptr->ind2) ) )
3710  {
3711  /* Go left */
3712  ptr=ptr->left;
3713  }
3714  else
3715  {
3716  /* Go right */
3717  ptr=ptr->right;
3718  }
3719  }
3720 
3721  /* If we found a match, assign the values */
3722  if (found)
3723  {
3724  *A = ptr->A;
3725  *B = ptr->B;
3726  *A14 = ptr->A14;
3727  *B14 = ptr->B14;
3728 
3729  return(TRUE);
3730  }
3731  else
3732  {
3733  return(FALSE);
3734  }
3735 }
struct indexed_vdw_pair * right
Definition: Parameters.h:168
const BigReal A
#define FALSE
Definition: common.h:116
IndexedVdwPair * vdw_pair_tree
Definition: Parameters.h:255
int Index
Definition: structures.h:26
struct indexed_vdw_pair * left
Definition: Parameters.h:169
const BigReal B
#define TRUE
Definition: common.h:117
void Parameters::get_vdw_params ( Real sigma,
Real epsilon,
Real sigma14,
Real epsilon14,
Index  index 
)
inline

Definition at line 497 of file Parameters.h.

References A, B, cbrt, vdw_val::epsilon, vdw_val::epsilon14, get_vdw_pair_params(), vdw_params::index, NAMD_die(), vdw_val::sigma, vdw_val::sigma14, and vdw_array.

Referenced by Molecule::print_atoms(), and ComputeNonbondedUtil::select().

499  {
500  if ( vdw_array ) {
505  } else {
506  // sigma and epsilon from A and B for each vdw type's interaction with itself
507  Real A; Real B; Real A14; Real B14;
508  if (get_vdw_pair_params(index, index, &A, &B, &A14, &B14) ) {
509  if (A && B) {
510  *sigma = sqrt(cbrt(A)) / sqrt(cbrt(B));
511  *epsilon = B*B / (4*A);
512  }
513  else {
514  *sigma = 0; *epsilon = 0;
515  }
516  if (A14 && B14) {
517  *sigma14 = sqrt(cbrt(A14)) / sqrt(cbrt(B14));
518  *epsilon14 = B14*B14 / (4*A14);
519  }
520  else {
521  *sigma14 = 0; *epsilon14 = 0;
522  }
523  }
524  else {NAMD_die ("Function get_vdw_params failed to derive Lennard-Jones sigma and epsilon from A and B values\n");}
525  }
526  }
Real epsilon
Definition: Parameters.h:151
const BigReal A
float Real
Definition: common.h:107
Real sigma14
Definition: Parameters.h:152
Real sigma14
Definition: Parameters.C:146
Real epsilon14
Definition: Parameters.C:147
Real sigma
Definition: Parameters.h:150
int get_vdw_pair_params(Index ind1, Index ind2, Real *, Real *, Real *, Real *)
Definition: Parameters.C:3680
void NAMD_die(const char *err_msg)
Definition: common.C:83
VdwValue * vdw_array
Definition: Parameters.h:249
Index index
Definition: Parameters.C:148
const BigReal B
Real epsilon14
Definition: Parameters.h:153
Real epsilon
Definition: Parameters.C:145
Real sigma
Definition: Parameters.C:144
#define cbrt(x)
Definition: Parameters.h:31
void Parameters::print_angle_params ( )

Definition at line 4865 of file Parameters.C.

References DebugM, and NumAngleParams.

4866 {
4867  DebugM(3,NumAngleParams << " ANGLE PARAMETERS\n"
4868  << "*****************************************" );
4869  traverse_angle_params(anglep);
4870 }
#define DebugM(x, y)
Definition: Debug.h:59
int NumAngleParams
Definition: Parameters.h:260
void Parameters::print_bond_params ( )

Definition at line 4847 of file Parameters.C.

References DebugM, and NumBondParams.

4848 {
4849  DebugM(3,NumBondParams << " BOND PARAMETERS\n" \
4850  << "*****************************************" \
4851  );
4852 
4853  traverse_bond_params(bondp);
4854 }
int NumBondParams
Definition: Parameters.h:259
#define DebugM(x, y)
Definition: Debug.h:59
void Parameters::print_crossterm_params ( )
void Parameters::print_dihedral_params ( )

Definition at line 4881 of file Parameters.C.

References DebugM, and NumDihedralParams.

4882 {
4883  DebugM(3,NumDihedralParams << " DIHEDRAL PARAMETERS\n" \
4884  << "*****************************************" );
4885 
4886  traverse_dihedral_params(dihedralp);
4887 }
#define DebugM(x, y)
Definition: Debug.h:59
int NumDihedralParams
Definition: Parameters.h:262
void Parameters::print_improper_params ( )

Definition at line 4898 of file Parameters.C.

References DebugM, and NumImproperParams.

4899 {
4900  DebugM(3,NumImproperParams << " IMPROPER PARAMETERS\n" \
4901  << "*****************************************" );
4902 
4903  traverse_improper_params(improperp);
4904 }
#define DebugM(x, y)
Definition: Debug.h:59
int NumImproperParams
Definition: Parameters.h:263
void Parameters::print_nbthole_pair_params ( )

Definition at line 4949 of file Parameters.C.

References DebugM, and NumNbtholePairParams.

4950 {
4951  DebugM(3,NumNbtholePairParams << " NBTHOLE PAIR PARAMETERS\n" \
4952  << "*****************************************" );
4953 
4954  traverse_nbthole_pair_params(nbthole_pairp);
4955 }
#define DebugM(x, y)
Definition: Debug.h:59
int NumNbtholePairParams
Definition: Parameters.h:272
void Parameters::print_param_summary ( )

Definition at line 4966 of file Parameters.C.

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

Referenced by NamdState::loadStructure().

4967 {
4968  iout << iINFO << "SUMMARY OF PARAMETERS:\n"
4969  << iINFO << NumBondParams << " BONDS\n"
4970  << iINFO << NumAngleParams << " ANGLES\n" << endi;
4971  if (cosAngles) {
4972  iout << iINFO << " " << (NumAngleParams - NumCosAngles) << " HARMONIC\n"
4973  << iINFO << " " << NumCosAngles << " COSINE-BASED\n" << endi;
4974  }
4975  iout << iINFO << NumDihedralParams << " DIHEDRAL\n"
4976  << iINFO << NumImproperParams << " IMPROPER\n"
4977  << iINFO << NumCrosstermParams << " CROSSTERM\n"
4978  << iINFO << NumVdwParams << " VDW\n"
4979  << iINFO << NumVdwPairParams << " VDW_PAIRS\n"
4980  << iINFO << NumNbtholePairParams << " NBTHOLE_PAIRS\n" << endi;
4981 }
std::ostream & iINFO(std::ostream &s)
Definition: InfoStream.C:107
int NumBondParams
Definition: Parameters.h:259
int NumVdwParams
Definition: Parameters.h:268
int NumVdwPairParams
Definition: Parameters.h:271
#define iout
Definition: InfoStream.h:87
int NumAngleParams
Definition: Parameters.h:260
int NumDihedralParams
Definition: Parameters.h:262
int NumNbtholePairParams
Definition: Parameters.h:272
int NumImproperParams
Definition: Parameters.h:263
int NumCosAngles
Definition: Parameters.h:261
int NumCrosstermParams
Definition: Parameters.h:264
infostream & endi(infostream &s)
Definition: InfoStream.C:38
void Parameters::print_vdw_pair_params ( )

Definition at line 4932 of file Parameters.C.

References DebugM, and NumVdwPairParams.

4933 {
4934  DebugM(3,NumVdwPairParams << " vdW PAIR PARAMETERS\n" \
4935  << "*****************************************" );
4936 
4937  traverse_vdw_pair_params(vdw_pairp);
4938 }
#define DebugM(x, y)
Definition: Debug.h:59
int NumVdwPairParams
Definition: Parameters.h:271
void Parameters::print_vdw_params ( )

Definition at line 4915 of file Parameters.C.

References DebugM, and NumVdwParams.

4916 {
4917  DebugM(3,NumVdwParams << " vdW PARAMETERS\n" \
4918  << "*****************************************" );
4919 
4920  traverse_vdw_params(vdwp);
4921 }
int NumVdwParams
Definition: Parameters.h:268
#define DebugM(x, y)
Definition: Debug.h:59
void Parameters::read_charmm_parameter_file ( char *  fname)

Definition at line 539 of file Parameters.C.

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

Referenced by Parameters().

541 {
542  int par_type=0; // What type of parameter are we currently
543  // dealing with? (vide infra)
544  int skipline; // skip this line?
545  int skipall = 0; // skip rest of file;
546  char buffer[512]; // Buffer to store each line of the file
547  char first_word[512]; // First word of the current line
548  FILE *pfile; // File descriptor for the parameter file
549 
550  /* Check to make sure that we haven't previously been told */
551  /* that all the files were read */
552  if (AllFilesRead)
553  {
554  NAMD_die("Tried to read another parameter file after being told that all files were read!");
555  }
556 
557  /* Try and open the file */
558  if ( (pfile = fopen(fname, "r")) == NULL)
559  {
560  char err_msg[256];
561 
562  sprintf(err_msg, "UNABLE TO OPEN CHARMM PARAMETER FILE %s\n", fname);
563  NAMD_die(err_msg);
564  }
565 
566  /* Keep reading in lines until we hit the EOF */
567  while (NAMD_read_line(pfile, buffer) != -1)
568  {
569  /* Get the first word of the line */
570  NAMD_find_first_word(buffer, first_word);
571  skipline=0;
572 
573  /* First, screen out things that we ignore. */
574  /* blank lines, lines that start with '!' or '*', lines that */
575  /* start with "END". */
576  if (!NAMD_blank_string(buffer) &&
577  (strncmp(first_word, "!", 1) != 0) &&
578  (strncmp(first_word, "*", 1) != 0) &&
579  (strncasecmp(first_word, "END", 3) != 0))
580  {
581  if ( skipall ) {
582  iout << iWARN << "SKIPPING PART OF PARAMETER FILE AFTER RETURN STATEMENT\n" << endi;
583  break;
584  }
585  /* Now, determine the apropriate parameter type. */
586  if (strncasecmp(first_word, "bond", 4)==0)
587  {
588  par_type=1; skipline=1;
589  }
590  else if (strncasecmp(first_word, "angl", 4)==0)
591  {
592  par_type=2; skipline=1;
593  }
594  else if (strncasecmp(first_word, "thet", 4)==0)
595  {
596  par_type=2; skipline=1;
597  }
598  else if (strncasecmp(first_word, "dihe", 4)==0)
599  {
600  par_type=3; skipline=1;
601  }
602  else if (strncasecmp(first_word, "phi", 3)==0)
603  {
604  par_type=3; skipline=1;
605  }
606  else if (strncasecmp(first_word, "impr", 4)==0)
607  {
608  par_type=4; skipline=1;
609  }
610  else if (strncasecmp(first_word, "imph", 4)==0)
611  {
612  par_type=4; skipline=1;
613  }
614  else if (strncasecmp(first_word, "nbon", 4)==0)
615  {
616  par_type=5; skipline=1;
617  }
618  else if (strncasecmp(first_word, "nonb", 4)==0)
619  {
620  par_type=5; skipline=1;
621  }
622  else if (strncasecmp(first_word, "nbfi", 4)==0)
623  {
624  par_type=6; skipline=1;
625  }
626  else if (strncasecmp(first_word, "hbon", 4)==0)
627  {
628  par_type=7; skipline=1;
629  }
630  else if (strncasecmp(first_word, "cmap", 4)==0)
631  {
632  par_type=8; skipline=1;
633  }
634  else if (strncasecmp(first_word, "nbta", 4)==0)
635  {
636  par_type=9; skipline=1;
637  }
638  else if (strncasecmp(first_word, "thol", 4)==0)
639  {
640  par_type=10; skipline=1;
641  }
642  else if (strncasecmp(first_word, "atom", 4)==0)
643  {
644  par_type=11; skipline=1;
645  }
646  else if (strncasecmp(first_word, "ioformat", 8)==0)
647  {
648  skipline=1;
649  }
650  else if (strncasecmp(first_word, "read", 4)==0)
651  {
652  skip_stream_read(buffer, pfile); skipline=1;
653  }
654  else if (strncasecmp(first_word, "return", 4)==0)
655  {
656  skipall=1; skipline=1;
657  }
658  else if ((strncasecmp(first_word, "nbxm", 4) == 0) ||
659  (strncasecmp(first_word, "grou", 4) == 0) ||
660  (strncasecmp(first_word, "cdie", 4) == 0) ||
661  (strncasecmp(first_word, "shif", 4) == 0) ||
662  (strncasecmp(first_word, "vgro", 4) == 0) ||
663  (strncasecmp(first_word, "vdis", 4) == 0) ||
664  (strncasecmp(first_word, "vswi", 4) == 0) ||
665  (strncasecmp(first_word, "cutn", 4) == 0) ||
666  (strncasecmp(first_word, "ctof", 4) == 0) ||
667  (strncasecmp(first_word, "cton", 4) == 0) ||
668  (strncasecmp(first_word, "eps", 3) == 0) ||
669  (strncasecmp(first_word, "e14f", 4) == 0) ||
670  (strncasecmp(first_word, "wmin", 4) == 0) ||
671  (strncasecmp(first_word, "aexp", 4) == 0) ||
672  (strncasecmp(first_word, "rexp", 4) == 0) ||
673  (strncasecmp(first_word, "haex", 4) == 0) ||
674  (strncasecmp(first_word, "aaex", 4) == 0) ||
675  (strncasecmp(first_word, "noac", 4) == 0) ||
676  (strncasecmp(first_word, "hbno", 4) == 0) ||
677  (strncasecmp(first_word, "cuth", 4) == 0) ||
678  (strncasecmp(first_word, "ctof", 4) == 0) ||
679  (strncasecmp(first_word, "cton", 4) == 0) ||
680  (strncasecmp(first_word, "cuth", 4) == 0) ||
681  (strncasecmp(first_word, "ctof", 4) == 0) ||
682  (strncasecmp(first_word, "cton", 4) == 0) )
683  {
684  if ((par_type != 5) && (par_type != 6) && (par_type != 7) && (par_type != 9))
685  {
686  char err_msg[512];
687 
688  sprintf(err_msg, "ERROR IN CHARMM PARAMETER FILE %s\nLINE=*%s*",fname, buffer);
689  NAMD_die(err_msg);
690  }
691  else
692  {
693  skipline = 1;
694  }
695  }
696  else if (par_type == 0)
697  {
698  /* This is an unknown paramter. */
699  /* This is BAD */
700  char err_msg[512];
701 
702  sprintf(err_msg, "UNKNOWN PARAMETER IN CHARMM PARAMETER FILE %s\nLINE=*%s*",fname, buffer);
703  NAMD_die(err_msg);
704  }
705  }
706  else
707  {
708  skipline=1;
709  }
710 
711  if ( (par_type != 0) && (!skipline) )
712  {
713  /* Now, call the appropriate function based */
714  /* on the type of parameter we have */
715  /* I know, this should really be a switch ... */
716  if (par_type == 1)
717  {
718  add_bond_param(buffer, TRUE);
719  NumBondParams++;
720  }
721  else if (par_type == 2)
722  {
723  add_angle_param(buffer);
724  NumAngleParams++;
725  }
726  else if (par_type == 3)
727  {
728  add_dihedral_param(buffer, pfile);
730  }
731  else if (par_type == 4)
732  {
733  add_improper_param(buffer, pfile);
735  }
736  else if (par_type == 5)
737  {
738  add_vdw_param(buffer);
739  NumVdwParams++;
740  }
741  else if (par_type == 6)
742  {
743  add_vdw_pair_param(buffer);
744  NumVdwPairParams++;
745  }
746  else if (par_type == 7)
747  {
748  add_hb_pair_param(buffer);
749  }
750  else if (par_type == 8)
751  {
752  add_crossterm_param(buffer, pfile);
754  }
755  else if (par_type == 9)
756  {
757 
758  if (table_ener == NULL) {
759  continue;
760  }
761 
762  add_table_pair_param(buffer);
764  }
765  else if (par_type == 10)
766  {
767  add_nbthole_pair_param(buffer);
769  }
770  else if (par_type == 11)
771  {
772  if ( strncasecmp(first_word, "mass", 4) != 0 ) {
773  char err_msg[512];
774  sprintf(err_msg, "UNKNOWN PARAMETER IN CHARMM PARAMETER FILE %s\nLINE=*%s*",fname, buffer);
775  NAMD_die(err_msg);
776  }
777  }
778  else
779  {
780  /* This really should not occour! */
781  /* This is an internal error. */
782  /* This is VERY BAD */
783  char err_msg[512];
784 
785  sprintf(err_msg, "INTERNAL ERROR IN CHARMM PARAMETER FILE %s\nLINE=*%s*",fname, buffer);
786  NAMD_die(err_msg);
787  }
788  }
789  }
790 
791  /* Close the file */
792  fclose(pfile);
793 
794  return;
795 }
int NumTablePairParams
Definition: Parameters.h:273
int NumBondParams
Definition: Parameters.h:259
int NAMD_read_line(FILE *fd, char *buf, int bufsize)
Definition: strlib.C:38
int NumVdwParams
Definition: Parameters.h:268
int NumVdwPairParams
Definition: Parameters.h:271
std::ostream & iWARN(std::ostream &s)
Definition: InfoStream.C:108
#define iout
Definition: InfoStream.h:87
void NAMD_find_first_word(char *source, char *word)
Definition: strlib.C:258
int NumAngleParams
Definition: Parameters.h:260
int NumDihedralParams
Definition: Parameters.h:262
int NumNbtholePairParams
Definition: Parameters.h:272
int NAMD_blank_string(char *str)
Definition: strlib.C:222
int NumImproperParams
Definition: Parameters.h:263
int NumCrosstermParams
Definition: Parameters.h:264
void NAMD_die(const char *err_msg)
Definition: common.C:83
BigReal * table_ener
Definition: Parameters.h:254
infostream & endi(infostream &s)
Definition: InfoStream.C:38
#define TRUE
Definition: common.h:117
void Parameters::read_ener_table ( SimParameters simParams)

Definition at line 6728 of file Parameters.C.

References SimParameters::cutoff, endi(), iout, NAMD_die(), namdnearbyint, numenerentries, read_energy_type(), read_energy_type_bothcubspline(), table_ener, table_types, SimParameters::tableInterpType, SimParameters::tableMaxDist, tablenumtypes, SimParameters::tableNumTypes, SimParameters::tableSpacing, and SimParameters::tabulatedEnergiesFile.

Referenced by Parameters().

6728  {
6729  char* table_file = simParams->tabulatedEnergiesFile;
6730  char* interp_type = simParams->tableInterpType;
6731  FILE* enertable;
6732  enertable = fopen(table_file, "r");
6733 
6734  if (enertable == NULL) {
6735  NAMD_die("ERROR: Could not open energy table file!\n");
6736  }
6737 
6738  char tableline[121];
6739  char* newtypename;
6740  int numtypes;
6741  int atom1;
6742  int atom2;
6743  int distbin;
6744  int readcount;
6745  Real dist;
6746  Real energy;
6747  Real force;
6748  Real table_spacing;
6749  Real maxdist;
6750 
6751 /* First find the header line and set the variables we need */
6752  iout << "Beginning table read\n" << endi;
6753  while(fgets(tableline,120,enertable)) {
6754  if (strncmp(tableline,"#",1)==0) {
6755  continue;
6756  }
6757  readcount = sscanf(tableline, "%i %f %f", &numtypes, &table_spacing, &maxdist);
6758  if (readcount != 3) {
6759  NAMD_die("ERROR: Couldn't parse table header line\n");
6760  }
6761  break;
6762  }
6763 
6764  if (maxdist < simParams->cutoff) {
6765  NAMD_die("Tabulated energies must at least extend to the cutoff distance\n");
6766  }
6767 
6768  if (maxdist > simParams->cutoff) {
6769  iout << "Info: Reducing max table size to match nonbond cutoff\n";
6770  maxdist = ceil(simParams->cutoff);
6771  }
6772 
6773 /* Now allocate memory for the table; we know what we should be getting */
6774  numenerentries = 2 * numtypes * int (namdnearbyint(maxdist/table_spacing) + 1);
6775  //Set up a full energy lookup table from a file
6776  //Allocate the table; layout is atom1 x atom2 x distance energy force
6777  fprintf(stdout, "Table has %i entries\n",numenerentries);
6778  //iout << "Allocating original energy table\n" << endi;
6779  if ( table_ener ) delete [] table_ener;
6781  if ( table_types ) delete [] table_types;
6782  table_types = new char*[numtypes];
6783 
6784 /* Finally, start reading the table */
6785  int numtypesread = 0; //Number of types read so far
6786 
6787  while(fgets(tableline,120,enertable)) {
6788  if (strncmp(tableline,"#",1)==0) {
6789  continue;
6790  }
6791  if (strncmp(tableline,"TYPE",4)==0) {
6792  // Read a new type
6793  newtypename = new char[5];
6794  int readcount = sscanf(tableline, "%*s %s", newtypename);
6795  if (readcount != 1) {
6796  NAMD_die("Couldn't read interaction type from TYPE line\n");
6797  }
6798 // printf("Setting table_types[%i] to %s\n", numtypesread, newtypename);
6799  table_types[numtypesread] = newtypename;
6800 
6801  if (numtypesread == numtypes) {
6802  NAMD_die("Error: Number of types in table doesn't match table header\n");
6803  }
6804 
6805  // Read the current energy type with the proper interpolation
6806  if (!strncasecmp(interp_type, "linear", 6)) {
6807  if (read_energy_type(enertable, numtypesread, table_ener, table_spacing, maxdist) != 0) {
6808  char err_msg[512];
6809  sprintf(err_msg, "Failed to read table energy (linear) type %s\n", newtypename);
6810  NAMD_die(err_msg);
6811  }
6812  } else if (!strncasecmp(interp_type, "cubic", 5)) {
6813  if (read_energy_type_bothcubspline(enertable, numtypesread, table_ener, table_spacing, maxdist) != 0) {
6814  char err_msg[512];
6815  sprintf(err_msg, "Failed to read table energy (cubic) type %s\n", newtypename);
6816  NAMD_die(err_msg);
6817  }
6818  } else {
6819  NAMD_die("Table interpolation type must be linear or cubic\n");
6820  }
6821 
6822  numtypesread++;
6823  continue;
6824  }
6825  //char err_msg[512];
6826  //sprintf(err_msg, "Unrecognized line in energy table file: %s\n", tableline);
6827  //NAMD_die(err_msg);
6828  }
6829 
6830  // See if we got what we expected
6831  if (numtypesread != numtypes) {
6832  char err_msg[512];
6833  sprintf(err_msg, "ERROR: Expected %i tabulated energy types but got %i\n", numtypes, numtypesread);
6834  NAMD_die(err_msg);
6835  }
6836 
6837  // Move the necessary information into simParams
6838  simParams->tableNumTypes = numtypes;
6839  simParams->tableSpacing = table_spacing;
6840  simParams->tableMaxDist = maxdist;
6841  tablenumtypes = numtypes;
6842 
6843  /*
6844 xxxxxx
6845  int numtypes = simParams->tableNumTypes;
6846 
6847  //This parameter controls the table spacing
6848  BigReal table_spacing = 0.01;
6849  BigReal maxdist = 20.0;
6850  if (maxdist > simParams->cutoff) {
6851  iout << "Info: Reducing max table size to match nonbond cutoff\n";
6852  maxdist = ceil(simParams->cutoff);
6853  }
6854 
6855  numenerentries = (numtypes + 1) * numtypes * int (ceil(maxdist/table_spacing));
6856 // fprintf(stdout, "Table arithmetic: maxdist %f, table_spacing %f, numtypes %i, numentries %i\n", maxdist, table_spacing, numtypes, numenerentries);
6857  columnsize = 2 * namdnearbyint(maxdist/table_spacing);
6858  rowsize = numtypes * columnsize;
6859  //Set up a full energy lookup table from a file
6860  //Allocate the table; layout is atom1 x atom2 x distance energy force
6861  fprintf(stdout, "Table has %i entries\n",numenerentries);
6862  //iout << "Allocating original energy table\n" << endi;
6863  if ( table_ener ) delete [] table_ener;
6864  table_ener = new Real[numenerentries];
6865  //
6866  //Set sentinel values for uninitialized table energies
6867  for (int i=0 ; i<numenerentries ; i++) {
6868  table_ener[i] = 1337.0;
6869  }
6870  Real compval = 1337.0;
6871 
6872  // iout << "Entering table reading\n" << endi;
6873  //iout << "Finished allocating table\n" << endi;
6874 
6875  //Counter to make sure we read the right number of energy entries
6876  int numentries = 0;
6877 
6878  //Now, start reading from the file
6879  char* table_file = simParams->tabulatedEnergiesFile;
6880  FILE* enertable;
6881  enertable = fopen(table_file, "r");
6882 
6883  if (enertable == NULL) {
6884  NAMD_die("ERROR: Could not open energy table file!\n");
6885  }
6886 
6887  char tableline[121];
6888  int atom1;
6889  int atom2;
6890  int distbin;
6891  Real dist;
6892  Real energy;
6893  Real force;
6894 
6895  iout << "Beginning table read\n" << endi;
6896  //Read the table entries
6897  while(fgets(tableline,120,enertable)) {
6898 // IOut << "Processing line " << tableline << "\n" << endi;
6899  if (strncmp(tableline,"#",1)==0) {
6900  continue;
6901  }
6902 
6903 
6904  sscanf(tableline, "%i %i %f %f %f\n", &atom1, &atom2, &dist, &energy, &force);
6905  distbin = int(namdnearbyint(dist/table_spacing));
6906 // fprintf(stdout, "Working on atoms %i and %i at distance %f\n",atom1,atom2,dist);
6907  if ((atom2 > atom1) || (distbin > int(namdnearbyint(maxdist/table_spacing)))) {
6908 // fprintf(stdout,"Rejected\n");
6909 // fprintf(stdout, "Error: Throwing out energy line beyond bounds\n");
6910  // fprintf(stdout, "Bad line: Atom 1: %i Atom 2: %i Distance Bin: %i Max Distance Bin: %i \n", atom1, atom2, distbin, int(namdnearbyint(maxdist/table_spacing)));
6911  } else {
6912  //The magic formula for the number of columns from previous rows
6913  //in the triangular matrix is (2ni+i-i**2)/2
6914  //Here n is the number of types, and i is atom2
6915 // fprintf(stdout, "Input: atom1 %f atom2 %f columnsize %f \n", float(atom1), float(atom2), float(columnsize));
6916 // 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);
6917  int eneraddress = columnsize*((2*numtypes*atom2 + atom2 - atom2*atom2)/2) + columnsize*(atom1-atom2) + 2*distbin - 2;
6918  int forceaddress = eneraddress + 1;
6919 // 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]);
6920  fflush(stdout);
6921 // fprintf(stdout, "Checking for dupes: Looking at: %f %f \n", table_ener[eneraddress], table_ener[forceaddress]);
6922  if ((table_ener[eneraddress] == compval && table_ener[forceaddress] == compval)) {
6923  numentries++;
6924  table_ener[eneraddress] = energy;
6925  table_ener[forceaddress] = force;
6926 // 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]);
6927  //table_ener[rowsize*atom2 + columnsize*atom1 + 2*distbin] = energy;
6928  //table_ener[rowsize*atom2 + columnsize*atom1 + 2*distbin + 1] = force;
6929 // fflush(stdout);
6930  } else {
6931 // fprintf(stdout,"Rejecting duplicate entry\n");
6932  }
6933  }
6934  // iout << "Done with line\n"<< endi;
6935  }
6936 
6937  // int should = numtypes * numtypes * (maxdist/table_spacing);
6938  // cout << "Entries: " << numentries << " should be " << should << "\n" << endi;
6939 // int exptypes = ceil((numtypes+1) * numtypes * (maxdist/table_spacing));
6940 //fprintf(stdout,"numtypes: %i maxdist: %f table_spacing: %f exptypes: %i\n",numtypes,maxdist,table_spacing);
6941  if (numentries != int(numenerentries/2)) {
6942  fprintf(stdout,"ERROR: Expected %i entries but got %i\n",numenerentries/2, numentries);
6943  NAMD_die("Number of energy table entries did not match expected value\n");
6944  }
6945  // iout << "Done with table\n"<< endi;
6946  fprintf(stdout, "Numenerentries: %i\n",numenerentries/2);
6947  */
6948 } /* END of function read_ener_table */
int numenerentries
Definition: Parameters.h:251
#define namdnearbyint(x)
Definition: common.h:74
static char ** table_types
Definition: Parameters.C:39
float Real
Definition: common.h:107
#define iout
Definition: InfoStream.h:87
int read_energy_type(FILE *, const int, BigReal *, const float, const float)
Definition: Parameters.C:7446
char tabulatedEnergiesFile[128]
void NAMD_die(const char *err_msg)
Definition: common.C:83
BigReal tableMaxDist
BigReal * table_ener
Definition: Parameters.h:254
int read_energy_type_bothcubspline(FILE *, const int, BigReal *, const float, const float)
Definition: Parameters.C:6971
int tablenumtypes
Definition: Parameters.h:258
infostream & endi(infostream &s)
Definition: InfoStream.C:38
char tableInterpType[128]
double BigReal
Definition: common.h:112
int Parameters::read_energy_type ( FILE *  enertable,
const int  typeindex,
BigReal table_ener,
const float  table_spacing,
const float  maxdist 
)

Definition at line 7446 of file Parameters.C.

References NAMD_die(), and namdnearbyint.

Referenced by read_ener_table().

7446  {
7447 
7448  char tableline[120];
7449 
7450  /* Last values read from table */
7451  BigReal readdist;
7452  BigReal readener;
7453  BigReal readforce;
7454  BigReal lastdist;
7455  BigReal lastener;
7456  BigReal lastforce;
7457  readdist = -1.0;
7458  readener = 0.0;
7459  readforce = 0.0;
7460 
7461  /* Keep track of where in the table we are */
7462  float currdist;
7463  int distbin;
7464  currdist = -1.0;
7465  distbin = -1;
7466 
7467  while(fgets(tableline,120,enertable) && distbin <= (int) (namdnearbyint(maxdist / table_spacing) + 1)) {
7468  printf("At distance %f + %f vs. %f\n", currdist, table_spacing, maxdist);
7469  if (strncmp(tableline,"#",1)==0) {
7470  continue;
7471  }
7472  if (strncmp(tableline,"TYPE",4)==0) {
7473  break;
7474  }
7475 
7476  // Read an energy line from the table
7477  lastdist = readdist;
7478  lastener = readener;
7479  lastforce = readforce;
7480  int readcount = sscanf(tableline, "%lf %lf %lf", &readdist, &readener, &readforce);
7481  if (distbin == -1) {
7482  if (readdist != 0.0) {
7483  NAMD_die("ERROR: Energy/force tables must start at d=0.0\n");
7484  } else {
7485  distbin = 0;
7486  continue;
7487  }
7488  }
7489  // printf("Read an energy line: %f %f %f\n", readdist, readener, readforce);
7490  if (readcount != 3) {
7491  char err_msg[512];
7492  sprintf(err_msg, "ERROR: Failed to parse table line %s!\n", tableline);
7493  NAMD_die(err_msg);
7494  }
7495 
7496  //Sanity check the current entry
7497  if (readdist < lastdist) {
7498  NAMD_die("ERROR: Encountered badly ordered entries in energy table!\n");
7499  }
7500 
7501  currdist = lastdist;
7502 
7503  while (currdist <= readdist && distbin <= (int) (namdnearbyint(maxdist / table_spacing))) {
7504  distbin = (int) (namdnearbyint(currdist / table_spacing));
7505  int table_loc = 2 * (distbin + (typeindex * (1 + (int) namdnearbyint(maxdist / table_spacing))));
7506  printf("Doing interpolation for energy between %f %f and %f %f: Dist %f\n", readener, readdist, lastener, lastdist, currdist);
7507  table_ener[table_loc] = interp_lin(readener, lastener, readdist, lastdist, currdist);
7508  table_ener[table_loc + 1] = interp_lin(readforce, lastforce, readdist, lastdist, currdist);
7509  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);
7510  currdist += table_spacing;
7511  distbin++;
7512  }
7513  }
7514 
7515  // Go back one line, since we should now be into the next TYPE block
7516  fseek(enertable, -1 * (long) strlen(tableline), SEEK_CUR);
7517 
7518  // Clean up and make sure everything worked ok
7519  distbin--;
7520  printf("Testing: %i vs %i (from %f / %f)\n", distbin, (int) (namdnearbyint(maxdist / table_spacing)), maxdist, table_spacing);
7521  if (distbin != (int) (namdnearbyint(maxdist / table_spacing))) return 1;
7522  return 0;
7523 }
#define namdnearbyint(x)
Definition: common.h:74
void NAMD_die(const char *err_msg)
Definition: common.C:83
BigReal * table_ener
Definition: Parameters.h:254
double BigReal
Definition: common.h:112
int Parameters::read_energy_type_bothcubspline ( FILE *  enertable,
const int  typeindex,
BigReal table_ener,
const float  table_spacing,
const float  maxdist 
)

Definition at line 6971 of file Parameters.C.

References INDEX, NAMD_die(), and namdnearbyint.

Referenced by read_ener_table().

6971  {
6972 
6973  char tableline[120];
6974  int i,j;
6975 
6976  /* Last values read from table */
6977  BigReal readdist;
6978  BigReal readener;
6979  BigReal readforce;
6980  BigReal spacing;
6981 // BigReal readforce;
6982  BigReal lastdist;
6983 // BigReal lastener;
6984 // BigReal lastforce;
6985 // readdist = -1.0;
6986 // readener = 0.0;
6987 // readforce = 0.0;
6988 
6989  /* Create arrays for holding the input values */
6990  std::vector<BigReal> dists;
6991  std::vector<BigReal> enervalues;
6992  std::vector<BigReal> forcevalues;
6993  int numentries = 0;
6994 
6995 
6996  /* Keep track of where in the table we are */
6997  BigReal currdist;
6998  int distbin;
6999  currdist = 0.0;
7000  lastdist = -1.0;
7001  distbin = 0;
7002 
7003  // Read all of the values first -- we'll interpolate later
7004  while(fgets(tableline,120,enertable) && distbin <= (int) (namdnearbyint(maxdist / table_spacing) + 1)) {
7005  if (strncmp(tableline,"#",1)==0) {
7006  continue;
7007  }
7008  if (strncmp(tableline,"TYPE",4)==0) {
7009  fseek(enertable, -1 * (long) strlen(tableline), SEEK_CUR);
7010  break;
7011  }
7012 
7013  // Read an energy line from the table
7014  int readcount = sscanf(tableline, "%lf %lf %lf", &readdist, &readener, &readforce);
7015 
7016  //printf("Read an energy line: %g %g %g\n", readdist, readener, readforce);
7017  if (readcount != 3) {
7018  char err_msg[512];
7019  sprintf(err_msg, "ERROR: Failed to parse table line %s!\n", tableline);
7020  NAMD_die(err_msg);
7021  }
7022 
7023  //Sanity check the current entry
7024  if (readdist < lastdist) {
7025  NAMD_die("ERROR: Encountered badly ordered entries in energy table!\n");
7026  }
7027 
7028  lastdist = readdist;
7029  dists.push_back(readdist);
7030  enervalues.push_back(readener);
7031  forcevalues.push_back(readforce);
7032  numentries++;
7033  }
7034 
7035  // Check the spacing and make sure it is uniform
7036  if (dists[0] != 0.0) {
7037  NAMD_die("ERROR: First data point for energy table must be at r=0\n");
7038  }
7039  spacing = dists[1] - dists[0];
7040  for (i=2; i<(numentries - 1); i++) {
7041  BigReal myspacing;
7042  myspacing = dists[i] - dists[i-1];
7043  if (fabs(myspacing - spacing) > 0.00001) {
7044  printf("Bad spacing in table: %f should be %f (between distances %f and %f)\n", myspacing, spacing, dists[i-1], dists[i]);
7045  NAMD_die("ERROR: Nonuniform table encountered on cubic interpolation. Use a uniformly spaced table or switch to linear interpolation.\n");
7046  }
7047  }
7048 
7049 // Perform cubic spline interpolation to get the energies and forces
7050 
7051  /* allocate spline coefficient matrix */
7052  // xe and xf / be and bf for energies and forces, respectively
7053  double* m = new double[numentries*numentries];
7054  double* xe = new double[numentries];
7055  double* xf = new double[numentries];
7056  double* be = new double[numentries];
7057  double* bf = new double[numentries];
7058 
7059  be[0] = 3 * (enervalues[1] - enervalues[0]);
7060  for (i=1; i < (numentries - 1); i++) {
7061 // printf("Control point %i at %f\n", i, enervalues[i]);
7062  be[i] = 3 * (enervalues[i+1] - enervalues[i-1]);
7063 // printf("be is %f\n", be[i]);
7064  }
7065  be[numentries - 1] = 3 * (enervalues[numentries - 1] - enervalues[numentries - 2]);
7066 
7067  bf[0] = 3 * (forcevalues[1] - forcevalues[0]);
7068  for (i=1; i < (numentries - 1); i++) {
7069 // printf("Control point %i at %f\n", i, forcevalues[i]);
7070  bf[i] = 3 * (forcevalues[i+1] - forcevalues[i-1]);
7071 // printf("bf is %f\n", bf[i]);
7072  }
7073  bf[numentries - 1] = 3 * (forcevalues[numentries - 1] - forcevalues[numentries - 2]);
7074 
7075  memset(m,0,numentries*numentries*sizeof(double));
7076 
7077  /* initialize spline coefficient matrix */
7078  m[0] = 2;
7079  for (i = 1; i < numentries; i++) {
7080  m[INDEX(numentries,i-1,i)] = 1;
7081  m[INDEX(numentries,i,i-1)] = 1;
7082  m[INDEX(numentries,i,i)] = 4;
7083  }
7084  m[INDEX(numentries,numentries-1,numentries-1)] = 2;
7085 
7086  /* Now we need to solve the equation M X = b for X */
7087  // Do this simultaneously for energy and force -- ONLY because we use the same matrix
7088 
7089  //Put m in upper triangular form and apply corresponding operations to b
7090  for (i=0; i<numentries; i++) {
7091  // zero the ith column in all rows below i
7092  const BigReal baseval = m[INDEX(numentries,i,i)];
7093  for (j=i+1; j<numentries; j++) {
7094  const BigReal myval = m[INDEX(numentries,j,i)];
7095  const BigReal factor = myval / baseval;
7096 
7097  for (int k=i; k<numentries; k++) {
7098  const BigReal subval = m[INDEX(numentries,i,k)];
7099  m[INDEX(numentries,j,k)] -= (factor * subval);
7100  }
7101 
7102  be[j] -= (factor * be[i]);
7103  bf[j] -= (factor * bf[i]);
7104 
7105  }
7106  }
7107 
7108  //Now work our way diagonally up from the bottom right to assign values to X
7109  for (i=numentries-1; i>=0; i--) {
7110 
7111  //Subtract the effects of previous columns
7112  for (j=i+1; j<numentries; j++) {
7113  be[i] -= ( m[INDEX(numentries,i,j)] * xe[j] );
7114  bf[i] -= ( m[INDEX(numentries,i,j)] * xf[j] );
7115  }
7116 
7117  xe[i] = be[i] / m[INDEX(numentries,i,i)];
7118  xf[i] = bf[i] / m[INDEX(numentries,i,i)];
7119 
7120  }
7121 
7122  // We now have the coefficient information we need to make the table
7123  // Now interpolate on each interval we want
7124 
7125  distbin = 0;
7126  int entriesperseg = (int) ceil(spacing / table_spacing);
7127  int table_prefix = 2 * typeindex * (int) (namdnearbyint(maxdist / table_spacing) + 1);
7128 
7129  for (i=0; i<numentries-1; i++) {
7130  BigReal Ae,Be,Ce,De;
7131  BigReal Af,Bf,Cf,Df;
7132  currdist = dists[i];
7133 
7134 // printf("Interpolating on interval %i\n", i);
7135 
7136  // Set up the coefficients for this section
7137  Ae = enervalues[i];
7138  Be = xe[i];
7139  Ce = 3 * (enervalues[i+1] - enervalues[i]) - (2 * xe[i]) - xe[i+1];
7140  De = 2 * (enervalues[i] - enervalues[i+1]) + xe[i] + xe[i+1];
7141 
7142  Af = forcevalues[i];
7143  Bf = xf[i];
7144  Cf = 3 * (forcevalues[i+1] - forcevalues[i]) - (2 * xf[i]) - xf[i+1];
7145  Df = 2 * (forcevalues[i] - forcevalues[i+1]) + xf[i] + xf[i+1];
7146 
7147  // Go over the region of interest and fill in the table
7148  for (j=0; j<entriesperseg; j++) {
7149  const BigReal mydist = currdist + (j * table_spacing);
7150  const BigReal mydistfrac = (float) j / (entriesperseg - 1);
7151  distbin = (int) namdnearbyint(mydist / table_spacing);
7152  if (distbin >= (int) namdnearbyint(maxdist / table_spacing)) break;
7153  BigReal energy;
7154  BigReal force;
7155 
7156  energy = Ae + (Be * mydistfrac) + (Ce * mydistfrac * mydistfrac) + (De * mydistfrac * mydistfrac * mydistfrac);
7157  force = Af + (Bf * mydistfrac) + (Cf * mydistfrac * mydistfrac) + (Df * mydistfrac * mydistfrac * mydistfrac);
7158 
7159 // 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);
7160  table_ener[table_prefix + 2 * distbin] = energy;
7161  table_ener[table_prefix + 2 * distbin + 1] = force;
7162  distbin++;
7163  }
7164  if (currdist >= maxdist) break;
7165  }
7166 
7167  //The procedure above leaves out the last entry -- add it explicitly
7168  distbin = (int) namdnearbyint(maxdist / table_spacing);
7169 // 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));
7170  table_ener[table_prefix + 2 * distbin] = enervalues[numentries - 1];
7171  table_ener[table_prefix + 2 * distbin + 1] = 0.0;
7172  distbin++;
7173 
7174 
7175  // Clean up and make sure everything worked ok
7176  delete m;
7177  delete xe;
7178  delete xf;
7179  delete be;
7180  delete bf;
7181  distbin--;
7182  printf("Testing: %i vs %i (from %f / %f)\n", distbin, (int) (namdnearbyint(maxdist / table_spacing)), maxdist, table_spacing);
7183  if (distbin != (int) (namdnearbyint(maxdist / table_spacing))) return 1;
7184  return 0;
7185 } /* end read_energy_type_bothcubspline */
#define namdnearbyint(x)
Definition: common.h:74
#define INDEX(ncols, i, j)
Definition: Parameters.C:35
void NAMD_die(const char *err_msg)
Definition: common.C:83
BigReal * table_ener
Definition: Parameters.h:254
double BigReal
Definition: common.h:112
int Parameters::read_energy_type_cubspline ( FILE *  enertable,
const int  typeindex,
BigReal table_ener,
const float  table_spacing,
const float  maxdist 
)

Definition at line 7207 of file Parameters.C.

References A, B, INDEX, NAMD_die(), namdnearbyint, and x.

7207  {
7208 
7209  char tableline[120];
7210  int i,j;
7211 
7212  /* Last values read from table */
7213  BigReal readdist;
7214  BigReal readener;
7215  BigReal spacing;
7216 // BigReal readforce;
7217  BigReal lastdist;
7218 // BigReal lastener;
7219 // BigReal lastforce;
7220 // readdist = -1.0;
7221 // readener = 0.0;
7222 // readforce = 0.0;
7223 
7224  /* Create arrays for holding the input values */
7225  std::vector<BigReal> dists;
7226  std::vector<BigReal> enervalues;
7227  int numentries = 0;
7228 
7229 
7230  /* Keep track of where in the table we are */
7231  BigReal currdist;
7232  int distbin;
7233  currdist = 0.0;
7234  lastdist = -1.0;
7235  distbin = 0;
7236 
7237  // Read all of the values first -- we'll interpolate later
7238  while(fgets(tableline,120,enertable) && distbin <= (int) (namdnearbyint(maxdist / table_spacing) + 1)) {
7239  if (strncmp(tableline,"#",1)==0) {
7240  continue;
7241  }
7242  if (strncmp(tableline,"TYPE",4)==0) {
7243  fseek(enertable, -1 * (long) strlen(tableline), SEEK_CUR);
7244  break;
7245  }
7246 
7247  // Read an energy line from the table
7248  int readcount = sscanf(tableline, "%lf %lf", &readdist, &readener);
7249 
7250  // printf("Read an energy line: %f %f %f\n", readdist, readener, readforce);
7251  if (readcount != 2) {
7252  char err_msg[512];
7253  sprintf(err_msg, "ERROR: Failed to parse table line %s!\n", tableline);
7254  NAMD_die(err_msg);
7255  }
7256 
7257  //Sanity check the current entry
7258  if (readdist < lastdist) {
7259  NAMD_die("ERROR: Encountered badly ordered entries in energy table!\n");
7260  }
7261 
7262  lastdist = readdist;
7263  dists.push_back(readdist);
7264  enervalues.push_back(readener);
7265  numentries++;
7266  }
7267 
7268  // Check the spacing and make sure it is uniform
7269  if (dists[0] != 0.0) {
7270  NAMD_die("ERROR: First data point for energy table must be at r=0\n");
7271  }
7272  spacing = dists[1] - dists[0];
7273  for (i=2; i<(numentries - 1); i++) {
7274  BigReal myspacing;
7275  myspacing = dists[i] - dists[i-1];
7276  if (fabs(myspacing - spacing) > 0.00001) {
7277  printf("Bad spacing in table: %f should be %f (between distances %f and %f)\n", myspacing, spacing, dists[i-1], dists[i]);
7278  NAMD_die("ERROR: Nonuniform table encountered on cubic interpolation. Use a uniformly spaced table or switch to linear interpolation.\n");
7279  }
7280  }
7281 
7282 // Perform cubic spline interpolation to get the energies and forces
7283 
7284  /* allocate spline coefficient matrix */
7285  double* m = new double[numentries*numentries];
7286  double* x = new double[numentries];
7287  double* b = new double[numentries];
7288 
7289  b[0] = 3 * (enervalues[1] - enervalues[0]);
7290  for (i=1; i < (numentries - 1); i++) {
7291  printf("Control point %i at %f\n", i, enervalues[i]);
7292  b[i] = 3 * (enervalues[i+1] - enervalues[i-1]);
7293  printf("b is %f\n", b[i]);
7294  }
7295  b[numentries - 1] = 3 * (enervalues[numentries - 1] - enervalues[numentries - 2]);
7296 
7297  memset(m,0,numentries*numentries*sizeof(double));
7298 
7299  /* initialize spline coefficient matrix */
7300  m[0] = 2;
7301  for (i = 1; i < numentries; i++) {
7302  m[INDEX(numentries,i-1,i)] = 1;
7303  m[INDEX(numentries,i,i-1)] = 1;
7304  m[INDEX(numentries,i,i)] = 4;
7305  }
7306  m[INDEX(numentries,numentries-1,numentries-1)] = 2;
7307 
7308  /* Now we need to solve the equation M X = b for X */
7309 
7310  printf("Solving the matrix equation: \n");
7311 
7312  for (i=0; i<numentries; i++) {
7313  printf(" ( ");
7314  for (j=0; j<numentries; j++) {
7315  printf(" %6.3f,", m[INDEX(numentries, i, j)]);
7316  }
7317  printf(" ) ( D%-3i ) = ( %6.3f )\n", i, b[i]);
7318  }
7319 
7320  //Put m in upper triangular form and apply corresponding operations to b
7321  for (i=0; i<numentries; i++) {
7322  // zero the ith column in all rows below i
7323  const BigReal baseval = m[INDEX(numentries,i,i)];
7324  for (j=i+1; j<numentries; j++) {
7325  const BigReal myval = m[INDEX(numentries,j,i)];
7326  const BigReal factor = myval / baseval;
7327 
7328  for (int k=i; k<numentries; k++) {
7329  const BigReal subval = m[INDEX(numentries,i,k)];
7330  m[INDEX(numentries,j,k)] -= (factor * subval);
7331  }
7332 
7333  b[j] -= (factor * b[i]);
7334 
7335  }
7336  }
7337 
7338  printf(" In upper diagonal form, equation is:\n");
7339  for (i=0; i<numentries; i++) {
7340  printf(" ( ");
7341  for (j=0; j<numentries; j++) {
7342  printf(" %6.3f,", m[INDEX(numentries, i, j)]);
7343  }
7344  printf(" ) ( D%-3i ) = ( %6.3f )\n", i, b[i]);
7345  }
7346 
7347  //Now work our way diagonally up from the bottom right to assign values to X
7348  for (i=numentries-1; i>=0; i--) {
7349 
7350  //Subtract the effects of previous columns
7351  for (j=i+1; j<numentries; j++) {
7352  b[i] -= ( m[INDEX(numentries,i,j)] * x[j] );
7353  }
7354 
7355  x[i] = b[i] / m[INDEX(numentries,i,i)];
7356 
7357  }
7358 
7359  printf(" Solution vector is:\n\t(");
7360  for (i=0; i<numentries; i++) {
7361  printf(" %6.3f ", x[i]);
7362  }
7363  printf(" ) \n");
7364 
7365  // We now have the coefficient information we need to make the table
7366  // Now interpolate on each interval we want
7367 
7368  distbin = 0;
7369  int entriesperseg = (int) ceil(spacing / table_spacing);
7370  int table_prefix = 2 * typeindex * (int) (namdnearbyint(maxdist / table_spacing) + 1);
7371 
7372  for (i=0; i<numentries-1; i++) {
7373  BigReal A,B,C,D;
7374  currdist = dists[i];
7375 
7376  printf("Interpolating on interval %i\n", i);
7377 
7378  // Set up the coefficients for this section
7379  A = enervalues[i];
7380  B = x[i];
7381  C = 3 * (enervalues[i+1] - enervalues[i]) - (2 * x[i]) - x[i+1];
7382  D = 2 * (enervalues[i] - enervalues[i+1]) + x[i] + x[i+1];
7383 
7384  printf("Coefficients for this interval: %f %f %f %f\n", A, B, C, D);
7385 
7386  // Go over the region of interest and fill in the table
7387  for (j=0; j<entriesperseg; j++) {
7388  const BigReal mydist = currdist + (j * table_spacing);
7389  const BigReal mydistfrac = (float) j / (entriesperseg - 1);
7390  distbin = (int) namdnearbyint(mydist / table_spacing);
7391  if (distbin >= (int) namdnearbyint(maxdist / table_spacing)) break;
7392  BigReal energy;
7393  BigReal force;
7394 
7395  energy = A + (B * mydistfrac) + (C * mydistfrac * mydistfrac) + (D * mydistfrac * mydistfrac * mydistfrac);
7396  force = B + (2 * C * mydistfrac) + (3 * D * mydistfrac * mydistfrac);
7397  // Multiply force by 1 / (interval length)
7398  force *= (1.0 / spacing);
7399 
7400  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);
7401  table_ener[table_prefix + 2 * distbin] = energy;
7402  table_ener[table_prefix + 2 * distbin + 1] = force;
7403  distbin++;
7404  }
7405  if (currdist >= maxdist) break;
7406  }
7407 
7408  //The procedure above leaves out the last entry -- add it explicitly
7409  distbin = (int) namdnearbyint(maxdist / table_spacing);
7410  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));
7411  table_ener[table_prefix + 2 * distbin] = enervalues[numentries - 1];
7412  table_ener[table_prefix + 2 * distbin + 1] = 0.0;
7413  distbin++;
7414 
7415 
7416  // Clean up and make sure everything worked ok
7417  delete m;
7418  delete x;
7419  delete b;
7420  distbin--;
7421  printf("Testing: %i vs %i (from %f / %f)\n", distbin, (int) (namdnearbyint(maxdist / table_spacing)), maxdist, table_spacing);
7422  if (distbin != (int) (namdnearbyint(maxdist / table_spacing))) return 1;
7423  return 0;
7424 } /* end read_energy_type_cubspline */
#define namdnearbyint(x)
Definition: common.h:74
const BigReal A
#define INDEX(ncols, i, j)
Definition: Parameters.C:35
void NAMD_die(const char *err_msg)
Definition: common.C:83
BigReal * table_ener
Definition: Parameters.h:254
const BigReal B
gridSize x
double BigReal
Definition: common.h:112
void Parameters::read_parameter_file ( char *  fname)

Definition at line 404 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, table_ener, and TRUE.

Referenced by Parameters().

406 {
407  char buffer[512]; // Buffer to store each line of the file
408  char first_word[512]; // First word of the current line
409  FILE *pfile; // File descriptor for the parameter file
410 
411  /* Check to make sure that we haven't previously been told */
412  /* that all the files were read */
413  if (AllFilesRead)
414  {
415  NAMD_die("Tried to read another parameter file after being told that all files were read!");
416  }
417 
418  /* Try and open the file */
419  if ( (pfile = Fopen(fname, "r")) == NULL)
420  {
421  char err_msg[256];
422 
423  sprintf(err_msg, "UNABLE TO OPEN XPLOR PARAMETER FILE %s\n", fname);
424  NAMD_die(err_msg);
425  }
426 
427  /* Keep reading in lines until we hit the EOF */
428  while (NAMD_read_line(pfile, buffer) != -1)
429  {
430  /* Get the first word of the line */
431  NAMD_find_first_word(buffer, first_word);
432 
433  /* First, screen out things that we ignore, such as */
434  /* blank lines, lines that start with '!', lines that */
435  /* start with "REMARK", lines that start with set", */
436  /* and most of the HBOND parameters which include */
437  /* AEXP, REXP, HAEX, AAEX, but not the HBOND statement */
438  /* which is parsed. */
439  if ((buffer[0] != '!') &&
440  !NAMD_blank_string(buffer) &&
441  (strncasecmp(first_word, "REMARK", 6) != 0) &&
442  (strcasecmp(first_word, "set")!=0) &&
443  (strncasecmp(first_word, "AEXP", 4) != 0) &&
444  (strncasecmp(first_word, "REXP", 4) != 0) &&
445  (strncasecmp(first_word, "HAEX", 4) != 0) &&
446  (strncasecmp(first_word, "AAEX", 4) != 0) &&
447  (strncasecmp(first_word, "NBOND", 5) != 0) &&
448  (strncasecmp(first_word, "CUTNB", 5) != 0) &&
449  (strncasecmp(first_word, "END", 3) != 0) &&
450  (strncasecmp(first_word, "CTONN", 5) != 0) &&
451  (strncasecmp(first_word, "EPS", 3) != 0) &&
452  (strncasecmp(first_word, "VSWI", 4) != 0) &&
453  (strncasecmp(first_word, "NBXM", 4) != 0) &&
454  (strncasecmp(first_word, "INHI", 4) != 0) )
455  {
456  /* Now, call the appropriate function based */
457  /* on the type of parameter we have */
458  if (strncasecmp(first_word, "bond", 4)==0)
459  {
460  add_bond_param(buffer, TRUE);
461  NumBondParams++;
462  }
463  else if (strncasecmp(first_word, "angl", 4)==0)
464  {
465  add_angle_param(buffer);
466  NumAngleParams++;
467  }
468  else if (strncasecmp(first_word, "dihe", 4)==0)
469  {
470  add_dihedral_param(buffer, pfile);
472  }
473  else if (strncasecmp(first_word, "impr", 4)==0)
474  {
475  add_improper_param(buffer, pfile);
477  }
478  else if (strncasecmp(first_word, "nonb", 4)==0)
479  {
480  add_vdw_param(buffer);
481  NumVdwParams++;
482  }
483  else if (strncasecmp(first_word, "nbfi", 4)==0)
484  {
485  add_vdw_pair_param(buffer);
486  NumVdwPairParams++;
487  }
488  else if (strncasecmp(first_word, "nbta", 4)==0)
489  {
490 
491  if (table_ener == NULL) {
492  continue;
493  }
494 
495  add_table_pair_param(buffer);
497  }
498  else if (strncasecmp(first_word, "hbon", 4)==0)
499  {
500  add_hb_pair_param(buffer);
501  }
502  else
503  {
504  /* This is an unknown paramter. */
505  /* This is BAD */
506  char err_msg[512];
507 
508  sprintf(err_msg, "UNKNOWN PARAMETER IN XPLOR PARAMETER FILE %s\nLINE=*%s*",
509  fname, buffer);
510  NAMD_die(err_msg);
511  }
512  }
513  }
514 
515  /* Close the file */
516  Fclose(pfile);
517 
518  return;
519 }
int NumTablePairParams
Definition: Parameters.h:273
int NumBondParams
Definition: Parameters.h:259
int NAMD_read_line(FILE *fd, char *buf, int bufsize)
Definition: strlib.C:38
int NumVdwParams
Definition: Parameters.h:268
int NumVdwPairParams
Definition: Parameters.h:271
void NAMD_find_first_word(char *source, char *word)
Definition: strlib.C:258
int NumAngleParams
Definition: Parameters.h:260
int NumDihedralParams
Definition: Parameters.h:262
int NAMD_blank_string(char *str)
Definition: strlib.C:222
int NumImproperParams
Definition: Parameters.h:263
FILE * Fopen(const char *filename, const char *mode)
Definition: common.C:265
void NAMD_die(const char *err_msg)
Definition: common.C:83
int Fclose(FILE *fout)
Definition: common.C:359
BigReal * table_ener
Definition: Parameters.h:254
#define TRUE
Definition: common.h:117
void Parameters::read_parm ( Ambertoppar amber_data,
BigReal  vdw14 
)

Definition at line 6366 of file Parameters.C.

References indexed_vdw_pair::A, indexed_vdw_pair::A14, angle_array, indexed_vdw_pair::B, indexed_vdw_pair::B14, bond_array, parm::Cn1, parm::Cn2, parm::Cno, parm::data_read, four_body_consts::delta, dihedral_array, parm::HB12, parm::HB6, indexed_vdw_pair::ind1, indexed_vdw_pair::ind2, iout, iWARN(), BondValue::k, AngleValue::k, four_body_consts::k, AngleValue::k_ub, indexed_vdw_pair::left, MAX_MULTIPLICITY, DihedralValue::multiplicity, four_body_consts::n, NAMD_die(), AngleValue::normal, parm::Nptra, parm::Ntypes, parm::Numang, NumAngleParams, parm::Numbnd, NumBondParams, NumDihedralParams, NumVdwPairParams, NumVdwParamsAssigned, parm::Phase, parm::Pk, parm::Pn, AngleValue::r_ub, parm::Req, indexed_vdw_pair::right, parm::Rk, parm::Teq, AngleValue::theta0, parm::Tk, DihedralValue::values, values, vdw_pair_tree, and BondValue::x0.

6367 {
6368  int i,j,ico,numtype,mul;
6369  IndexedVdwPair *new_node;
6370 
6371  if (!amber_data->data_read)
6372  NAMD_die("No data read from parm file yet!");
6373 
6374  // Copy bond parameters
6375  NumBondParams = amber_data->Numbnd;
6376  if (NumBondParams)
6378  if (bond_array == NULL)
6379  NAMD_die("memory allocation of bond_array failed!");
6380  }
6381  for (i=0; i<NumBondParams; ++i)
6382  { bond_array[i].k = amber_data->Rk[i];
6383  bond_array[i].x0 = amber_data->Req[i];
6384  }
6385 
6386  // Copy angle parameters
6387  NumAngleParams = amber_data->Numang;
6388  if (NumAngleParams)
6390  if (angle_array == NULL)
6391  NAMD_die("memory allocation of angle_array failed!");
6392  }
6393  for (i=0; i<NumAngleParams; ++i)
6394  { angle_array[i].k = amber_data->Tk[i];
6395  angle_array[i].theta0 = amber_data->Teq[i];
6396  // Amber has no k_ub and r_ub for angle parameters, so they're set to 0
6397  angle_array[i].k_ub = angle_array[i].r_ub = 0;
6398  // All angles are harmonic
6399  angle_array[i].normal = 1;
6400  }
6401 
6402  // Copy dihedral parameters
6403  // Note: If the periodicity is negative, it means the following
6404  // entry is another term in a multitermed dihedral, and the
6405  // absolute value is the true periodicity; in this case the
6406  // following entry in "dihedral_array" should be left empty,
6407  // NOT be skipped, in order to keep the next dihedral's index
6408  // unchanged.
6409  NumDihedralParams = amber_data->Nptra;
6410  if (NumDihedralParams)
6411  { dihedral_array = new DihedralValue[amber_data->Nptra];
6412  if (dihedral_array == NULL)
6413  NAMD_die("memory allocation of dihedral_array failed!");
6414  }
6415  mul = 0;
6416  for (i=0; i<NumDihedralParams; ++i)
6417  { dihedral_array[i-mul].values[mul].k = amber_data->Pk[i];
6418  dihedral_array[i-mul].values[mul].n = int(fabs(amber_data->Pn[i])+0.5);
6419  if (dihedral_array[i-mul].values[mul].n == 0)
6420  { char err_msg[128];
6421  sprintf(err_msg, "The periodicity of dihedral # %d is zero!", i+1);
6422  NAMD_die(err_msg);
6423  }
6424  dihedral_array[i-mul].values[mul].delta = amber_data->Phase[i];
6425  // If the periodicity is positive, it means the following
6426  // entry is a new dihedral term.
6427  if (amber_data->Pn[i] > 0)
6428  { dihedral_array[i-mul].multiplicity = mul+1;
6429  mul = 0;
6430  }
6431  else if (++mul >= MAX_MULTIPLICITY)
6432  { char err_msg[181];
6433  sprintf(err_msg, "Multiple dihedral with multiplicity of %d greater than max of %d",
6434  mul+1, MAX_MULTIPLICITY);
6435  NAMD_die(err_msg);
6436  }
6437  }
6438  if (mul > 0)
6439  NAMD_die("Negative periodicity in the last dihedral entry!");
6440 
6441  // Copy VDW parameters: AMBER explicitly gives VDW parameters between every
6442  // 2 atom types, so the data are copied to vdw_pair_tree
6443  // In AMBER, all 1-4 VDW interactions are divided by factor vdw14
6444  NumVdwParamsAssigned = numtype = amber_data->Ntypes; // Number of atom types
6445  NumVdwPairParams = numtype * (numtype+1) / 2;
6446  for (i=0; i<numtype; ++i)
6447  for (j=i; j<numtype; ++j)
6448  { new_node = new IndexedVdwPair;
6449  if (new_node == NULL)
6450  NAMD_die("memory allocation of vdw_pair_tree failed!");
6451  new_node->ind1 = i;
6452  new_node->ind2 = j;
6453  new_node->left = new_node->right = NULL;
6454  // ico is the index of interaction between atom type i and j into
6455  // the parameter arrays. If it's positive, the interaction is
6456  // 6-12 VDW; otherwise it's 10-12 H-bond interaction. NAMD doesn't
6457  // have 10-12 term, so if ico is negative, then the 10-12
6458  // coefficients must be 0, otherwise die.
6459  ico = amber_data->Cno[numtype*i+j];
6460  if (ico>0)
6461  { new_node->A = amber_data->Cn1[ico-1];
6462  new_node->A14 = new_node->A / vdw14;
6463  new_node->B = amber_data->Cn2[ico-1];
6464  new_node->B14 = new_node->B / vdw14;
6465  }
6466  else if (amber_data->HB12[abs(ico)-1]==0.0 && amber_data->HB6[abs(ico)-1]==0.0)
6467  { new_node->A = new_node->A14 = new_node->B = new_node->B14 = 0.0;
6468  iout << iWARN << "Encounter 10-12 H-bond term\n";
6469  }
6470  else
6471  NAMD_die("Encounter non-zero 10-12 H-bond term!");
6472  // Add the node to the binary tree
6473  vdw_pair_tree = add_to_indexed_vdw_pairs(new_node, vdw_pair_tree);
6474  }
6475 }
_REAL * Pk
Definition: parm.h:24
_REAL * Teq
Definition: parm.h:24
int NumBondParams
Definition: Parameters.h:259
_REAL * Phase
Definition: parm.h:24
struct indexed_vdw_pair * right
Definition: Parameters.h:168
Real r_ub
Definition: Parameters.h:94
int NumVdwPairParams
Definition: Parameters.h:271
std::ostream & iWARN(std::ostream &s)
Definition: InfoStream.C:108
DihedralValue * dihedral_array
Definition: Parameters.h:243
_REAL * Tk
Definition: parm.h:24
#define iout
Definition: InfoStream.h:87
int NumAngleParams
Definition: Parameters.h:260
IndexedVdwPair * vdw_pair_tree
Definition: Parameters.h:255
int NumDihedralParams
Definition: Parameters.h:262
int Numbnd
Definition: parm.h:17
FourBodyConsts values[MAX_MULTIPLICITY]
Definition: Parameters.h:108
AngleValue * angle_array
Definition: Parameters.h:242
_REAL * HB6
Definition: parm.h:24
_REAL * HB12
Definition: parm.h:24
struct indexed_vdw_pair * left
Definition: Parameters.h:169
#define MAX_MULTIPLICITY
Definition: Parameters.h:68
int Ntypes
Definition: parm.h:17
int normal
Definition: Parameters.h:95
_REAL * Cn2
Definition: parm.h:24
struct indexed_vdw_pair IndexedVdwPair
void NAMD_die(const char *err_msg)
Definition: common.C:83
int NumVdwParamsAssigned
Definition: Parameters.h:270
_REAL * Rk
Definition: parm.h:24
int data_read
Definition: parm.h:34
_REAL * Cn1
Definition: parm.h:24
Real x0
Definition: Parameters.h:85
BondValue * bond_array
Definition: Parameters.h:241
int Numang
Definition: parm.h:17
_REAL * Req
Definition: parm.h:24
Real theta0
Definition: Parameters.h:92
valT values[SORTTILELISTSKERNEL_ITEMS_PER_THREAD]
Real k_ub
Definition: Parameters.h:93
int Nptra
Definition: parm.h:17
Real k
Definition: Parameters.h:84
int * Cno
Definition: parm.h:27
_REAL * Pn
Definition: parm.h:24
void Parameters::receive_Parameters ( MIStream msg)

Definition at line 5424 of file Parameters.C.

References indexed_vdw_pair::A, indexed_vdw_pair::A14, nbthole_pair_value::alphai, nbthole_pair_value::alphaj, angle_array, indexed_vdw_pair::B, indexed_vdw_pair::B14, bond_array, crossterm_array, four_body_consts::delta, dihedral_array, CrosstermValue::dim, MIStream::get(), gromacsPair_array, improper_array, indexed_vdw_pair::ind1, nbthole_pair_value::ind1, indexed_table_pair::ind1, indexed_vdw_pair::ind2, nbthole_pair_value::ind2, indexed_table_pair::ind2, BondValue::k, AngleValue::k, four_body_consts::k, indexed_table_pair::K, AngleValue::k_ub, indexed_vdw_pair::left, indexed_table_pair::left, MAX_ATOMTYPE_CHARS, MAX_MULTIPLICITY, DihedralValue::multiplicity, ImproperValue::multiplicity, four_body_consts::n, NAMD_die(), nbthole_array, AngleValue::normal, NumAngleParams, NumBondParams, NumCrosstermParams, NumDihedralParams, numenerentries, NumGromacsPairParams, NumImproperParams, NumNbtholePairParams, NumTablePairParams, NumVdwPairParams, NumVdwParams, NumVdwParamsAssigned, GromacsPairValue::pairC12, GromacsPairValue::pairC6, AngleValue::r_ub, indexed_vdw_pair::right, indexed_table_pair::right, tab_pair_tree, table_ener, AngleValue::theta0, nbthole_pair_value::tholeij, TRUE, DihedralValue::values, ImproperValue::values, vdw_array, vdw_pair_tree, and BondValue::x0.

Referenced by Node::resendMolecule().

5426 {
5427  int i, j; // Loop counters
5428  Real *a1, *a2, *a3, *a4; // Temporary arrays to get data from message in
5429  int *i1, *i2, *i3; // Temporary int array to get data from message in
5430  IndexedVdwPair *new_node; // New node for vdw pair param tree
5431  IndexedTablePair *new_tab_node;
5432  Real **kvals; // Force constant values for dihedrals and impropers
5433  int **nvals; // Periodicity values for dihedrals and impropers
5434  Real **deltavals; // Phase shift values for dihedrals and impropers
5435  BigReal *pairC6, *pairC12; // JLai
5436 
5437  // Get the bonded parameters
5438  msg->get(NumBondParams);
5439 
5440  if (NumBondParams)
5441  {
5443  a1 = new Real[NumBondParams];
5444  a2 = new Real[NumBondParams];
5445 
5446  if ( (bond_array == NULL) || (a1 == NULL) || (a2 == NULL) )
5447  {
5448  NAMD_die("memory allocation failed in Parameters::receive_Parameters");
5449  }
5450 
5451  msg->get(NumBondParams, a1);
5452  msg->get(NumBondParams, a2);
5453 
5454  for (i=0; i<NumBondParams; i++)
5455  {
5456  bond_array[i].k = a1[i];
5457  bond_array[i].x0 = a2[i];
5458  }
5459 
5460  delete [] a1;
5461  delete [] a2;
5462  }
5463 
5464  // Get the angle parameters
5465  msg->get(NumAngleParams);
5466 
5467  if (NumAngleParams)
5468  {
5470  a1 = new Real[NumAngleParams];
5471  a2 = new Real[NumAngleParams];
5472  a3 = new Real[NumAngleParams];
5473  a4 = new Real[NumAngleParams];
5474  i1 = new int[NumAngleParams];
5475 
5476  if ( (angle_array == NULL) || (a1 == NULL) || (a2 == NULL) ||
5477  (a3 == NULL) || (a4 == NULL) || (i1 == NULL))
5478  {
5479  NAMD_die("memory allocation failed in Parameters::receive_Parameters");
5480  }
5481 
5482  msg->get(NumAngleParams, a1);
5483  msg->get(NumAngleParams, a2);
5484  msg->get(NumAngleParams, a3);
5485  msg->get(NumAngleParams, a4);
5486  msg->get(NumAngleParams, i1);
5487 
5488  for (i=0; i<NumAngleParams; i++)
5489  {
5490  angle_array[i].k = a1[i];
5491  angle_array[i].theta0 = a2[i];
5492  angle_array[i].k_ub = a3[i];
5493  angle_array[i].r_ub = a4[i];
5494  angle_array[i].normal = i1[i];
5495  }
5496 
5497  delete [] a1;
5498  delete [] a2;
5499  delete [] a3;
5500  delete [] a4;
5501  delete [] i1;
5502  }
5503 
5504  // Get the dihedral parameters
5505  msg->get(NumDihedralParams);
5506 
5507  if (NumDihedralParams)
5508  {
5510 
5511  i1 = new int[NumDihedralParams];
5512  kvals = new Real *[MAX_MULTIPLICITY];
5513  nvals = new int *[MAX_MULTIPLICITY];
5514  deltavals = new Real *[MAX_MULTIPLICITY];
5515 
5516  if ( (i1 == NULL) || (kvals == NULL) || (nvals == NULL) ||
5517  (deltavals == NULL) || (dihedral_array == NULL) )
5518  {
5519  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5520  }
5521 
5522  for (i=0; i<MAX_MULTIPLICITY; i++)
5523  {
5524  kvals[i] = new Real[NumDihedralParams];
5525  nvals[i] = new int[NumDihedralParams];
5526  deltavals[i] = new Real[NumDihedralParams];
5527 
5528  if ( (kvals[i] == NULL) || (nvals[i] == NULL) || (deltavals[i] == NULL) )
5529  {
5530  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5531  }
5532  }
5533 
5534  msg->get(NumDihedralParams, i1);
5535 
5536  for (i=0; i<MAX_MULTIPLICITY; i++)
5537  {
5538  msg->get(NumDihedralParams, kvals[i]);
5539  msg->get(NumDihedralParams, nvals[i]);
5540  msg->get(NumDihedralParams, deltavals[i]);
5541  }
5542 
5543  for (i=0; i<NumDihedralParams; i++)
5544  {
5545  dihedral_array[i].multiplicity = i1[i];
5546 
5547  for (j=0; j<MAX_MULTIPLICITY; j++)
5548  {
5549  dihedral_array[i].values[j].k = kvals[j][i];
5550  dihedral_array[i].values[j].n = nvals[j][i];
5551  dihedral_array[i].values[j].delta = deltavals[j][i];
5552  }
5553  }
5554 
5555  for (i=0; i<MAX_MULTIPLICITY; i++)
5556  {
5557  delete [] kvals[i];
5558  delete [] nvals[i];
5559  delete [] deltavals[i];
5560  }
5561 
5562  delete [] i1;
5563  delete [] kvals;
5564  delete [] nvals;
5565  delete [] deltavals;
5566  }
5567 
5568  // Get the improper parameters
5569  msg->get(NumImproperParams);
5570 
5571  if (NumImproperParams)
5572  {
5574  i1 = new int[NumImproperParams];
5575  kvals = new Real *[MAX_MULTIPLICITY];
5576  nvals = new int *[MAX_MULTIPLICITY];
5577  deltavals = new Real *[MAX_MULTIPLICITY];
5578 
5579  if ( (i1 == NULL) || (kvals == NULL) || (nvals == NULL) ||
5580  (deltavals == NULL) || (improper_array==NULL) )
5581  {
5582  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5583  }
5584 
5585  for (i=0; i<MAX_MULTIPLICITY; i++)
5586  {
5587  kvals[i] = new Real[NumImproperParams];
5588  nvals[i] = new int[NumImproperParams];
5589  deltavals[i] = new Real[NumImproperParams];
5590 
5591  if ( (kvals[i] == NULL) || (nvals[i] == NULL) || (deltavals[i] == NULL) )
5592  {
5593  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5594  }
5595  }
5596 
5597  msg->get(NumImproperParams,i1);
5598 
5599  for (i=0; i<MAX_MULTIPLICITY; i++)
5600  {
5601  msg->get(NumImproperParams,kvals[i]);
5602  msg->get(NumImproperParams,nvals[i]);
5603  msg->get(NumImproperParams,deltavals[i]);
5604  }
5605 
5606  for (i=0; i<NumImproperParams; i++)
5607  {
5608  improper_array[i].multiplicity = i1[i];
5609 
5610  for (j=0; j<MAX_MULTIPLICITY; j++)
5611  {
5612  improper_array[i].values[j].k = kvals[j][i];
5613  improper_array[i].values[j].n = nvals[j][i];
5614  improper_array[i].values[j].delta = deltavals[j][i];
5615  }
5616  }
5617 
5618  for (i=0; i<MAX_MULTIPLICITY; i++)
5619  {
5620  delete [] kvals[i];
5621  delete [] nvals[i];
5622  delete [] deltavals[i];
5623  }
5624 
5625  delete [] i1;
5626  delete [] kvals;
5627  delete [] nvals;
5628  delete [] deltavals;
5629  }
5630 
5631  // Get the crossterm parameters
5632  msg->get(NumCrosstermParams);
5633 
5634  if (NumCrosstermParams)
5635  {
5637 
5638  for (i=0; i<NumCrosstermParams; ++i) {
5639  int nvals = CrosstermValue::dim * CrosstermValue::dim * 2 * 2;
5640  msg->get(nvals,&crossterm_array[i].c[0][0].d00);
5641  }
5642  }
5643 
5644  // Get GromacsPairs parameters
5645  // JLai
5646  msg->get(NumGromacsPairParams);
5647 
5649  {
5651  pairC6 = new BigReal[NumGromacsPairParams];
5652  pairC12 = new BigReal[NumGromacsPairParams];
5653 
5654  if ( (pairC6 == NULL) || (pairC12 == NULL) ) {
5655  NAMD_die("memory allocation failed in Parameters::receive_Parameters");
5656  }
5657 
5658  msg->get(NumGromacsPairParams,pairC6);
5659  msg->get(NumGromacsPairParams,pairC12);
5660 
5661  for (i=0; i<NumGromacsPairParams; ++i) {
5662  gromacsPair_array[i].pairC6 = pairC6[i];
5663  gromacsPair_array[i].pairC12 = pairC12[i];
5664  }
5665 
5666  delete [] pairC6;
5667  delete [] pairC12;
5668  }
5669  // JLai
5670 
5671  //Get the energy table
5672  msg->get(numenerentries);
5673  if (numenerentries > 0) {
5674  //fprintf(stdout, "Getting tables\n");
5675  //fprintf(infofile, "Trying to allocate table\n");
5677  //fprintf(infofile, "Finished table allocation\n");
5678  if (table_ener==NULL)
5679  {
5680  NAMD_die("memory allocation failed in Parameters::receive_Parameters");
5681  }
5682 
5683  msg->get(numenerentries, table_ener);
5684  }
5685 
5686  // Get the vdw parameters
5687  msg->get(NumVdwParams);
5688  msg->get(NumVdwParamsAssigned);
5689 
5690  if (NumVdwParams)
5691  {
5692  atomTypeNames = new char[NumVdwParams*(MAX_ATOMTYPE_CHARS+1)];
5694  a1 = new Real[NumVdwParams];
5695  a2 = new Real[NumVdwParams];
5696  a3 = new Real[NumVdwParams];
5697  a4 = new Real[NumVdwParams];
5698 
5699  if ( (vdw_array==NULL) || (a1==NULL) || (a2==NULL) || (a3==NULL)
5700  || (a4==NULL) || (atomTypeNames==NULL))
5701  {
5702  NAMD_die("memory allocation failed in Parameters::receive_Parameters");
5703  }
5704 
5705  msg->get(NumVdwParams * (MAX_ATOMTYPE_CHARS+1), atomTypeNames);
5706  msg->get(NumVdwParams, a1);
5707  msg->get(NumVdwParams, a2);
5708  msg->get(NumVdwParams, a3);
5709  msg->get(NumVdwParams, a4);
5710 
5711  for (i=0; i<NumVdwParams; i++)
5712  {
5713  vdw_array[i].sigma = a1[i];
5714  vdw_array[i].epsilon = a2[i];
5715  vdw_array[i].sigma14 = a3[i];
5716  vdw_array[i].epsilon14 = a4[i];
5717  }
5718 
5719  delete [] a1;
5720  delete [] a2;
5721  delete [] a3;
5722  delete [] a4;
5723  }
5724 
5725  // Get the vdw_pair_parameters
5726  msg->get(NumVdwPairParams);
5727 
5728  if (NumVdwPairParams)
5729  {
5730  a1 = new Real[NumVdwPairParams];
5731  a2 = new Real[NumVdwPairParams];
5732  a3 = new Real[NumVdwPairParams];
5733  a4 = new Real[NumVdwPairParams];
5734  i1 = new int[NumVdwPairParams];
5735  i2 = new int[NumVdwPairParams];
5736 
5737  if ( (a1 == NULL) || (a2 == NULL) || (a3 == NULL) || (a4 == NULL) ||
5738  (i1 == NULL) || (i2 == NULL) )
5739  {
5740  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5741  }
5742 
5743  msg->get(NumVdwPairParams, i1);
5744  msg->get(NumVdwPairParams, i2);
5745  msg->get(NumVdwPairParams, a1);
5746  msg->get(NumVdwPairParams, a2);
5747  msg->get(NumVdwPairParams, a3);
5748  msg->get(NumVdwPairParams, a4);
5749 
5750  for (i=0; i<NumVdwPairParams; i++)
5751  {
5752  new_node = (IndexedVdwPair *) malloc(sizeof(IndexedVdwPair));
5753 
5754  if (new_node == NULL)
5755  {
5756  NAMD_die("memory allocation failed in Parameters::receive_Parameters");
5757  }
5758 
5759  new_node->ind1 = i1[i];
5760  new_node->ind2 = i2[i];
5761  new_node->A = a1[i];
5762  new_node->A14 = a2[i];
5763  new_node->B = a3[i];
5764  new_node->B14 = a4[i];
5765  new_node->left = NULL;
5766  new_node->right = NULL;
5767 
5768  vdw_pair_tree = add_to_indexed_vdw_pairs(new_node, vdw_pair_tree);
5769  }
5770 
5771  delete [] i1;
5772  delete [] i2;
5773  delete [] a1;
5774  delete [] a2;
5775  delete [] a3;
5776  delete [] a4;
5777  }
5778 
5779  // Get the nbthole_pair_parameters
5780  msg->get(NumNbtholePairParams);
5781 
5783  {
5785  a1 = new Real[NumNbtholePairParams];
5786  a2 = new Real[NumNbtholePairParams];
5787  a3 = new Real[NumNbtholePairParams];
5788  i1 = new int[NumNbtholePairParams];
5789  i2 = new int[NumNbtholePairParams];
5790 
5791  if ( (nbthole_array == NULL) || (a1 == NULL) || (a2 == NULL) || (a3 == NULL)
5792  || (i1 == NULL) || (i2 == NULL) )
5793  {
5794  NAMD_die("memory allocation failed in Parameters::receive_Parameters");
5795  }
5796 
5797  msg->get(NumNbtholePairParams, i1);
5798  msg->get(NumNbtholePairParams, i2);
5799  msg->get(NumNbtholePairParams, a1);
5800  msg->get(NumNbtholePairParams, a2);
5801  msg->get(NumNbtholePairParams, a3);
5802 
5803  for (i=0; i<NumNbtholePairParams; i++)
5804  {
5805 
5806  nbthole_array[i].ind1 = i1[i];
5807  nbthole_array[i].ind2 = i2[i];
5808  nbthole_array[i].alphai = a1[i];
5809  nbthole_array[i].alphaj = a2[i];
5810  nbthole_array[i].tholeij = a3[i];
5811 
5812  }
5813 
5814  delete [] i1;
5815  delete [] i2;
5816  delete [] a1;
5817  delete [] a2;
5818  delete [] a3;
5819  }
5820 
5821  // Get the table_pair_parameters
5822  msg->get(NumTablePairParams);
5823 
5824  if (NumTablePairParams)
5825  {
5826  i1 = new int[NumTablePairParams];
5827  i2 = new int[NumTablePairParams];
5828  i3 = new int[NumTablePairParams];
5829 
5830  if ( (i3 == NULL) || (i1 == NULL) || (i2 == NULL) )
5831  {
5832  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5833  }
5834 
5835  msg->get(NumTablePairParams, i1);
5836  msg->get(NumTablePairParams, i2);
5837  msg->get(NumTablePairParams, i3);
5838 
5839  for (i=0; i<NumTablePairParams; i++)
5840  {
5841  new_tab_node = (IndexedTablePair *) malloc(sizeof(IndexedTablePair));
5842 
5843  if (new_tab_node == NULL)
5844  {
5845  NAMD_die("memory allocation failed in Parameters::receive_Parameters");
5846  }
5847 
5848 // printf("Adding new table pair with ind1 %i ind2 %i k %i\n", i1[i], i2[i],i3[i]);
5849  new_tab_node->ind1 = i1[i];
5850  new_tab_node->ind2 = i2[i];
5851  new_tab_node->K = i3[i];
5852  new_tab_node->left = NULL;
5853  new_tab_node->right = NULL;
5854 
5855  tab_pair_tree = add_to_indexed_table_pairs(new_tab_node, tab_pair_tree);
5856  }
5857 
5858  delete [] i1;
5859  delete [] i2;
5860  delete [] i3;
5861  }
5862 
5863  // receive the hydrogen bond parameters
5864  // hbondParams.receive_message(msg);
5865 
5866  AllFilesRead = TRUE;
5867 
5868  delete msg;
5869 }
int numenerentries
Definition: Parameters.h:251
int NumTablePairParams
Definition: Parameters.h:273
GromacsPairValue * gromacsPair_array
Definition: Parameters.h:247
int NumBondParams
Definition: Parameters.h:259
struct indexed_vdw_pair * right
Definition: Parameters.h:168
struct indexed_table_pair * right
Definition: Parameters.h:201
int NumVdwParams
Definition: Parameters.h:268
Real r_ub
Definition: Parameters.h:94
float Real
Definition: common.h:107
int NumVdwPairParams
Definition: Parameters.h:271
MIStream * get(char &data)
Definition: MStream.h:29
DihedralValue * dihedral_array
Definition: Parameters.h:243
int NumAngleParams
Definition: Parameters.h:260
CrosstermValue * crossterm_array
Definition: Parameters.h:245
IndexedVdwPair * vdw_pair_tree
Definition: Parameters.h:255
int NumDihedralParams
Definition: Parameters.h:262
FourBodyConsts values[MAX_MULTIPLICITY]
Definition: Parameters.h:108
AngleValue * angle_array
Definition: Parameters.h:242
struct indexed_vdw_pair * left
Definition: Parameters.h:169
int NumNbtholePairParams
Definition: Parameters.h:272
#define MAX_MULTIPLICITY
Definition: Parameters.h:68
int normal
Definition: Parameters.h:95
int NumImproperParams
Definition: Parameters.h:263
struct indexed_table_pair * left
Definition: Parameters.h:202
int NumCrosstermParams
Definition: Parameters.h:264
ImproperValue * improper_array
Definition: Parameters.h:244
void NAMD_die(const char *err_msg)
Definition: common.C:83
int NumVdwParamsAssigned
Definition: Parameters.h:270
NbtholePairValue * nbthole_array
Definition: Parameters.h:250
FourBodyConsts values[MAX_MULTIPLICITY]
Definition: Parameters.h:114
VdwValue * vdw_array
Definition: Parameters.h:249
BigReal * table_ener
Definition: Parameters.h:254
Real x0
Definition: Parameters.h:85
BondValue * bond_array
Definition: Parameters.h:241
int NumGromacsPairParams
Definition: Parameters.h:266
Real theta0
Definition: Parameters.h:92
IndexedTablePair * tab_pair_tree
Definition: Parameters.h:257
#define TRUE
Definition: common.h:117
Real k_ub
Definition: Parameters.h:93
Real k
Definition: Parameters.h:84
#define MAX_ATOMTYPE_CHARS
Definition: Parameters.h:71
double BigReal
Definition: common.h:112
void Parameters::send_Parameters ( MOStream msg)

Definition at line 5048 of file Parameters.C.

References nbthole_pair_value::alphai, nbthole_pair_value::alphaj, angle_array, bond_array, crossterm_array, four_body_consts::delta, dihedral_array, CrosstermValue::dim, MOStream::end(), vdw_val::epsilon, vdw_val::epsilon14, gromacsPair_array, improper_array, nbthole_pair_value::ind1, nbthole_pair_value::ind2, BondValue::k, AngleValue::k, four_body_consts::k, AngleValue::k_ub, MAX_ATOMTYPE_CHARS, MAX_MULTIPLICITY, DihedralValue::multiplicity, ImproperValue::multiplicity, four_body_consts::n, NAMD_die(), nbthole_array, nbthole_pair_tree, AngleValue::normal, NumAngleParams, NumBondParams, NumCrosstermParams, NumDihedralParams, numenerentries, NumGromacsPairParams, NumImproperParams, NumNbtholePairParams, NumTablePairParams, NumVdwPairParams, NumVdwParams, NumVdwParamsAssigned, GromacsPairValue::pairC12, GromacsPairValue::pairC6, MOStream::put(), AngleValue::r_ub, vdw_val::sigma, vdw_val::sigma14, tab_pair_tree, table_ener, AngleValue::theta0, nbthole_pair_value::tholeij, DihedralValue::values, ImproperValue::values, vdw_array, vdw_pair_tree, and BondValue::x0.

Referenced by Node::resendMolecule().

5049 {
5050  Real *a1, *a2, *a3, *a4; // Temporary arrays for sending messages
5051  int *i1, *i2, *i3; // Temporary int array
5052  int i, j; // Loop counters
5053  Real **kvals; // Force constant values for dihedrals and impropers
5054  int **nvals; // Periodicity values for dihedrals and impropers
5055  Real **deltavals; // Phase shift values for dihedrals and impropers
5056  BigReal *pairC6, *pairC12; // JLai
5057  /*MOStream *msg=comm_obj->newOutputStream(ALLBUTME, STATICPARAMSTAG, BUFSIZE);
5058  if ( msg == NULL )
5059  {
5060  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5061  }*/
5062 
5063  // Send the bond parameters
5064  msg->put(NumBondParams);
5065 
5066  if (NumBondParams)
5067  {
5068  a1 = new Real[NumBondParams];
5069  a2 = new Real[NumBondParams];
5070 
5071  if ( (a1 == NULL) || (a2 == NULL) )
5072  {
5073  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5074  }
5075 
5076  for (i=0; i<NumBondParams; i++)
5077  {
5078  a1[i] = bond_array[i].k;
5079  a2[i] = bond_array[i].x0;
5080  }
5081 
5082  msg->put(NumBondParams, a1)->put(NumBondParams, a2);
5083 
5084  delete [] a1;
5085  delete [] a2;
5086  }
5087 
5088  // Send the angle parameters
5089  msg->put(NumAngleParams);
5090 
5091  if (NumAngleParams)
5092  {
5093  a1 = new Real[NumAngleParams];
5094  a2 = new Real[NumAngleParams];
5095  a3 = new Real[NumAngleParams];
5096  a4 = new Real[NumAngleParams];
5097  i1 = new int[NumAngleParams];
5098 
5099  if ( (a1 == NULL) || (a2 == NULL) || (a3 == NULL) ||
5100  (a4 == NULL) || (i1 == NULL))
5101  {
5102  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5103  }
5104 
5105  for (i=0; i<NumAngleParams; i++)
5106  {
5107  a1[i] = angle_array[i].k;
5108  a2[i] = angle_array[i].theta0;
5109  a3[i] = angle_array[i].k_ub;
5110  a4[i] = angle_array[i].r_ub;
5111  i1[i] = angle_array[i].normal;
5112  }
5113 
5114  msg->put(NumAngleParams, a1)->put(NumAngleParams, a2);
5115  msg->put(NumAngleParams, a3)->put(NumAngleParams, a4);
5116  msg->put(NumAngleParams, i1);
5117 
5118  delete [] a1;
5119  delete [] a2;
5120  delete [] a3;
5121  delete [] a4;
5122  delete [] i1;
5123  }
5124 
5125  // Send the dihedral parameters
5126  msg->put(NumDihedralParams);
5127 
5128  if (NumDihedralParams)
5129  {
5130  i1 = new int[NumDihedralParams];
5131  kvals = new Real *[MAX_MULTIPLICITY];
5132  nvals = new int *[MAX_MULTIPLICITY];
5133  deltavals = new Real *[MAX_MULTIPLICITY];
5134 
5135  if ( (i1 == NULL) || (kvals == NULL) || (nvals == NULL) ||
5136  (deltavals == NULL) )
5137  {
5138  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5139  }
5140 
5141  for (i=0; i<MAX_MULTIPLICITY; i++)
5142  {
5143  kvals[i] = new Real[NumDihedralParams];
5144  nvals[i] = new int[NumDihedralParams];
5145  deltavals[i] = new Real[NumDihedralParams];
5146 
5147  if ( (kvals[i] == NULL) || (nvals[i] == NULL) || (deltavals[i] == NULL) )
5148  {
5149  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5150  }
5151  }
5152 
5153  for (i=0; i<NumDihedralParams; i++)
5154  {
5155  i1[i] = dihedral_array[i].multiplicity;
5156 
5157  for (j=0; j<MAX_MULTIPLICITY; j++)
5158  {
5159  kvals[j][i] = dihedral_array[i].values[j].k;
5160  nvals[j][i] = dihedral_array[i].values[j].n;
5161  deltavals[j][i] = dihedral_array[i].values[j].delta;
5162  }
5163  }
5164 
5165  msg->put(NumDihedralParams, i1);
5166 
5167  for (i=0; i<MAX_MULTIPLICITY; i++)
5168  {
5169  msg->put(NumDihedralParams, kvals[i]);
5170  msg->put(NumDihedralParams, nvals[i]);
5171  msg->put(NumDihedralParams, deltavals[i]);
5172 
5173  delete [] kvals[i];
5174  delete [] nvals[i];
5175  delete [] deltavals[i];
5176  }
5177 
5178  delete [] i1;
5179  delete [] kvals;
5180  delete [] nvals;
5181  delete [] deltavals;
5182  }
5183 
5184  // Send the improper parameters
5185  msg->put(NumImproperParams);
5186 
5187  if (NumImproperParams)
5188  {
5189  i1 = new int[NumImproperParams];
5190  kvals = new Real *[MAX_MULTIPLICITY];
5191  nvals = new int *[MAX_MULTIPLICITY];
5192  deltavals = new Real *[MAX_MULTIPLICITY];
5193 
5194  if ( (i1 == NULL) || (kvals == NULL) || (nvals == NULL) ||
5195  (deltavals == NULL) )
5196  {
5197  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5198  }
5199 
5200  for (i=0; i<MAX_MULTIPLICITY; i++)
5201  {
5202  kvals[i] = new Real[NumImproperParams];
5203  nvals[i] = new int[NumImproperParams];
5204  deltavals[i] = new Real[NumImproperParams];
5205 
5206  if ( (kvals[i] == NULL) || (nvals[i] == NULL) || (deltavals[i] == NULL) )
5207  {
5208  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5209  }
5210  }
5211 
5212  for (i=0; i<NumImproperParams; i++)
5213  {
5214  i1[i] = improper_array[i].multiplicity;
5215 
5216  for (j=0; j<MAX_MULTIPLICITY; j++)
5217  {
5218  kvals[j][i] = improper_array[i].values[j].k;
5219  nvals[j][i] = improper_array[i].values[j].n;
5220  deltavals[j][i] = improper_array[i].values[j].delta;
5221  }
5222  }
5223 
5224  msg->put(NumImproperParams, i1);
5225 
5226  for (i=0; i<MAX_MULTIPLICITY; i++)
5227  {
5228  msg->put(NumImproperParams, kvals[i]);
5229  msg->put(NumImproperParams, nvals[i]);
5230  msg->put(NumImproperParams, deltavals[i]);
5231 
5232  delete [] kvals[i];
5233  delete [] nvals[i];
5234  delete [] deltavals[i];
5235  }
5236 
5237  delete [] i1;
5238  delete [] kvals;
5239  delete [] nvals;
5240  delete [] deltavals;
5241  }
5242 
5243  // Send the crossterm parameters
5244  msg->put(NumCrosstermParams);
5245 
5246  if (NumCrosstermParams)
5247  {
5248  for (i=0; i<NumCrosstermParams; ++i) {
5249  int nvals = CrosstermValue::dim * CrosstermValue::dim * 2 * 2;
5250  msg->put(nvals,&crossterm_array[i].c[0][0].d00);
5251  }
5252  }
5253  // Send the GromacsPairs parameters
5254  // JLai
5255  msg->put(NumGromacsPairParams);
5256 
5258  {
5259  pairC6 = new BigReal[NumGromacsPairParams];
5260  pairC12 = new BigReal[NumGromacsPairParams];
5261  if ( (pairC6 == NULL) || (pairC12 == NULL) ) {
5262  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5263  }
5264 
5265  for (i=0; i<NumGromacsPairParams; i++) {
5266  pairC6[i] = gromacsPair_array[i].pairC6;
5267  pairC12[i] = gromacsPair_array[i].pairC12;
5268  }
5269 
5270  msg->put(NumGromacsPairParams,pairC6);
5271  msg->put(NumGromacsPairParams,pairC12);
5272 
5273  delete [] pairC6;
5274  delete [] pairC12;
5275  }
5276  // End of JLai
5277 
5278  //
5279  //Send the energy table parameters
5280  msg->put(numenerentries);
5281 
5282  if (numenerentries) {
5283  /*
5284  b1 = new Real[numenerentries];
5285  if (b1 == NULL)
5286  {
5287  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5288  }
5289  */
5290 
5291  msg->put(numenerentries, table_ener);
5292  }
5293 
5294  // Send the vdw parameters
5295  msg->put(NumVdwParams);
5296  msg->put(NumVdwParamsAssigned);
5297 
5298  if (NumVdwParams)
5299  {
5300  a1 = new Real[NumVdwParams];
5301  a2 = new Real[NumVdwParams];
5302  a3 = new Real[NumVdwParams];
5303  a4 = new Real[NumVdwParams];
5304 
5305  if ( (a1 == NULL) || (a2 == NULL) || (a3 == NULL) )
5306  {
5307  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5308  }
5309 
5310  for (i=0; i<NumVdwParams; i++)
5311  {
5312  a1[i] = vdw_array[i].sigma;
5313  a2[i] = vdw_array[i].epsilon;
5314  a3[i] = vdw_array[i].sigma14;
5315  a4[i] = vdw_array[i].epsilon14;
5316  }
5317 
5318  msg->put(NumVdwParams * (MAX_ATOMTYPE_CHARS+1), atomTypeNames);
5319  msg->put(NumVdwParams, a1);
5320  msg->put(NumVdwParams, a2);
5321  msg->put(NumVdwParams, a3);
5322  msg->put(NumVdwParams, a4);
5323 
5324  delete [] a1;
5325  delete [] a2;
5326  delete [] a3;
5327  delete [] a4;
5328  }
5329 
5330  // Send the vdw pair parameters
5331  msg->put(NumVdwPairParams);
5332 
5333  if (NumVdwPairParams)
5334  {
5335  a1 = new Real[NumVdwPairParams];
5336  a2 = new Real[NumVdwPairParams];
5337  a3 = new Real[NumVdwPairParams];
5338  a4 = new Real[NumVdwPairParams];
5339  i1 = new int[NumVdwPairParams];
5340  i2 = new int[NumVdwPairParams];
5341 
5342  if ( (a1 == NULL) || (a2 == NULL) || (a3 == NULL) || (a4 == NULL) ||
5343  (i1 == NULL) || (i2 == NULL) )
5344  {
5345  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5346  }
5347 
5348  vdw_pair_to_arrays(i1, i2, a1, a2, a3, a4, 0, vdw_pair_tree);
5349 
5350  msg->put(NumVdwPairParams, i1)->put(NumVdwPairParams, i2);
5351  msg->put(NumVdwPairParams, a1);
5352  msg->put(NumVdwPairParams, a2)->put(NumVdwPairParams, a3);
5353  msg->put(NumVdwPairParams, a4);
5354  }
5355 
5356  // Send the nbthole pair parameters
5357  msg->put(NumNbtholePairParams);
5358 
5360  {
5361  a1 = new Real[NumNbtholePairParams];
5362  a2 = new Real[NumNbtholePairParams];
5363  a3 = new Real[NumNbtholePairParams];
5364  i1 = new int[NumNbtholePairParams];
5365  i2 = new int[NumNbtholePairParams];
5366 
5367  if ( (a1 == NULL) || (a2 == NULL) || (a3 == NULL) || (i1 == NULL) || (i2 == NULL) )
5368  {
5369  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5370  }
5371 
5372  nbthole_pair_to_arrays(i1, i2, a1, a2, a3, 0, nbthole_pair_tree);
5373 
5374  for (i=0; i<NumNbtholePairParams; i++)
5375  {
5376  nbthole_array[i].ind1 = i1[i];
5377  nbthole_array[i].ind2 = i2[i];
5378  nbthole_array[i].alphai = a1[i];
5379  nbthole_array[i].alphaj = a2[i];
5380  nbthole_array[i].tholeij = a3[i];
5381  }
5382 
5383  msg->put(NumNbtholePairParams, i1)->put(NumNbtholePairParams, i2);
5384  msg->put(NumNbtholePairParams, a1);
5385  msg->put(NumNbtholePairParams, a2)->put(NumNbtholePairParams, a3);
5386  }
5387 
5388  // Send the table pair parameters
5389  //printf("Pairs: %i\n", NumTablePairParams);
5390  msg->put(NumTablePairParams);
5391 
5392  if (NumTablePairParams)
5393  {
5394  i1 = new int[NumTablePairParams];
5395  i2 = new int[NumTablePairParams];
5396  i3 = new int[NumTablePairParams];
5397 
5398  if ( (i3 == NULL) || (i1 == NULL) || (i2 == NULL) )
5399  {
5400  NAMD_die("memory allocation failed in Parameters::send_Parameters");
5401  }
5402 
5403  table_pair_to_arrays(i1, i2, i3, 0, tab_pair_tree);
5404 
5405  msg->put(NumTablePairParams, i1)->put(NumTablePairParams, i2);
5406  msg->put(NumTablePairParams, i3);
5407  }
5408 
5409  // send the hydrogen bond parameters
5410  // hbondParams.create_message(msg);
5411  msg->end();
5412  delete msg;
5413 }
int numenerentries
Definition: Parameters.h:251
int NumTablePairParams
Definition: Parameters.h:273
GromacsPairValue * gromacsPair_array
Definition: Parameters.h:247
int NumBondParams
Definition: Parameters.h:259
void end(void)
Definition: MStream.C:176
Real epsilon
Definition: Parameters.h:151
int NumVdwParams
Definition: Parameters.h:268
Real r_ub
Definition: Parameters.h:94
float Real
Definition: common.h:107
int NumVdwPairParams
Definition: Parameters.h:271
Real sigma14
Definition: Parameters.h:152
DihedralValue * dihedral_array
Definition: Parameters.h:243
int NumAngleParams
Definition: Parameters.h:260
CrosstermValue * crossterm_array
Definition: Parameters.h:245
IndexedVdwPair * vdw_pair_tree
Definition: Parameters.h:255
int NumDihedralParams
Definition: Parameters.h:262
IndexedNbtholePair * nbthole_pair_tree
Definition: Parameters.h:256
FourBodyConsts values[MAX_MULTIPLICITY]
Definition: Parameters.h:108
AngleValue * angle_array
Definition: Parameters.h:242
int NumNbtholePairParams
Definition: Parameters.h:272
#define MAX_MULTIPLICITY
Definition: Parameters.h:68
int normal
Definition: Parameters.h:95
Real sigma
Definition: Parameters.h:150
int NumImproperParams
Definition: Parameters.h:263
int NumCrosstermParams
Definition: Parameters.h:264
ImproperValue * improper_array
Definition: Parameters.h:244
void NAMD_die(const char *err_msg)
Definition: common.C:83
int NumVdwParamsAssigned
Definition: Parameters.h:270
NbtholePairValue * nbthole_array
Definition: Parameters.h:250
FourBodyConsts values[MAX_MULTIPLICITY]
Definition: Parameters.h:114
VdwValue * vdw_array
Definition: Parameters.h:249
BigReal * table_ener
Definition: Parameters.h:254
Real x0
Definition: Parameters.h:85
BondValue * bond_array
Definition: Parameters.h:241
int NumGromacsPairParams
Definition: Parameters.h:266
Real epsilon14
Definition: Parameters.h:153
MOStream * put(char data)
Definition: MStream.h:112
Real theta0
Definition: Parameters.h:92
IndexedTablePair * tab_pair_tree
Definition: Parameters.h:257
Real k_ub
Definition: Parameters.h:93
Real k
Definition: Parameters.h:84
#define MAX_ATOMTYPE_CHARS
Definition: Parameters.h:71
double BigReal
Definition: common.h:112

Member Data Documentation

AngleValue* Parameters::angle_array
BondValue* Parameters::bond_array
int Parameters::columnsize

Definition at line 253 of file Parameters.h.

Referenced by ComputeNonbondedUtil::select().

CrosstermValue* Parameters::crossterm_array
DihedralValue* Parameters::dihedral_array
GromacsPairValue* Parameters::gromacsPair_array
ImproperValue* Parameters::improper_array
NbtholePairValue* Parameters::nbthole_array

Definition at line 250 of file Parameters.h.

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

IndexedNbtholePair* Parameters::nbthole_pair_tree

Definition at line 256 of file Parameters.h.

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

int Parameters::NumAngleParams
int Parameters::NumBondParams
int Parameters::NumCosAngles

Definition at line 261 of file Parameters.h.

Referenced by print_param_summary().

int Parameters::NumCrosstermParams
int Parameters::NumDihedralParams
int Parameters::numenerentries

Definition at line 251 of file Parameters.h.

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

int Parameters::NumGromacsPairParams

Definition at line 266 of file Parameters.h.

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

int Parameters::NumImproperParams
int Parameters::NumNbtholePairParams
int Parameters::NumTablePairParams
int Parameters::NumTableParams

Definition at line 269 of file Parameters.h.

int Parameters::NumVdwPairParams
int Parameters::NumVdwParams
int Parameters::NumVdwParamsAssigned
int Parameters::rowsize

Definition at line 252 of file Parameters.h.

Referenced by ComputeNonbondedUtil::select().

IndexedTablePair* Parameters::tab_pair_tree
BigReal* Parameters::table_ener
int Parameters::tablenumtypes

Definition at line 258 of file Parameters.h.

Referenced by get_int_table_type(), and read_ener_table().

VdwValue* Parameters::vdw_array
IndexedVdwPair* Parameters::vdw_pair_tree

The documentation for this class was generated from the following files: