ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/WATER.cpp
Revision: 989
Committed: Tue Jan 27 19:38:27 2004 UTC (20 years, 5 months ago) by gezelter
File size: 27948 byte(s)
Log Message:
More BASS changes to do new rigidBody scheme
a copy of WATER.cpp from this morning

File Contents

# User Rev Content
1 chrisfen 976 #include <stdlib.h>
2     #include <stdio.h>
3     #include <string.h>
4    
5     #include <iostream>
6     using namespace std;
7    
8     #ifdef IS_MPI
9     #include <mpi.h>
10     #endif //is_mpi
11    
12     #include "ForceFields.hpp"
13     #include "SRI.hpp"
14     #include "simError.h"
15    
16     #include "fortranWrappers.hpp"
17    
18     #ifdef IS_MPI
19     #include "mpiForceField.h"
20     #endif // is_mpi
21    
22    
23    
24     namespace WATER_NS{
25    
26 gezelter 989 // Declare the structures that will be passed by the parser and MPI
27 chrisfen 976
28     typedef struct{
29     char name[15];
30     double mass;
31     double epslon;
32     double sigma;
33     double charge;
34     int isDirectional;
35     int isLJ;
36     int isCharge;
37     int ident;
38     int last; // 0 -> default
39     // 1 -> in MPI: tells nodes to stop listening
40     } atomStruct;
41    
42     typedef struct{
43     char name[15];
44 gezelter 989 double Ixx;
45     double Iyy;
46     double Izz;
47 chrisfen 976 double dipole;
48     double w0;
49     double v0;
50     double v0p;
51     double rl;
52     double ru;
53     double rlp;
54     double rup;
55     int isDipole;
56     int isSticky;
57     int last; // 0 -> default
58     // 1 -> in MPI: tells nodes to stop listening
59     } directionalStruct;
60    
61     int parseAtom( char *lineBuffer, int lineNum, atomStruct &info );
62     int parseDirectional( char *lineBuffer, int lineNum, directionalStruct &info );
63    
64    
65     #ifdef IS_MPI
66    
67     MPI_Datatype mpiAtomStructType;
68     MPI_Datatype mpiDirectionalStructType;
69    
70     #endif
71    
72     class LinkedAtomType {
73     public:
74     LinkedAtomType(){
75     next = NULL;
76     name[0] = '\0';
77     }
78     ~LinkedAtomType(){ if( next != NULL ) delete next; }
79    
80     LinkedAtomType* find(char* key){
81     if( !strcmp(name, key) ) return this;
82     if( next != NULL ) return next->find(key);
83     return NULL;
84     }
85    
86    
87     void add( atomStruct &info ){
88    
89     // check for duplicates
90    
91     if( !strcmp( info.name, name ) ){
92     sprintf( painCave.errMsg,
93     "Duplicate WATER atom type \"%s\" found in "
94     "the WATER param file./n",
95     name );
96     painCave.isFatal = 1;
97     simError();
98     }
99    
100     if( next != NULL ) next->add(info);
101     else{
102     next = new LinkedAtomType();
103     strcpy(next->name, info.name);
104     next->isDirectional = info.isDirectional;
105     next->isLJ = info.isLJ;
106     next->isCharge = info.isCharge;
107     next->mass = info.mass;
108     next->epslon = info.epslon;
109     next->sigma = info.sigma;
110     next->charge = info.charge;
111     next->ident = info.ident;
112     }
113     }
114    
115     #ifdef IS_MPI
116    
117     void duplicate( atomStruct &info ){
118     strcpy(info.name, name);
119     info.isDirectional = isDirectional;
120     info.isLJ = isLJ;
121     info.isCharge = isCharge;
122     info.mass = mass;
123     info.epslon = epslon;
124     info.sigma = sigma;
125     info.charge = charge;
126     info.ident = ident;
127     info.last = 0;
128     }
129    
130     #endif
131    
132     char name[15];
133     int isDirectional;
134     int isLJ;
135     int isCharge;
136     double mass;
137     double epslon;
138     double sigma;
139     double charge;
140     int ident;
141     LinkedAtomType* next;
142     };
143    
144     class LinkedDirectionalType {
145     public:
146     LinkedDirectionalType(){
147     next = NULL;
148     name[0] = '\0';
149     }
150     ~LinkedDirectionalType(){ if( next != NULL ) delete next; }
151    
152     LinkedDirectionalType* find(char* key){
153     if( !strcmp(name, key) ) return this;
154     if( next != NULL ) return next->find(key);
155     return NULL;
156     }
157    
158    
159     void add( directionalStruct &info ){
160    
161     // check for duplicates
162    
163     if( !strcmp( info.name, name ) ){
164     sprintf( painCave.errMsg,
165     "Duplicate WATER directional type \"%s\" found in "
166     "the WATER param file./n",
167     name );
168     painCave.isFatal = 1;
169     simError();
170     }
171    
172     if( next != NULL ) next->add(info);
173     else{
174     next = new LinkedDirectionalType();
175     strcpy(next->name, info.name);
176     next->isDipole = info.isDipole;
177     next->isSticky = info.isSticky;
178 gezelter 989 next->Ixx = info.Ixx;
179     next->Iyy = info.Iyy;
180     next->Izz = info.Izz;
181     next->dipole = info.dipole;
182 chrisfen 976 next->w0 = info.w0;
183     next->v0 = info.v0;
184     next->v0p = info.v0p;
185     next->rl = info.rl;
186     next->ru = info.ru;
187     next->rlp = info.rlp;
188     next->rup = info.rup;
189     }
190     }
191    
192     #ifdef IS_MPI
193    
194     void duplicate( directionalStruct &info ){
195     strcpy(info.name, name);
196     info.isDipole = isDipole;
197     info.isSticky = isSticky;
198 gezelter 989 info.Ixx = Ixx;
199     info.Iyy = Iyy;
200     info.Izz = Izz;
201 chrisfen 976 info.dipole = dipole;
202     info.w0 = w0;
203     info.v0 = v0;
204     info.v0p = v0p;
205     info.rl = rl;
206     info.ru = ru;
207     info.rlp = rlp;
208     info.rup = rup;
209     info.last = 0;
210     }
211    
212     #endif
213    
214     char name[15];
215     int isDipole;
216     int isSticky;
217 gezelter 989 double Ixx;
218     double Iyy;
219     double Izz;
220 chrisfen 976 double dipole;
221     double w0;
222     double v0;
223     double v0p;
224     double rl;
225     double ru;
226     double rlp;
227     double rup;
228     LinkedDirectionalType* next;
229     };
230    
231     LinkedAtomType* headAtomType;
232     LinkedAtomType* currentAtomType;
233     LinkedDirectionalType* headDirectionalType;
234     LinkedDirectionalType* currentDirectionalType;
235 gezelter 989 } // namespace
236 chrisfen 976
237     using namespace WATER_NS;
238    
239     //****************************************************************
240     // begins the actual forcefield stuff.
241     //****************************************************************
242    
243    
244     WATER::WATER(){
245    
246     char fileName[200];
247     char* ffPath_env = "FORCE_PARAM_PATH";
248     char* ffPath;
249     char temp[200];
250    
251     headAtomType = NULL;
252     currentAtomType = NULL;
253     headDirectionalType = NULL;
254     currentDirectionalType = NULL;
255    
256     // do the funtion wrapping
257     wrapMeFF( this );
258    
259     #ifdef IS_MPI
260     int i;
261    
262     // **********************************************************************
263     // Init the atomStruct mpi type
264    
265     atomStruct atomProto; // mpiPrototype
266     int atomBC[3] = {15,4,5}; // block counts
267     MPI_Aint atomDspls[3]; // displacements
268     MPI_Datatype atomMbrTypes[3]; // member mpi types
269    
270     MPI_Address(&atomProto.name, &atomDspls[0]);
271     MPI_Address(&atomProto.mass, &atomDspls[1]);
272 gezelter 989 MPI_Address(&atomProto.isDirectional, &atomDspls[2]);
273 chrisfen 976
274     atomMbrTypes[0] = MPI_CHAR;
275     atomMbrTypes[1] = MPI_DOUBLE;
276     atomMbrTypes[2] = MPI_INT;
277    
278     for (i=2; i >= 0; i--) atomDspls[i] -= atomDspls[0];
279    
280     MPI_Type_struct(3, atomBC, atomDspls, atomMbrTypes, &mpiAtomStructType);
281     MPI_Type_commit(&mpiAtomStructType);
282    
283     // ***********************************************************************
284    
285     // **********************************************************************
286     // Init the directionalStruct mpi type
287    
288     directionalStruct directionalProto; // mpiPrototype
289 gezelter 989 int directionalBC[3] = {15,11,3}; // block counts
290 chrisfen 976 MPI_Aint directionalDspls[3]; // displacements
291     MPI_Datatype directionalMbrTypes[3]; // member mpi types
292    
293     MPI_Address(&directionalProto.name, &directionalDspls[0]);
294 gezelter 989 MPI_Address(&directionalProto.Ixx, &directionalDspls[1]);
295     MPI_Address(&directionalProto.isDipole, &directionalDspls[2]);
296 chrisfen 976
297     directionalMbrTypes[0] = MPI_CHAR;
298     directionalMbrTypes[1] = MPI_DOUBLE;
299     directionalMbrTypes[2] = MPI_INT;
300    
301     for (i=2; i >= 0; i--) directionalDspls[i] -= directionalDspls[0];
302    
303     MPI_Type_struct(3, directionalBC, directionalDspls, directionalMbrTypes,
304     &mpiDirectionalStructType);
305     MPI_Type_commit(&mpiDirectionalStructType);
306    
307     // ***********************************************************************
308    
309     if( worldRank == 0 ){
310     #endif
311    
312     // generate the force file name
313    
314     strcpy( fileName, "WATER.frc" );
315     // fprintf( stderr,"Trying to open %s\n", fileName );
316    
317     // attempt to open the file in the current directory first.
318    
319     frcFile = fopen( fileName, "r" );
320    
321     if( frcFile == NULL ){
322    
323     // next see if the force path enviorment variable is set
324    
325     ffPath = getenv( ffPath_env );
326     if( ffPath == NULL ) {
327     STR_DEFINE(ffPath, FRC_PATH );
328     }
329    
330    
331     strcpy( temp, ffPath );
332     strcat( temp, "/" );
333     strcat( temp, fileName );
334     strcpy( fileName, temp );
335    
336     frcFile = fopen( fileName, "r" );
337    
338     if( frcFile == NULL ){
339    
340     sprintf( painCave.errMsg,
341     "Error opening the force field parameter file: %s\n"
342     "Have you tried setting the FORCE_PARAM_PATH environment "
343     "variable?\n",
344     fileName );
345     painCave.isFatal = 1;
346     simError();
347     }
348     }
349    
350     #ifdef IS_MPI
351     }
352    
353     sprintf( checkPointMsg, "WATER file opened sucessfully." );
354     MPIcheckPoint();
355    
356     #endif // is_mpi
357     }
358    
359    
360     WATER::~WATER(){
361    
362     if( headAtomType != NULL ) delete headAtomType;
363     if( headDirectionalType != NULL ) delete headDirectionalType;
364    
365     #ifdef IS_MPI
366     if( worldRank == 0 ){
367     #endif // is_mpi
368    
369     fclose( frcFile );
370    
371     #ifdef IS_MPI
372     }
373     #endif // is_mpi
374     }
375    
376     void WATER::cleanMe( void ){
377    
378     #ifdef IS_MPI
379    
380     // keep the linked list in the mpi version
381    
382     #else // is_mpi
383    
384     // delete the linked list in the single processor version
385    
386     if( headAtomType != NULL ) delete headAtomType;
387     if( headDirectionalType != NULL ) delete headDirectionalType;
388    
389     #endif // is_mpi
390     }
391    
392    
393     void WATER::initForceField( int ljMixRule ){
394    
395     initFortran( ljMixRule, entry_plug->useReactionField );
396     }
397    
398    
399     void WATER::readParams( void ){
400    
401     int identNum;
402    
403     atomStruct atomInfo;
404     directionalStruct directionalInfo;
405 gezelter 989 fpos_t *atomPos;
406 chrisfen 976
407     atomInfo.last = 1; // initialize last to have the last set.
408     directionalInfo.last = 1; // if things go well, last will be set to 0
409    
410     atomPos = new fpos_t;
411     bigSigma = 0.0;
412    
413     #ifdef IS_MPI
414     if( worldRank == 0 ){
415     #endif
416    
417     // read in the atom types.
418    
419     headAtomType = new LinkedAtomType;
420 gezelter 989 headDirectionalType = new LinkedDirectionalType;
421 chrisfen 976
422     fastForward( "AtomTypes", "initializeAtoms" );
423    
424     // we are now at the AtomTypes section.
425    
426     eof_test = fgets( readLine, sizeof(readLine), frcFile );
427     lineNum++;
428    
429    
430     // read a line, and start parsing out the atom types
431    
432     if( eof_test == NULL ){
433     sprintf( painCave.errMsg,
434     "Error in reading Atoms from force file at line %d.\n",
435     lineNum );
436     painCave.isFatal = 1;
437     simError();
438     }
439    
440     identNum = 1;
441     // stop reading at end of file, or at next section
442     while( readLine[0] != '#' && eof_test != NULL ){
443    
444     // toss comment lines
445     if( readLine[0] != '!' ){
446    
447     // the parser returns 0 if the line was blank
448     if( parseAtom( readLine, lineNum, atomInfo ) ){
449     atomInfo.ident = identNum;
450     headAtomType->add( atomInfo );
451     if ( atomInfo.isDirectional ) {
452    
453     // if the atom is directional, skip to the directional section
454     // and parse directional info.
455     fgetpos( frcFile, atomPos );
456     sectionSearch( "DirectionalTypes", atomInfo.name,
457     "initializeDirectional" );
458     parseDirectional( readLine, lineNum, directionalInfo );
459 gezelter 989 headDirectionalType->add( directionalInfo );
460 chrisfen 976
461     // return to the AtomTypes section
462     fsetpos( frcFile, atomPos );
463     }
464     identNum++;
465     }
466     }
467     eof_test = fgets( readLine, sizeof(readLine), frcFile );
468     lineNum++;
469     }
470    
471     #ifdef IS_MPI
472    
473     // send out the linked list to all the other processes
474    
475     sprintf( checkPointMsg,
476     "WATER atom structures read successfully." );
477     MPIcheckPoint();
478     currentAtomType = headAtomType->next; //skip the first element who is a place holder.
479     while( currentAtomType != NULL ){
480     currentAtomType->duplicate( atomInfo );
481    
482     sendFrcStruct( &atomInfo, mpiAtomStructType );
483    
484     sprintf( checkPointMsg,
485     "successfully sent WATER force type: \"%s\"\n",
486     atomInfo.name );
487     MPIcheckPoint();
488    
489     currentAtomType = currentAtomType->next;
490     }
491 gezelter 989
492 chrisfen 976 atomInfo.last = 1;
493     sendFrcStruct( &atomInfo, mpiAtomStructType );
494    
495     if ( atomInfo.isDirectional ){
496     // send out the linked list to all the other processes
497    
498     sprintf( checkPointMsg,
499     "WATER directional structures read successfully." );
500     MPIcheckPoint();
501    
502     currentDirectionalType = headDirectionalType->next;
503     while( currentDirectionalType != NULL ){
504     currentDirectionalType->duplicate( directionalInfo );
505     sendFrcStruct( &directionalInfo, mpiDirectionalStructType );
506     currentDirectionalType = currentDirectionalType->next;
507     }
508     directionalInfo.last = 1;
509     sendFrcStruct( &directionalInfo, mpiDirectionalStructType );
510     }
511     }
512 gezelter 989
513 chrisfen 976 else{
514    
515     // listen for node 0 to send out the force params
516    
517     MPIcheckPoint();
518 gezelter 989
519 chrisfen 976 headAtomType = new LinkedAtomType;
520     receiveFrcStruct( &atomInfo, mpiAtomStructType );
521    
522     while( !atomInfo.last ){
523    
524     headAtomType->add( atomInfo );
525    
526     MPIcheckPoint();
527    
528     receiveFrcStruct( &atomInfo, mpiAtomStructType );
529     }
530    
531     if ( atomInfo.isDirectional ) {
532     // listen for node 0 to send out the force params
533    
534     MPIcheckPoint();
535    
536     headDirectionalType = new LinkedDirectionalType;
537     receiveFrcStruct( &directionalInfo, mpiDirectionalStructType );
538     while( !directionalInfo.last ){
539    
540     headDirectionalType->add( directionalInfo );
541     receiveFrcStruct( &directionalInfo, mpiDirectionalStructType );
542     }
543    
544     sprintf( checkPointMsg,
545     "WATER directional structures broadcast successfully." );
546     MPIcheckPoint();
547     }
548     }
549    
550     #endif // is_mpi
551    
552     // call new A_types in fortran
553    
554     int isError;
555    
556     // dummy variables
557     int isGB = 0;
558     int isEAM = 0;
559 gezelter 989 int notDipole = 0;
560     int notSSD = 0;
561     double noDipMoment = 0.0;
562    
563    
564     currentAtomType = headAtomType->next;
565     currentDirectionalType = headDirectionalType->next;
566    
567 chrisfen 976 while( currentAtomType != NULL ){
568    
569     if( currentAtomType->isLJ ) entry_plug->useLJ = 1;
570 gezelter 989 if( currentAtomType->isCharge ) entry_plug->useCharges = 1;
571 chrisfen 976 if( currentAtomType->isDirectional ){
572 gezelter 989 // only directional atoms can have dipoles or be sticky
573     if ( currentDirectionalType->isDipole ) entry_plug->useDipoles = 1;
574     if ( currentDirectionalType->isSticky ) {
575 chrisfen 976 entry_plug->useSticky = 1;
576 gezelter 989 set_sticky_params( &(currentDirectionalType->w0),
577     &(currentDirectionalType->v0),
578     &(currentDirectionalType->v0p),
579     &(currentDirectionalType->rl),
580     &(currentDirectionalType->ru),
581     &(currentDirectionalType->rlp),
582     &(currentDirectionalType->rup));
583 chrisfen 976 }
584 gezelter 989 if( currentAtomType->name[0] != '\0' ){
585     isError = 0;
586     makeAtype( &(currentAtomType->ident),
587     &(currentAtomType->isLJ),
588     &(currentDirectionalType->isSticky),
589     &(currentDirectionalType->isDipole),
590     &isGB,
591     &isEAM,
592     &(currentAtomType->isCharge),
593     &(currentAtomType->epslon),
594     &(currentAtomType->sigma),
595     &(currentAtomType->charge),
596     &(currentDirectionalType->dipole),
597     &isError );
598     if( isError ){
599     sprintf( painCave.errMsg,
600     "Error initializing the \"%s\" atom type in fortran\n",
601     currentAtomType->name );
602     painCave.isFatal = 1;
603     simError();
604     }
605     }
606     currentDirectionalType->next;
607 chrisfen 976 }
608 gezelter 989
609     else {
610     // use all dummy variables if this is not a directional atom
611     if( currentAtomType->name[0] != '\0' ){
612     isError = 0;
613     makeAtype( &(currentAtomType->ident),
614     &(currentAtomType->isLJ),
615     &notSSD,
616     &notDipole,
617     &isGB,
618     &isEAM,
619     &(currentAtomType->isCharge),
620     &(currentAtomType->epslon),
621     &(currentAtomType->sigma),
622     &(currentAtomType->charge),
623     &noDipMoment,
624     &isError );
625     if( isError ){
626     sprintf( painCave.errMsg,
627     "Error initializing the \"%s\" atom type in fortran\n",
628     currentAtomType->name );
629     painCave.isFatal = 1;
630     simError();
631     }
632 chrisfen 976 }
633     }
634     currentAtomType = currentAtomType->next;
635     }
636    
637     #ifdef IS_MPI
638     sprintf( checkPointMsg,
639     "WATER atom structures successfully sent to fortran\n" );
640     MPIcheckPoint();
641     #endif // is_mpi
642    
643     }
644    
645    
646     void WATER::initializeAtoms( int nAtoms, Atom** the_atoms ){
647    
648 gezelter 989 int i,j;
649 chrisfen 976
650     // initialize the atoms
651 gezelter 989 DirectionalAtom* dAtom;
652     double inertialMat[3][3];
653 chrisfen 976
654     for( i=0; i<nAtoms; i++ ){
655 gezelter 989 fprintf(stderr, "flag 1\n");
656 chrisfen 976 currentAtomType = headAtomType->find( the_atoms[i]->getType() );
657 gezelter 989 fprintf(stderr, "%s is the type\n", the_atoms[i]->getType());
658 chrisfen 976 if( currentAtomType == NULL ){
659     sprintf( painCave.errMsg,
660     "AtomType error, %s not found in force file.\n",
661     the_atoms[i]->getType() );
662     painCave.isFatal = 1;
663     simError();
664     }
665 gezelter 989 fprintf(stderr, "flag 2\n");
666     if( currentAtomType->isLJ ) the_atoms[i]->setLJ();
667     if( currentAtomType->isCharge ) the_atoms[i]->setCharged();
668 chrisfen 976 the_atoms[i]->setMass( currentAtomType->mass );
669     the_atoms[i]->setIdent( currentAtomType->ident );
670 gezelter 989 fprintf(stderr, "flag 3\n");
671     if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
672 chrisfen 976
673 gezelter 989 if( currentAtomType->isDirectional ){
674     fprintf(stderr, "flag 4\n");
675     currentDirectionalType =
676     headDirectionalType->find( the_atoms[i]->getType() );
677     fprintf(stderr, "%s is the type\n", the_atoms[i]->getType());
678     if( currentDirectionalType == NULL ){
679     sprintf( painCave.errMsg,
680     "DirectionalType error, %s not found in force file.\n",
681     the_atoms[i]->getType() );
682     painCave.isFatal = 1;
683     simError();
684     }
685 chrisfen 976
686 gezelter 989 // zero out the moments of inertia matrix
687     for( i=0; i<3; i++ )
688     for( j=0; j<3; j++ )
689     inertialMat[i][j] = 0.0;
690 chrisfen 976
691 gezelter 989 // load the force file moments of inertia
692     inertialMat[0][0] = currentDirectionalType->Ixx;
693     inertialMat[1][1] = currentDirectionalType->Iyy;
694     inertialMat[2][2] = currentDirectionalType->Izz;
695     fprintf(stderr, "Let's try pointing to isDirectional\n");
696     fprintf(stderr, "%i what is this\n",the_atoms[i]->isDirectional());
697     dAtom = (DirectionalAtom *) the_atoms[i];
698     fprintf(stderr, "%i is isDipole\n", currentDirectionalType->isDipole);
699     dAtom->setHasDipole( currentDirectionalType->isDipole );
700     dAtom->setMu( currentDirectionalType->dipole );
701     fprintf(stderr, "flag 5\n");
702     dAtom->setMu( currentDirectionalType->dipole );
703     fprintf(stderr,"flag 6\n");
704     // if it's sticky then it's an SSD type
705     dAtom->setSSD( currentDirectionalType->isSticky );
706     dAtom->setJx( 0.0 );
707     dAtom->setJy( 0.0 );
708     dAtom->setJz( 0.0 );
709     dAtom->setI( inertialMat );
710     fprintf(stderr, "flag 7\n");
711     entry_plug->n_dipoles++;
712     fprintf(stderr, "flag 8\n");
713     }
714     else{
715     sprintf( painCave.errMsg,
716     "WATER error: Atom \"%s\" is directional, yet no standard"
717     " orientation was specifed in the BASS file.\n",
718     currentAtomType->name );
719 chrisfen 976 painCave.isFatal = 1;
720     simError();
721     }
722     }
723     }
724    
725     void WATER::initializeBonds( int nBonds, Bond** BondArray,
726     bond_pair* the_bonds ){
727    
728     if( nBonds ){
729     sprintf( painCave.errMsg,
730     "WATER does not support bonds.\n" );
731     painCave.isFatal = 1;
732     simError();
733     }
734     }
735    
736     void WATER::initializeBends( int nBends, Bend** bendArray,
737     bend_set* the_bends ){
738    
739     if( nBends ){
740     sprintf( painCave.errMsg,
741     "WATER does not support bends.\n" );
742     painCave.isFatal = 1;
743     simError();
744     }
745     }
746    
747     void WATER::initializeTorsions( int nTorsions, Torsion** torsionArray,
748     torsion_set* the_torsions ){
749    
750     if( nTorsions ){
751     sprintf( painCave.errMsg,
752     "WATER does not support torsions.\n" );
753     painCave.isFatal = 1;
754     simError();
755     }
756     }
757    
758     void WATER::fastForward( char* stopText, char* searchOwner ){
759    
760     int foundText = 0;
761     char* the_token;
762    
763     rewind( frcFile );
764     lineNum = 0;
765    
766     eof_test = fgets( readLine, sizeof(readLine), frcFile );
767     lineNum++;
768     if( eof_test == NULL ){
769     sprintf( painCave.errMsg, "Error fast forwarding force file for %s: "
770     " file is empty.\n",
771     searchOwner );
772     painCave.isFatal = 1;
773     simError();
774     }
775    
776    
777     while( !foundText ){
778     while( eof_test != NULL && readLine[0] != '#' ){
779     eof_test = fgets( readLine, sizeof(readLine), frcFile );
780     lineNum++;
781     }
782     if( eof_test == NULL ){
783     sprintf( painCave.errMsg,
784     "Error fast forwarding force file for %s at "
785     "line %d: file ended unexpectedly.\n",
786     searchOwner,
787     lineNum );
788     painCave.isFatal = 1;
789     simError();
790     }
791    
792     the_token = strtok( readLine, " ,;\t#\n" );
793     foundText = !strcmp( stopText, the_token );
794    
795     if( !foundText ){
796     eof_test = fgets( readLine, sizeof(readLine), frcFile );
797     lineNum++;
798    
799     if( eof_test == NULL ){
800     sprintf( painCave.errMsg,
801     "Error fast forwarding force file for %s at "
802     "line %d: file ended unexpectedly.\n",
803     searchOwner,
804     lineNum );
805     painCave.isFatal = 1;
806     simError();
807     }
808     }
809     }
810     }
811    
812     void WATER::sectionSearch( char* secHead, char* stopText, char* searchOwner ){
813    
814     int foundSection = 0;
815     int foundText = 0;
816     char* the_token;
817 gezelter 989 fpos_t *tempPos;
818 chrisfen 976
819     rewind( frcFile );
820     lineNum = 0;
821 gezelter 989 tempPos = new fpos_t;
822 chrisfen 976
823     eof_test = fgets( readLine, sizeof(readLine), frcFile );
824     lineNum++;
825     if( eof_test == NULL ){
826     sprintf( painCave.errMsg, "Error fast forwarding force file for %s: "
827     " file is empty.\n",
828     searchOwner );
829     painCave.isFatal = 1;
830     simError();
831     }
832    
833     while( !foundSection ){
834     while( eof_test != NULL && readLine[0] != '#' ){
835     eof_test = fgets( readLine, sizeof(readLine), frcFile );
836     lineNum++;
837     }
838     if( eof_test == NULL ){
839     sprintf( painCave.errMsg,
840     "Error fast forwarding force file for %s at "
841     "line %d: file ended unexpectedly.\n",
842     searchOwner,
843     lineNum );
844     painCave.isFatal = 1;
845     simError();
846     }
847     the_token = strtok( readLine, " ,;\t#\n" );
848     foundSection = !strcmp( secHead, the_token );
849    
850     if( !foundSection ){
851     eof_test = fgets( readLine, sizeof(readLine), frcFile );
852     lineNum++;
853    
854     if( eof_test == NULL ){
855     sprintf( painCave.errMsg,
856     "Error section searching force file for %s at "
857     "line %d: file ended unexpectedly.\n",
858     searchOwner,
859     lineNum );
860     painCave.isFatal = 1;
861     simError();
862     }
863     }
864 gezelter 989 }
865 chrisfen 976
866     while( !foundText ){
867 gezelter 989
868     fgetpos( frcFile, tempPos );
869     eof_test = fgets( readLine, sizeof(readLine), frcFile );
870     lineNum++;
871    
872     if( eof_test == NULL ){
873     sprintf( painCave.errMsg,
874     "Error fast forwarding force file for %s at "
875     "line %d: file ended unexpectedly.\n",
876     searchOwner,
877     lineNum );
878     painCave.isFatal = 1;
879 chrisfen 976 simError();
880 gezelter 989 }
881 chrisfen 976
882     the_token = strtok( readLine, " ,;\t#\n" );
883 gezelter 989 if( the_token != NULL ){
884     foundText = !strcmp( stopText, the_token );
885     }
886 chrisfen 976 }
887 gezelter 989 fsetpos( frcFile, tempPos );
888     eof_test = fgets( readLine, sizeof(readLine), frcFile );
889     lineNum++;
890 chrisfen 976 }
891    
892    
893     int WATER_NS::parseAtom( char *lineBuffer, int lineNum, atomStruct &info ){
894    
895     char* the_token;
896    
897     the_token = strtok( lineBuffer, " \n\t,;" );
898     if( the_token != NULL ){
899    
900     strcpy( info.name, the_token );
901    
902     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
903     sprintf( painCave.errMsg,
904     "Error parseing AtomTypes: line %d\n", lineNum );
905     painCave.isFatal = 1;
906     simError();
907     }
908    
909     info.isDirectional = atoi( the_token );
910    
911     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
912     sprintf( painCave.errMsg,
913     "Error parseing AtomTypes: line %d\n", lineNum );
914     painCave.isFatal = 1;
915     simError();
916     }
917    
918     info.isLJ = atoi( the_token );
919    
920     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
921     sprintf( painCave.errMsg,
922     "Error parseing AtomTypes: line %d\n", lineNum );
923     painCave.isFatal = 1;
924     simError();
925     }
926    
927     info.isCharge = atoi( the_token );
928    
929     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
930     sprintf( painCave.errMsg,
931     "Error parseing AtomTypes: line %d\n", lineNum );
932     painCave.isFatal = 1;
933     simError();
934     }
935    
936     info.mass = atof( the_token );
937    
938     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
939     sprintf( painCave.errMsg,
940     "Error parseing AtomTypes: line %d\n", lineNum );
941     painCave.isFatal = 1;
942     simError();
943     }
944    
945     info.epslon = atof( the_token );
946    
947     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
948     sprintf( painCave.errMsg,
949     "Error parseing AtomTypes: line %d\n", lineNum );
950     painCave.isFatal = 1;
951     simError();
952     }
953    
954     info.sigma = atof( the_token );
955    
956     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
957     sprintf( painCave.errMsg,
958     "Error parseing AtomTypes: line %d\n", lineNum );
959     painCave.isFatal = 1;
960     simError();
961     }
962    
963     info.charge = atof( the_token );
964    
965     return 1;
966     }
967     else return 0;
968     }
969    
970     int WATER_NS::parseDirectional( char *lineBuffer, int lineNum, directionalStruct &info ){
971    
972     char* the_token;
973    
974     the_token = strtok( lineBuffer, " \n\t,;" );
975     if( the_token != NULL ){
976    
977     strcpy( info.name, the_token );
978    
979     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
980     sprintf( painCave.errMsg,
981     "Error parseing DirectionalTypes: line %d\n", lineNum );
982     painCave.isFatal = 1;
983     simError();
984     }
985    
986     info.isDipole = atoi( the_token );
987    
988     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
989     sprintf( painCave.errMsg,
990     "Error parseing DirectionalTypes: line %d\n", lineNum );
991     painCave.isFatal = 1;
992     simError();
993     }
994    
995     info.isSticky = atoi( the_token );
996    
997     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
998     sprintf( painCave.errMsg,
999     "Error parseing DirectionalTypes: line %d\n", lineNum );
1000     painCave.isFatal = 1;
1001     simError();
1002     }
1003    
1004 gezelter 989 info.Ixx = atof( the_token );
1005 chrisfen 976
1006     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1007     sprintf( painCave.errMsg,
1008     "Error parseing DirectionalTypes: line %d\n", lineNum );
1009     painCave.isFatal = 1;
1010     simError();
1011     }
1012    
1013 gezelter 989 info.Iyy = atof( the_token );
1014 chrisfen 976
1015     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1016     sprintf( painCave.errMsg,
1017     "Error parseing DirectionalTypes: line %d\n", lineNum );
1018     painCave.isFatal = 1;
1019     simError();
1020     }
1021    
1022 gezelter 989 info.Izz = atof( the_token );
1023 chrisfen 976
1024     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1025     sprintf( painCave.errMsg,
1026     "Error parseing DirectionalTypes: line %d\n", lineNum );
1027     painCave.isFatal = 1;
1028     simError();
1029     }
1030    
1031    
1032     info.dipole = atof( the_token );
1033    
1034     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1035     sprintf( painCave.errMsg,
1036     "Error parseing DirectionalTypes: line %d\n", lineNum );
1037     painCave.isFatal = 1;
1038     simError();
1039     }
1040    
1041     info.w0 = atof( the_token );
1042    
1043     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1044     sprintf( painCave.errMsg,
1045     "Error parseing DirectionalTypes: line %d\n", lineNum );
1046     painCave.isFatal = 1;
1047     simError();
1048     }
1049    
1050     info.v0 = atof( the_token );
1051    
1052     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1053     sprintf( painCave.errMsg,
1054     "Error parseing DirectionalTypes: line %d\n", lineNum );
1055     painCave.isFatal = 1;
1056     simError();
1057     }
1058    
1059     info.v0p = atof( the_token );
1060    
1061     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1062     sprintf( painCave.errMsg,
1063     "Error parseing DirectionalTypes: line %d\n", lineNum );
1064     painCave.isFatal = 1;
1065     simError();
1066     }
1067    
1068     info.rl = atof( the_token );
1069    
1070     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1071     sprintf( painCave.errMsg,
1072     "Error parseing DirectionalTypes: line %d\n", lineNum );
1073     painCave.isFatal = 1;
1074     simError();
1075     }
1076    
1077     info.ru = atof( the_token );
1078    
1079     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1080     sprintf( painCave.errMsg,
1081     "Error parseing DirectionalTypes: line %d\n", lineNum );
1082     painCave.isFatal = 1;
1083     simError();
1084     }
1085    
1086     info.rlp = atof( the_token );
1087    
1088     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
1089     sprintf( painCave.errMsg,
1090     "Error parseing DirectionalTypes: line %d\n", lineNum );
1091     painCave.isFatal = 1;
1092     simError();
1093     }
1094    
1095     info.rup = atof( the_token );
1096    
1097     return 1;
1098     }
1099     else return 0;
1100     }