ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/TraPPE_ExFF.cpp
Revision: 976
Committed: Thu Jan 22 17:34:20 2004 UTC (20 years, 5 months ago) by chrisfen
File size: 40994 byte(s)
Log Message: 
Corrected spelling in several directories, and stated WATER.cpp

File Contents

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