ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/DumpReader.cpp
Revision: 1203
Committed: Thu May 27 18:59:17 2004 UTC (20 years, 1 month ago) by gezelter
File size: 17964 byte(s)
Log Message:
Cutoff group changes under MPI

File Contents

# Content
1 #define _LARGEFILE_SOURCE64
2 #define _FILE_OFFSET_BITS 64
3
4 #include <sys/types.h>
5 #include <sys/stat.h>
6
7 #include <iostream>
8 #include <math.h>
9
10 #include <stdio.h>
11 #include <stdlib.h>
12 #include <string.h>
13
14
15 #include "ReadWrite.hpp"
16 #include "simError.h"
17
18 #ifdef IS_MPI
19 #include <mpi.h>
20 #include "mpiSimulation.hpp"
21 #define TAKE_THIS_TAG_CHAR 0
22 #define TAKE_THIS_TAG_INT 1
23 #endif // is_mpi
24
25
26 DumpReader :: DumpReader(const char *in_name ){
27
28 isScanned = false;
29
30 #ifdef IS_MPI
31 if (worldRank == 0) {
32 #endif
33
34 inFile = fopen(in_name, "r");
35 if(inFile == NULL){
36 sprintf(painCave.errMsg,
37 "Cannot open file: %s\n", in_name);
38 painCave.isFatal = 1;
39 simError();
40 }
41
42 inFileName = in_name;
43 #ifdef IS_MPI
44 }
45 strcpy( checkPointMsg, "Dump file opened for reading successfully." );
46 MPIcheckPoint();
47 #endif
48 return;
49 }
50
51 DumpReader :: ~DumpReader( ){
52 #ifdef IS_MPI
53 if (worldRank == 0) {
54 #endif
55 vector<fpos_t*>::iterator i;
56
57 int error;
58 error = fclose( inFile );
59 if( error ){
60 sprintf( painCave.errMsg,
61 "Error closing %s\n", inFileName.c_str());
62 simError();
63 }
64
65 for(i = framePos.begin(); i != framePos.end(); ++i)
66 delete *i;
67 framePos.clear();
68
69 #ifdef IS_MPI
70 }
71 strcpy( checkPointMsg, "Dump file closed successfully." );
72 MPIcheckPoint();
73 #endif
74
75 return;
76 }
77
78 int DumpReader::getNframes( void ){
79
80 if( !isScanned )
81 scanFile();
82 return framePos.size();
83 }
84
85 void DumpReader::scanFile( void ){
86
87 int vectorSize;
88 int i, j, k;
89 int lineNum = 0;
90 char readBuffer[2000];
91 char* foo;
92 fpos_t *currPos;
93 double time;
94
95
96
97 #ifdef IS_MPI
98 if( worldRank == 0 ){
99 #endif // is_mpi
100
101 rewind( inFile );
102
103 currPos = new fpos_t;
104 fgetpos( inFile, currPos );
105 fgets( readBuffer, sizeof( readBuffer ), inFile );
106 lineNum++;
107 if( feof( inFile ) ){
108 sprintf( painCave.errMsg,
109 "File \"%s\" ended unexpectedly at line %d\n",
110 inFileName.c_str(),
111 lineNum );
112 painCave.isFatal = 1;
113 simError();
114 }
115
116 while( !feof( inFile ) ){
117
118 framePos.push_back(currPos);
119
120 i = atoi(readBuffer);
121
122 fgets( readBuffer, sizeof( readBuffer ), inFile );
123 lineNum++;
124 if( feof( inFile ) ){
125 sprintf( painCave.errMsg,
126 "File \"%s\" ended unexpectedly at line %d\n",
127 inFileName.c_str(),
128 lineNum );
129 painCave.isFatal = 1;
130 simError();
131 }
132
133 for(j=0; j<i; j++){
134
135 fgets( readBuffer, sizeof( readBuffer ), inFile );
136 lineNum++;
137 if( feof( inFile ) ){
138 sprintf( painCave.errMsg,
139 "File \"%s\" ended unexpectedly at line %d,"
140 " with atom %d\n",
141 inFileName.c_str(),
142 lineNum,
143 j );
144 painCave.isFatal = 1;
145 simError();
146 }
147
148 }
149
150 currPos = new fpos_t;
151 fgetpos( inFile, currPos );
152 fgets( readBuffer, sizeof( readBuffer ), inFile );
153 lineNum++;
154 }
155
156 delete currPos;
157 rewind( inFile );
158
159 isScanned = true;
160
161 #ifdef IS_MPI
162 }
163 strcpy( checkPointMsg, "Successfully scanned DumpFile\n" );
164 MPIcheckPoint();
165 #endif // is_mpi
166 }
167
168 void DumpReader :: readFrame( SimInfo* the_simnfo, int whichFrame){
169
170 simnfo = the_simnfo;
171
172 this->readSet( whichFrame );
173 }
174
175
176
177 void DumpReader :: readSet( int whichFrame ){
178
179 int i, j;
180
181 #ifdef IS_MPI
182 int done, which_node, which_atom; // loop counter
183 #endif //is_mpi
184
185 const int BUFFERSIZE = 2000; // size of the read buffer
186 int nTotObjs; // the number of atoms
187 char read_buffer[BUFFERSIZE]; //the line buffer for reading
188
189 char *eof_test; // ptr to see when we reach the end of the file
190 char *parseErr;
191
192 vector<StuntDouble*> integrableObjects;
193
194
195 #ifndef IS_MPI
196
197 fsetpos(inFile, framePos[whichFrame]);
198 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
199 if( eof_test == NULL ){
200 sprintf( painCave.errMsg,
201 "DumpReader error: error reading 1st line of \"%s\"\n",
202 inFileName.c_str() );
203 painCave.isFatal = 1;
204 simError();
205 }
206
207 nTotObjs = atoi( read_buffer );
208
209 if( nTotObjs != simnfo->getTotIntegrableObjects() ){
210 sprintf( painCave.errMsg,
211 "DumpReader error. %s n_atoms, %d, "
212 "does not match the BASS file's n_atoms, %d.\n",
213 inFileName.c_str(), nTotObjs, simnfo->getTotIntegrableObjects());
214 painCave.isFatal = 1;
215 simError();
216 }
217
218 //read the box mat from the comment line
219
220 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
221 if(eof_test == NULL){
222 sprintf( painCave.errMsg,
223 "error in reading commment in %s\n", inFileName.c_str());
224 painCave.isFatal = 1;
225 simError();
226 }
227
228 parseErr = parseCommentLine( read_buffer, simnfo);
229 if( parseErr != NULL ){
230 strcpy( painCave.errMsg, parseErr );
231 painCave.isFatal = 1;
232 simError();
233 }
234
235 //parse dump lines
236
237 for( i=0; i < simnfo->n_mol; i++){
238
239 integrableObjects = (simnfo->molecules[i]).getIntegrableObjects();
240
241 for(j = 0; j < integrableObjects.size(); j++){
242
243 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
244 if(eof_test == NULL){
245 sprintf(painCave.errMsg,
246 "error in reading file %s\n"
247 "natoms = %d; index = %d\n"
248 "error reading the line from the file.\n",
249 inFileName.c_str(), nTotObjs, i );
250 painCave.isFatal = 1;
251 simError();
252 }
253
254 parseErr = parseDumpLine( read_buffer, integrableObjects[j]);
255 if( parseErr != NULL ){
256 strcpy( painCave.errMsg, parseErr );
257 painCave.isFatal = 1;
258 simError();
259 }
260 }
261 }
262
263 // MPI Section of code..........
264 #else //IS_MPI
265
266 // first thing first, suspend fatalities.
267 painCave.isEventLoop = 1;
268
269 int myStatus; // 1 = wakeup & success; 0 = error; -1 = AllDone
270 int haveError;
271
272 MPI_Status istatus;
273 int *MolToProcMap = mpiSim->getMolToProcMap();
274 int localIndex;
275 int nCurObj;
276 int nitems;
277
278 nTotObjs = simnfo->getTotIntegrableObjects();
279 haveError = 0;
280 if (worldRank == 0) {
281 fsetpos(inFile, framePos[whichFrame]);
282
283 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
284 if( eof_test == NULL ){
285 sprintf( painCave.errMsg,
286 "Error reading 1st line of %s \n ",inFileName.c_str());
287 haveError = 1;
288 simError();
289 }
290
291 nitems = atoi( read_buffer );
292
293 // Check to see that the number of integrable objects in the intial configuration file is the
294 // same as declared in simBass.
295
296 if( nTotObjs != nitems){
297 sprintf( painCave.errMsg,
298 "DumpReadererror. %s n_atoms, %d, "
299 "does not match the BASS file's n_atoms, %d.\n",
300 inFileName.c_str(), nTotObjs, simnfo->getTotIntegrableObjects());
301 haveError= 1;
302 simError();
303 }
304
305 //read the boxMat from the comment line
306
307 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
308 if(eof_test == NULL){
309 sprintf( painCave.errMsg,
310 "error in reading commment in %s\n", inFileName.c_str());
311 haveError = 1;
312 simError();
313 }
314
315 //Every single processor will parse the comment line by itself
316 //By using this way, we might lose some efficiency, but if we want to add
317 //more parameters into comment line, we only need to modify function
318 //parseCommentLine
319
320 MPI_Bcast(read_buffer, BUFFERSIZE, MPI_CHAR, 0, MPI_COMM_WORLD);
321
322 parseErr = parseCommentLine( read_buffer, simnfo);
323
324 if( parseErr != NULL ){
325 strcpy( painCave.errMsg, parseErr );
326 haveError = 1;
327 simError();
328 }
329
330 for (i=0 ; i < mpiSim->getNMolGlobal(); i++) {
331 which_node = MolToProcMap[i];
332 if(which_node == 0){
333 //molecules belong to master node
334
335 localIndex = mpiSim->getGlobalToLocalMol(i);
336
337 if(localIndex == -1) {
338 strcpy(painCave.errMsg, "Molecule not found on node 0!");
339 haveError = 1;
340 simError();
341 }
342
343 integrableObjects = (simnfo->molecules[localIndex]).getIntegrableObjects();
344 for(j=0; j < integrableObjects.size(); j++){
345
346 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
347 if(eof_test == NULL){
348 sprintf(painCave.errMsg,
349 "error in reading file %s\n"
350 "natoms = %d; index = %d\n"
351 "error reading the line from the file.\n",
352 inFileName.c_str(), nTotObjs, i );
353 haveError= 1;
354 simError();
355 }
356
357 if(haveError) nodeZeroError();
358
359 parseDumpLine(read_buffer, integrableObjects[j]);
360
361 }
362
363
364 }
365 else{
366 //molecule belongs to slave nodes
367
368 MPI_Recv(&nCurObj, 1, MPI_INT, which_node,
369 TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus);
370
371 for(j=0; j < nCurObj; j++){
372
373 eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
374 if(eof_test == NULL){
375 sprintf(painCave.errMsg,
376 "error in reading file %s\n"
377 "natoms = %d; index = %d\n"
378 "error reading the line from the file.\n",
379 inFileName.c_str(), nTotObjs, i );
380 haveError= 1;
381 simError();
382 }
383
384 if(haveError) nodeZeroError();
385
386 MPI_Send(read_buffer, BUFFERSIZE, MPI_CHAR, which_node,
387 TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD);
388
389 }
390
391 }
392
393 }
394
395 }
396 else{
397 //actions taken at slave nodes
398 MPI_Bcast(read_buffer, BUFFERSIZE, MPI_CHAR, 0, MPI_COMM_WORLD);
399
400 parseErr = parseCommentLine( read_buffer, simnfo);
401
402 if( parseErr != NULL ){
403 strcpy( painCave.errMsg, parseErr );
404 haveError = 1;
405 simError();
406 }
407
408 for (i=0 ; i < mpiSim->getNMolGlobal(); i++) {
409 which_node = MolToProcMap[i];
410
411 if(which_node == worldRank){
412 //molecule with global index i belongs to this processor
413
414 localIndex = mpiSim->getGlobalToLocalMol(i);
415
416 if(localIndex == -1) {
417 sprintf(painCave.errMsg, "Molecule not found on node %d\n", worldRank);
418 haveError = 1;
419 simError();
420 }
421
422 integrableObjects = (simnfo->molecules[localIndex]).getIntegrableObjects();
423
424 nCurObj = integrableObjects.size();
425
426 MPI_Send(&nCurObj, 1, MPI_INT, 0,
427 TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
428
429 for(j = 0; j < integrableObjects.size(); j++){
430
431 MPI_Recv(read_buffer, BUFFERSIZE, MPI_CHAR, 0,
432 TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD, &istatus);
433
434 parseErr = parseDumpLine(read_buffer, integrableObjects[j]);
435
436 if( parseErr != NULL ){
437 strcpy( painCave.errMsg, parseErr );
438 simError();
439 }
440
441 }
442
443 }
444
445 }
446
447 }
448
449 #endif
450 }
451
452 char* DumpReader::parseDumpLine(char* readLine, StuntDouble* sd){
453
454 char *foo; // the pointer to the current string token
455
456 double pos[3]; // position place holders
457 double vel[3]; // velocity placeholders
458 double q[4]; // the quaternions
459 double ji[3]; // angular velocity placeholders;
460 double qSqr, qLength; // needed to normalize the quaternion vector.
461
462
463 // set the string tokenizer
464
465 foo = strtok(readLine, " ,;\t");
466
467 // check the atom name to the current atom
468
469 if( strcmp( foo, sd->getType() ) ){
470 sprintf( painCave.errMsg,
471 "DumpReader error. Does not"
472 " match the BASS atom %s.\n",
473 sd->getType() );
474 return strdup( painCave.errMsg );
475 }
476
477 // get the positions
478
479 foo = strtok(NULL, " ,;\t");
480 if(foo == NULL){
481 sprintf( painCave.errMsg,
482 "error in reading postition x from %s\n",
483 inFileName.c_str());
484 return strdup( painCave.errMsg );
485 }
486 pos[0] = atof( foo );
487
488 foo = strtok(NULL, " ,;\t");
489 if(foo == NULL){
490 sprintf( painCave.errMsg,
491 "error in reading postition y from %s\n",
492 inFileName.c_str());
493 return strdup( painCave.errMsg );
494 }
495 pos[1] = atof( foo );
496
497 foo = strtok(NULL, " ,;\t");
498 if(foo == NULL){
499 sprintf( painCave.errMsg,
500 "error in reading postition z from %s\n",
501 inFileName.c_str());
502 return strdup( painCave.errMsg );
503 }
504 pos[2] = atof( foo );
505
506
507 // get the velocities
508
509 foo = strtok(NULL, " ,;\t");
510 if(foo == NULL){
511 sprintf( painCave.errMsg,
512 "error in reading velocity x from %s\n",
513 inFileName.c_str() );
514 return strdup( painCave.errMsg );
515 }
516 vel[0] = atof( foo );
517
518 foo = strtok(NULL, " ,;\t");
519 if(foo == NULL){
520 sprintf( painCave.errMsg,
521 "error in reading velocity x from %s\n",
522 inFileName.c_str() );
523 return strdup( painCave.errMsg );
524 }
525 vel[1] = atof( foo );
526
527 foo = strtok(NULL, " ,;\t");
528 if(foo == NULL){
529 sprintf( painCave.errMsg,
530 "error in reading velocity x from %s\n",
531 inFileName.c_str() );
532 return strdup( painCave.errMsg );
533 }
534 vel[2] = atof( foo );
535
536
537 // add the positions and velocities to the atom
538
539 sd->setPos( pos );
540 sd->setVel( vel );
541
542 if (!sd->isDirectional())
543 return NULL;
544
545 // get the quaternions
546
547 if( sd->isDirectional() ){
548
549 foo = strtok(NULL, " ,;\t");
550 if(foo == NULL){
551 sprintf( painCave.errMsg,
552 "error in reading velocity x from %s\n",
553 inFileName.c_str() );
554 return strdup( painCave.errMsg );
555 }
556 q[0] = atof( foo );
557
558 foo = strtok(NULL, " ,;\t");
559 if(foo == NULL){
560 sprintf( painCave.errMsg,
561 "error in reading velocity x from %s\n",
562 inFileName.c_str() );
563 return strdup( painCave.errMsg );
564 }
565 q[1] = atof( foo );
566
567 foo = strtok(NULL, " ,;\t");
568 if(foo == NULL){
569 sprintf( painCave.errMsg,
570 "error in reading velocity x from %s\n",
571 inFileName.c_str() );
572 return strdup( painCave.errMsg );
573 }
574 q[2] = atof( foo );
575
576 foo = strtok(NULL, " ,;\t");
577 if(foo == NULL){
578 sprintf( painCave.errMsg,
579 "error in reading velocity x from %s\n",
580 inFileName.c_str() );
581 return strdup( painCave.errMsg );
582 }
583 q[3] = atof( foo );
584
585 // get the angular velocities
586
587 foo = strtok(NULL, " ,;\t");
588 if(foo == NULL){
589 sprintf( painCave.errMsg,
590 "error in reading velocity x from %s\n",
591 inFileName.c_str() );
592 return strdup( painCave.errMsg );
593 }
594 ji[0] = atof( foo );
595
596 foo = strtok(NULL, " ,;\t");
597 if(foo == NULL){
598 sprintf( painCave.errMsg,
599 "error in reading velocity x from %s\n",
600 inFileName.c_str() );
601 return strdup( painCave.errMsg );
602 }
603 ji[1] = atof(foo );
604
605 foo = strtok(NULL, " ,;\t");
606 if(foo == NULL){
607 sprintf( painCave.errMsg,
608 "error in reading velocity x from %s\n",
609 inFileName.c_str() );
610 return strdup( painCave.errMsg );
611 }
612 ji[2] = atof( foo );
613
614
615 // check that the quaternion vector is normalized
616
617 qSqr = (q[0] * q[0]) + (q[1] * q[1]) + (q[2] * q[2]) + (q[3] * q[3]);
618
619 qLength = sqrt( qSqr );
620 q[0] = q[0] / qLength;
621 q[1] = q[1] / qLength;
622 q[2] = q[2] / qLength;
623 q[3] = q[3] / qLength;
624
625 // add quaternion and angular velocities
626
627 sd->setQ( q );
628 sd->setJ( ji );
629 }
630
631
632
633 return NULL;
634 }
635
636
637 char* DumpReader::parseCommentLine(char* readLine, SimInfo* entry_plug){
638
639 double currTime;
640 double boxMat[9];
641 double theBoxMat3[3][3];
642 double chi;
643 double integralOfChidt;
644 double eta[9];
645
646 char *foo; // the pointer to the current string token
647
648 // set the string tokenizer
649
650 foo = strtok(readLine, " ,;\t");
651 // set the timeToken.
652
653 if(foo == NULL){
654 sprintf( painCave.errMsg,
655 "error in reading Time from %s\n",
656 inFileName.c_str() );
657 return strdup( painCave.errMsg );
658 }
659
660 currTime = atof( foo );
661 entry_plug->setTime( currTime );
662
663 //get H-Matrix
664
665 for(int i = 0 ; i < 9; i++){
666 foo = strtok(NULL, " ,;\t");
667 if(foo == NULL){
668 sprintf( painCave.errMsg,
669 "error in reading H[%d] from %s\n", i, inFileName.c_str() );
670 return strdup( painCave.errMsg );
671 }
672 boxMat[i] = atof( foo );
673 }
674
675 for(int i=0;i<3;i++)
676 for(int j=0;j<3;j++) theBoxMat3[i][j] = boxMat[3*j+i];
677
678 //set H-Matrix
679 entry_plug->setBoxM( theBoxMat3 );
680
681 //get chi and integralOfChidt, they should appear by pair
682
683 if( entry_plug->useInitXSstate ){
684 foo = strtok(NULL, " ,;\t\n");
685 if(foo != NULL){
686 chi = atof(foo);
687
688 foo = strtok(NULL, " ,;\t\n");
689 if(foo == NULL){
690 sprintf( painCave.errMsg,
691 "chi and integralOfChidt should appear by pair in %s\n", inFileName.c_str() );
692 return strdup( painCave.errMsg );
693 }
694 integralOfChidt = atof( foo );
695
696 //push chi and integralOfChidt into SimInfo::properties which can be
697 //retrieved by integrator later
698 DoubleData* chiValue = new DoubleData();
699 chiValue->setID(CHIVALUE_ID);
700 chiValue->setData(chi);
701 entry_plug->addProperty(chiValue);
702
703 DoubleData* integralOfChidtValue = new DoubleData();
704 integralOfChidtValue->setID(INTEGRALOFCHIDT_ID);
705 integralOfChidtValue->setData(integralOfChidt);
706 entry_plug->addProperty(integralOfChidtValue);
707
708 }
709 else
710 return NULL;
711
712 //get eta
713 foo = strtok(NULL, " ,;\t\n");
714 if(foo != NULL ){
715
716 for(int i = 0 ; i < 9; i++){
717
718 if(foo == NULL){
719 sprintf( painCave.errMsg,
720 "error in reading eta[%d] from %s\n", i, inFileName.c_str() );
721 return strdup( painCave.errMsg );
722 }
723 eta[i] = atof( foo );
724 foo = strtok(NULL, " ,;\t\n");
725 }
726 }
727 else
728 return NULL;
729
730 //push eta into SimInfo::properties which can be
731 //retrieved by integrator later
732 //entry_plug->setBoxM( theBoxMat3 );
733 DoubleArrayData* etaValue = new DoubleArrayData();
734 etaValue->setID(ETAVALUE_ID);
735 etaValue->setData(eta, 9);
736 entry_plug->addProperty(etaValue);
737 }
738
739 return NULL;
740 }
741
742 #ifdef IS_MPI
743 void DumpReader::nodeZeroError( void ){
744 int j, myStatus;
745
746 myStatus = 0;
747 for (j = 0; j < mpiSim->getNProcessors(); j++) {
748 MPI_Send( &myStatus, 1, MPI_INT, j,
749 TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
750 }
751
752
753 MPI_Finalize();
754 exit (0);
755
756 }
757
758 void DumpReader::anonymousNodeDie( void ){
759
760 MPI_Finalize();
761 exit (0);
762 }
763 #endif