ViewVC Help
View File | Revision Log | Show Annotations | View Changeset | Root Listing
root/group/trunk/OOPSE/libmdtools/InitializeFromFile.cpp
Revision: 436
Committed: Fri Mar 28 21:45:03 2003 UTC (21 years, 3 months ago) by chuckv
File size: 13245 byte(s)
Log Message:
Bug fixes in read-write routines.

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 <mpi++.h>
17 #include "mpiSimulation.hpp"
18 #define TAKE_THIS_TAG_CHAR 0
19 #define TAKE_THIS_TAG_INT 1
20
21 void nodeZeroError( void );
22 void anonymousNodeDie( void );
23
24 #endif // is_mpi
25
26 InitializeFromFile :: InitializeFromFile( char *in_name ){
27 #ifdef IS_MPI
28 if (worldRank == 0) {
29 #endif
30
31 c_in_file = fopen(in_name, "r");
32 if(c_in_file == NULL){
33 sprintf(painCave.errMsg,
34 "Cannot open file: %s\n", in_name);
35 painCave.isFatal = 1;
36 simError();
37 }
38
39 strcpy( c_in_name, in_name);
40 #ifdef IS_MPI
41 }
42 strcpy( checkPointMsg, "Infile opened for reading successfully." );
43 MPIcheckPoint();
44 #endif
45 return;
46 }
47
48 InitializeFromFile :: ~InitializeFromFile( ){
49 #ifdef IS_MPI
50 if (worldRank == 0) {
51 #endif
52 int error;
53 error = fclose( c_in_file );
54 if( error ){
55 sprintf( painCave.errMsg,
56 "Error closing %s\n", c_in_name );
57 simError();
58 }
59 #ifdef IS_MPI
60 }
61 strcpy( checkPointMsg, "Infile closed successfully." );
62 MPIcheckPoint();
63 #endif
64
65 return;
66 }
67
68
69 void InitializeFromFile :: read_xyz( SimInfo* the_entry_plug ){
70
71 int i, j, done, which_node, which_atom; // loop counter
72
73 const int BUFFERSIZE = 2000; // size of the read buffer
74 int n_atoms; // the number of atoms
75 char read_buffer[BUFFERSIZE]; //the line buffer for reading
76 #ifdef IS_MPI
77 char send_buffer[BUFFERSIZE];
78 #endif
79
80 char *eof_test; // ptr to see when we reach the end of the file
81 char *parseErr;
82 int procIndex;
83
84 entry_plug = the_entry_plug;
85
86
87 #ifndef IS_MPI
88 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
89 if( eof_test == NULL ){
90 sprintf( painCave.errMsg,
91 "InitializeFromFile error: error reading 1st line of \"%s\"\n",
92 c_in_name );
93 painCave.isFatal = 1;
94 simError();
95 }
96
97 n_atoms = atoi( read_buffer );
98
99 Atom **atoms = entry_plug->atoms;
100 DirectionalAtom* dAtom;
101
102 if( n_atoms != entry_plug->n_atoms ){
103 sprintf( painCave.errMsg,
104 "Initialize from File error. %s n_atoms, %d, "
105 "does not match the BASS file's n_atoms, %d.\n",
106 c_in_name, n_atoms, entry_plug->n_atoms );
107 painCave.isFatal = 1;
108 simError();
109 }
110
111 //read and toss the comment line
112
113 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
114 if(eof_test == NULL){
115 sprintf( painCave.errMsg,
116 "error in reading commment in %s\n", c_in_name);
117 painCave.isFatal = 1;
118 simError();
119 }
120
121 for( i=0; i < n_atoms; i++){
122
123 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
124 if(eof_test == NULL){
125 sprintf(painCave.errMsg,
126 "error in reading file %s\n"
127 "natoms = %d; index = %d\n"
128 "error reading the line from the file.\n",
129 c_in_name, n_atoms, i );
130 painCave.isFatal = 1;
131 simError();
132 }
133
134
135 parseErr = parseDumpLine( read_buffer, i );
136 if( parseErr != NULL ){
137 strcpy( painCave.errMsg, parseErr );
138 painCave.isFatal = 1;
139 simError();
140 }
141 }
142
143
144 // MPI Section of code..........
145 #else //IS_MPI
146
147 // first thing first, suspend fatalities.
148 painCave.isEventLoop = 1;
149
150 int myStatus; // 1 = wakeup & success; 0 = error; -1 = AllDone
151 int haveError;
152
153 MPI::Status istatus;
154 int *AtomToProcMap = mpiSim->getAtomToProcMap();
155
156
157 haveError = 0;
158 if (worldRank == 0) {
159
160 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
161 if( eof_test == NULL ){
162 sprintf( painCave.errMsg,
163 "Error reading 1st line of %d \n ",c_in_name);
164 haveError = 1;
165 simError();
166 }
167
168 n_atoms = atoi( read_buffer );
169
170 Atom **atoms = entry_plug->atoms;
171 DirectionalAtom* dAtom;
172
173 // Check to see that the number of atoms in the intial configuration file is the
174 // same as declared in simBass.
175
176 if( n_atoms != mpiSim->getTotAtoms() ){
177 sprintf( painCave.errMsg,
178 "Initialize from File error. %s n_atoms, %d, "
179 "does not match the BASS file's n_atoms, %d.\n",
180 c_in_name, n_atoms, entry_plug->n_atoms );
181 haveError= 1;
182 simError();
183 }
184
185 //read and toss the comment line
186
187 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
188 if(eof_test == NULL){
189 sprintf( painCave.errMsg,
190 "error in reading commment in %s\n", c_in_name);
191 haveError= 1;
192 simError();
193 }
194
195 if(haveError) nodeZeroError();
196
197 for (i=0 ; i < mpiSim->getTotAtoms(); i++) {
198
199 eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
200 if(eof_test == NULL){
201 sprintf(painCave.errMsg,
202 "error in reading file %s\n"
203 "natoms = %d; index = %d\n"
204 "error reading the line from the file.\n",
205 c_in_name, n_atoms, i );
206 haveError= 1;
207 simError();
208 }
209
210 if(haveError) nodeZeroError();
211
212 // Get the Node number which wants this atom:
213 which_node = AtomToProcMap[i];
214 if (which_node == 0) {
215 parseErr = parseDumpLine( read_buffer, i );
216 if( parseErr != NULL ){
217 strcpy( painCave.errMsg, parseErr );
218 haveError = 1;
219 simError();
220 }
221 if(haveError) nodeZeroError();
222 }
223
224 else {
225
226 myStatus = 1;
227 MPI::COMM_WORLD.Send(&myStatus, 1, MPI_INT, which_node,
228 TAKE_THIS_TAG_INT);
229 MPI::COMM_WORLD.Send(read_buffer, BUFFERSIZE, MPI_CHAR, which_node,
230 TAKE_THIS_TAG_CHAR);
231 MPI::COMM_WORLD.Send(&i, 1, MPI_INT, which_node, TAKE_THIS_TAG_INT);
232 MPI::COMM_WORLD.Recv(&myStatus, 1, MPI_INT, which_node, TAKE_THIS_TAG_INT, istatus);
233
234 if(!myStatus) nodeZeroError();
235 }
236 }
237 myStatus = -1;
238 for (j = 0; j < mpiSim->getNumberProcessors(); j++) {
239 MPI::COMM_WORLD.Send( &myStatus, 1, MPI_INT, j,
240 TAKE_THIS_TAG_INT);
241 }
242
243 } else {
244
245 done = 0;
246 while (!done) {
247
248 MPI::COMM_WORLD.Recv(&myStatus, 1, MPI_INT, 0,
249 TAKE_THIS_TAG_INT, istatus);
250
251 if(!myStatus) anonymousNodeDie();
252
253 if(myStatus < 0) break;
254
255 MPI::COMM_WORLD.Recv(read_buffer, BUFFERSIZE, MPI_CHAR, 0,
256 TAKE_THIS_TAG_CHAR, istatus);
257 MPI::COMM_WORLD.Recv(&which_atom, 1, MPI_INT, 0,
258 TAKE_THIS_TAG_INT, istatus);
259
260 myStatus = 1;
261 parseErr = parseDumpLine( read_buffer, which_atom );
262 if( parseErr != NULL ){
263 strcpy( painCave.errMsg, parseErr );
264 myStatus = 0;;
265 simError();
266 }
267
268 MPI::COMM_WORLD.Send( &myStatus, 1, MPI_INT, 0,
269 TAKE_THIS_TAG_INT);
270
271 }
272 }
273
274 // last thing last, enable fatalities.
275 painCave.isEventLoop = 0;
276
277 #endif
278 }
279
280 char* InitializeFromFile::parseDumpLine(char* readLine, int globalIndex){
281
282 char *foo; // the pointer to the current string token
283
284 double rx, ry, rz; // position place holders
285 double vx, vy, vz; // velocity placeholders
286 double q[4]; // the quaternions
287 double jx, jy, jz; // angular velocity placeholders;
288 double qSqr, qLength; // needed to normalize the quaternion vector.
289
290 Atom **atoms = entry_plug->atoms;
291 DirectionalAtom* dAtom;
292
293 int j, n_atoms, atomIndex;
294
295 #ifdef IS_MPI
296 n_atoms = mpiSim->getTotAtoms();
297 atomIndex=-1;
298 for (j=0; j < mpiSim->getMyNlocal(); j++) {
299 if (atoms[j]->getGlobalIndex() == globalIndex) atomIndex = j;
300 }
301 if (atomIndex == -1) {
302 sprintf( painCave.errMsg,
303 "Initialize from file error. Atom at index %d "
304 "in file %s does not exist on processor %d .\n",
305 globalIndex, c_in_name, mpiSim->getMyNode() );
306 return strdup( painCave.errMsg );
307 }
308 #else
309 n_atoms = entry_plug->n_atoms;
310 atomIndex = globalIndex;
311 #endif // is_mpi
312
313 // set the string tokenizer
314
315 foo = strtok(readLine, " ,;\t");
316
317 // check the atom name to the current atom
318
319 if( strcmp( foo, atoms[atomIndex]->getType() ) ){
320 sprintf( painCave.errMsg,
321 "Initialize from file error. Atom %s at index %d "
322 "in file %s does not"
323 " match the BASS atom %s.\n",
324 foo, atomIndex, c_in_name, atoms[atomIndex]->getType() );
325 return strdup( painCave.errMsg );
326 }
327
328 // get the positions
329
330 foo = strtok(NULL, " ,;\t");
331 if(foo == NULL){
332 sprintf( painCave.errMsg,
333 "error in reading postition x from %s\n"
334 "natoms = %d, index = %d\n",
335 c_in_name, n_atoms, atomIndex );
336 return strdup( painCave.errMsg );
337 }
338 rx = atof( foo );
339
340 foo = strtok(NULL, " ,;\t");
341 if(foo == NULL){
342 sprintf( painCave.errMsg,
343 "error in reading postition y from %s\n"
344 "natoms = %d, index = %d\n",
345 c_in_name, n_atoms, atomIndex );
346 return strdup( painCave.errMsg );
347 }
348 ry = atof( foo );
349
350 foo = strtok(NULL, " ,;\t");
351 if(foo == NULL){
352 sprintf( painCave.errMsg,
353 "error in reading postition z from %s\n"
354 "natoms = %d, index = %d\n",
355 c_in_name, n_atoms, atomIndex );
356 return strdup( painCave.errMsg );
357 }
358 rz = atof( foo );
359
360
361 // get the velocities
362
363 foo = strtok(NULL, " ,;\t");
364 if(foo == NULL){
365 sprintf( painCave.errMsg,
366 "error in reading velocity x from %s\n"
367 "natoms = %d, index = %d\n",
368 c_in_name, n_atoms, atomIndex );
369 return strdup( painCave.errMsg );
370 }
371 vx = atof( foo );
372
373 foo = strtok(NULL, " ,;\t");
374 if(foo == NULL){
375 sprintf( painCave.errMsg,
376 "error in reading velocity y from %s\n"
377 "natoms = %d, index = %d\n",
378 c_in_name, n_atoms, atomIndex );
379 return strdup( painCave.errMsg );
380 }
381 vy = atof( foo );
382
383 foo = strtok(NULL, " ,;\t");
384 if(foo == NULL){
385 sprintf( painCave.errMsg,
386 "error in reading velocity z from %s\n"
387 "natoms = %d, index = %d\n",
388 c_in_name, n_atoms, atomIndex );
389 return strdup( painCave.errMsg );
390 }
391 vz = atof( foo );
392
393
394 // get the quaternions
395
396 if( atoms[atomIndex]->isDirectional() ){
397
398 foo = strtok(NULL, " ,;\t");
399 if(foo == NULL){
400 sprintf(painCave.errMsg,
401 "error in reading quaternion 0 from %s\n"
402 "natoms = %d, index = %d\n",
403 c_in_name, n_atoms, atomIndex );
404 return strdup( painCave.errMsg );
405 }
406 q[0] = atof( foo );
407
408 foo = strtok(NULL, " ,;\t");
409 if(foo == NULL){
410 sprintf( painCave.errMsg,
411 "error in reading quaternion 1 from %s\n"
412 "natoms = %d, index = %d\n",
413 c_in_name, n_atoms, atomIndex );
414 return strdup( painCave.errMsg );
415 }
416 q[1] = atof( foo );
417
418 foo = strtok(NULL, " ,;\t");
419 if(foo == NULL){
420 sprintf( painCave.errMsg,
421 "error in reading quaternion 2 from %s\n"
422 "natoms = %d, index = %d\n",
423 c_in_name, n_atoms, atomIndex );
424 return strdup( painCave.errMsg );
425 }
426 q[2] = atof( foo );
427
428 foo = strtok(NULL, " ,;\t");
429 if(foo == NULL){
430 sprintf( painCave.errMsg,
431 "error in reading quaternion 3 from %s\n"
432 "natoms = %d, index = %d\n",
433 c_in_name, n_atoms, atomIndex );
434 return strdup( painCave.errMsg );
435 }
436 q[3] = atof( foo );
437
438 // get the angular velocities
439
440 foo = strtok(NULL, " ,;\t");
441 if(foo == NULL){
442 sprintf( painCave.errMsg,
443 "error in reading angular momentum jx from %s\n"
444 "natoms = %d, index = %d\n",
445 c_in_name, n_atoms, atomIndex );
446 return strdup( painCave.errMsg );
447 }
448 jx = atof( foo );
449
450 foo = strtok(NULL, " ,;\t");
451 if(foo == NULL){
452 sprintf( painCave.errMsg,
453 "error in reading angular momentum jy from %s\n"
454 "natoms = %d, index = %d\n",
455 c_in_name, n_atoms, atomIndex );
456 return strdup( painCave.errMsg );
457 }
458 jy = atof(foo );
459
460 foo = strtok(NULL, " ,;\t");
461 if(foo == NULL){
462 sprintf( painCave.errMsg,
463 "error in reading angular momentum jz from %s\n"
464 "natoms = %d, index = %d\n",
465 c_in_name, n_atoms, atomIndex );
466 return strdup( painCave.errMsg );
467 }
468 jz = atof( foo );
469
470 dAtom = ( DirectionalAtom* )atoms[atomIndex];
471
472 // check that the quaternion vector is normalized
473
474 qSqr = (q[0] * q[0]) + (q[1] * q[1]) + (q[2] * q[2]) + (q[3] * q[3]);
475
476 qLength = sqrt( qSqr );
477 q[0] = q[0] / qLength;
478 q[1] = q[1] / qLength;
479 q[2] = q[2] / qLength;
480 q[3] = q[3] / qLength;
481
482 dAtom->setQ( q );
483
484 // add the angular velocities
485
486 dAtom->setJx( jx );
487 dAtom->setJy( jy );
488 dAtom->setJz( jz );
489 }
490
491 // add the positions and velocities to the atom
492
493 atoms[atomIndex]->setX( rx );
494 atoms[atomIndex]->setY( ry );
495 atoms[atomIndex]->setZ( rz );
496
497 atoms[atomIndex]->set_vx( vx );
498 atoms[atomIndex]->set_vy( vy );
499 atoms[atomIndex]->set_vz( vz );
500
501 return NULL;
502 }
503
504
505 #ifdef IS_MPI
506
507 // a couple of functions to let us escape the read loop
508
509 void nodeZeroError( void ){
510 int j, myStatus;
511
512 myStatus = 0;
513 for (j = 0; j < mpiSim->getNumberProcessors(); j++) {
514 MPI::COMM_WORLD.Send( &myStatus, 1, MPI_INT, j,
515 TAKE_THIS_TAG_INT);
516 }
517
518
519 MPI_Finalize();
520 exit (0);
521
522 }
523
524 void anonymousNodeDie( void ){
525
526 MPI_Finalize();
527 exit (0);
528 }
529
530 #endif //is_mpi