ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/TraPPE_ExFF.cpp
Revision: 442
Committed: Wed Apr 2 15:01:06 2003 UTC (21 years, 3 months ago) by mmeineke
File size: 40590 byte(s)
Log Message: 
Fixed a bug where MPI was not getting the proper atomIdents.

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