OOPSE/libmdtools/DumpReader.cpp

#include <iostream>
#include <cmath>

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <unistd.h>
#include <sys/types.h>
#include <sys/stat.h>

#include "ReadWrite.hpp"
#include "simError.h"

#ifdef IS_MPI
#include <mpi.h>
#include "mpiSimulation.hpp"
#define TAKE_THIS_TAG_CHAR 0
#define TAKE_THIS_TAG_INT 1

namespace dumpRead{
  void nodeZeroError( void );
  void anonymousNodeDie( void );
}

using namespace dumpRead;

#endif // is_mpi

DumpReader :: DumpReader( char *in_name ){
#ifdef IS_MPI
  if (worldRank == 0) {
#endif

  c_in_file = fopen(in_name, "r");
  if(c_in_file == NULL){
    sprintf(painCave.errMsg,
            "Cannot open file: %s\n", in_name);
    painCave.isFatal = 1;
    simError();
  }
  
  strcpy( c_in_name, in_name);
#ifdef IS_MPI
  }
  strcpy( checkPointMsg, "Dump file opened for reading successfully." );
  MPIcheckPoint();
#endif
  return;  
}

DumpReader :: ~DumpReader( ){
#ifdef IS_MPI
  if (worldRank == 0) {
#endif
  int error;
  error = fclose( c_in_file );
  if( error ){
    sprintf( painCave.errMsg,
             "Error closing %s\n", c_in_name );
    simError();
  }
#ifdef IS_MPI
  }
  strcpy( checkPointMsg, "Dump file closed successfully." );
  MPIcheckPoint();
#endif

  return;
}


