ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE_old/src/mdtools/libmdCode/TraPPE_ExFF.cpp
Revision: 296
Committed: Thu Mar 6 20:05:39 2003 UTC (21 years, 4 months ago) by mmeineke
File size: 41532 byte(s)
Log Message:
dded more static arrays to the Atom class;

File Contents

# Content
1 #include <cstdlib>
2 #include <cstdio>
3 #include <cstring>
4
5 #include <iostream>
6 using namespace std;
7
8 #include "ForceFields.hpp"
9 #include "SRI.hpp"
10 #include "simError.h"
11
12 #include <fortranWrappers.hpp>
13
14 #ifdef IS_MPI
15 #include "mpiForceField.h"
16 #endif // is_mpi
17
18 namespace TPE { // restrict the access of the folowing to this file only.
19
20
21 // Declare the structures that will be passed by MPI
22
23 typedef struct{
24 char name[15];
25 double mass;
26 double epslon;
27 double sigma;
28 double dipole;
29 double w0;
30 double v0;
31 int isSSD;
32 int isDipole;
33 int ident;
34 int last; // 0 -> default
35 // 1 -> tells nodes to stop listening
36 } atomStruct;
37
38
39 typedef struct{
40 char nameA[15];
41 char nameB[15];
42 char type[30];
43 double d0;
44 int last; // 0 -> default
45 // 1 -> tells nodes to stop listening
46 } bondStruct;
47
48
49 typedef struct{
50 char nameA[15];
51 char nameB[15];
52 char nameC[15];
53 char type[30];
54 double k1, k2, k3, t0;
55 int last; // 0 -> default
56 // 1 -> tells nodes to stop listening
57 } bendStruct;
58
59
60 typedef struct{
61 char nameA[15];
62 char nameB[15];
63 char nameC[15];
64 char nameD[15];
65 char type[30];
66 double k1, k2, k3, k4;
67 int last; // 0 -> default
68 // 1 -> tells nodes to stop listening
69 } torsionStruct;
70
71
72 int parseAtom( char *lineBuffer, int lineNum, atomStruct &info );
73 int parseBond( char *lineBuffer, int lineNum, bondStruct &info );
74 int parseBend( char *lineBuffer, int lineNum, bendStruct &info );
75 int parseTorsion( char *lineBuffer, int lineNum, torsionStruct &info );
76
77
78 #ifdef IS_MPI
79
80 MPI_Datatype mpiAtomStructType;
81 MPI_Datatype mpiBondStructType;
82 MPI_Datatype mpiBendStructType;
83 MPI_Datatype mpiTorsionStructType;
84
85 #endif
86
87 } // namespace
88
89 using namespace TPE;
90
91
92 //****************************************************************
93 // begins the actual forcefield stuff.
94 //****************************************************************
95
96
97 TraPPE_ExFF::TraPPE_ExFF(){
98
99 char fileName[200];
100 char* ffPath_env = "FORCE_PARAM_PATH";
101 char* ffPath;
102 char temp[200];
103 char errMsg[1000];
104
105 // do the funtion wrapping
106 wrapMeFF( this );
107
108
109 #ifdef IS_MPI
110 int i;
111
112 // **********************************************************************
113 // Init the atomStruct mpi type
114
115 atomStruct atomProto; // mpiPrototype
116 int atomBC[3] = {15,6,4}; // block counts
117 MPI_Aint atomDspls[3]; // displacements
118 MPI_Datatype atomMbrTypes[3]; // member mpi types
119
120 MPI_Address(&atomProto.name, &atomDspls[0]);
121 MPI_Address(&atomProto.mass, &atomDspls[1]);
122 MPI_Address(&atomProto.isSSD, &atomDspls[2]);
123
124 atomMbrTypes[0] = MPI_CHAR;
125 atomMbrTypes[1] = MPI_DOUBLE;
126 atomMbrTypes[2] = MPI_INT;
127
128 for (i=2; i >= 0; i--) atomDspls[i] -= atomDspls[0];
129
130 MPI_Type_struct(3, atomBC, atomDspls, atomMbrTypes, &mpiAtomStructType);
131 MPI_Type_commit(&mpiAtomStructType);
132
133
134 // **********************************************************************
135 // Init the bondStruct mpi type
136
137 bondStruct bondProto; // mpiPrototype
138 int bondBC[3] = {60,1,1}; // block counts
139 MPI_Aint bondDspls[3]; // displacements
140 MPI_Datatype bondMbrTypes[3]; // member mpi types
141
142 MPI_Address(&bondProto.nameA, &bondDspls[0]);
143 MPI_Address(&bondProto.d0, &bondDspls[1]);
144 MPI_Address(&bondProto.last, &bondDspls[2]);
145
146 bondMbrTypes[0] = MPI_CHAR;
147 bondMbrTypes[1] = MPI_DOUBLE;
148 bondMbrTypes[2] = MPI_INT;
149
150 for (i=2; i >= 0; i--) bondDspls[i] -= bondDspls[0];
151
152 MPI_Type_struct(3, bondBC, bondDspls, bondMbrTypes, &mpiBondStructType);
153 MPI_Type_commit(&mpiBondStructType);
154
155
156 // **********************************************************************
157 // Init the bendStruct mpi type
158
159 bendStruct bendProto; // mpiPrototype
160 int bendBC[3] = {75,4,1}; // block counts
161 MPI_Aint bendDspls[3]; // displacements
162 MPI_Datatype bendMbrTypes[3]; // member mpi types
163
164 MPI_Address(&bendProto.nameA, &bendDspls[0]);
165 MPI_Address(&bendProto.k1, &bendDspls[1]);
166 MPI_Address(&bendProto.last, &bendDspls[2]);
167
168 bendMbrTypes[0] = MPI_CHAR;
169 bendMbrTypes[1] = MPI_DOUBLE;
170 bendMbrTypes[2] = MPI_INT;
171
172 for (i=2; i >= 0; i--) bendDspls[i] -= bendDspls[0];
173
174 MPI_Type_struct(3, bendBC, bendDspls, bendMbrTypes, &mpiBendStructType);
175 MPI_Type_commit(&mpiBendStructType);
176
177
178 // **********************************************************************
179 // Init the torsionStruct mpi type
180
181 torsionStruct torsionProto; // mpiPrototype
182 int torsionBC[3] = {90,4,1}; // block counts
183 MPI_Aint torsionDspls[3]; // displacements
184 MPI_Datatype torsionMbrTypes[3]; // member mpi types
185
186 MPI_Address(&torsionProto.nameA, &torsionDspls[0]);
187 MPI_Address(&torsionProto.k1, &torsionDspls[1]);
188 MPI_Address(&torsionProto.last, &torsionDspls[2]);
189
190 torsionMbrTypes[0] = MPI_CHAR;
191 torsionMbrTypes[1] = MPI_DOUBLE;
192 torsionMbrTypes[2] = MPI_INT;
193
194 for (i=2; i >= 0; i--) torsionDspls[i] -= torsionDspls[0];
195
196 MPI_Type_struct(3, torsionBC, torsionDspls, torsionMbrTypes,
197 &mpiTorsionStructType);
198 MPI_Type_commit(&mpiTorsionStructType);
199
200 // ***********************************************************************
201
202 if( worldRank == 0 ){
203 #endif
204
205 // generate the force file name
206
207 strcpy( fileName, "TraPPE_Ex.frc" );
208 // fprintf( stderr,"Trying to open %s\n", fileName );
209
210 // attempt to open the file in the current directory first.
211
212 frcFile = fopen( fileName, "r" );
213
214 if( frcFile == NULL ){
215
216 // next see if the force path enviorment variable is set
217
218 ffPath = getenv( ffPath_env );
219 if( ffPath == NULL ) {
220 sprintf( painCave.errMsg,
221 "Error opening the force field parameter file: %s\n"
222 "Have you tried setting the FORCE_PARAM_PATH environment "
223 "vairable?\n",
224 fileName );
225 painCave.isFatal = 1;
226 simError();
227 }
228
229
230 strcpy( temp, ffPath );
231 strcat( temp, "/" );
232 strcat( temp, fileName );
233 strcpy( fileName, temp );
234
235 frcFile = fopen( fileName, "r" );
236
237 if( frcFile == NULL ){
238
239 sprintf( painCave.errMsg,
240 "Error opening the force field parameter file: %s\n"
241 "Have you tried setting the FORCE_PARAM_PATH environment "
242 "vairable?\n",
243 fileName );
244 painCave.isFatal = 1;
245 simError();
246 }
247 }
248
249 #ifdef IS_MPI
250 }
251
252 sprintf( checkPointMsg, "TraPPE_ExFF file opened sucessfully." );
253 MPIcheckPoint();
254
255 #endif // is_mpi
256 }
257
258
259 TraPPE_ExFF::~TraPPE_ExFF(){
260
261 #ifdef IS_MPI
262 if( worldRank == 0 ){
263 #endif // is_mpi
264
265 fclose( frcFile );
266
267 #ifdef IS_MPI
268 }
269 #endif // is_mpi
270 }
271
272 void TraPPE_ExFF::doForces( int calcPot ){
273
274 int i, isError;
275 double* frc;
276 double* pos;
277 double* tau;
278 short int passedCalcPot = (short int)calcPot;
279
280 // forces are zeroed here, before any are acumulated.
281 // NOTE: do not rezero the forces in Fortran.
282
283 for(i=0; i<entry_plug->n_atoms; i++){
284 entry_plug->atoms[i]->zeroForces();
285 }
286
287 frc = Atom::getFrcArray();
288 pos = Atom::getPosArray();
289 tau = entry_plug->tau;
290
291 isError = 0;
292 fortranForceLoop( pos, frc, &(entry_plug->lrPot), tau,
293 &passedCalcPot, &isError );
294
295
296 if( isError ){
297 sprintf( painCave.errMsg,
298 "Error returned from the fortran force calculation.\n" );
299 painCave.isFatal = 1;
300 simError();
301 }
302
303 #ifdef IS_MPI
304 sprintf( checkPointMsg,
305 "successfully returned from the force calculation.\n" );
306 MPIcheckPoint();
307 #endif // is_mpi
308
309 }
310
311 void TraPPE_ExFF::initFortran( void ){
312
313 int nLocal = entry_plug->n_atoms;
314 int *ident;
315 int isError;
316 int i;
317
318 ident = new int[nLocal];
319
320 for(i=0; i<nLocal; i++){
321 ident[i] = entry_plug->atoms[i]->getIdent();
322 }
323
324 isError = 0;
325 initfortran( &nLocal, ident, &isError );
326
327 if(isError){
328 sprintf( painCave.errMsg,
329 "TraPPE_ExFF error: There was an error initializing the component list in fortran.\n" );
330 painCave.isFatal = 1;
331 simError();
332 }
333
334
335 #ifdef IS_MPI
336 sprintf( checkPointMsg, "TraPPE_ExFF successfully initialized the fortran component list.\n" );
337 MPIcheckPoint();
338 #endif // is_mpi
339
340 delete[] ident;
341
342 }
343
344
345
346
347 void TraPPE_ExFF::initializeAtoms( void ){
348
349 class LinkedType {
350 public:
351 LinkedType(){
352 next = NULL;
353 name[0] = '\0';
354 }
355 ~LinkedType(){ if( next != NULL ) delete next; }
356
357 LinkedType* find(char* key){
358 if( !strcmp(name, key) ) return this;
359 if( next != NULL ) return next->find(key);
360 return NULL;
361 }
362
363 void add( atomStruct &info ){
364
365 // check for duplicates
366
367 if( !strcmp( info.name, name ) ){
368 sprintf( painCave.errMsg,
369 "Duplicate TraPPE_Ex atom type \"%s\" found in "
370 "the TraPPE_ExFF param file./n",
371 name );
372 painCave.isFatal = 1;
373 simError();
374 }
375
376 if( next != NULL ) next->add(info);
377 else{
378 next = new LinkedType();
379 strcpy(next->name, info.name);
380 next->isDipole = info.isDipole;
381 next->isSSD = info.isSSD;
382 next->mass = info.mass;
383 next->epslon = info.epslon;
384 next->sigma = info.sigma;
385 next->dipole = info.dipole;
386 next->w0 = info.w0;
387 next->v0 = info.v0;
388 next->ident = info.ident;
389 }
390 }
391
392 #ifdef IS_MPI
393
394 void duplicate( atomStruct &info ){
395 strcpy(info.name, name);
396 info.isDipole = isDipole;
397 info.isSSD = isSSD;
398 info.mass = mass;
399 info.epslon = epslon;
400 info.sigma = sigma;
401 info.dipole = dipole;
402 info.w0 = w0;
403 info.v0 = v0;
404 info.last = 0;
405 }
406
407
408 #endif
409
410 char name[15];
411 int isDipole;
412 int isSSD;
413 double mass;
414 double epslon;
415 double sigma;
416 double dipole;
417 double w0;
418 double v0;
419 int ident;
420 LinkedType* next;
421 };
422
423 LinkedType* headAtomType;
424 LinkedType* currentAtomType;
425 atomStruct info;
426 info.last = 1; // initialize last to have the last set.
427 // if things go well, last will be set to 0
428
429
430
431 int i;
432 int identNum;
433
434 Atom** the_atoms;
435 int nAtoms;
436 the_atoms = entry_plug->atoms;
437 nAtoms = entry_plug->n_atoms;
438
439
440 //////////////////////////////////////////////////
441 // a quick water fix
442
443 double waterI[3][3];
444 waterI[0][0] = 1.76958347772500;
445 waterI[0][1] = 0.0;
446 waterI[0][2] = 0.0;
447
448 waterI[1][0] = 0.0;
449 waterI[1][1] = 0.614537057924513;
450 waterI[1][2] = 0.0;
451
452 waterI[2][0] = 0.0;
453 waterI[2][1] = 0.0;
454 waterI[2][2] = 1.15504641980049;
455
456
457 double headI[3][3];
458 headI[0][0] = 1125;
459 headI[0][1] = 0.0;
460 headI[0][2] = 0.0;
461
462 headI[1][0] = 0.0;
463 headI[1][1] = 1125;
464 headI[1][2] = 0.0;
465
466 headI[2][0] = 0.0;
467 headI[2][1] = 0.0;
468 headI[2][2] = 250;
469
470
471
472 //////////////////////////////////////////////////
473
474
475 #ifdef IS_MPI
476 if( worldRank == 0 ){
477 #endif
478
479 // read in the atom types.
480
481 headAtomType = new LinkedType;
482
483 fastForward( "AtomTypes", "initializeAtoms" );
484
485 // we are now at the AtomTypes section.
486
487 eof_test = fgets( readLine, sizeof(readLine), frcFile );
488 lineNum++;
489
490
491 // read a line, and start parseing out the atom types
492
493 if( eof_test == NULL ){
494 sprintf( painCave.errMsg,
495 "Error in reading Atoms from force file at line %d.\n",
496 lineNum );
497 painCave.isFatal = 1;
498 simError();
499 }
500
501 identNum = 1;
502 // stop reading at end of file, or at next section
503 while( readLine[0] != '#' && eof_test != NULL ){
504
505 // toss comment lines
506 if( readLine[0] != '!' ){
507
508 // the parser returns 0 if the line was blank
509 if( parseAtom( readLine, lineNum, info ) ){
510 info.ident = identNum;
511 headAtomType->add( info );;
512 identNum++;
513 }
514 }
515 eof_test = fgets( readLine, sizeof(readLine), frcFile );
516 lineNum++;
517 }
518
519 #ifdef IS_MPI
520
521 // send out the linked list to all the other processes
522
523 sprintf( checkPointMsg,
524 "TraPPE_ExFF atom structures read successfully." );
525 MPIcheckPoint();
526
527 currentAtomType = headAtomType->next; //skip the first element who is a place holder.
528 while( currentAtomType != NULL ){
529 currentAtomType->duplicate( info );
530
531
532
533 sendFrcStruct( &info, mpiAtomStructType );
534
535 sprintf( checkPointMsg,
536 "successfully sent TraPPE_Ex force type: \"%s\"\n",
537 info.name );
538 MPIcheckPoint();
539
540 currentAtomType = currentAtomType->next;
541 }
542 info.last = 1;
543 sendFrcStruct( &info, mpiAtomStructType );
544
545 }
546
547 else{
548
549 // listen for node 0 to send out the force params
550
551 MPIcheckPoint();
552
553 headAtomType = new LinkedType;
554 recieveFrcStruct( &info, mpiAtomStructType );
555
556 while( !info.last ){
557
558
559
560 headAtomType->add( info );
561
562 MPIcheckPoint();
563
564 recieveFrcStruct( &info, mpiAtomStructType );
565 }
566 }
567 #endif // is_mpi
568
569 // call new A_types in fortran
570
571 int isError;
572
573 // dummy variables
574
575 int isGB = 0;
576 int isLJ = 1;
577
578
579 currentAtomType = headAtomType;
580 while( currentAtomType != NULL ){
581
582 if( currentAtomType->name[0] != '\0' ){
583 isError = 0;
584 newTPEtype( &(currentAtomType->ident),
585 &(currentAtomType->mass),
586 &(currentAtomType->epslon),
587 &(currentAtomType->sigma),
588 &isLJ,
589 &(currentAtomType->isSSD),
590 &(currentAtomType->isDipole),
591 &isGB,
592 &(currentAtomType->w0),
593 &(currentAtomType->v0),
594 &(currentAtomType->dipole),
595 &isError );
596 if( isError ){
597 sprintf( painCave.errMsg,
598 "Error initializing the \"%s\" atom type in fortran\n",
599 currentAtomType->name );
600 painCave.isFatal = 1;
601 simError();
602 }
603 }
604 currentAtomType = currentAtomType->next;
605 }
606
607 #ifdef IS_MPI
608 sprintf( checkPointMsg,
609 "TraPPE_ExFF atom structures successfully sent to fortran\n" );
610 MPIcheckPoint();
611 #endif // is_mpi
612
613
614 // initialize the atoms
615
616 double bigSigma = 0.0;
617 DirectionalAtom* dAtom;
618
619 for( i=0; i<nAtoms; i++ ){
620
621 currentAtomType = headAtomType->find( the_atoms[i]->getType() );
622 if( currentAtomType == NULL ){
623 sprintf( painCave.errMsg,
624 "AtomType error, %s not found in force file.\n",
625 the_atoms[i]->getType() );
626 painCave.isFatal = 1;
627 simError();
628 }
629
630 the_atoms[i]->setMass( currentAtomType->mass );
631 the_atoms[i]->setEpslon( currentAtomType->epslon );
632 the_atoms[i]->setSigma( currentAtomType->sigma );
633 the_atoms[i]->setLJ();
634
635 if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
636
637 if( currentAtomType->isDipole ){
638 if( the_atoms[i]->isDirectional() ){
639
640 dAtom = (DirectionalAtom *) the_atoms[i];
641 dAtom->setMu( currentAtomType->dipole );
642 dAtom->setHasDipole( 1 );
643 dAtom->setJx( 0.0 );
644 dAtom->setJy( 0.0 );
645 dAtom->setJz( 0.0 );
646
647 if(!strcmp("SSD",the_atoms[i]->getType())){
648 dAtom->setI( waterI );
649 dAtom->setSSD( 1 );
650 }
651 else if(!strcmp("HEAD",the_atoms[i]->getType())){
652 dAtom->setI( headI );
653 dAtom->setSSD( 0 );
654 }
655 else{
656 sprintf(painCave.errMsg,
657 "AtmType error, %s does not have a moment of inertia set.\n",
658 the_atoms[i]->getType() );
659 painCave.isFatal = 1;
660 simError();
661 }
662 entry_plug->n_dipoles++;
663 }
664 else{
665
666 sprintf( painCave.errMsg,
667 "TraPPE_ExFF error: Atom \"%s\" is a dipole, yet no standard"
668 " orientation was specifed in the BASS file.\n",
669 currentAtomType->name );
670 painCave.isFatal = 1;
671 simError();
672 }
673 }
674 else{
675 if( the_atoms[i]->isDirectional() ){
676 sprintf( painCave.errMsg,
677 "TraPPE_ExFF error: Atom \"%s\" was given a standard"
678 "orientation in the BASS file, yet it is not a dipole.\n",
679 currentAtomType->name);
680 painCave.isFatal = 1;
681 simError();
682 }
683 }
684 }
685
686 #ifdef IS_MPI
687 double tempBig = bigSigma;
688 MPI::COMM_WORLD.Allreduce( &tempBig, &bigSigma, 1, MPI_DOUBLE, MPI_MAX );
689 #endif //is_mpi
690
691 //calc rCut and rList
692
693 entry_plug->rCut = 2.5 * bigSigma;
694
695 if(entry_plug->rCut > (entry_plug->box_x / 2.0))
696 entry_plug->rCut = entry_plug->box_x / 2.0;
697
698 if(entry_plug->rCut > (entry_plug->box_y / 2.0))
699 entry_plug->rCut = entry_plug->box_y / 2.0;
700
701 if(entry_plug->rCut > (entry_plug->box_z / 2.0))
702 entry_plug->rCut = entry_plug->box_z / 2.0;
703
704 entry_plug->rList = entry_plug->rCut + 1.0;
705
706
707 // clean up the memory
708
709 delete headAtomType;
710
711 #ifdef IS_MPI
712 sprintf( checkPointMsg, "TraPPE_Ex atoms initialized succesfully" );
713 MPIcheckPoint();
714 #endif // is_mpi
715
716 initFortran();
717 entry_plug->refreshSim();
718
719 }
720
721 void TraPPE_ExFF::initializeBonds( bond_pair* the_bonds ){
722
723 class LinkedType {
724 public:
725 LinkedType(){
726 next = NULL;
727 nameA[0] = '\0';
728 nameB[0] = '\0';
729 type[0] = '\0';
730 }
731 ~LinkedType(){ if( next != NULL ) delete next; }
732
733 LinkedType* find(char* key1, char* key2){
734 if( !strcmp(nameA, key1 ) && !strcmp( nameB, key2 ) ) return this;
735 if( !strcmp(nameA, key2 ) && !strcmp( nameB, key1 ) ) return this;
736 if( next != NULL ) return next->find(key1, key2);
737 return NULL;
738 }
739
740
741 void add( bondStruct &info ){
742
743 // check for duplicates
744 int dup = 0;
745
746 if( !strcmp(nameA, info.nameA ) && !strcmp( nameB, info.nameB ) ) dup = 1;
747 if( !strcmp(nameA, info.nameB ) && !strcmp( nameB, info.nameA ) ) dup = 1;
748
749 if(dup){
750 sprintf( painCave.errMsg,
751 "Duplicate TraPPE_Ex bond type \"%s - %s\" found in "
752 "the TraPPE_ExFF param file./n",
753 nameA, nameB );
754 painCave.isFatal = 1;
755 simError();
756 }
757
758
759 if( next != NULL ) next->add(info);
760 else{
761 next = new LinkedType();
762 strcpy(next->nameA, info.nameA);
763 strcpy(next->nameB, info.nameB);
764 strcpy(next->type, info.type);
765 next->d0 = info.d0;
766 }
767 }
768
769 #ifdef IS_MPI
770 void duplicate( bondStruct &info ){
771 strcpy(info.nameA, nameA);
772 strcpy(info.nameB, nameB);
773 strcpy(info.type, type);
774 info.d0 = d0;
775 info.last = 0;
776 }
777
778
779 #endif
780
781 char nameA[15];
782 char nameB[15];
783 char type[30];
784 double d0;
785
786 LinkedType* next;
787 };
788
789
790
791 LinkedType* headBondType;
792 LinkedType* currentBondType;
793 bondStruct info;
794 info.last = 1; // initialize last to have the last set.
795 // if things go well, last will be set to 0
796
797 SRI **the_sris;
798 Atom** the_atoms;
799 int nBonds;
800 the_sris = entry_plug->sr_interactions;
801 the_atoms = entry_plug->atoms;
802 nBonds = entry_plug->n_bonds;
803
804 int i, a, b;
805 char* atomA;
806 char* atomB;
807
808 #ifdef IS_MPI
809 if( worldRank == 0 ){
810 #endif
811
812 // read in the bond types.
813
814 headBondType = new LinkedType;
815
816 fastForward( "BondTypes", "initializeBonds" );
817
818 // we are now at the bondTypes section
819
820 eof_test = fgets( readLine, sizeof(readLine), frcFile );
821 lineNum++;
822
823
824 // read a line, and start parseing out the atom types
825
826 if( eof_test == NULL ){
827 sprintf( painCave.errMsg,
828 "Error in reading bonds from force file at line %d.\n",
829 lineNum );
830 painCave.isFatal = 1;
831 simError();
832 }
833
834 // stop reading at end of file, or at next section
835 while( readLine[0] != '#' && eof_test != NULL ){
836
837 // toss comment lines
838 if( readLine[0] != '!' ){
839
840 // the parser returns 0 if the line was blank
841 if( parseBond( readLine, lineNum, info ) ){
842 headBondType->add( info );
843 }
844 }
845 eof_test = fgets( readLine, sizeof(readLine), frcFile );
846 lineNum++;
847 }
848
849 #ifdef IS_MPI
850
851 // send out the linked list to all the other processes
852
853 sprintf( checkPointMsg,
854 "TraPPE_Ex bond structures read successfully." );
855 MPIcheckPoint();
856
857 currentBondType = headBondType;
858 while( currentBondType != NULL ){
859 currentBondType->duplicate( info );
860 sendFrcStruct( &info, mpiBondStructType );
861 currentBondType = currentBondType->next;
862 }
863 info.last = 1;
864 sendFrcStruct( &info, mpiBondStructType );
865
866 }
867
868 else{
869
870 // listen for node 0 to send out the force params
871
872 MPIcheckPoint();
873
874 headBondType = new LinkedType;
875 recieveFrcStruct( &info, mpiBondStructType );
876 while( !info.last ){
877
878 headBondType->add( info );
879 recieveFrcStruct( &info, mpiBondStructType );
880 }
881 }
882 #endif // is_mpi
883
884
885 // initialize the Bonds
886
887
888 for( i=0; i<nBonds; i++ ){
889
890 a = the_bonds[i].a;
891 b = the_bonds[i].b;
892
893 atomA = the_atoms[a]->getType();
894 atomB = the_atoms[b]->getType();
895 currentBondType = headBondType->find( atomA, atomB );
896 if( currentBondType == NULL ){
897 sprintf( painCave.errMsg,
898 "BondType error, %s - %s not found in force file.\n",
899 atomA, atomB );
900 painCave.isFatal = 1;
901 simError();
902 }
903
904 if( !strcmp( currentBondType->type, "fixed" ) ){
905
906 the_sris[i] = new ConstrainedBond( *the_atoms[a],
907 *the_atoms[b],
908 currentBondType->d0 );
909 entry_plug->n_constraints++;
910 }
911 }
912
913
914 // clean up the memory
915
916 delete headBondType;
917
918 #ifdef IS_MPI
919 sprintf( checkPointMsg, "TraPPE_Ex bonds initialized succesfully" );
920 MPIcheckPoint();
921 #endif // is_mpi
922
923 }
924
925 void TraPPE_ExFF::initializeBends( bend_set* the_bends ){
926
927 class LinkedType {
928 public:
929 LinkedType(){
930 next = NULL;
931 nameA[0] = '\0';
932 nameB[0] = '\0';
933 nameC[0] = '\0';
934 type[0] = '\0';
935 }
936 ~LinkedType(){ if( next != NULL ) delete next; }
937
938 LinkedType* find( char* key1, char* key2, char* key3 ){
939 if( !strcmp( nameA, key1 ) && !strcmp( nameB, key2 )
940 && !strcmp( nameC, key3 ) ) return this;
941 if( !strcmp( nameA, key3 ) && !strcmp( nameB, key2 )
942 && !strcmp( nameC, key1 ) ) return this;
943 if( next != NULL ) return next->find(key1, key2, key3);
944 return NULL;
945 }
946
947 void add( bendStruct &info ){
948
949 // check for duplicates
950 int dup = 0;
951
952 if( !strcmp( nameA, info.nameA ) && !strcmp( nameB, info.nameB )
953 && !strcmp( nameC, info.nameC ) ) dup = 1;
954 if( !strcmp( nameA, info.nameC ) && !strcmp( nameB, info.nameB )
955 && !strcmp( nameC, info.nameA ) ) dup = 1;
956
957 if(dup){
958 sprintf( painCave.errMsg,
959 "Duplicate TraPPE_Ex bend type \"%s - %s - %s\" found in "
960 "the TraPPE_ExFF param file./n",
961 nameA, nameB, nameC );
962 painCave.isFatal = 1;
963 simError();
964 }
965
966 if( next != NULL ) next->add(info);
967 else{
968 next = new LinkedType();
969 strcpy(next->nameA, info.nameA);
970 strcpy(next->nameB, info.nameB);
971 strcpy(next->nameC, info.nameC);
972 strcpy(next->type, info.type);
973 next->k1 = info.k1;
974 next->k2 = info.k2;
975 next->k3 = info.k3;
976 next->t0 = info.t0;
977 }
978 }
979
980 #ifdef IS_MPI
981
982 void duplicate( bendStruct &info ){
983 strcpy(info.nameA, nameA);
984 strcpy(info.nameB, nameB);
985 strcpy(info.nameC, nameC);
986 strcpy(info.type, type);
987 info.k1 = k1;
988 info.k2 = k2;
989 info.k3 = k3;
990 info.t0 = t0;
991 info.last = 0;
992 }
993
994 #endif // is_mpi
995
996 char nameA[15];
997 char nameB[15];
998 char nameC[15];
999 char type[30];
1000 double k1, k2, k3, t0;
1001
1002 LinkedType* next;
1003 };
1004
1005 LinkedType* headBendType;
1006 LinkedType* currentBendType;
1007 bendStruct info;
1008 info.last = 1; // initialize last to have the last set.
1009 // if things go well, last will be set to 0
1010
1011 QuadraticBend* qBend;
1012 SRI **the_sris;
1013 Atom** the_atoms;
1014 int nBends;
1015 the_sris = entry_plug->sr_interactions;
1016 the_atoms = entry_plug->atoms;
1017 nBends = entry_plug->n_bends;
1018
1019 int i, a, b, c;
1020 char* atomA;
1021 char* atomB;
1022 char* atomC;
1023
1024
1025 #ifdef IS_MPI
1026 if( worldRank == 0 ){
1027 #endif
1028
1029 // read in the bend types.
1030
1031 headBendType = new LinkedType;
1032
1033 fastForward( "BendTypes", "initializeBends" );
1034
1035 // we are now at the bendTypes section
1036
1037 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1038 lineNum++;
1039
1040 // read a line, and start parseing out the bend types
1041
1042 if( eof_test == NULL ){
1043 sprintf( painCave.errMsg,
1044 "Error in reading bends from force file at line %d.\n",
1045 lineNum );
1046 painCave.isFatal = 1;
1047 simError();
1048 }
1049
1050 // stop reading at end of file, or at next section
1051 while( readLine[0] != '#' && eof_test != NULL ){
1052
1053 // toss comment lines
1054 if( readLine[0] != '!' ){
1055
1056 // the parser returns 0 if the line was blank
1057 if( parseBend( readLine, lineNum, info ) ){
1058 headBendType->add( info );
1059 }
1060 }
1061 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1062 lineNum++;
1063 }
1064
1065 #ifdef IS_MPI
1066
1067 // send out the linked list to all the other processes
1068
1069 sprintf( checkPointMsg,
1070 "TraPPE_Ex bend structures read successfully." );
1071 MPIcheckPoint();
1072
1073 currentBendType = headBendType;
1074 while( currentBendType != NULL ){
1075 currentBendType->duplicate( info );
1076 sendFrcStruct( &info, mpiBendStructType );
1077 currentBendType = currentBendType->next;
1078 }
1079 info.last = 1;
1080 sendFrcStruct( &info, mpiBendStructType );
1081
1082 }
1083
1084 else{
1085
1086 // listen for node 0 to send out the force params
1087
1088 MPIcheckPoint();
1089
1090 headBendType = new LinkedType;
1091 recieveFrcStruct( &info, mpiBendStructType );
1092 while( !info.last ){
1093
1094 headBendType->add( info );
1095 recieveFrcStruct( &info, mpiBendStructType );
1096 }
1097 }
1098 #endif // is_mpi
1099
1100 // initialize the Bends
1101
1102 int index;
1103
1104 for( i=0; i<nBends; i++ ){
1105
1106 a = the_bends[i].a;
1107 b = the_bends[i].b;
1108 c = the_bends[i].c;
1109
1110 atomA = the_atoms[a]->getType();
1111 atomB = the_atoms[b]->getType();
1112 atomC = the_atoms[c]->getType();
1113 currentBendType = headBendType->find( atomA, atomB, atomC );
1114 if( currentBendType == NULL ){
1115 sprintf( painCave.errMsg, "BendType error, %s - %s - %s not found"
1116 " in force file.\n",
1117 atomA, atomB, atomC );
1118 painCave.isFatal = 1;
1119 simError();
1120 }
1121
1122 if( !strcmp( currentBendType->type, "quadratic" ) ){
1123
1124 index = i + entry_plug->n_bonds;
1125 qBend = new QuadraticBend( *the_atoms[a],
1126 *the_atoms[b],
1127 *the_atoms[c] );
1128 qBend->setConstants( currentBendType->k1,
1129 currentBendType->k2,
1130 currentBendType->k3,
1131 currentBendType->t0 );
1132 the_sris[index] = qBend;
1133 }
1134 }
1135
1136
1137 // clean up the memory
1138
1139 delete headBendType;
1140
1141 #ifdef IS_MPI
1142 sprintf( checkPointMsg, "TraPPE_Ex bends initialized succesfully" );
1143 MPIcheckPoint();
1144 #endif // is_mpi
1145
1146 }
1147
1148 void TraPPE_ExFF::initializeTorsions( torsion_set* the_torsions ){
1149
1150 class LinkedType {
1151 public:
1152 LinkedType(){
1153 next = NULL;
1154 nameA[0] = '\0';
1155 nameB[0] = '\0';
1156 nameC[0] = '\0';
1157 type[0] = '\0';
1158 }
1159 ~LinkedType(){ if( next != NULL ) delete next; }
1160
1161 LinkedType* find( char* key1, char* key2, char* key3, char* key4 ){
1162 if( !strcmp( nameA, key1 ) && !strcmp( nameB, key2 ) &&
1163 !strcmp( nameC, key3 ) && !strcmp( nameD, key4 ) ) return this;
1164
1165 if( !strcmp( nameA, key4 ) && !strcmp( nameB, key3 ) &&
1166 !strcmp( nameC, key2 ) && !strcmp( nameD, key1 ) ) return this;
1167
1168 if( next != NULL ) return next->find(key1, key2, key3, key4);
1169 return NULL;
1170 }
1171
1172 void add( torsionStruct &info ){
1173
1174 // check for duplicates
1175 int dup = 0;
1176
1177 if( !strcmp( nameA, info.nameA ) && !strcmp( nameB, info.nameB ) &&
1178 !strcmp( nameC, info.nameC ) && !strcmp( nameD, info.nameD ) ) dup = 1;
1179
1180 if( !strcmp( nameA, info.nameD ) && !strcmp( nameB, info.nameC ) &&
1181 !strcmp( nameC, info.nameB ) && !strcmp( nameD, info.nameA ) ) dup = 1;
1182
1183 if(dup){
1184 sprintf( painCave.errMsg,
1185 "Duplicate TraPPE_Ex torsion type \"%s - %s - %s - %s\" found in "
1186 "the TraPPE_ExFF param file./n", nameA, nameB, nameC, nameD );
1187 painCave.isFatal = 1;
1188 simError();
1189 }
1190
1191 if( next != NULL ) next->add(info);
1192 else{
1193 next = new LinkedType();
1194 strcpy(next->nameA, info.nameA);
1195 strcpy(next->nameB, info.nameB);
1196 strcpy(next->nameC, info.nameC);
1197 strcpy(next->type, info.type);
1198 next->k1 = info.k1;
1199 next->k2 = info.k2;
1200 next->k3 = info.k3;
1201 next->k4 = info.k4;
1202 }
1203 }
1204
1205 #ifdef IS_MPI
1206
1207 void duplicate( torsionStruct &info ){
1208 strcpy(info.nameA, nameA);
1209 strcpy(info.nameB, nameB);
1210 strcpy(info.nameC, nameC);
1211 strcpy(info.nameD, nameD);
1212 strcpy(info.type, type);
1213 info.k1 = k1;
1214 info.k2 = k2;
1215 info.k3 = k3;
1216 info.k4 = k4;
1217 info.last = 0;
1218 }
1219
1220 #endif
1221
1222 char nameA[15];
1223 char nameB[15];
1224 char nameC[15];
1225 char nameD[15];
1226 char type[30];
1227 double k1, k2, k3, k4;
1228
1229 LinkedType* next;
1230 };
1231
1232 LinkedType* headTorsionType;
1233 LinkedType* currentTorsionType;
1234 torsionStruct info;
1235 info.last = 1; // initialize last to have the last set.
1236 // if things go well, last will be set to 0
1237
1238 int i, a, b, c, d, index;
1239 char* atomA;
1240 char* atomB;
1241 char* atomC;
1242 char* atomD;
1243 CubicTorsion* cTors;
1244
1245 SRI **the_sris;
1246 Atom** the_atoms;
1247 int nTorsions;
1248 the_sris = entry_plug->sr_interactions;
1249 the_atoms = entry_plug->atoms;
1250 nTorsions = entry_plug->n_torsions;
1251
1252 #ifdef IS_MPI
1253 if( worldRank == 0 ){
1254 #endif
1255
1256 // read in the torsion types.
1257
1258 headTorsionType = new LinkedType;
1259
1260 fastForward( "TorsionTypes", "initializeTorsions" );
1261
1262 // we are now at the torsionTypes section
1263
1264 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1265 lineNum++;
1266
1267
1268 // read a line, and start parseing out the atom types
1269
1270 if( eof_test == NULL ){
1271 sprintf( painCave.errMsg,
1272 "Error in reading torsions from force file at line %d.\n",
1273 lineNum );
1274 painCave.isFatal = 1;
1275 simError();
1276 }
1277
1278 // stop reading at end of file, or at next section
1279 while( readLine[0] != '#' && eof_test != NULL ){
1280
1281 // toss comment lines
1282 if( readLine[0] != '!' ){
1283
1284 // the parser returns 0 if the line was blank
1285 if( parseTorsion( readLine, lineNum, info ) ){
1286 headTorsionType->add( info );
1287 }
1288 }
1289 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1290 lineNum++;
1291 }
1292
1293 #ifdef IS_MPI
1294
1295 // send out the linked list to all the other processes
1296
1297 sprintf( checkPointMsg,
1298 "TraPPE_Ex torsion structures read successfully." );
1299 MPIcheckPoint();
1300
1301 currentTorsionType = headTorsionType;
1302 while( currentTorsionType != NULL ){
1303 currentTorsionType->duplicate( info );
1304 sendFrcStruct( &info, mpiTorsionStructType );
1305 currentTorsionType = currentTorsionType->next;
1306 }
1307 info.last = 1;
1308 sendFrcStruct( &info, mpiTorsionStructType );
1309
1310 }
1311
1312 else{
1313
1314 // listen for node 0 to send out the force params
1315
1316 MPIcheckPoint();
1317
1318 headTorsionType = new LinkedType;
1319 recieveFrcStruct( &info, mpiTorsionStructType );
1320 while( !info.last ){
1321
1322 headTorsionType->add( info );
1323 recieveFrcStruct( &info, mpiTorsionStructType );
1324 }
1325 }
1326 #endif // is_mpi
1327
1328 // initialize the Torsions
1329
1330 for( i=0; i<nTorsions; i++ ){
1331
1332 a = the_torsions[i].a;
1333 b = the_torsions[i].b;
1334 c = the_torsions[i].c;
1335 d = the_torsions[i].d;
1336
1337 atomA = the_atoms[a]->getType();
1338 atomB = the_atoms[b]->getType();
1339 atomC = the_atoms[c]->getType();
1340 atomD = the_atoms[d]->getType();
1341 currentTorsionType = headTorsionType->find( atomA, atomB, atomC, atomD );
1342 if( currentTorsionType == NULL ){
1343 sprintf( painCave.errMsg,
1344 "TorsionType error, %s - %s - %s - %s not found"
1345 " in force file.\n",
1346 atomA, atomB, atomC, atomD );
1347 painCave.isFatal = 1;
1348 simError();
1349 }
1350
1351 if( !strcmp( currentTorsionType->type, "cubic" ) ){
1352 index = i + entry_plug->n_bonds + entry_plug->n_bends;
1353
1354 cTors = new CubicTorsion( *the_atoms[a], *the_atoms[b],
1355 *the_atoms[c], *the_atoms[d] );
1356 cTors->setConstants( currentTorsionType->k1, currentTorsionType->k2,
1357 currentTorsionType->k3, currentTorsionType->k4 );
1358 the_sris[index] = cTors;
1359 }
1360 }
1361
1362
1363 // clean up the memory
1364
1365 delete headTorsionType;
1366
1367 #ifdef IS_MPI
1368 sprintf( checkPointMsg, "TraPPE_Ex torsions initialized succesfully" );
1369 MPIcheckPoint();
1370 #endif // is_mpi
1371
1372 }
1373
1374 void TraPPE_ExFF::fastForward( char* stopText, char* searchOwner ){
1375
1376 int foundText = 0;
1377 char* the_token;
1378
1379 rewind( frcFile );
1380 lineNum = 0;
1381
1382 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1383 lineNum++;
1384 if( eof_test == NULL ){
1385 sprintf( painCave.errMsg, "Error fast forwarding force file for %s: "
1386 " file is empty.\n",
1387 searchOwner );
1388 painCave.isFatal = 1;
1389 simError();
1390 }
1391
1392
1393 while( !foundText ){
1394 while( eof_test != NULL && readLine[0] != '#' ){
1395 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1396 lineNum++;
1397 }
1398 if( eof_test == NULL ){
1399 sprintf( painCave.errMsg,
1400 "Error fast forwarding force file for %s at "
1401 "line %d: file ended unexpectedly.\n",
1402 searchOwner,
1403 lineNum );
1404 painCave.isFatal = 1;
1405 simError();
1406 }
1407
1408 the_token = strtok( readLine, " ,;\t#\n" );
1409 foundText = !strcmp( stopText, the_token );
1410
1411 if( !foundText ){
1412 eof_test = fgets( readLine, sizeof(readLine), frcFile );
1413 lineNum++;
1414
1415 if( eof_test == NULL ){
1416 sprintf( painCave.errMsg,
1417 "Error fast forwarding force file for %s at "
1418 "line %d: file ended unexpectedly.\n",
1419 searchOwner,
1420 lineNum );
1421 painCave.isFatal = 1;
1422 simError();
1423 }
1424 }
1425 }
1426 }
1427
1428
1429 int parseAtom( char *lineBuffer, int lineNum, atomStruct &info ){
1430
1431 char* the_token;
1432
1433 the_token = strtok( lineBuffer, " \n\t,;" );
1434 if( the_token != NULL ){
1435
1436 strcpy( info.name, the_token );
1437
1438 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1439 sprintf( painCave.errMsg,
1440 "Error parseing AtomTypes: line %d\n", lineNum );
1441 painCave.isFatal = 1;
1442 simError();
1443 }
1444
1445 info.isDipole = atoi( the_token );
1446
1447 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1448 sprintf( painCave.errMsg,
1449 "Error parseing AtomTypes: line %d\n", lineNum );
1450 painCave.isFatal = 1;
1451 simError();
1452 }
1453
1454 info.isSSD = atoi( the_token );
1455
1456 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1457 sprintf( painCave.errMsg,
1458 "Error parseing AtomTypes: line %d\n", lineNum );
1459 painCave.isFatal = 1;
1460 simError();
1461 }
1462
1463 info.mass = atof( the_token );
1464
1465 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1466 sprintf( painCave.errMsg,
1467 "Error parseing AtomTypes: line %d\n", lineNum );
1468 painCave.isFatal = 1;
1469 simError();
1470 }
1471
1472 info.epslon = atof( the_token );
1473
1474 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1475 sprintf( painCave.errMsg,
1476 "Error parseing AtomTypes: line %d\n", lineNum );
1477 painCave.isFatal = 1;
1478 simError();
1479 }
1480
1481 info.sigma = atof( the_token );
1482
1483 if( info.isDipole ){
1484
1485 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1486 sprintf( painCave.errMsg,
1487 "Error parseing AtomTypes: line %d\n", lineNum );
1488 painCave.isFatal = 1;
1489 simError();
1490 }
1491
1492 info.dipole = atof( the_token );
1493 }
1494 else info.dipole = 0.0;
1495
1496 if( info.isSSD ){
1497
1498 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1499 sprintf( painCave.errMsg,
1500 "Error parseing AtomTypes: line %d\n", lineNum );
1501 painCave.isFatal = 1;
1502 simError();
1503 }
1504
1505 info.w0 = atof( the_token );
1506
1507 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1508 sprintf( painCave.errMsg,
1509 "Error parseing AtomTypes: line %d\n", lineNum );
1510 painCave.isFatal = 1;
1511 simError();
1512 }
1513
1514 info.v0 = atof( the_token );
1515 }
1516 else info.v0 = info.w0 = 0.0;
1517
1518 return 1;
1519 }
1520 else return 0;
1521 }
1522
1523 int parseBond( char *lineBuffer, int lineNum, bondStruct &info ){
1524
1525 char* the_token;
1526
1527 the_token = strtok( lineBuffer, " \n\t,;" );
1528 if( the_token != NULL ){
1529
1530 strcpy( info.nameA, the_token );
1531
1532 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1533 sprintf( painCave.errMsg,
1534 "Error parseing BondTypes: line %d\n", lineNum );
1535 painCave.isFatal = 1;
1536 simError();
1537 }
1538
1539 strcpy( info.nameB, the_token );
1540
1541 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1542 sprintf( painCave.errMsg,
1543 "Error parseing BondTypes: line %d\n", lineNum );
1544 painCave.isFatal = 1;
1545 simError();
1546 }
1547
1548 strcpy( info.type, the_token );
1549
1550 if( !strcmp( info.type, "fixed" ) ){
1551 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1552 sprintf( painCave.errMsg,
1553 "Error parseing BondTypes: line %d\n", lineNum );
1554 painCave.isFatal = 1;
1555 simError();
1556 }
1557
1558 info.d0 = atof( the_token );
1559 }
1560 else{
1561 sprintf( painCave.errMsg,
1562 "Unknown TraPPE_Ex bond type \"%s\" at line %d\n",
1563 info.type,
1564 lineNum );
1565 painCave.isFatal = 1;
1566 simError();
1567 }
1568
1569 return 1;
1570 }
1571 else return 0;
1572 }
1573
1574
1575 int parseBend( char *lineBuffer, int lineNum, bendStruct &info ){
1576
1577 char* the_token;
1578
1579 the_token = strtok( lineBuffer, " \n\t,;" );
1580 if( the_token != NULL ){
1581
1582 strcpy( info.nameA, the_token );
1583
1584 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1585 sprintf( painCave.errMsg,
1586 "Error parseing BondTypes: line %d\n", lineNum );
1587 painCave.isFatal = 1;
1588 simError();
1589 }
1590
1591 strcpy( info.nameB, the_token );
1592
1593 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1594 sprintf( painCave.errMsg,
1595 "Error parseing BondTypes: line %d\n", lineNum );
1596 painCave.isFatal = 1;
1597 simError();
1598 }
1599
1600 strcpy( info.nameC, the_token );
1601
1602 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1603 sprintf( painCave.errMsg,
1604 "Error parseing BondTypes: line %d\n", lineNum );
1605 painCave.isFatal = 1;
1606 simError();
1607 }
1608
1609 strcpy( info.type, the_token );
1610
1611 if( !strcmp( info.type, "quadratic" ) ){
1612 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1613 sprintf( painCave.errMsg,
1614 "Error parseing BendTypes: line %d\n", lineNum );
1615 painCave.isFatal = 1;
1616 simError();
1617 }
1618
1619 info.k1 = atof( the_token );
1620
1621 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1622 sprintf( painCave.errMsg,
1623 "Error parseing BendTypes: line %d\n", lineNum );
1624 painCave.isFatal = 1;
1625 simError();
1626 }
1627
1628 info.k2 = atof( the_token );
1629
1630 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1631 sprintf( painCave.errMsg,
1632 "Error parseing BendTypes: line %d\n", lineNum );
1633 painCave.isFatal = 1;
1634 simError();
1635 }
1636
1637 info.k3 = atof( the_token );
1638
1639 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1640 sprintf( painCave.errMsg,
1641 "Error parseing BendTypes: line %d\n", lineNum );
1642 painCave.isFatal = 1;
1643 simError();
1644 }
1645
1646 info.t0 = atof( the_token );
1647 }
1648
1649 else{
1650 sprintf( painCave.errMsg,
1651 "Unknown TraPPE_Ex bend type \"%s\" at line %d\n",
1652 info.type,
1653 lineNum );
1654 painCave.isFatal = 1;
1655 simError();
1656 }
1657
1658 return 1;
1659 }
1660 else return 0;
1661 }
1662
1663 int parseTorsion( char *lineBuffer, int lineNum, torsionStruct &info ){
1664
1665 char* the_token;
1666
1667 the_token = strtok( lineBuffer, " \n\t,;" );
1668 if( the_token != NULL ){
1669
1670 strcpy( info.nameA, the_token );
1671
1672 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1673 sprintf( painCave.errMsg,
1674 "Error parseing TorsionTypes: line %d\n", lineNum );
1675 painCave.isFatal = 1;
1676 simError();
1677 }
1678
1679 strcpy( info.nameB, the_token );
1680
1681 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1682 sprintf( painCave.errMsg,
1683 "Error parseing TorsionTypes: line %d\n", lineNum );
1684 painCave.isFatal = 1;
1685 simError();
1686 }
1687
1688 strcpy( info.nameC, the_token );
1689
1690 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1691 sprintf( painCave.errMsg,
1692 "Error parseing TorsionTypes: line %d\n", lineNum );
1693 painCave.isFatal = 1;
1694 simError();
1695 }
1696
1697 strcpy( info.nameD, the_token );
1698
1699 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1700 sprintf( painCave.errMsg,
1701 "Error parseing TorsionTypes: line %d\n", lineNum );
1702 painCave.isFatal = 1;
1703 simError();
1704 }
1705
1706 strcpy( info.type, the_token );
1707
1708 if( !strcmp( info.type, "cubic" ) ){
1709 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1710 sprintf( painCave.errMsg,
1711 "Error parseing TorsionTypes: line %d\n", lineNum );
1712 painCave.isFatal = 1;
1713 simError();
1714 }
1715
1716 info.k1 = atof( the_token );
1717
1718 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1719 sprintf( painCave.errMsg,
1720 "Error parseing TorsionTypes: line %d\n", lineNum );
1721 painCave.isFatal = 1;
1722 simError();
1723 }
1724
1725 info.k2 = atof( the_token );
1726
1727 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1728 sprintf( painCave.errMsg,
1729 "Error parseing TorsionTypes: line %d\n", lineNum );
1730 painCave.isFatal = 1;
1731 simError();
1732 }
1733
1734 info.k3 = atof( the_token );
1735
1736 if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1737 sprintf( painCave.errMsg,
1738 "Error parseing TorsionTypes: line %d\n", lineNum );
1739 painCave.isFatal = 1;
1740 simError();
1741 }
1742
1743 info.k4 = atof( the_token );
1744
1745 }
1746
1747 else{
1748 sprintf( painCave.errMsg,
1749 "Unknown TraPPE_Ex torsion type \"%s\" at line %d\n",
1750 info.type,
1751 lineNum );
1752 painCave.isFatal = 1;
1753 simError();
1754 }
1755
1756 return 1;
1757 }
1758
1759 else return 0;
1760 }