ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/TraPPE_ExFF.cpp
Revision: 493
Committed: Mon Apr 14 19:03:38 2003 UTC (21 years, 2 months ago) by mmeineke
File size: 40525 byte(s)
Log Message: 
added Ghost bends to the TraPPE_Ex forceField

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 gezelter 463 if(currentAtomType->isSSD) {
765     entry_plug->useSticky = 1;
766     set_sticky_params( &(currentAtomType->w0), &(currentAtomType->v0));
767     }
768 mmeineke 377
769     if( currentAtomType->name[0] != '\0' ){
770     isError = 0;
771     makeAtype( &(currentAtomType->ident),
772     &isLJ,
773     &(currentAtomType->isSSD),
774     &(currentAtomType->isDipole),
775     &isGB,
776     &(currentAtomType->epslon),
777     &(currentAtomType->sigma),
778     &(currentAtomType->dipole),
779     &isError );
780     if( isError ){
781     sprintf( painCave.errMsg,
782     "Error initializing the \"%s\" atom type in fortran\n",
783     currentAtomType->name );
784     painCave.isFatal = 1;
785     simError();
786     }
787     }
788     currentAtomType = currentAtomType->next;
789     }
790    
791     #ifdef IS_MPI
792     sprintf( checkPointMsg,
793     "TraPPE_ExFF atom structures successfully sent to fortran\n" );
794     MPIcheckPoint();
795     #endif // is_mpi
796    
797    
798 mmeineke 420
799     // read in the bonds
800    
801 mmeineke 421 #ifdef IS_MPI
802     if( worldRank == 0 ){
803     #endif
804    
805     // read in the bond types.
806    
807     headBondType = new LinkedBondType;
808    
809     fastForward( "BondTypes", "initializeBonds" );
810 mmeineke 420
811 mmeineke 421 // we are now at the bondTypes section
812 mmeineke 420
813 mmeineke 421 eof_test = fgets( readLine, sizeof(readLine), frcFile );
814     lineNum++;
815    
816    
817     // read a line, and start parseing out the atom types
818 mmeineke 420
819 mmeineke 421 if( eof_test == NULL ){
820     sprintf( painCave.errMsg,
821     "Error in reading bonds from force file at line %d.\n",
822     lineNum );
823     painCave.isFatal = 1;
824     simError();
825     }
826    
827     // stop reading at end of file, or at next section
828     while( readLine[0] != '#' && eof_test != NULL ){
829 mmeineke 420
830 mmeineke 421 // toss comment lines
831     if( readLine[0] != '!' ){
832    
833     // the parser returns 0 if the line was blank
834     if( parseBond( readLine, lineNum, bondInfo ) ){
835     headBondType->add( bondInfo );
836     }
837     }
838     eof_test = fgets( readLine, sizeof(readLine), frcFile );
839     lineNum++;
840     }
841    
842     #ifdef IS_MPI
843    
844     // send out the linked list to all the other processes
845    
846     sprintf( checkPointMsg,
847     "TraPPE_Ex bond structures read successfully." );
848     MPIcheckPoint();
849    
850 chuckv 438 currentBondType = headBondType->next;
851 mmeineke 421 while( currentBondType != NULL ){
852     currentBondType->duplicate( bondInfo );
853     sendFrcStruct( &bondInfo, mpiBondStructType );
854     currentBondType = currentBondType->next;
855     }
856     bondInfo.last = 1;
857     sendFrcStruct( &bondInfo, mpiBondStructType );
858    
859     }
860    
861     else{
862    
863     // listen for node 0 to send out the force params
864    
865 chuckv 438 MPIcheckPoint();
866 mmeineke 421
867     headBondType = new LinkedBondType;
868     recieveFrcStruct( &bondInfo, mpiBondStructType );
869     while( !bondInfo.last ){
870    
871     headBondType->add( bondInfo );
872     recieveFrcStruct( &bondInfo, mpiBondStructType );
873     }
874     }
875     #endif // is_mpi
876 mmeineke 420
877    
878 mmeineke 421 // read in the bends
879 mmeineke 420
880 mmeineke 421 #ifdef IS_MPI
881     if( worldRank == 0 ){
882     #endif
883    
884     // read in the bend types.
885    
886     headBendType = new LinkedBendType;
887    
888     fastForward( "BendTypes", "initializeBends" );
889    
890     // we are now at the bendTypes section
891    
892     eof_test = fgets( readLine, sizeof(readLine), frcFile );
893     lineNum++;
894    
895     // read a line, and start parseing out the bend types
896    
897     if( eof_test == NULL ){
898     sprintf( painCave.errMsg,
899     "Error in reading bends from force file at line %d.\n",
900     lineNum );
901     painCave.isFatal = 1;
902     simError();
903     }
904    
905     // stop reading at end of file, or at next section
906     while( readLine[0] != '#' && eof_test != NULL ){
907    
908     // toss comment lines
909     if( readLine[0] != '!' ){
910    
911     // the parser returns 0 if the line was blank
912     if( parseBend( readLine, lineNum, bendInfo ) ){
913     headBendType->add( bendInfo );
914     }
915     }
916     eof_test = fgets( readLine, sizeof(readLine), frcFile );
917     lineNum++;
918     }
919    
920     #ifdef IS_MPI
921    
922     // send out the linked list to all the other processes
923    
924     sprintf( checkPointMsg,
925     "TraPPE_Ex bend structures read successfully." );
926     MPIcheckPoint();
927    
928 chuckv 438 currentBendType = headBendType->next;
929 mmeineke 421 while( currentBendType != NULL ){
930     currentBendType->duplicate( bendInfo );
931     sendFrcStruct( &bendInfo, mpiBendStructType );
932     currentBendType = currentBendType->next;
933     }
934     bendInfo.last = 1;
935     sendFrcStruct( &bendInfo, mpiBendStructType );
936    
937     }
938    
939     else{
940    
941     // listen for node 0 to send out the force params
942    
943     MPIcheckPoint();
944    
945     headBendType = new LinkedBendType;
946     recieveFrcStruct( &bendInfo, mpiBendStructType );
947     while( !bendInfo.last ){
948    
949     headBendType->add( bendInfo );
950     recieveFrcStruct( &bendInfo, mpiBendStructType );
951     }
952     }
953     #endif // is_mpi
954    
955    
956     // read in the torsions
957    
958     #ifdef IS_MPI
959     if( worldRank == 0 ){
960     #endif
961    
962     // read in the torsion types.
963    
964     headTorsionType = new LinkedTorsionType;
965    
966     fastForward( "TorsionTypes", "initializeTorsions" );
967    
968     // we are now at the torsionTypes section
969    
970     eof_test = fgets( readLine, sizeof(readLine), frcFile );
971     lineNum++;
972    
973    
974     // read a line, and start parseing out the atom types
975    
976     if( eof_test == NULL ){
977     sprintf( painCave.errMsg,
978     "Error in reading torsions from force file at line %d.\n",
979     lineNum );
980     painCave.isFatal = 1;
981     simError();
982     }
983    
984     // stop reading at end of file, or at next section
985     while( readLine[0] != '#' && eof_test != NULL ){
986    
987     // toss comment lines
988     if( readLine[0] != '!' ){
989    
990     // the parser returns 0 if the line was blank
991     if( parseTorsion( readLine, lineNum, torsionInfo ) ){
992     headTorsionType->add( torsionInfo );
993    
994     }
995     }
996     eof_test = fgets( readLine, sizeof(readLine), frcFile );
997     lineNum++;
998     }
999    
1000     #ifdef IS_MPI
1001    
1002     // send out the linked list to all the other processes
1003    
1004     sprintf( checkPointMsg,
1005     "TraPPE_Ex torsion structures read successfully." );
1006     MPIcheckPoint();
1007    
1008 chuckv 438 currentTorsionType = headTorsionType->next;
1009 mmeineke 421 while( currentTorsionType != NULL ){
1010     currentTorsionType->duplicate( torsionInfo );
1011     sendFrcStruct( &torsionInfo, mpiTorsionStructType );
1012     currentTorsionType = currentTorsionType->next;
1013     }
1014     torsionInfo.last = 1;
1015     sendFrcStruct( &torsionInfo, mpiTorsionStructType );
1016    
1017     }
1018    
1019     else{
1020    
1021     // listen for node 0 to send out the force params
1022    
1023     MPIcheckPoint();
1024    
1025     headTorsionType = new LinkedTorsionType;
1026     recieveFrcStruct( &torsionInfo, mpiTorsionStructType );
1027     while( !torsionInfo.last ){
1028    
1029     headTorsionType->add( torsionInfo );
1030     recieveFrcStruct( &torsionInfo, mpiTorsionStructType );
1031     }
1032     }
1033     #endif // is_mpi
1034    
1035 mmeineke 433 entry_plug->useLJ = 1;
1036 mmeineke 421 }
1037    
1038    
1039    
1040     void TraPPE_ExFF::initializeAtoms( int nAtoms, Atom** the_atoms ){
1041 mmeineke 420
1042    
1043     //////////////////////////////////////////////////
1044     // a quick water fix
1045    
1046     double waterI[3][3];
1047     waterI[0][0] = 1.76958347772500;
1048     waterI[0][1] = 0.0;
1049     waterI[0][2] = 0.0;
1050    
1051     waterI[1][0] = 0.0;
1052     waterI[1][1] = 0.614537057924513;
1053     waterI[1][2] = 0.0;
1054    
1055     waterI[2][0] = 0.0;
1056     waterI[2][1] = 0.0;
1057     waterI[2][2] = 1.15504641980049;
1058    
1059    
1060     double headI[3][3];
1061     headI[0][0] = 1125;
1062     headI[0][1] = 0.0;
1063     headI[0][2] = 0.0;
1064    
1065     headI[1][0] = 0.0;
1066     headI[1][1] = 1125;
1067     headI[1][2] = 0.0;
1068    
1069     headI[2][0] = 0.0;
1070     headI[2][1] = 0.0;
1071     headI[2][2] = 250;
1072    
1073     //////////////////////////////////////////////////
1074    
1075    
1076 mmeineke 377 // initialize the atoms
1077    
1078     DirectionalAtom* dAtom;
1079    
1080 mmeineke 428 for(int i=0; i<nAtoms; i++ ){
1081 gezelter 394
1082 mmeineke 377 currentAtomType = headAtomType->find( the_atoms[i]->getType() );
1083     if( currentAtomType == NULL ){
1084     sprintf( painCave.errMsg,
1085     "AtomType error, %s not found in force file.\n",
1086     the_atoms[i]->getType() );
1087     painCave.isFatal = 1;
1088     simError();
1089     }
1090    
1091     the_atoms[i]->setMass( currentAtomType->mass );
1092     the_atoms[i]->setEpslon( currentAtomType->epslon );
1093     the_atoms[i]->setSigma( currentAtomType->sigma );
1094     the_atoms[i]->setIdent( currentAtomType->ident );
1095     the_atoms[i]->setLJ();
1096    
1097     if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
1098    
1099     if( currentAtomType->isDipole ){
1100     if( the_atoms[i]->isDirectional() ){
1101    
1102     dAtom = (DirectionalAtom *) the_atoms[i];
1103     dAtom->setMu( currentAtomType->dipole );
1104     dAtom->setHasDipole( 1 );
1105     dAtom->setJx( 0.0 );
1106     dAtom->setJy( 0.0 );
1107     dAtom->setJz( 0.0 );
1108    
1109     if(!strcmp("SSD",the_atoms[i]->getType())){
1110     dAtom->setI( waterI );
1111     dAtom->setSSD( 1 );
1112     }
1113     else if(!strcmp("HEAD",the_atoms[i]->getType())){
1114     dAtom->setI( headI );
1115     dAtom->setSSD( 0 );
1116     }
1117     else{
1118     sprintf(painCave.errMsg,
1119     "AtmType error, %s does not have a moment of inertia set.\n",
1120     the_atoms[i]->getType() );
1121     painCave.isFatal = 1;
1122     simError();
1123     }
1124     entry_plug->n_dipoles++;
1125     }
1126     else{
1127    
1128     sprintf( painCave.errMsg,
1129     "TraPPE_ExFF error: Atom \"%s\" is a dipole, yet no standard"
1130     " orientation was specifed in the BASS file.\n",
1131     currentAtomType->name );
1132     painCave.isFatal = 1;
1133     simError();
1134     }
1135     }
1136     else{
1137     if( the_atoms[i]->isDirectional() ){
1138     sprintf( painCave.errMsg,
1139     "TraPPE_ExFF error: Atom \"%s\" was given a standard"
1140     "orientation in the BASS file, yet it is not a dipole.\n",
1141     currentAtomType->name);
1142     painCave.isFatal = 1;
1143     simError();
1144     }
1145     }
1146     }
1147     }
1148    
1149 mmeineke 421 void TraPPE_ExFF::initializeBonds( int nBonds, Bond** bondArray,
1150     bond_pair* the_bonds ){
1151     int i,a,b;
1152     char* atomA;
1153     char* atomB;
1154 mmeineke 377
1155     Atom** the_atoms;
1156     the_atoms = entry_plug->atoms;
1157    
1158    
1159     // initialize the Bonds
1160    
1161     for( i=0; i<nBonds; i++ ){
1162    
1163     a = the_bonds[i].a;
1164     b = the_bonds[i].b;
1165    
1166     atomA = the_atoms[a]->getType();
1167     atomB = the_atoms[b]->getType();
1168     currentBondType = headBondType->find( atomA, atomB );
1169     if( currentBondType == NULL ){
1170     sprintf( painCave.errMsg,
1171     "BondType error, %s - %s not found in force file.\n",
1172     atomA, atomB );
1173     painCave.isFatal = 1;
1174     simError();
1175     }
1176    
1177     if( !strcmp( currentBondType->type, "fixed" ) ){
1178    
1179 mmeineke 421 bondArray[i] = new ConstrainedBond( *the_atoms[a],
1180     *the_atoms[b],
1181     currentBondType->d0 );
1182 mmeineke 377 entry_plug->n_constraints++;
1183     }
1184     }
1185     }
1186    
1187 mmeineke 421 void TraPPE_ExFF::initializeBends( int nBends, Bend** bendArray,
1188     bend_set* the_bends ){
1189 mmeineke 377
1190     QuadraticBend* qBend;
1191     GhostBend* gBend;
1192     Atom** the_atoms;
1193     the_atoms = entry_plug->atoms;
1194 mmeineke 421
1195 mmeineke 377 int i, a, b, c;
1196     char* atomA;
1197     char* atomB;
1198     char* atomC;
1199    
1200     // initialize the Bends
1201    
1202     for( i=0; i<nBends; i++ ){
1203    
1204     a = the_bends[i].a;
1205     b = the_bends[i].b;
1206     c = the_bends[i].c;
1207    
1208     atomA = the_atoms[a]->getType();
1209     atomB = the_atoms[b]->getType();
1210    
1211     if( the_bends[i].isGhost ) atomC = "GHOST";
1212     else atomC = the_atoms[c]->getType();
1213    
1214     currentBendType = headBendType->find( atomA, atomB, atomC );
1215     if( currentBendType == NULL ){
1216     sprintf( painCave.errMsg, "BendType error, %s - %s - %s not found"
1217     " in force file.\n",
1218     atomA, atomB, atomC );
1219     painCave.isFatal = 1;
1220     simError();
1221     }
1222    
1223     if( !strcmp( currentBendType->type, "quadratic" ) ){
1224    
1225     if( the_bends[i].isGhost){
1226    
1227     if( the_bends[i].ghost == b ){
1228     // do nothing
1229     }
1230     else if( the_bends[i].ghost == a ){
1231     c = a;
1232     a = b;
1233 mmeineke 493 b = c;
1234 mmeineke 377 }
1235     else{
1236     sprintf( painCave.errMsg,
1237     "BendType error, %s - %s - %s,\n"
1238     " --> central atom is not "
1239     "correctly identified with the "
1240     "\"ghostVectorSource = \" tag.\n",
1241     atomA, atomB, atomC );
1242     painCave.isFatal = 1;
1243     simError();
1244     }
1245    
1246     gBend = new GhostBend( *the_atoms[a],
1247     *the_atoms[b] );
1248     gBend->setConstants( currentBendType->k1,
1249     currentBendType->k2,
1250     currentBendType->k3,
1251     currentBendType->t0 );
1252 mmeineke 421 bendArray[i] = gBend;
1253 mmeineke 377 }
1254     else{
1255     qBend = new QuadraticBend( *the_atoms[a],
1256     *the_atoms[b],
1257     *the_atoms[c] );
1258     qBend->setConstants( currentBendType->k1,
1259     currentBendType->k2,
1260     currentBendType->k3,
1261     currentBendType->t0 );
1262 mmeineke 421 bendArray[i] = qBend;
1263 mmeineke 377 }
1264     }
1265     }
1266     }
1267    
1268 mmeineke 421 void TraPPE_ExFF::initializeTorsions( int nTorsions, Torsion** torsionArray,
1269     torsion_set* the_torsions ){
1270 mmeineke 377
1271 mmeineke 421 int i, a, b, c, d;
1272 mmeineke 377 char* atomA;
1273     char* atomB;
1274     char* atomC;
1275     char* atomD;
1276 mmeineke 421
1277 mmeineke 377 CubicTorsion* cTors;
1278     Atom** the_atoms;
1279     the_atoms = entry_plug->atoms;
1280    
1281     // initialize the Torsions
1282    
1283     for( i=0; i<nTorsions; i++ ){
1284    
1285     a = the_torsions[i].a;
1286     b = the_torsions[i].b;
1287     c = the_torsions[i].c;
1288     d = the_torsions[i].d;
1289    
1290     atomA = the_atoms[a]->getType();
1291     atomB = the_atoms[b]->getType();
1292     atomC = the_atoms[c]->getType();
1293     atomD = the_atoms[d]->getType();
1294     currentTorsionType = headTorsionType->find( atomA, atomB, atomC, atomD );
1295     if( currentTorsionType == NULL ){
1296     sprintf( painCave.errMsg,
1297     "TorsionType error, %s - %s - %s - %s not found"
1298     " in force file.\n",
1299     atomA, atomB, atomC, atomD );
1300     painCave.isFatal = 1;
1301     simError();
1302     }
1303    
1304     if( !strcmp( currentTorsionType->type, "cubic" ) ){
1305    
1306     cTors = new CubicTorsion( *the_atoms[a], *the_atoms[b],
1307     *the_atoms[c], *the_atoms[d] );
1308     cTors->setConstants( currentTorsionType->k1, currentTorsionType->k2,
1309     currentTorsionType->k3, currentTorsionType->k4 );
1310 mmeineke 421 torsionArray[i] = cTors;
1311 mmeineke 377 }
1312     }
1313     }
1314    
1315     void TraPPE_ExFF::fastForward( char* stopText, char* searchOwner ){
1316    
1317     int foundText = 0;
1318     char* the_token;
1319    
1320     rewind( frcFile );
1321     lineNum = 0;
1322    
1323     eof_test = fgets( readLine, sizeof(readLine), frcFile );
1324     lineNum++;
1325     if( eof_test == NULL ){
1326     sprintf( painCave.errMsg, "Error fast forwarding force file for %s: "
1327     " file is empty.\n",
1328     searchOwner );
1329     painCave.isFatal = 1;
1330     simError();
1331     }
1332    
1333    
1334     while( !foundText ){
1335     while( eof_test != NULL && readLine[0] != '#' ){
1336     eof_test = fgets( readLine, sizeof(readLine), frcFile );
1337     lineNum++;
1338     }
1339     if( eof_test == NULL ){
1340     sprintf( painCave.errMsg,
1341     "Error fast forwarding force file for %s at "
1342     "line %d: file ended unexpectedly.\n",
1343     searchOwner,
1344     lineNum );
1345     painCave.isFatal = 1;
1346     simError();
1347     }
1348    
1349     the_token = strtok( readLine, " ,;\t#\n" );
1350     foundText = !strcmp( stopText, the_token );
1351    
1352     if( !foundText ){
1353     eof_test = fgets( readLine, sizeof(readLine), frcFile );
1354     lineNum++;
1355    
1356     if( eof_test == NULL ){
1357     sprintf( painCave.errMsg,
1358     "Error fast forwarding force file for %s at "
1359     "line %d: file ended unexpectedly.\n",
1360     searchOwner,
1361     lineNum );
1362     painCave.isFatal = 1;
1363     simError();
1364     }
1365     }
1366     }
1367     }
1368    
1369    
1370     int TPE::parseAtom( char *lineBuffer, int lineNum, atomStruct &info ){
1371    
1372     char* the_token;
1373    
1374     the_token = strtok( lineBuffer, " \n\t,;" );
1375     if( the_token != NULL ){
1376    
1377     strcpy( info.name, the_token );
1378    
1379     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1380     sprintf( painCave.errMsg,
1381     "Error parseing AtomTypes: line %d\n", lineNum );
1382     painCave.isFatal = 1;
1383     simError();
1384     }
1385    
1386     info.isDipole = atoi( the_token );
1387    
1388     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1389     sprintf( painCave.errMsg,
1390     "Error parseing AtomTypes: line %d\n", lineNum );
1391     painCave.isFatal = 1;
1392     simError();
1393     }
1394    
1395     info.isSSD = atoi( the_token );
1396    
1397     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1398     sprintf( painCave.errMsg,
1399     "Error parseing AtomTypes: line %d\n", lineNum );
1400     painCave.isFatal = 1;
1401     simError();
1402     }
1403    
1404     info.mass = atof( the_token );
1405    
1406     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1407     sprintf( painCave.errMsg,
1408     "Error parseing AtomTypes: line %d\n", lineNum );
1409     painCave.isFatal = 1;
1410     simError();
1411     }
1412    
1413     info.epslon = atof( the_token );
1414    
1415     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1416     sprintf( painCave.errMsg,
1417     "Error parseing AtomTypes: line %d\n", lineNum );
1418     painCave.isFatal = 1;
1419     simError();
1420     }
1421    
1422     info.sigma = atof( the_token );
1423    
1424     if( info.isDipole ){
1425    
1426     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1427     sprintf( painCave.errMsg,
1428     "Error parseing AtomTypes: line %d\n", lineNum );
1429     painCave.isFatal = 1;
1430     simError();
1431     }
1432    
1433     info.dipole = atof( the_token );
1434     }
1435     else info.dipole = 0.0;
1436    
1437     if( info.isSSD ){
1438    
1439     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1440     sprintf( painCave.errMsg,
1441     "Error parseing AtomTypes: line %d\n", lineNum );
1442     painCave.isFatal = 1;
1443     simError();
1444     }
1445    
1446     info.w0 = atof( the_token );
1447    
1448     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1449     sprintf( painCave.errMsg,
1450     "Error parseing AtomTypes: line %d\n", lineNum );
1451     painCave.isFatal = 1;
1452     simError();
1453     }
1454    
1455     info.v0 = atof( the_token );
1456     }
1457     else info.v0 = info.w0 = 0.0;
1458    
1459     return 1;
1460     }
1461     else return 0;
1462     }
1463    
1464     int TPE::parseBond( char *lineBuffer, int lineNum, bondStruct &info ){
1465    
1466     char* the_token;
1467    
1468     the_token = strtok( lineBuffer, " \n\t,;" );
1469     if( the_token != NULL ){
1470    
1471     strcpy( info.nameA, the_token );
1472    
1473     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1474     sprintf( painCave.errMsg,
1475     "Error parseing BondTypes: line %d\n", lineNum );
1476     painCave.isFatal = 1;
1477     simError();
1478     }
1479    
1480     strcpy( info.nameB, the_token );
1481    
1482     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1483     sprintf( painCave.errMsg,
1484     "Error parseing BondTypes: line %d\n", lineNum );
1485     painCave.isFatal = 1;
1486     simError();
1487     }
1488    
1489     strcpy( info.type, the_token );
1490    
1491     if( !strcmp( info.type, "fixed" ) ){
1492     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1493     sprintf( painCave.errMsg,
1494     "Error parseing BondTypes: line %d\n", lineNum );
1495     painCave.isFatal = 1;
1496     simError();
1497     }
1498    
1499     info.d0 = atof( the_token );
1500     }
1501     else{
1502     sprintf( painCave.errMsg,
1503     "Unknown TraPPE_Ex bond type \"%s\" at line %d\n",
1504     info.type,
1505     lineNum );
1506     painCave.isFatal = 1;
1507     simError();
1508     }
1509    
1510     return 1;
1511     }
1512     else return 0;
1513     }
1514    
1515    
1516     int TPE::parseBend( char *lineBuffer, int lineNum, bendStruct &info ){
1517    
1518     char* the_token;
1519    
1520     the_token = strtok( lineBuffer, " \n\t,;" );
1521     if( the_token != NULL ){
1522    
1523     strcpy( info.nameA, the_token );
1524    
1525     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1526     sprintf( painCave.errMsg,
1527 chuckv 388 "Error parseing BendTypes: line %d\n", lineNum );
1528 mmeineke 377 painCave.isFatal = 1;
1529     simError();
1530     }
1531    
1532     strcpy( info.nameB, the_token );
1533    
1534     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1535     sprintf( painCave.errMsg,
1536 chuckv 388 "Error parseing BendTypes: line %d\n", lineNum );
1537 mmeineke 377 painCave.isFatal = 1;
1538     simError();
1539     }
1540    
1541     strcpy( info.nameC, the_token );
1542    
1543     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1544     sprintf( painCave.errMsg,
1545 chuckv 388 "Error parseing BendTypes: line %d\n", lineNum );
1546 mmeineke 377 painCave.isFatal = 1;
1547     simError();
1548     }
1549    
1550     strcpy( info.type, the_token );
1551    
1552     if( !strcmp( info.type, "quadratic" ) ){
1553     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1554     sprintf( painCave.errMsg,
1555     "Error parseing BendTypes: line %d\n", lineNum );
1556     painCave.isFatal = 1;
1557     simError();
1558     }
1559    
1560     info.k1 = atof( the_token );
1561    
1562     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1563     sprintf( painCave.errMsg,
1564     "Error parseing BendTypes: line %d\n", lineNum );
1565     painCave.isFatal = 1;
1566     simError();
1567     }
1568    
1569     info.k2 = atof( the_token );
1570    
1571     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1572     sprintf( painCave.errMsg,
1573     "Error parseing BendTypes: line %d\n", lineNum );
1574     painCave.isFatal = 1;
1575     simError();
1576     }
1577    
1578     info.k3 = atof( the_token );
1579    
1580     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1581     sprintf( painCave.errMsg,
1582     "Error parseing BendTypes: line %d\n", lineNum );
1583     painCave.isFatal = 1;
1584     simError();
1585     }
1586    
1587     info.t0 = atof( the_token );
1588     }
1589    
1590     else{
1591     sprintf( painCave.errMsg,
1592     "Unknown TraPPE_Ex bend type \"%s\" at line %d\n",
1593     info.type,
1594     lineNum );
1595     painCave.isFatal = 1;
1596     simError();
1597     }
1598    
1599     return 1;
1600     }
1601     else return 0;
1602     }
1603    
1604     int TPE::parseTorsion( char *lineBuffer, int lineNum, torsionStruct &info ){
1605    
1606     char* the_token;
1607    
1608     the_token = strtok( lineBuffer, " \n\t,;" );
1609     if( the_token != NULL ){
1610    
1611     strcpy( info.nameA, the_token );
1612    
1613     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1614     sprintf( painCave.errMsg,
1615     "Error parseing TorsionTypes: line %d\n", lineNum );
1616     painCave.isFatal = 1;
1617     simError();
1618     }
1619    
1620     strcpy( info.nameB, the_token );
1621    
1622     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1623     sprintf( painCave.errMsg,
1624     "Error parseing TorsionTypes: line %d\n", lineNum );
1625     painCave.isFatal = 1;
1626     simError();
1627     }
1628    
1629     strcpy( info.nameC, the_token );
1630    
1631     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1632     sprintf( painCave.errMsg,
1633     "Error parseing TorsionTypes: line %d\n", lineNum );
1634     painCave.isFatal = 1;
1635     simError();
1636     }
1637    
1638     strcpy( info.nameD, the_token );
1639    
1640     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1641     sprintf( painCave.errMsg,
1642     "Error parseing TorsionTypes: line %d\n", lineNum );
1643     painCave.isFatal = 1;
1644     simError();
1645     }
1646    
1647     strcpy( info.type, the_token );
1648    
1649     if( !strcmp( info.type, "cubic" ) ){
1650     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1651     sprintf( painCave.errMsg,
1652     "Error parseing TorsionTypes: line %d\n", lineNum );
1653     painCave.isFatal = 1;
1654     simError();
1655     }
1656    
1657     info.k1 = atof( the_token );
1658    
1659     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1660     sprintf( painCave.errMsg,
1661     "Error parseing TorsionTypes: line %d\n", lineNum );
1662     painCave.isFatal = 1;
1663     simError();
1664     }
1665    
1666     info.k2 = atof( the_token );
1667    
1668     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1669     sprintf( painCave.errMsg,
1670     "Error parseing TorsionTypes: line %d\n", lineNum );
1671     painCave.isFatal = 1;
1672     simError();
1673     }
1674    
1675     info.k3 = atof( the_token );
1676    
1677     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1678     sprintf( painCave.errMsg,
1679     "Error parseing TorsionTypes: line %d\n", lineNum );
1680     painCave.isFatal = 1;
1681     simError();
1682     }
1683    
1684     info.k4 = atof( the_token );
1685    
1686     }
1687    
1688     else{
1689     sprintf( painCave.errMsg,
1690     "Unknown TraPPE_Ex torsion type \"%s\" at line %d\n",
1691     info.type,
1692     lineNum );
1693     painCave.isFatal = 1;
1694     simError();
1695     }
1696    
1697     return 1;
1698     }
1699    
1700     else return 0;
1701     }