void DumpReader :: getFrame( SimInfo* the_simnfo, int whichFrame){

  int i, j, done, which_node, which_atom; // loop counter

  const int BUFFERSIZE = 2000; // size of the read buffer
  int n_atoms; // the number of atoms
  char read_buffer[BUFFERSIZE]; //the line buffer for reading 
#ifdef IS_MPI
  char send_buffer[BUFFERSIZE];
#endif

  char *eof_test; // ptr to see when we reach the end of the file 
  char *parseErr;
  int procIndex;
  double boxMat[9];
  double theBoxMat3[3][3];

  simnfo = the_simnfo;

#ifndef IS_MPI
  eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
  if( eof_test == NULL ){
    sprintf( painCave.errMsg,
             "DumpReader error: error reading 1st line of \"%s\"\n",
             c_in_name );
    painCave.isFatal = 1;
    simError();
  }
  
  n_atoms = atoi( read_buffer );

  Atom **atoms = simnfo->atoms;
  DirectionalAtom* dAtom;

  if( n_atoms != simnfo->n_atoms ){
    sprintf( painCave.errMsg,
             "DumpReader error. %s n_atoms, %d, "
             "does not match the BASS file's n_atoms, %d.\n",
             c_in_name, n_atoms, simnfo->n_atoms );
    painCave.isFatal = 1;
    simError();
  }
  
  //read the box mat from the comment line 
  
  eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
  if(eof_test == NULL){
    sprintf( painCave.errMsg,
             "error in reading commment in %s\n", c_in_name);
    painCave.isFatal = 1;
    simError();
  }

  parseErr = parseCommentLine( read_buffer, time, boxMat );
  if( parseErr != NULL ){
    strcpy( painCave.errMsg, parseErr );
    painCave.isFatal = 1;
    simError();
  }

  simnfo->setTime( time );
  
  for(i=0;i<3;i++)
    for(j=0;j<3;j++) theBoxMat3[i][j] = boxMat[3*j+i];
  

  simnfo->setBoxM( theBoxMat3 );
  

  for( i=0; i < n_atoms; i++){
    
    eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
    if(eof_test == NULL){
      sprintf(painCave.errMsg,
              "error in reading file %s\n"
              "natoms  = %d; index = %d\n"
              "error reading the line from the file.\n",
              c_in_name, n_atoms, i );
      painCave.isFatal = 1;
      simError();
    }

    
    parseErr = parseDumpLine( read_buffer, i );
    if( parseErr != NULL ){
      strcpy( painCave.errMsg, parseErr );
      painCave.isFatal = 1;
      simError();
    }
  }


  // MPI Section of code..........
#else //IS_MPI

  // first thing first, suspend fatalities.
  painCave.isEventLoop = 1;

  int myStatus; // 1 = wakeup & success; 0 = error; -1 = AllDone
  int haveError;
  
  MPI_Status istatus;
  int *AtomToProcMap = mpiSim->getAtomToProcMap();

  
  haveError = 0;
  if (worldRank == 0) {

    eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
    if( eof_test == NULL ){
      sprintf( painCave.errMsg,
               "Error reading 1st line of %d \n ",c_in_name);
      haveError = 1;
      simError();
    }
    
    n_atoms = atoi( read_buffer );
    
    Atom **atoms = simnfo->atoms;
    DirectionalAtom* dAtom;
 
    // Check to see that the number of atoms in the intial configuration file is the
    // same as declared in simBass.
    
    if( n_atoms != mpiSim->getTotAtoms() ){
      sprintf( painCave.errMsg,
               "Initialize from File error. %s n_atoms, %d, "
               "does not match the BASS file's n_atoms, %d.\n",
               c_in_name, n_atoms, simnfo->n_atoms );
      haveError= 1;
      simError();
    }
    
    //read the time and boxMat from the comment line 
    
    eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
    if(eof_test == NULL){
      sprintf( painCave.errMsg,
               "error in reading commment in %s\n", c_in_name);
      haveError = 1;
      simError();
    }
    
    parseErr = parseCommentLine( read_buffer, time, boxMat );
    if( parseErr != NULL ){
      strcpy( painCave.errMsg, parseErr );
      haveError = 1;
      simError();
    }

    MPI_Bcast(time, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD );
    
    MPI_Bcast(boxMat, 9, MPI_DOUBLE, 0, MPI_COMM_WORLD );
    
    if(haveError) nodeZeroError();
    
    for (i=0 ; i < mpiSim->getTotAtoms(); i++) {
      
      eof_test = fgets(read_buffer, sizeof(read_buffer), c_in_file);
      if(eof_test == NULL){
        sprintf(painCave.errMsg,
                "error in reading file %s\n"
                "natoms  = %d; index = %d\n"
                "error reading the line from the file.\n",
                c_in_name, n_atoms, i );
        haveError= 1;
        simError();
      }
    
      if(haveError) nodeZeroError();

      // Get the Node number which wants this atom:
      which_node = AtomToProcMap[i];     
      if (which_node == 0) {
        parseErr = parseDumpLine( read_buffer, i );
        if( parseErr != NULL ){
          strcpy( painCave.errMsg, parseErr );
          haveError = 1;
          simError();
        }    
        if(haveError) nodeZeroError();
      } 
      
      else {
       
        myStatus = 1;
        MPI_Send(&myStatus, 1, MPI_INT, which_node, 
                 TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
        MPI_Send(read_buffer, BUFFERSIZE, MPI_CHAR, which_node, 
                 TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD);
        MPI_Send(&i, 1, MPI_INT, which_node, TAKE_THIS_TAG_INT,
                 MPI_COMM_WORLD);
        MPI_Recv(&myStatus, 1, MPI_INT, which_node, TAKE_THIS_TAG_INT, 
                 MPI_COMM_WORLD, &istatus);
        
        if(!myStatus) nodeZeroError();
      }
    }
    myStatus = -1;
    for (j = 0; j < mpiSim->getNumberProcessors(); j++) {      
      MPI_Send( &myStatus, 1, MPI_INT, j, 
                TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
    }
    
  } else {
    
    MPI_Bcast(time, 1, MPI_DOUBLE, 0, MPI_COMM_WORLD);
    MPI_Bcast(boxMat, 9, MPI_DOUBLE, 0, MPI_COMM_WORLD);

    done = 0;
    while (!done) {

      MPI_Recv(&myStatus, 1, MPI_INT, 0, 
               TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus);
      
      if(!myStatus) anonymousNodeDie();
      
      if(myStatus < 0) break;

      MPI_Recv(read_buffer, BUFFERSIZE, MPI_CHAR, 0, 
               TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD, &istatus);
      MPI_Recv(&which_atom, 1, MPI_INT, 0, 
               TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus);
      
      myStatus = 1;
      parseErr = parseDumpLine( read_buffer, which_atom );
      if( parseErr != NULL ){
        strcpy( painCave.errMsg, parseErr );
        myStatus = 0;;
        simError();
      }
      
      MPI_Send( &myStatus, 1, MPI_INT, 0, 
                TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
      
    }
  }
     
 // last  thing last, enable  fatalities.
  painCave.isEventLoop = 0;

  simnfo->setTime( time );

  for(i=0;i<3;i++)
    for(j=0;j<3;j++) theBoxMat3[i][j] = boxMat[3*j+i];

  simnfo->setBoxM( theBoxMat3 );

   
#endif
}

char* DumpReader::parseDumpLine(char* readLine, int globalIndex){

  char *foo; // the pointer to the current string token 
  
  double rx, ry, rz; // position place holders
  double vx, vy, vz; // velocity placeholders
  double q[4]; // the quaternions
  double jx, jy, jz; // angular velocity placeholders;
  double qSqr, qLength; // needed to normalize the quaternion vector.
  
  Atom **atoms = simnfo->atoms;
  DirectionalAtom* dAtom;
  
  int j, n_atoms, atomIndex;

#ifdef IS_MPI
  n_atoms = mpiSim->getTotAtoms();
  atomIndex=-1;        
  for (j=0; j < mpiSim->getMyNlocal(); j++) {
    if (atoms[j]->getGlobalIndex() == globalIndex) atomIndex = j;
  }
  if (atomIndex == -1) {
    sprintf( painCave.errMsg,
             "Initialize from file error. Atom at index %d "
             "in file %s does not exist on processor %d .\n",
             globalIndex, c_in_name, mpiSim->getMyNode() );
    return strdup( painCave.errMsg );
  }  
#else
  n_atoms = simnfo->n_atoms;
  atomIndex = globalIndex;
#endif // is_mpi

  // set the string tokenizer
  
  foo = strtok(readLine, " ,;\t");
  
  // check the atom name to the current atom
  
  if( strcmp( foo, atoms[atomIndex]->getType() ) ){
    sprintf( painCave.errMsg,
             "Initialize from file error. Atom %s at index %d "
             "in file %s does not"
             " match the BASS atom %s.\n",
             foo, atomIndex, c_in_name, atoms[atomIndex]->getType() );
    return strdup( painCave.errMsg );
  }
    
  // get the positions

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading postition x from %s\n"
             "natoms  = %d, index = %d\n",
             c_in_name, n_atoms, atomIndex );
    return strdup( painCave.errMsg );
  }
  rx = atof( foo );
   
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading postition y from %s\n"
             "natoms  = %d, index = %d\n",
             c_in_name, n_atoms, atomIndex );
    return strdup( painCave.errMsg );
  }
  ry = atof( foo );
     
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading postition z from %s\n"
             "natoms  = %d, index = %d\n",
             c_in_name, n_atoms, atomIndex );
    return strdup( painCave.errMsg );
  }
  rz = atof( foo );    


  // get the velocities

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading velocity x from %s\n"
             "natoms  = %d, index = %d\n",
             c_in_name, n_atoms, atomIndex );
    return strdup( painCave.errMsg );
  }
  vx = atof( foo );
    
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading velocity y from %s\n"
             "natoms  = %d, index = %d\n",
             c_in_name, n_atoms, atomIndex );
    return strdup( painCave.errMsg );
  }
  vy = atof( foo );
    
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading velocity z from %s\n"
             "natoms  = %d, index = %d\n",
             c_in_name, n_atoms, atomIndex );
    return strdup( painCave.errMsg );
  }
  vz = atof( foo );
    
    
  // get the quaternions
    
  if( atoms[atomIndex]->isDirectional() ){
      
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf(painCave.errMsg,
              "error in reading quaternion 0 from %s\n"
              "natoms  = %d, index = %d\n",
              c_in_name, n_atoms, atomIndex );
      return strdup( painCave.errMsg );
    }
    q[0] = atof( foo );
      
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading quaternion 1 from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, atomIndex );
      return strdup( painCave.errMsg );
    }
    q[1] = atof( foo );
      
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading quaternion 2 from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, atomIndex );
      return strdup( painCave.errMsg );
    }
    q[2] = atof( foo );
      
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading quaternion 3 from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, atomIndex );
      return strdup( painCave.errMsg );
    }
    q[3] = atof( foo );
      
    // get the angular velocities
      
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading angular momentum jx from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, atomIndex );
      return strdup( painCave.errMsg );
    }
    jx = atof( foo );
      
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading angular momentum jy from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, atomIndex );
      return strdup( painCave.errMsg );
    }
    jy = atof(foo );
      
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading angular momentum jz from %s\n"
               "natoms  = %d, index = %d\n",
               c_in_name, n_atoms, atomIndex );
      return strdup( painCave.errMsg );
    }
    jz = atof( foo );
      
    dAtom = ( DirectionalAtom* )atoms[atomIndex];

    // check that the quaternion vector is normalized

    qSqr = (q[0] * q[0]) + (q[1] * q[1]) + (q[2] * q[2]) + (q[3] * q[3]);
      
    qLength = sqrt( qSqr );
    q[0] = q[0] / qLength;
    q[1] = q[1] / qLength;
    q[2] = q[2] / qLength;
    q[3] = q[3] / qLength;
      
    dAtom->setQ( q );
      
    // add the angular velocities

    dAtom->setJx( jx );
    dAtom->setJy( jy );
    dAtom->setJz( jz );
  }
    
  // add the positions and velocities to the atom
    
  atoms[atomIndex]->setX( rx );
  atoms[atomIndex]->setY( ry );
  atoms[atomIndex]->setZ( rz );
    
  atoms[atomIndex]->set_vx( vx );
  atoms[atomIndex]->set_vy( vy );
  atoms[atomIndex]->set_vz( vz );

  return NULL;
}


