ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/InitializeFromFile.cpp
Revision: 689
Committed: Tue Aug 12 19:56:49 2003 UTC (20 years, 10 months ago) by tim
File size: 16584 byte(s)
Log Message:
debugging globals

File Contents

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