ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/LJFF.cpp
Revision: 787
Committed: Thu Sep 25 19:27:15 2003 UTC (20 years, 9 months ago) by mmeineke
File size: 11833 byte(s)
Log Message:
cleaned things with gcc -Wall and g++ -Wall

File Contents

# User Rev Content
1 mmeineke 559 #include <cstdlib>
2     #include <cstdio>
3     #include <cstring>
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 LJ_NS{
25    
26     // Declare the structures that will be passed by the parser and MPI
27    
28     typedef struct{
29     char name[15];
30     double mass;
31     double epslon;
32     double sigma;
33     int ident;
34     int last; // 0 -> default
35     // 1 -> in MPI: tells nodes to stop listening
36     } atomStruct;
37    
38     int parseAtom( char *lineBuffer, int lineNum, atomStruct &info );
39    
40     #ifdef IS_MPI
41    
42     MPI_Datatype mpiAtomStructType;
43    
44     #endif
45    
46     class LinkedAtomType {
47     public:
48     LinkedAtomType(){
49     next = NULL;
50     name[0] = '\0';
51     }
52     ~LinkedAtomType(){ if( next != NULL ) delete next; }
53    
54     LinkedAtomType* find(char* key){
55     if( !strcmp(name, key) ) return this;
56     if( next != NULL ) return next->find(key);
57     return NULL;
58     }
59    
60    
61     void add( atomStruct &info ){
62    
63     // check for duplicates
64    
65     if( !strcmp( info.name, name ) ){
66     sprintf( painCave.errMsg,
67     "Duplicate LJ atom type \"%s\" found in "
68 mmeineke 561 "the LJFF param file./n",
69 mmeineke 559 name );
70     painCave.isFatal = 1;
71     simError();
72     }
73    
74     if( next != NULL ) next->add(info);
75     else{
76     next = new LinkedAtomType();
77     strcpy(next->name, info.name);
78     next->mass = info.mass;
79     next->epslon = info.epslon;
80     next->sigma = info.sigma;
81     next->ident = info.ident;
82     }
83     }
84    
85    
86     #ifdef IS_MPI
87    
88     void duplicate( atomStruct &info ){
89     strcpy(info.name, name);
90     info.mass = mass;
91     info.epslon = epslon;
92     info.sigma = sigma;
93     info.ident = ident;
94     info.last = 0;
95     }
96    
97    
98     #endif
99    
100     char name[15];
101     double mass;
102     double epslon;
103     double sigma;
104     int ident;
105     LinkedAtomType* next;
106     };
107    
108     LinkedAtomType* headAtomType;
109     LinkedAtomType* currentAtomType;
110    
111     }
112    
113     using namespace LJ_NS;
114    
115     //****************************************************************
116     // begins the actual forcefield stuff.
117     //****************************************************************
118    
119    
120 mmeineke 561 LJFF::LJFF(){
121 mmeineke 559
122     char fileName[200];
123     char* ffPath_env = "FORCE_PARAM_PATH";
124     char* ffPath;
125     char temp[200];
126    
127     headAtomType = NULL;
128     currentAtomType = NULL;
129    
130     // do the funtion wrapping
131     wrapMeFF( this );
132    
133     #ifdef IS_MPI
134     int i;
135    
136     // **********************************************************************
137     // Init the atomStruct mpi type
138    
139     atomStruct atomProto; // mpiPrototype
140     int atomBC[3] = {15,3,2}; // block counts
141     MPI_Aint atomDspls[3]; // displacements
142     MPI_Datatype atomMbrTypes[3]; // member mpi types
143    
144     MPI_Address(&atomProto.name, &atomDspls[0]);
145     MPI_Address(&atomProto.mass, &atomDspls[1]);
146     MPI_Address(&atomProto.ident, &atomDspls[2]);
147    
148     atomMbrTypes[0] = MPI_CHAR;
149     atomMbrTypes[1] = MPI_DOUBLE;
150     atomMbrTypes[2] = MPI_INT;
151    
152     for (i=2; i >= 0; i--) atomDspls[i] -= atomDspls[0];
153    
154     MPI_Type_struct(3, atomBC, atomDspls, atomMbrTypes, &mpiAtomStructType);
155     MPI_Type_commit(&mpiAtomStructType);
156    
157     // ***********************************************************************
158    
159     if( worldRank == 0 ){
160     #endif
161    
162     // generate the force file name
163    
164 mmeineke 561 strcpy( fileName, "LJFF.frc" );
165     fprintf( stderr,"Trying to open %s\n", fileName );
166 mmeineke 559
167     // attempt to open the file in the current directory first.
168    
169     frcFile = fopen( fileName, "r" );
170    
171     if( frcFile == NULL ){
172    
173     // next see if the force path enviorment variable is set
174    
175     ffPath = getenv( ffPath_env );
176     if( ffPath == NULL ) {
177     STR_DEFINE(ffPath, FRC_PATH );
178     }
179    
180    
181     strcpy( temp, ffPath );
182     strcat( temp, "/" );
183     strcat( temp, fileName );
184     strcpy( fileName, temp );
185    
186     frcFile = fopen( fileName, "r" );
187    
188     if( frcFile == NULL ){
189    
190     sprintf( painCave.errMsg,
191     "Error opening the force field parameter file: %s\n"
192     "Have you tried setting the FORCE_PARAM_PATH environment "
193     "vairable?\n",
194     fileName );
195     painCave.isFatal = 1;
196     simError();
197     }
198     }
199    
200     #ifdef IS_MPI
201     }
202    
203 mmeineke 561 sprintf( checkPointMsg, "LJFF file opened sucessfully." );
204 mmeineke 559 MPIcheckPoint();
205    
206     #endif // is_mpi
207     }
208    
209    
210 mmeineke 561 LJFF::~LJFF(){
211 mmeineke 559
212     if( headAtomType != NULL ) delete headAtomType;
213    
214     #ifdef IS_MPI
215     if( worldRank == 0 ){
216     #endif // is_mpi
217    
218     fclose( frcFile );
219    
220     #ifdef IS_MPI
221     }
222     #endif // is_mpi
223     }
224    
225 mmeineke 561 void LJFF::initForceField( int ljMixRule ){
226 mmeineke 559
227     initFortran( ljMixRule, 0 );
228     }
229    
230 mmeineke 561 void LJFF::cleanMe( void ){
231 mmeineke 559
232     #ifdef IS_MPI
233    
234     // keep the linked list in the mpi version
235    
236     #else // is_mpi
237    
238     // delete the linked list in the single processor version
239    
240     if( headAtomType != NULL ) delete headAtomType;
241    
242     #endif // is_mpi
243     }
244    
245 mmeineke 561 void LJFF::readParams( void ){
246 mmeineke 559
247     atomStruct info;
248     info.last = 1; // initialize last to have the last set.
249     // if things go well, last will be set to 0
250    
251     int identNum;
252    
253    
254     bigSigma = 0.0;
255     #ifdef IS_MPI
256     if( worldRank == 0 ){
257     #endif
258    
259     // read in the atom types.
260    
261     headAtomType = new LinkedAtomType;
262    
263     fastForward( "AtomTypes", "initializeAtoms" );
264    
265     // we are now at the AtomTypes section.
266    
267     eof_test = fgets( readLine, sizeof(readLine), frcFile );
268     lineNum++;
269    
270    
271     // read a line, and start parseing out the atom types
272    
273     if( eof_test == NULL ){
274     sprintf( painCave.errMsg,
275     "Error in reading Atoms from force file at line %d.\n",
276     lineNum );
277     painCave.isFatal = 1;
278     simError();
279     }
280    
281     identNum = 1;
282     // stop reading at end of file, or at next section
283     while( readLine[0] != '#' && eof_test != NULL ){
284    
285     // toss comment lines
286     if( readLine[0] != '!' ){
287    
288     // the parser returns 0 if the line was blank
289     if( parseAtom( readLine, lineNum, info ) ){
290     info.ident = identNum;
291     headAtomType->add( info );;
292     identNum++;
293     }
294     }
295     eof_test = fgets( readLine, sizeof(readLine), frcFile );
296     lineNum++;
297     }
298    
299     #ifdef IS_MPI
300    
301     // send out the linked list to all the other processes
302    
303     sprintf( checkPointMsg,
304 mmeineke 561 "LJFF atom structures read successfully." );
305 mmeineke 559 MPIcheckPoint();
306    
307     currentAtomType = headAtomType->next; //skip the first element who is a place holder.
308     while( currentAtomType != NULL ){
309     currentAtomType->duplicate( info );
310    
311    
312    
313     sendFrcStruct( &info, mpiAtomStructType );
314    
315     sprintf( checkPointMsg,
316     "successfully sent lJ force type: \"%s\"\n",
317     info.name );
318     MPIcheckPoint();
319    
320     currentAtomType = currentAtomType->next;
321     }
322     info.last = 1;
323     sendFrcStruct( &info, mpiAtomStructType );
324    
325     }
326    
327     else{
328    
329     // listen for node 0 to send out the force params
330    
331     MPIcheckPoint();
332    
333     headAtomType = new LinkedAtomType;
334     recieveFrcStruct( &info, mpiAtomStructType );
335    
336     while( !info.last ){
337    
338    
339    
340     headAtomType->add( info );
341    
342     MPIcheckPoint();
343    
344     recieveFrcStruct( &info, mpiAtomStructType );
345     }
346     }
347     #endif // is_mpi
348    
349     // call new A_types in fortran
350    
351     int isError;
352    
353     // dummy variables
354     int isLJ = 1;
355     int isDipole = 0;
356     int isSSD = 0;
357     int isGB = 0;
358 chuckv 631 int isEAM = 0;
359 mmeineke 559 double dipole = 0.0;
360    
361     currentAtomType = headAtomType;
362     while( currentAtomType != NULL ){
363    
364     if( currentAtomType->name[0] != '\0' ){
365     isError = 0;
366     makeAtype( &(currentAtomType->ident),
367     &isLJ,
368     &isSSD,
369     &isDipole,
370     &isGB,
371 chuckv 631 &isEAM,
372 mmeineke 559 &(currentAtomType->epslon),
373     &(currentAtomType->sigma),
374     &dipole,
375     &isError );
376     if( isError ){
377     sprintf( painCave.errMsg,
378     "Error initializing the \"%s\" atom type in fortran\n",
379     currentAtomType->name );
380     painCave.isFatal = 1;
381     simError();
382     }
383     }
384     currentAtomType = currentAtomType->next;
385     }
386    
387     entry_plug->useLJ = 1;
388    
389     #ifdef IS_MPI
390     sprintf( checkPointMsg,
391 mmeineke 561 "LJFF atom structures successfully sent to fortran\n" );
392 mmeineke 559 MPIcheckPoint();
393     #endif // is_mpi
394    
395     }
396    
397    
398 mmeineke 561 void LJFF::initializeAtoms( int nAtoms, Atom** the_atoms ){
399 mmeineke 559
400     int i;
401    
402     // initialize the atoms
403    
404    
405     for( i=0; i<nAtoms; i++ ){
406    
407     currentAtomType = headAtomType->find( the_atoms[i]->getType() );
408     if( currentAtomType == NULL ){
409     sprintf( painCave.errMsg,
410     "AtomType error, %s not found in force file.\n",
411     the_atoms[i]->getType() );
412     painCave.isFatal = 1;
413     simError();
414     }
415    
416     the_atoms[i]->setMass( currentAtomType->mass );
417     the_atoms[i]->setEpslon( currentAtomType->epslon );
418     the_atoms[i]->setSigma( currentAtomType->sigma );
419     the_atoms[i]->setIdent( currentAtomType->ident );
420     the_atoms[i]->setLJ();
421    
422     if( bigSigma < currentAtomType->sigma ) bigSigma = currentAtomType->sigma;
423     }
424     }
425    
426 mmeineke 561 void LJFF::initializeBonds( int nBonds, Bond** BondArray,
427 mmeineke 559 bond_pair* the_bonds ){
428    
429     if( nBonds ){
430     sprintf( painCave.errMsg,
431 mmeineke 561 "LJFF does not support bonds.\n" );
432 mmeineke 559 painCave.isFatal = 1;
433     simError();
434     }
435     }
436    
437 mmeineke 561 void LJFF::initializeBends( int nBends, Bend** bendArray,
438 mmeineke 559 bend_set* the_bends ){
439    
440     if( nBends ){
441     sprintf( painCave.errMsg,
442 mmeineke 561 "LJFF does not support bends.\n" );
443 mmeineke 559 painCave.isFatal = 1;
444     simError();
445     }
446     }
447    
448 mmeineke 561 void LJFF::initializeTorsions( int nTorsions, Torsion** torsionArray,
449 mmeineke 559 torsion_set* the_torsions ){
450    
451     if( nTorsions ){
452     sprintf( painCave.errMsg,
453 mmeineke 561 "LJFF does not support torsions.\n" );
454 mmeineke 559 painCave.isFatal = 1;
455     simError();
456     }
457     }
458    
459 mmeineke 561 void LJFF::fastForward( char* stopText, char* searchOwner ){
460 mmeineke 559
461     int foundText = 0;
462     char* the_token;
463    
464     rewind( frcFile );
465     lineNum = 0;
466    
467     eof_test = fgets( readLine, sizeof(readLine), frcFile );
468     lineNum++;
469     if( eof_test == NULL ){
470     sprintf( painCave.errMsg, "Error fast forwarding force file for %s: "
471     " file is empty.\n",
472     searchOwner );
473     painCave.isFatal = 1;
474     simError();
475     }
476    
477    
478     while( !foundText ){
479     while( eof_test != NULL && readLine[0] != '#' ){
480     eof_test = fgets( readLine, sizeof(readLine), frcFile );
481     lineNum++;
482     }
483     if( eof_test == NULL ){
484     sprintf( painCave.errMsg,
485     "Error fast forwarding force file for %s at "
486     "line %d: file ended unexpectedly.\n",
487     searchOwner,
488     lineNum );
489     painCave.isFatal = 1;
490     simError();
491     }
492    
493     the_token = strtok( readLine, " ,;\t#\n" );
494     foundText = !strcmp( stopText, the_token );
495    
496     if( !foundText ){
497     eof_test = fgets( readLine, sizeof(readLine), frcFile );
498     lineNum++;
499    
500     if( eof_test == NULL ){
501     sprintf( painCave.errMsg,
502     "Error fast forwarding force file for %s at "
503     "line %d: file ended unexpectedly.\n",
504     searchOwner,
505     lineNum );
506     painCave.isFatal = 1;
507     simError();
508     }
509     }
510     }
511     }
512    
513    
514    
515     int LJ_NS::parseAtom( char *lineBuffer, int lineNum, atomStruct &info ){
516    
517     char* the_token;
518    
519     the_token = strtok( lineBuffer, " \n\t,;" );
520     if( the_token != NULL ){
521    
522     strcpy( info.name, the_token );
523    
524     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
525     sprintf( painCave.errMsg,
526     "Error parseing AtomTypes: line %d\n", lineNum );
527     painCave.isFatal = 1;
528     simError();
529     }
530    
531     info.mass = atof( the_token );
532    
533     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
534     sprintf( painCave.errMsg,
535     "Error parseing AtomTypes: line %d\n", lineNum );
536     painCave.isFatal = 1;
537     simError();
538     }
539    
540     info.epslon = atof( the_token );
541    
542     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
543     sprintf( painCave.errMsg,
544     "Error parseing AtomTypes: line %d\n", lineNum );
545     painCave.isFatal = 1;
546     simError();
547     }
548    
549     info.sigma = atof( the_token );
550    
551     return 1;
552     }
553     else return 0;
554     }