ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/LJFF.cpp
Revision: 829
Committed: Tue Oct 28 16:03:37 2003 UTC (20 years, 8 months ago) by gezelter
File size: 11836 byte(s)
Log Message:
replace c++ header stuff with more portable c header stuff
Also, mod file fixes and portability changes
Some fortran changes will need to be reversed.

File Contents

# Content
1 #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 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 "the LJFF param file./n",
69 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 LJFF::LJFF(){
121
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 strcpy( fileName, "LJFF.frc" );
165 fprintf( stderr,"Trying to open %s\n", fileName );
166
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 sprintf( checkPointMsg, "LJFF file opened sucessfully." );
204 MPIcheckPoint();
205
206 #endif // is_mpi
207 }
208
209
210 LJFF::~LJFF(){
211
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 void LJFF::initForceField( int ljMixRule ){
226
227 initFortran( ljMixRule, 0 );
228 }
229
230 void LJFF::cleanMe( void ){
231
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 void LJFF::readParams( void ){
246
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 "LJFF atom structures read successfully." );
305 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 int isEAM = 0;
359 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 &isEAM,
372 &(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 "LJFF atom structures successfully sent to fortran\n" );
392 MPIcheckPoint();
393 #endif // is_mpi
394
395 }
396
397
398 void LJFF::initializeAtoms( int nAtoms, Atom** the_atoms ){
399
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 void LJFF::initializeBonds( int nBonds, Bond** BondArray,
427 bond_pair* the_bonds ){
428
429 if( nBonds ){
430 sprintf( painCave.errMsg,
431 "LJFF does not support bonds.\n" );
432 painCave.isFatal = 1;
433 simError();
434 }
435 }
436
437 void LJFF::initializeBends( int nBends, Bend** bendArray,
438 bend_set* the_bends ){
439
440 if( nBends ){
441 sprintf( painCave.errMsg,
442 "LJFF does not support bends.\n" );
443 painCave.isFatal = 1;
444 simError();
445 }
446 }
447
448 void LJFF::initializeTorsions( int nTorsions, Torsion** torsionArray,
449 torsion_set* the_torsions ){
450
451 if( nTorsions ){
452 sprintf( painCave.errMsg,
453 "LJFF does not support torsions.\n" );
454 painCave.isFatal = 1;
455 simError();
456 }
457 }
458
459 void LJFF::fastForward( char* stopText, char* searchOwner ){
460
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 }