char* DumpReader::parseCommentLine(char* readLine, double time, 
                                   double boxMat[9]){

  char *foo; // the pointer to the current string token 
  int j;

  // set the string tokenizer
  
  foo = strtok(readLine, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading time from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  time = atof( foo );

  // get the Hx vector

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Hx[0] from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  boxMat[0] = atof( foo );
   
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Hx[1] from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  boxMat[1] = atof( foo );
     
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Hx[2] from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  boxMat[2] = atof( foo );    

  // get the Hy vector

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Hy[0] from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  boxMat[3] = atof( foo );
   
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Hy[1] from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  boxMat[4] = atof( foo );
     
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Hy[2] from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  boxMat[5] = atof( foo );    

  // get the Hz vector

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Hz[0] from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  boxMat[6] = atof( foo );
   
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Hz[1] from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  boxMat[7] = atof( foo );
     
  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Hz[2] from %s\n",
             c_in_name );
    return strdup( painCave.errMsg );
  }
  boxMat[8] = atof( foo );    

  return NULL;
}


#ifdef IS_MPI

// a couple of functions to let us escape the read loop

void initFile::nodeZeroError( void ){
  int j, myStatus;
  
  myStatus = 0;
  for (j = 0; j < mpiSim->getNumberProcessors(); j++) {      
    MPI_Send( &myStatus, 1, MPI_INT, j, 
              TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
  }  
  

  MPI_Finalize();
  exit (0);
  
}

void initFile::anonymousNodeDie( void ){

  MPI_Finalize();
  exit (0);
}

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