ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/TraPPE_ExFF.cpp
Revision: 451
Committed: Thu Apr 3 22:01:29 2003 UTC (21 years, 3 months ago) by mmeineke
File size: 40578 byte(s)
Log Message: 
fixed a possible call before initialize bug.

File Contents

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