ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/mdtools/interface_implementation/LJ_FF.cpp
Revision: 224
Committed: Wed Jan 8 21:53:53 2003 UTC (21 years, 8 months ago) by mmeineke
File size: 10253 byte(s)
Log Message:
fixed up the LJ_FF implementation. It is now the standard for building a force routine.

File Contents

# User Rev Content
1 chuckv 215 #include <cstdlib>
2     #include <cstdio>
3     #include <cstring>
4    
5     #include <iostream>
6     using namespace std;
7    
8     #include "ForceFields.hpp"
9     #include "SRI.hpp"
10     #include "simError.h"
11    
12    
13 mmeineke 224 // Declare the structures that will be passed by the parser and MPI
14 chuckv 215
15     typedef struct{
16     char name[15];
17     double mass;
18     double epslon;
19     double sigma;
20 mmeineke 224 int ident;
21 chuckv 215 int last; // 0 -> default
22 mmeineke 224 // 1 -> in MPI: tells nodes to stop listening
23 chuckv 215 } atomStruct;
24 mmeineke 224
25     int parseAtomLJ( char *lineBuffer, int lineNum, atomStruct &info );
26    
27     #ifdef IS_MPI
28    
29     #include "mpiForceField.h"
30    
31 chuckv 215 MPI_Datatype mpiAtomStructType;
32    
33     #endif
34    
35    
36     LJ_FF::LJ_FF(){
37    
38     char fileName[200];
39     char* ffPath_env = "FORCE_PARAM_PATH";
40     char* ffPath;
41     char temp[200];
42     char errMsg[1000];
43    
44     #ifdef IS_MPI
45     int i;
46    
47     // **********************************************************************
48     // Init the atomStruct mpi type
49    
50     atomStruct atomProto; // mpiPrototype
51 mmeineke 224 int atomBC[3] = {15,3,2}; // block counts
52 chuckv 215 MPI_Aint atomDspls[3]; // displacements
53     MPI_Datatype atomMbrTypes[3]; // member mpi types
54    
55     MPI_Address(&atomProto.name, &atomDspls[0]);
56     MPI_Address(&atomProto.mass, &atomDspls[1]);
57 mmeineke 224 MPI_Address(&atomProto.ident, &atomDspls[2]);
58 chuckv 215
59     atomMbrTypes[0] = MPI_CHAR;
60     atomMbrTypes[1] = MPI_DOUBLE;
61     atomMbrTypes[2] = MPI_INT;
62    
63     for (i=2; i >= 0; i--) atomDspls[i] -= atomDspls[0];
64    
65     MPI_Type_struct(3, atomBC, atomDspls, atomMbrTypes, &mpiAtomStructType);
66     MPI_Type_commit(&mpiAtomStructType);
67    
68     // ***********************************************************************
69    
70     if( worldRank == 0 ){
71     #endif
72    
73     // generate the force file name
74    
75     strcpy( fileName, "LJ_FF.frc" );
76     // fprintf( stderr,"Trying to open %s\n", fileName );
77    
78     // attempt to open the file in the current directory first.
79    
80     frcFile = fopen( fileName, "r" );
81    
82     if( frcFile == NULL ){
83    
84     // next see if the force path enviorment variable is set
85    
86     ffPath = getenv( ffPath_env );
87     if( ffPath == NULL ) {
88     sprintf( painCave.errMsg,
89     "Error opening the force field parameter file: %s\n"
90     "Have you tried setting the FORCE_PARAM_PATH environment "
91     "vairable?\n",
92     fileName );
93     painCave.isFatal = 1;
94     simError();
95     }
96    
97    
98     strcpy( temp, ffPath );
99     strcat( temp, "/" );
100     strcat( temp, fileName );
101     strcpy( fileName, temp );
102    
103     frcFile = fopen( fileName, "r" );
104    
105     if( frcFile == NULL ){
106    
107     sprintf( painCave.errMsg,
108     "Error opening the force field parameter file: %s\n"
109     "Have you tried setting the FORCE_PARAM_PATH environment "
110     "vairable?\n",
111     fileName );
112     painCave.isFatal = 1;
113     simError();
114     }
115     }
116    
117     #ifdef IS_MPI
118     }
119    
120     sprintf( checkPointMsg, "LJ_FF file opened sucessfully." );
121     MPIcheckPoint();
122    
123     #endif // is_mpi
124     }
125    
126    
127     LJ_FF::~LJ_FF(){
128    
129     #ifdef IS_MPI
130     if( worldRank == 0 ){
131     #endif // is_mpi
132    
133     fclose( frcFile );
134    
135     #ifdef IS_MPI
136     }
137     #endif // is_mpi
138     }
139    
140     void LJ_FF::initializeAtoms( void ){
141    
142     class LinkedType {
143     public:
144     LinkedType(){
145     next = NULL;
146     name[0] = '\0';
147     }
148     ~LinkedType(){ if( next != NULL ) delete next; }
149    
150     LinkedType* find(char* key){
151     if( !strcmp(name, key) ) return this;
152     if( next != NULL ) return next->find(key);
153     return NULL;
154     }
155    
156 mmeineke 224
157 chuckv 215 void add( atomStruct &info ){
158 mmeineke 224
159     // check for duplicates
160    
161     if( !strcmp( info.name, name ) ){
162     sprintf( simError.painCave,
163     "Duplicate LJ atom type \"%s\" found in "
164     "the LJ_FF param file./n",
165     name );
166     painCave.isFatal = 1;
167     simError();
168     }
169    
170 chuckv 215 if( next != NULL ) next->add(info);
171     else{
172     next = new LinkedType();
173     strcpy(next->name, info.name);
174     next->mass = info.mass;
175     next->epslon = info.epslon;
176     next->sigma = info.sigma;
177 mmeineke 224 next->ident = info.ident;
178 chuckv 215 }
179     }
180    
181 mmeineke 224
182     #ifdef IS_MPI
183    
184 chuckv 215 void duplicate( atomStruct &info ){
185     strcpy(info.name, name);
186     info.mass = mass;
187     info.epslon = epslon;
188     info.sigma = sigma;
189 mmeineke 224 info.ident = ident;
190 chuckv 215 info.last = 0;
191     }
192    
193    
194     #endif
195    
196     char name[15];
197     double mass;
198     double epslon;
199     double sigma;
200 mmeineke 224 int ident;
201 chuckv 215 LinkedType* next;
202     };
203    
204     LinkedType* headAtomType;
205     LinkedType* currentAtomType;
206     atomStruct info;
207     info.last = 1; // initialize last to have the last set.
208     // if things go well, last will be set to 0
209    
210     int i;
211 mmeineke 224 int identNum;
212 chuckv 215
213     Atom** the_atoms;
214     int nAtoms;
215     the_atoms = entry_plug->atoms;
216     nAtoms = entry_plug->n_atoms;
217    
218    
219     #ifdef IS_MPI
220     if( worldRank == 0 ){
221     #endif
222    
223     // read in the atom types.
224 mmeineke 224
225 chuckv 215 headAtomType = new LinkedType;
226    
227 mmeineke 224 fastFoward( "AtomTypes", "initializeAtoms" );
228    
229 chuckv 215 // we are now at the AtomTypes section.
230    
231     eof_test = fgets( readLine, sizeof(readLine), frcFile );
232     lineNum++;
233    
234 mmeineke 224
235     // read a line, and start parseing out the atom types
236    
237 chuckv 215 if( eof_test == NULL ){
238     sprintf( painCave.errMsg,
239     "Error in reading Atoms from force file at line %d.\n",
240     lineNum );
241     painCave.isFatal = 1;
242     simError();
243     }
244    
245 mmeineke 224 identNum = 1;
246     // stop reading at end of file, or at next section
247 chuckv 215 while( readLine[0] != '#' && eof_test != NULL ){
248 mmeineke 224
249     // toss comment lines
250 chuckv 215 if( readLine[0] != '!' ){
251    
252 mmeineke 224 // the parser returns 0 if the line was blank
253     if( parseAtomLJ( readLine, lineNum, info ) ){
254     info.ident = identNum;
255     headAtomType->add( info );;
256     identNum++;
257 chuckv 215 }
258     }
259     eof_test = fgets( readLine, sizeof(readLine), frcFile );
260     lineNum++;
261     }
262    
263     #ifdef IS_MPI
264    
265     // send out the linked list to all the other processes
266    
267     sprintf( checkPointMsg,
268     "LJ_FF atom structures read successfully." );
269     MPIcheckPoint();
270    
271     currentAtomType = headAtomType;
272     while( currentAtomType != NULL ){
273     currentAtomType->duplicate( info );
274     sendFrcStruct( &info, mpiAtomStructType );
275     currentAtomType = currentAtomType->next;
276     }
277     info.last = 1;
278     sendFrcStruct( &info, mpiAtomStructType );
279    
280     }
281    
282     else{
283    
284     // listen for node 0 to send out the force params
285    
286     MPIcheckPoint();
287    
288     headAtomType = new LinkedType;
289     recieveFrcStruct( &info, mpiAtomStructType );
290     while( !info.last ){
291    
292     headAtomType->add( info );
293     recieveFrcStruct( &info, mpiAtomStructType );
294     }
295     }
296     #endif // is_mpi
297    
298    
299     // initialize the atoms
300    
301     Atom* thisAtom;
302    
303     for( i=0; i<nAtoms; i++ ){
304    
305     currentAtomType = headAtomType->find( the_atoms[i]->getType() );
306     if( currentAtomType == NULL ){
307     sprintf( painCave.errMsg,
308     "AtomType error, %s not found in force file.\n",
309     the_atoms[i]->getType() );
310     painCave.isFatal = 1;
311     simError();
312     }
313    
314     the_atoms[i]->setMass( currentAtomType->mass );
315     the_atoms[i]->setEpslon( currentAtomType->epslon );
316     the_atoms[i]->setSigma( currentAtomType->sigma );
317 mmeineke 224 the_atoms[i]->setIdent( currentAtomType->ident );
318 chuckv 215 the_atoms[i]->setLJ();
319     }
320    
321    
322     // clean up the memory
323    
324     delete headAtomType;
325    
326     #ifdef IS_MPI
327     sprintf( checkPointMsg, "LJ_FF atoms initialized succesfully" );
328     MPIcheckPoint();
329     #endif // is_mpi
330    
331     }
332    
333     void LJ_FF::initializeBonds( bond_pair* the_bonds ){
334    
335     if( entry_plug->n_bonds ){
336     sprintf( painCave.errMsg,
337     "LJ_FF does not support bonds.\n" );
338     painCave.isFatal = 1;
339     simError();
340     }
341     #ifdef IS_MPI
342     MPIcheckPoint();
343     #endif // is_mpi
344    
345     }
346    
347     void LJ_FF::initializeBends( bend_set* the_bends ){
348    
349     if( entry_plug->n_bends ){
350     sprintf( painCave.errMsg,
351     "LJ_FF does not support bends.\n" );
352     painCave.isFatal = 1;
353     simError();
354     }
355     #ifdef IS_MPI
356     MPIcheckPoint();
357     #endif // is_mpi
358    
359     }
360    
361     void LJ_FF::initializeTorsions( torsion_set* the_torsions ){
362    
363     if( entry_plug->n_torsions ){
364     sprintf( painCave.errMsg,
365     "LJ_FF does not support torsions.\n" );
366     painCave.isFatal = 1;
367     simError();
368     }
369     #ifdef IS_MPI
370     MPIcheckPoint();
371     #endif // is_mpi
372    
373     }
374 mmeineke 224
375    
376     void LJ_FF::fastForward( char* stopText, char* searchOwner ){
377    
378     int foundText = 0;
379     char* the_token;
380    
381     rewind( frcFile );
382     lineNum = 0;
383    
384     eof_test = fgets( readLine, sizeof(readLine), frcFile );
385     lineNum++;
386     if( eof_test == NULL ){
387     sprintf( painCave.errMsg, "Error fast forwarding force file for %s: "
388     " file is empty.\n",
389     searchOwner );
390     painCave.isFatal = 1;
391     simError();
392     }
393    
394    
395     while( !foundText ){
396     while( eof_test != NULL && readLine[0] != '#' ){
397     eof_test = fgets( readLine, sizeof(readLine), frcFile );
398     lineNum++;
399     }
400     if( eof_test == NULL ){
401     sprintf( painCave.errMsg,
402     "Error fast forwarding force file for %s at "
403     "line %d: file ended unexpectedly.\n",
404     searchOwner,
405     lineNum );
406     painCave.isFatal = 1;
407     simError();
408     }
409    
410     the_token = strtok( readLine, " ,;\t#\n" );
411     foundText = !strcmp( stopText, the_token );
412    
413     if( !foundText ){
414     eof_test = fgets( readLine, sizeof(readLine), frcFile );
415     lineNum++;
416    
417     if( eof_test == NULL ){
418     sprintf( painCave.errMsg,
419     "Error fast forwarding force file for %s at "
420     "line %d: file ended unexpectedly.\n",
421     searchOwner,
422     lineNum );
423     painCave.isFatal = 1;
424     simError();
425     }
426     }
427     }
428     }
429    
430    
431    
432     int parseAtomLJ( char *lineBuffer, int lineNum, atomStruct &info ){
433    
434     char* the_token;
435    
436     the_token = strtok( lineBuffer, " \n\t,;" );
437     if( the_token != NULL ){
438    
439     strcpy( info.name, the_token );
440    
441     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
442     sprintf( painCave.errMsg,
443     "Error parseing AtomTypes: line %d\n", lineNum );
444     painCave.isFatal = 1;
445     simError();
446     }
447    
448     info.mass = atof( the_token );
449    
450     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
451     sprintf( painCave.errMsg,
452     "Error parseing AtomTypes: line %d\n", lineNum );
453     painCave.isFatal = 1;
454     simError();
455     }
456    
457     info.epslon = atof( the_token );
458    
459     if( ( the_token = strtok( NULL, " \n\t,;" ) ) == NULL ){
460     sprintf( painCave.errMsg,
461     "Error parseing AtomTypes: line %d\n", lineNum );
462     painCave.isFatal = 1;
463     simError();
464     }
465    
466     info.sigma = atof( the_token );
467    
468     return 1;
469     }
470     else return 0;
471     }