ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/TraPPE_ExFF.cpp
Revision: 438
Committed: Mon Mar 31 21:50:59 2003 UTC (21 years, 3 months ago) by chuckv
File size: 40399 byte(s)
Log Message:
Fixes in MPI force calc and in Trappe_Ex parsing.

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