ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/DumpReader.cpp
Revision: 643
Committed: Mon Jul 21 21:27:40 2003 UTC (20 years, 11 months ago) by mmeineke
File size: 16624 byte(s)
Log Message:
some initial changes to Dumpwriter and friends to accomadate random file access

File Contents

# Content
1 #define _FILE_OFFSET_BITS 64
2
3 #include <iostream>
4 #include <cmath>
5
6 #include <stdio.h>
7 #include <stdlib.h>
8 #include <string.h>
9 #include <unistd.h>
10 #include <sys/types.h>
11 #include <sys/stat.h>
12
13 #include "ReadWrite.hpp"
14 #include "simError.h"
15
16 #ifdef IS_MPI
17 #include <mpi.h>
18 #include "mpiSimulation.hpp"
19 #define TAKE_THIS_TAG_CHAR 0
20 #define TAKE_THIS_TAG_INT 1
21
22 namespace dumpRead{
23 void nodeZeroError( void );
24 void anonymousNodeDie( void );
25 }
26
27 using namespace dumpRead;
28
29 #endif // is_mpi
30
31 DumpReader :: DumpReader( char *in_name ){
32 #ifdef IS_MPI
33 if (worldRank == 0) {
34 #endif
35
36 c_in_file = fopen(in_name, "r");
37 if(c_in_file == NULL){
38 sprintf(painCave.errMsg,
39 "Cannot open file: %s\n", in_name);
40 painCave.isFatal = 1;
41 simError();
42 }
43
44 strcpy( c_in_name, in_name);
45 #ifdef IS_MPI
46 }
47 strcpy( checkPointMsg, "Dump file opened for reading successfully." );
48 MPIcheckPoint();
49 #endif
50 return;
51 }
52
53 DumpReader :: ~DumpReader( ){
54 #ifdef IS_MPI
55 if (worldRank == 0) {
56 #endif
57 int error;
58 error = fclose( c_in_file );
59 if( error ){
60 sprintf( painCave.errMsg,
61 "Error closing %s\n", c_in_name );
62 simError();
63 }
64 #ifdef IS_MPI
65 }
66 strcpy( checkPointMsg, "Dump file closed successfully." );
67 MPIcheckPoint();
68 #endif
69
70 return;
71 }
72
73
74
75 void DumpReader :: getFrame( SimInfo* the_simnfo, int whichFrame){
76
77 int i, j, done, which_node, which_atom; // loop counter
78
79 const int BUFFERSIZE = 2000; // size of the read buffer
80 int n_atoms; // the number of atoms
81 char read_buffer[BUFFERSIZE]; //the line buffer for reading
82 #ifdef IS_MPI
83 char send_buffer[BUFFERSIZE];
84 #endif
85
86 char *eof_test; // ptr to see when we reach the end of the file
87 char *parseErr;
88 int procIndex;
89 double boxMat[9];
90 double theBoxMat3[3][3];
91
92 simnfo = the_simnfo;
93
94 #ifndef IS_MPI
95 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
96 if( eof_test == NULL ){
97 sprintf( painCave.errMsg,
98 "DumpReader error: error reading 1st line of \"%s\"\n",
99 c_in_name );
100 painCave.isFatal = 1;
101 simError();
102 }
103
104 n_atoms = atoi( read_buffer );
105
106 Atom **atoms = simnfo->atoms;
107 DirectionalAtom* dAtom;
108
109 if( n_atoms != simnfo->n_atoms ){
110 sprintf( painCave.errMsg,
111 "DumpReader error. %s n_atoms, %d, "
112 "does not match the BASS file's n_atoms, %d.\n",
113 c_in_name, n_atoms, simnfo->n_atoms );
114 painCave.isFatal = 1;
115 simError();
116 }
117
118 //read the box mat from the comment line
119
120 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
121 if(eof_test == NULL){
122 sprintf( painCave.errMsg,
123 "error in reading commment in %s\n", c_in_name);
124 painCave.isFatal = 1;
125 simError();
126 }
127
128 parseErr = parseCommentLine( read_buffer, time, boxMat );
129 if( parseErr != NULL ){
130 strcpy( painCave.errMsg, parseErr );
131 painCave.isFatal = 1;
132 simError();
133 }
134
135 simnfo->setTime( time );
136
137 for(i=0;i<3;i++)
138 for(j=0;j<3;j++) theBoxMat3[i][j] = boxMat[3*j+i];
139
140
141 simnfo->setBoxM( theBoxMat3 );
142
143
144 for( i=0; i < n_atoms; i++){
145
146 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
147 if(eof_test == NULL){
148 sprintf(painCave.errMsg,
149 "error in reading file %s\n"
150 "natoms = %d; index = %d\n"
151 "error reading the line from the file.\n",
152 c_in_name, n_atoms, i );
153 painCave.isFatal = 1;
154 simError();
155 }
156
157
158 parseErr = parseDumpLine( read_buffer, i );
159 if( parseErr != NULL ){
160 strcpy( painCave.errMsg, parseErr );
161 painCave.isFatal = 1;
162 simError();
163 }
164 }
165
166
167 // MPI Section of code..........
168 #else //IS_MPI
169
170 // first thing first, suspend fatalities.
171 painCave.isEventLoop = 1;
172
173 int myStatus; // 1 = wakeup & success; 0 = error; -1 = AllDone
174 int haveError;
175
176 MPI_Status istatus;
177 int *AtomToProcMap = mpiSim->getAtomToProcMap();
178
179
180 haveError = 0;
181 if (worldRank == 0) {
182
183 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
184 if( eof_test == NULL ){
185 sprintf( painCave.errMsg,
186 "Error reading 1st line of %d \n ",c_in_name);
187 haveError = 1;
188 simError();
189 }
190
191 n_atoms = atoi( read_buffer );
192
193 Atom **atoms = simnfo->atoms;
194 DirectionalAtom* dAtom;
195
196 // Check to see that the number of atoms in the intial configuration file is the
197 // same as declared in simBass.
198
199 if( n_atoms != mpiSim->getTotAtoms() ){
200 sprintf( painCave.errMsg,
201 "Initialize from File error. %s n_atoms, %d, "
202 "does not match the BASS file's n_atoms, %d.\n",
203 c_in_name, n_atoms, simnfo->n_atoms );
204 haveError= 1;
205 simError();
206 }
207
208 //read the time and boxMat from the comment line
209
210 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
211 if(eof_test == NULL){
212 sprintf( painCave.errMsg,
213 "error in reading commment in %s\n", c_in_name);
214 haveError = 1;
215 simError();
216 }
217
218 parseErr = parseCommentLine( read_buffer, time, boxMat );
219 if( parseErr != NULL ){
220 strcpy( painCave.errMsg, parseErr );
221 haveError = 1;
222 simError();
223 }
224
225 MPI_Bcast(time, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD );
226
227 MPI_Bcast(boxMat, 9, MPI_DOUBLE, 0, MPI_COMM_WORLD );
228
229 if(haveError) nodeZeroError();
230
231 for (i=0 ; i < mpiSim->getTotAtoms(); i++) {
232
233 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
234 if(eof_test == NULL){
235 sprintf(painCave.errMsg,
236 "error in reading file %s\n"
237 "natoms = %d; index = %d\n"
238 "error reading the line from the file.\n",
239 c_in_name, n_atoms, i );
240 haveError= 1;
241 simError();
242 }
243
244 if(haveError) nodeZeroError();
245
246 // Get the Node number which wants this atom:
247 which_node = AtomToProcMap[i];
248 if (which_node == 0) {
249 parseErr = parseDumpLine( read_buffer, i );
250 if( parseErr != NULL ){
251 strcpy( painCave.errMsg, parseErr );
252 haveError = 1;
253 simError();
254 }
255 if(haveError) nodeZeroError();
256 }
257
258 else {
259
260 myStatus = 1;
261 MPI_Send(&myStatus, 1, MPI_INT, which_node,
262 TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
263 MPI_Send(read_buffer, BUFFERSIZE, MPI_CHAR, which_node,
264 TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD);
265 MPI_Send(&i, 1, MPI_INT, which_node, TAKE_THIS_TAG_INT,
266 MPI_COMM_WORLD);
267 MPI_Recv(&myStatus, 1, MPI_INT, which_node, TAKE_THIS_TAG_INT,
268 MPI_COMM_WORLD, &istatus);
269
270 if(!myStatus) nodeZeroError();
271 }
272 }
273 myStatus = -1;
274 for (j = 0; j < mpiSim->getNumberProcessors(); j++) {
275 MPI_Send( &myStatus, 1, MPI_INT, j,
276 TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
277 }
278
279 } else {
280
281 MPI_Bcast(time, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
282 MPI_Bcast(boxMat, 9, MPI_DOUBLE, 0, MPI_COMM_WORLD);
283
284 done = 0;
285 while (!done) {
286
287 MPI_Recv(&myStatus, 1, MPI_INT, 0,
288 TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus);
289
290 if(!myStatus) anonymousNodeDie();
291
292 if(myStatus < 0) break;
293
294 MPI_Recv(read_buffer, BUFFERSIZE, MPI_CHAR, 0,
295 TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD, &istatus);
296 MPI_Recv(&which_atom, 1, MPI_INT, 0,
297 TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus);
298
299 myStatus = 1;
300 parseErr = parseDumpLine( read_buffer, which_atom );
301 if( parseErr != NULL ){
302 strcpy( painCave.errMsg, parseErr );
303 myStatus = 0;;
304 simError();
305 }
306
307 MPI_Send( &myStatus, 1, MPI_INT, 0,
308 TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
309
310 }
311 }
312
313 // last thing last, enable fatalities.
314 painCave.isEventLoop = 0;
315
316 simnfo->setTime( time );
317
318 for(i=0;i<3;i++)
319 for(j=0;j<3;j++) theBoxMat3[i][j] = boxMat[3*j+i];
320
321 simnfo->setBoxM( theBoxMat3 );
322
323
324 #endif
325 }
326
327 char* DumpReader::parseDumpLine(char* readLine, int globalIndex){
328
329 char *foo; // the pointer to the current string token
330
331 double rx, ry, rz; // position place holders
332 double vx, vy, vz; // velocity placeholders
333 double q[4]; // the quaternions
334 double jx, jy, jz; // angular velocity placeholders;
335 double qSqr, qLength; // needed to normalize the quaternion vector.
336
337 Atom **atoms = simnfo->atoms;
338 DirectionalAtom* dAtom;
339
340 int j, n_atoms, atomIndex;
341
342 #ifdef IS_MPI
343 n_atoms = mpiSim->getTotAtoms();
344 atomIndex=-1;
345 for (j=0; j < mpiSim->getMyNlocal(); j++) {
346 if (atoms[j]->getGlobalIndex() == globalIndex) atomIndex = j;
347 }
348 if (atomIndex == -1) {
349 sprintf( painCave.errMsg,
350 "Initialize from file error. Atom at index %d "
351 "in file %s does not exist on processor %d .\n",
352 globalIndex, c_in_name, mpiSim->getMyNode() );
353 return strdup( painCave.errMsg );
354 }
355 #else
356 n_atoms = simnfo->n_atoms;
357 atomIndex = globalIndex;
358 #endif // is_mpi
359
360 // set the string tokenizer
361
362 foo = strtok(readLine, " ,;\t");
363
364 // check the atom name to the current atom
365
366 if( strcmp( foo, atoms[atomIndex]->getType() ) ){
367 sprintf( painCave.errMsg,
368 "Initialize from file error. Atom %s at index %d "
369 "in file %s does not"
370 " match the BASS atom %s.\n",
371 foo, atomIndex, c_in_name, atoms[atomIndex]->getType() );
372 return strdup( painCave.errMsg );
373 }
374
375 // get the positions
376
377 foo = strtok(NULL, " ,;\t");
378 if(foo == NULL){
379 sprintf( painCave.errMsg,
380 "error in reading postition x from %s\n"
381 "natoms = %d, index = %d\n",
382 c_in_name, n_atoms, atomIndex );
383 return strdup( painCave.errMsg );
384 }
385 rx = atof( foo );
386
387 foo = strtok(NULL, " ,;\t");
388 if(foo == NULL){
389 sprintf( painCave.errMsg,
390 "error in reading postition y from %s\n"
391 "natoms = %d, index = %d\n",
392 c_in_name, n_atoms, atomIndex );
393 return strdup( painCave.errMsg );
394 }
395 ry = atof( foo );
396
397 foo = strtok(NULL, " ,;\t");
398 if(foo == NULL){
399 sprintf( painCave.errMsg,
400 "error in reading postition z from %s\n"
401 "natoms = %d, index = %d\n",
402 c_in_name, n_atoms, atomIndex );
403 return strdup( painCave.errMsg );
404 }
405 rz = atof( foo );
406
407
408 // get the velocities
409
410 foo = strtok(NULL, " ,;\t");
411 if(foo == NULL){
412 sprintf( painCave.errMsg,
413 "error in reading velocity x from %s\n"
414 "natoms = %d, index = %d\n",
415 c_in_name, n_atoms, atomIndex );
416 return strdup( painCave.errMsg );
417 }
418 vx = atof( foo );
419
420 foo = strtok(NULL, " ,;\t");
421 if(foo == NULL){
422 sprintf( painCave.errMsg,
423 "error in reading velocity y from %s\n"
424 "natoms = %d, index = %d\n",
425 c_in_name, n_atoms, atomIndex );
426 return strdup( painCave.errMsg );
427 }
428 vy = atof( foo );
429
430 foo = strtok(NULL, " ,;\t");
431 if(foo == NULL){
432 sprintf( painCave.errMsg,
433 "error in reading velocity z from %s\n"
434 "natoms = %d, index = %d\n",
435 c_in_name, n_atoms, atomIndex );
436 return strdup( painCave.errMsg );
437 }
438 vz = atof( foo );
439
440
441 // get the quaternions
442
443 if( atoms[atomIndex]->isDirectional() ){
444
445 foo = strtok(NULL, " ,;\t");
446 if(foo == NULL){
447 sprintf(painCave.errMsg,
448 "error in reading quaternion 0 from %s\n"
449 "natoms = %d, index = %d\n",
450 c_in_name, n_atoms, atomIndex );
451 return strdup( painCave.errMsg );
452 }
453 q[0] = atof( foo );
454
455 foo = strtok(NULL, " ,;\t");
456 if(foo == NULL){
457 sprintf( painCave.errMsg,
458 "error in reading quaternion 1 from %s\n"
459 "natoms = %d, index = %d\n",
460 c_in_name, n_atoms, atomIndex );
461 return strdup( painCave.errMsg );
462 }
463 q[1] = atof( foo );
464
465 foo = strtok(NULL, " ,;\t");
466 if(foo == NULL){
467 sprintf( painCave.errMsg,
468 "error in reading quaternion 2 from %s\n"
469 "natoms = %d, index = %d\n",
470 c_in_name, n_atoms, atomIndex );
471 return strdup( painCave.errMsg );
472 }
473 q[2] = atof( foo );
474
475 foo = strtok(NULL, " ,;\t");
476 if(foo == NULL){
477 sprintf( painCave.errMsg,
478 "error in reading quaternion 3 from %s\n"
479 "natoms = %d, index = %d\n",
480 c_in_name, n_atoms, atomIndex );
481 return strdup( painCave.errMsg );
482 }
483 q[3] = atof( foo );
484
485 // get the angular velocities
486
487 foo = strtok(NULL, " ,;\t");
488 if(foo == NULL){
489 sprintf( painCave.errMsg,
490 "error in reading angular momentum jx from %s\n"
491 "natoms = %d, index = %d\n",
492 c_in_name, n_atoms, atomIndex );
493 return strdup( painCave.errMsg );
494 }
495 jx = atof( foo );
496
497 foo = strtok(NULL, " ,;\t");
498 if(foo == NULL){
499 sprintf( painCave.errMsg,
500 "error in reading angular momentum jy from %s\n"
501 "natoms = %d, index = %d\n",
502 c_in_name, n_atoms, atomIndex );
503 return strdup( painCave.errMsg );
504 }
505 jy = atof(foo );
506
507 foo = strtok(NULL, " ,;\t");
508 if(foo == NULL){
509 sprintf( painCave.errMsg,
510 "error in reading angular momentum jz from %s\n"
511 "natoms = %d, index = %d\n",
512 c_in_name, n_atoms, atomIndex );
513 return strdup( painCave.errMsg );
514 }
515 jz = atof( foo );
516
517 dAtom = ( DirectionalAtom* )atoms[atomIndex];
518
519 // check that the quaternion vector is normalized
520
521 qSqr = (q[0] * q[0]) + (q[1] * q[1]) + (q[2] * q[2]) + (q[3] * q[3]);
522
523 qLength = sqrt( qSqr );
524 q[0] = q[0] / qLength;
525 q[1] = q[1] / qLength;
526 q[2] = q[2] / qLength;
527 q[3] = q[3] / qLength;
528
529 dAtom->setQ( q );
530
531 // add the angular velocities
532
533 dAtom->setJx( jx );
534 dAtom->setJy( jy );
535 dAtom->setJz( jz );
536 }
537
538 // add the positions and velocities to the atom
539
540 atoms[atomIndex]->setX( rx );
541 atoms[atomIndex]->setY( ry );
542 atoms[atomIndex]->setZ( rz );
543
544 atoms[atomIndex]->set_vx( vx );
545 atoms[atomIndex]->set_vy( vy );
546 atoms[atomIndex]->set_vz( vz );
547
548 return NULL;
549 }
550
551
552 char* DumpReader::parseCommentLine(char* readLine, double time,
553 double boxMat[9]){
554
555 char *foo; // the pointer to the current string token
556 int j;
557
558 // set the string tokenizer
559
560 foo = strtok(readLine, " ,;\t");
561 if(foo == NULL){
562 sprintf( painCave.errMsg,
563 "error in reading time from %s\n",
564 c_in_name );
565 return strdup( painCave.errMsg );
566 }
567 time = atof( foo );
568
569 // get the Hx vector
570
571 foo = strtok(NULL, " ,;\t");
572 if(foo == NULL){
573 sprintf( painCave.errMsg,
574 "error in reading Hx[0] from %s\n",
575 c_in_name );
576 return strdup( painCave.errMsg );
577 }
578 boxMat[0] = atof( foo );
579
580 foo = strtok(NULL, " ,;\t");
581 if(foo == NULL){
582 sprintf( painCave.errMsg,
583 "error in reading Hx[1] from %s\n",
584 c_in_name );
585 return strdup( painCave.errMsg );
586 }
587 boxMat[1] = atof( foo );
588
589 foo = strtok(NULL, " ,;\t");
590 if(foo == NULL){
591 sprintf( painCave.errMsg,
592 "error in reading Hx[2] from %s\n",
593 c_in_name );
594 return strdup( painCave.errMsg );
595 }
596 boxMat[2] = atof( foo );
597
598 // get the Hy vector
599
600 foo = strtok(NULL, " ,;\t");
601 if(foo == NULL){
602 sprintf( painCave.errMsg,
603 "error in reading Hy[0] from %s\n",
604 c_in_name );
605 return strdup( painCave.errMsg );
606 }
607 boxMat[3] = atof( foo );
608
609 foo = strtok(NULL, " ,;\t");
610 if(foo == NULL){
611 sprintf( painCave.errMsg,
612 "error in reading Hy[1] from %s\n",
613 c_in_name );
614 return strdup( painCave.errMsg );
615 }
616 boxMat[4] = atof( foo );
617
618 foo = strtok(NULL, " ,;\t");
619 if(foo == NULL){
620 sprintf( painCave.errMsg,
621 "error in reading Hy[2] from %s\n",
622 c_in_name );
623 return strdup( painCave.errMsg );
624 }
625 boxMat[5] = atof( foo );
626
627 // get the Hz vector
628
629 foo = strtok(NULL, " ,;\t");
630 if(foo == NULL){
631 sprintf( painCave.errMsg,
632 "error in reading Hz[0] from %s\n",
633 c_in_name );
634 return strdup( painCave.errMsg );
635 }
636 boxMat[6] = atof( foo );
637
638 foo = strtok(NULL, " ,;\t");
639 if(foo == NULL){
640 sprintf( painCave.errMsg,
641 "error in reading Hz[1] from %s\n",
642 c_in_name );
643 return strdup( painCave.errMsg );
644 }
645 boxMat[7] = atof( foo );
646
647 foo = strtok(NULL, " ,;\t");
648 if(foo == NULL){
649 sprintf( painCave.errMsg,
650 "error in reading Hz[2] from %s\n",
651 c_in_name );
652 return strdup( painCave.errMsg );
653 }
654 boxMat[8] = atof( foo );
655
656 return NULL;
657 }
658
659
660 #ifdef IS_MPI
661
662 // a couple of functions to let us escape the read loop
663
664 void initFile::nodeZeroError( void ){
665 int j, myStatus;
666
667 myStatus = 0;
668 for (j = 0; j < mpiSim->getNumberProcessors(); j++) {
669 MPI_Send( &myStatus, 1, MPI_INT, j,
670 TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
671 }
672
673
674 MPI_Finalize();
675 exit (0);
676
677 }
678
679 void initFile::anonymousNodeDie( void ){
680
681 MPI_Finalize();
682 exit (0);
683 }
684
685 #endif //is_mpi