ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/TraPPE_ExFF.cpp
Revision: 433
Committed: Fri Mar 28 15:28:53 2003 UTC (21 years, 3 months ago) by mmeineke
File size: 40878 byte(s)
Log Message: 
fixed long range interactions in Trappe

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     currentBondType = headBondType;
846     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     MPIcheckPoint();
861    
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     currentBendType = headBendType;
924     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     currentTorsionType = headTorsionType;
1004     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     #ifdef IS_MPI
1145     sprintf( checkPointMsg, "TraPPE_Ex atoms initialized succesfully" );
1146     MPIcheckPoint();
1147     #endif // is_mpi
1148    
1149     }
1150    
1151 mmeineke 421 void TraPPE_ExFF::initializeBonds( int nBonds, Bond** bondArray,
1152     bond_pair* the_bonds ){
1153     int i,a,b;
1154     char* atomA;
1155     char* atomB;
1156 mmeineke 377
1157     Atom** the_atoms;
1158     the_atoms = entry_plug->atoms;
1159    
1160    
1161     // initialize the Bonds
1162    
1163     for( i=0; i<nBonds; i++ ){
1164    
1165     a = the_bonds[i].a;
1166     b = the_bonds[i].b;
1167    
1168     atomA = the_atoms[a]->getType();
1169     atomB = the_atoms[b]->getType();
1170     currentBondType = headBondType->find( atomA, atomB );
1171     if( currentBondType == NULL ){
1172     sprintf( painCave.errMsg,
1173     "BondType error, %s - %s not found in force file.\n",
1174     atomA, atomB );
1175     painCave.isFatal = 1;
1176     simError();
1177     }
1178    
1179     if( !strcmp( currentBondType->type, "fixed" ) ){
1180    
1181 mmeineke 421 bondArray[i] = new ConstrainedBond( *the_atoms[a],
1182     *the_atoms[b],
1183     currentBondType->d0 );
1184 mmeineke 377 entry_plug->n_constraints++;
1185     }
1186     }
1187    
1188     #ifdef IS_MPI
1189     sprintf( checkPointMsg, "TraPPE_Ex bonds initialized succesfully" );
1190     MPIcheckPoint();
1191     #endif // is_mpi
1192    
1193     }
1194    
1195 mmeineke 421 void TraPPE_ExFF::initializeBends( int nBends, Bend** bendArray,
1196     bend_set* the_bends ){
1197 mmeineke 377
1198     QuadraticBend* qBend;
1199     GhostBend* gBend;
1200     Atom** the_atoms;
1201     the_atoms = entry_plug->atoms;
1202 mmeineke 421
1203 mmeineke 377 int i, a, b, c;
1204     char* atomA;
1205     char* atomB;
1206     char* atomC;
1207    
1208     // initialize the Bends
1209    
1210     for( i=0; i<nBends; i++ ){
1211    
1212     a = the_bends[i].a;
1213     b = the_bends[i].b;
1214     c = the_bends[i].c;
1215    
1216     atomA = the_atoms[a]->getType();
1217     atomB = the_atoms[b]->getType();
1218    
1219     if( the_bends[i].isGhost ) atomC = "GHOST";
1220     else atomC = the_atoms[c]->getType();
1221    
1222     currentBendType = headBendType->find( atomA, atomB, atomC );
1223     if( currentBendType == NULL ){
1224     sprintf( painCave.errMsg, "BendType error, %s - %s - %s not found"
1225     " in force file.\n",
1226     atomA, atomB, atomC );
1227     painCave.isFatal = 1;
1228     simError();
1229     }
1230    
1231     if( !strcmp( currentBendType->type, "quadratic" ) ){
1232    
1233     if( the_bends[i].isGhost){
1234    
1235     if( the_bends[i].ghost == b ){
1236     // do nothing
1237     }
1238     else if( the_bends[i].ghost == a ){
1239     c = a;
1240     a = b;
1241     b = a;
1242     }
1243     else{
1244     sprintf( painCave.errMsg,
1245     "BendType error, %s - %s - %s,\n"
1246     " --> central atom is not "
1247     "correctly identified with the "
1248     "\"ghostVectorSource = \" tag.\n",
1249     atomA, atomB, atomC );
1250     painCave.isFatal = 1;
1251     simError();
1252     }
1253    
1254     gBend = new GhostBend( *the_atoms[a],
1255     *the_atoms[b] );
1256     gBend->setConstants( currentBendType->k1,
1257     currentBendType->k2,
1258     currentBendType->k3,
1259     currentBendType->t0 );
1260 mmeineke 421 bendArray[i] = gBend;
1261 mmeineke 377 }
1262     else{
1263     qBend = new QuadraticBend( *the_atoms[a],
1264     *the_atoms[b],
1265     *the_atoms[c] );
1266     qBend->setConstants( currentBendType->k1,
1267     currentBendType->k2,
1268     currentBendType->k3,
1269     currentBendType->t0 );
1270 mmeineke 421 bendArray[i] = qBend;
1271 mmeineke 377 }
1272     }
1273     }
1274    
1275     #ifdef IS_MPI
1276     sprintf( checkPointMsg, "TraPPE_Ex bends initialized succesfully" );
1277     MPIcheckPoint();
1278     #endif // is_mpi
1279    
1280     }
1281    
1282 mmeineke 421 void TraPPE_ExFF::initializeTorsions( int nTorsions, Torsion** torsionArray,
1283     torsion_set* the_torsions ){
1284 mmeineke 377
1285 mmeineke 421 int i, a, b, c, d;
1286 mmeineke 377 char* atomA;
1287     char* atomB;
1288     char* atomC;
1289     char* atomD;
1290 mmeineke 421
1291 mmeineke 377 CubicTorsion* cTors;
1292     Atom** the_atoms;
1293     the_atoms = entry_plug->atoms;
1294    
1295     // initialize the Torsions
1296    
1297     for( i=0; i<nTorsions; i++ ){
1298    
1299     a = the_torsions[i].a;
1300     b = the_torsions[i].b;
1301     c = the_torsions[i].c;
1302     d = the_torsions[i].d;
1303    
1304     atomA = the_atoms[a]->getType();
1305     atomB = the_atoms[b]->getType();
1306     atomC = the_atoms[c]->getType();
1307     atomD = the_atoms[d]->getType();
1308     currentTorsionType = headTorsionType->find( atomA, atomB, atomC, atomD );
1309     if( currentTorsionType == NULL ){
1310     sprintf( painCave.errMsg,
1311     "TorsionType error, %s - %s - %s - %s not found"
1312     " in force file.\n",
1313     atomA, atomB, atomC, atomD );
1314     painCave.isFatal = 1;
1315     simError();
1316     }
1317    
1318     if( !strcmp( currentTorsionType->type, "cubic" ) ){
1319    
1320     cTors = new CubicTorsion( *the_atoms[a], *the_atoms[b],
1321     *the_atoms[c], *the_atoms[d] );
1322     cTors->setConstants( currentTorsionType->k1, currentTorsionType->k2,
1323     currentTorsionType->k3, currentTorsionType->k4 );
1324 mmeineke 421 torsionArray[i] = cTors;
1325 mmeineke 377 }
1326     }
1327    
1328     #ifdef IS_MPI
1329     sprintf( checkPointMsg, "TraPPE_Ex torsions initialized succesfully" );
1330     MPIcheckPoint();
1331     #endif // is_mpi
1332    
1333     }
1334    
1335     void TraPPE_ExFF::fastForward( char* stopText, char* searchOwner ){
1336    
1337     int foundText = 0;
1338     char* the_token;
1339    
1340     rewind( frcFile );
1341     lineNum = 0;
1342    
1343     eof_test = fgets( readLine, sizeof(readLine), frcFile );
1344     lineNum++;
1345     if( eof_test == NULL ){
1346     sprintf( painCave.errMsg, "Error fast forwarding force file for %s: "
1347     " file is empty.\n",
1348     searchOwner );
1349     painCave.isFatal = 1;
1350     simError();
1351     }
1352    
1353    
1354     while( !foundText ){
1355     while( eof_test != NULL && readLine[0] != '#' ){
1356     eof_test = fgets( readLine, sizeof(readLine), frcFile );
1357     lineNum++;
1358     }
1359     if( eof_test == NULL ){
1360     sprintf( painCave.errMsg,
1361     "Error fast forwarding force file for %s at "
1362     "line %d: file ended unexpectedly.\n",
1363     searchOwner,
1364     lineNum );
1365     painCave.isFatal = 1;
1366     simError();
1367     }
1368    
1369     the_token = strtok( readLine, " ,;\t#\n" );
1370     foundText = !strcmp( stopText, the_token );
1371    
1372     if( !foundText ){
1373     eof_test = fgets( readLine, sizeof(readLine), frcFile );
1374     lineNum++;
1375    
1376     if( eof_test == NULL ){
1377     sprintf( painCave.errMsg,
1378     "Error fast forwarding force file for %s at "
1379     "line %d: file ended unexpectedly.\n",
1380     searchOwner,
1381     lineNum );
1382     painCave.isFatal = 1;
1383     simError();
1384     }
1385     }
1386     }
1387     }
1388    
1389    
1390     int TPE::parseAtom( char *lineBuffer, int lineNum, atomStruct &info ){
1391    
1392     char* the_token;
1393    
1394     the_token = strtok( lineBuffer, " \n\t,;" );
1395     if( the_token != NULL ){
1396    
1397     strcpy( info.name, the_token );
1398    
1399     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1400     sprintf( painCave.errMsg,
1401     "Error parseing AtomTypes: line %d\n", lineNum );
1402     painCave.isFatal = 1;
1403     simError();
1404     }
1405    
1406     info.isDipole = atoi( the_token );
1407    
1408     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1409     sprintf( painCave.errMsg,
1410     "Error parseing AtomTypes: line %d\n", lineNum );
1411     painCave.isFatal = 1;
1412     simError();
1413     }
1414    
1415     info.isSSD = atoi( the_token );
1416    
1417     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1418     sprintf( painCave.errMsg,
1419     "Error parseing AtomTypes: line %d\n", lineNum );
1420     painCave.isFatal = 1;
1421     simError();
1422     }
1423    
1424     info.mass = atof( the_token );
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.epslon = atof( the_token );
1434    
1435     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1436     sprintf( painCave.errMsg,
1437     "Error parseing AtomTypes: line %d\n", lineNum );
1438     painCave.isFatal = 1;
1439     simError();
1440     }
1441    
1442     info.sigma = atof( the_token );
1443    
1444     if( info.isDipole ){
1445    
1446     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1447     sprintf( painCave.errMsg,
1448     "Error parseing AtomTypes: line %d\n", lineNum );
1449     painCave.isFatal = 1;
1450     simError();
1451     }
1452    
1453     info.dipole = atof( the_token );
1454     }
1455     else info.dipole = 0.0;
1456    
1457     if( info.isSSD ){
1458    
1459     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1460     sprintf( painCave.errMsg,
1461     "Error parseing AtomTypes: line %d\n", lineNum );
1462     painCave.isFatal = 1;
1463     simError();
1464     }
1465    
1466     info.w0 = atof( the_token );
1467    
1468     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1469     sprintf( painCave.errMsg,
1470     "Error parseing AtomTypes: line %d\n", lineNum );
1471     painCave.isFatal = 1;
1472     simError();
1473     }
1474    
1475     info.v0 = atof( the_token );
1476     }
1477     else info.v0 = info.w0 = 0.0;
1478    
1479     return 1;
1480     }
1481     else return 0;
1482     }
1483    
1484     int TPE::parseBond( char *lineBuffer, int lineNum, bondStruct &info ){
1485    
1486     char* the_token;
1487    
1488     the_token = strtok( lineBuffer, " \n\t,;" );
1489     if( the_token != NULL ){
1490    
1491     strcpy( info.nameA, the_token );
1492    
1493     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1494     sprintf( painCave.errMsg,
1495     "Error parseing BondTypes: line %d\n", lineNum );
1496     painCave.isFatal = 1;
1497     simError();
1498     }
1499    
1500     strcpy( info.nameB, the_token );
1501    
1502     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1503     sprintf( painCave.errMsg,
1504     "Error parseing BondTypes: line %d\n", lineNum );
1505     painCave.isFatal = 1;
1506     simError();
1507     }
1508    
1509     strcpy( info.type, the_token );
1510    
1511     if( !strcmp( info.type, "fixed" ) ){
1512     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1513     sprintf( painCave.errMsg,
1514     "Error parseing BondTypes: line %d\n", lineNum );
1515     painCave.isFatal = 1;
1516     simError();
1517     }
1518    
1519     info.d0 = atof( the_token );
1520     }
1521     else{
1522     sprintf( painCave.errMsg,
1523     "Unknown TraPPE_Ex bond type \"%s\" at line %d\n",
1524     info.type,
1525     lineNum );
1526     painCave.isFatal = 1;
1527     simError();
1528     }
1529    
1530     return 1;
1531     }
1532     else return 0;
1533     }
1534    
1535    
1536     int TPE::parseBend( char *lineBuffer, int lineNum, bendStruct &info ){
1537    
1538     char* the_token;
1539    
1540     the_token = strtok( lineBuffer, " \n\t,;" );
1541     if( the_token != NULL ){
1542    
1543     strcpy( info.nameA, the_token );
1544    
1545     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1546     sprintf( painCave.errMsg,
1547 chuckv 388 "Error parseing BendTypes: line %d\n", lineNum );
1548 mmeineke 377 painCave.isFatal = 1;
1549     simError();
1550     }
1551    
1552     strcpy( info.nameB, the_token );
1553    
1554     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1555     sprintf( painCave.errMsg,
1556 chuckv 388 "Error parseing BendTypes: line %d\n", lineNum );
1557 mmeineke 377 painCave.isFatal = 1;
1558     simError();
1559     }
1560    
1561     strcpy( info.nameC, the_token );
1562    
1563     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1564     sprintf( painCave.errMsg,
1565 chuckv 388 "Error parseing BendTypes: line %d\n", lineNum );
1566 mmeineke 377 painCave.isFatal = 1;
1567     simError();
1568     }
1569    
1570     strcpy( info.type, the_token );
1571    
1572     if( !strcmp( info.type, "quadratic" ) ){
1573     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1574     sprintf( painCave.errMsg,
1575     "Error parseing BendTypes: line %d\n", lineNum );
1576     painCave.isFatal = 1;
1577     simError();
1578     }
1579    
1580     info.k1 = atof( the_token );
1581    
1582     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1583     sprintf( painCave.errMsg,
1584     "Error parseing BendTypes: line %d\n", lineNum );
1585     painCave.isFatal = 1;
1586     simError();
1587     }
1588    
1589     info.k2 = atof( the_token );
1590    
1591     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1592     sprintf( painCave.errMsg,
1593     "Error parseing BendTypes: line %d\n", lineNum );
1594     painCave.isFatal = 1;
1595     simError();
1596     }
1597    
1598     info.k3 = atof( the_token );
1599    
1600     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1601     sprintf( painCave.errMsg,
1602     "Error parseing BendTypes: line %d\n", lineNum );
1603     painCave.isFatal = 1;
1604     simError();
1605     }
1606    
1607     info.t0 = atof( the_token );
1608     }
1609    
1610     else{
1611     sprintf( painCave.errMsg,
1612     "Unknown TraPPE_Ex bend type \"%s\" at line %d\n",
1613     info.type,
1614     lineNum );
1615     painCave.isFatal = 1;
1616     simError();
1617     }
1618    
1619     return 1;
1620     }
1621     else return 0;
1622     }
1623    
1624     int TPE::parseTorsion( char *lineBuffer, int lineNum, torsionStruct &info ){
1625    
1626     char* the_token;
1627    
1628     the_token = strtok( lineBuffer, " \n\t,;" );
1629     if( the_token != NULL ){
1630    
1631     strcpy( info.nameA, the_token );
1632    
1633     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1634     sprintf( painCave.errMsg,
1635     "Error parseing TorsionTypes: line %d\n", lineNum );
1636     painCave.isFatal = 1;
1637     simError();
1638     }
1639    
1640     strcpy( info.nameB, the_token );
1641    
1642     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1643     sprintf( painCave.errMsg,
1644     "Error parseing TorsionTypes: line %d\n", lineNum );
1645     painCave.isFatal = 1;
1646     simError();
1647     }
1648    
1649     strcpy( info.nameC, the_token );
1650    
1651     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1652     sprintf( painCave.errMsg,
1653     "Error parseing TorsionTypes: line %d\n", lineNum );
1654     painCave.isFatal = 1;
1655     simError();
1656     }
1657    
1658     strcpy( info.nameD, the_token );
1659    
1660     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1661     sprintf( painCave.errMsg,
1662     "Error parseing TorsionTypes: line %d\n", lineNum );
1663     painCave.isFatal = 1;
1664     simError();
1665     }
1666    
1667     strcpy( info.type, the_token );
1668    
1669     if( !strcmp( info.type, "cubic" ) ){
1670     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1671     sprintf( painCave.errMsg,
1672     "Error parseing TorsionTypes: line %d\n", lineNum );
1673     painCave.isFatal = 1;
1674     simError();
1675     }
1676    
1677     info.k1 = atof( the_token );
1678    
1679     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1680     sprintf( painCave.errMsg,
1681     "Error parseing TorsionTypes: line %d\n", lineNum );
1682     painCave.isFatal = 1;
1683     simError();
1684     }
1685    
1686     info.k2 = atof( the_token );
1687    
1688     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1689     sprintf( painCave.errMsg,
1690     "Error parseing TorsionTypes: line %d\n", lineNum );
1691     painCave.isFatal = 1;
1692     simError();
1693     }
1694    
1695     info.k3 = atof( the_token );
1696    
1697     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1698     sprintf( painCave.errMsg,
1699     "Error parseing TorsionTypes: line %d\n", lineNum );
1700     painCave.isFatal = 1;
1701     simError();
1702     }
1703    
1704     info.k4 = atof( the_token );
1705    
1706     }
1707    
1708     else{
1709     sprintf( painCave.errMsg,
1710     "Unknown TraPPE_Ex torsion type \"%s\" at line %d\n",
1711     info.type,
1712     lineNum );
1713     painCave.isFatal = 1;
1714     simError();
1715     }
1716    
1717     return 1;
1718     }
1719    
1720     else return 0;
1721     }