OOPSE/libmdtools/DumpReader.cpp

#define _FILE_OFFSET_BITS 64

#include <iostream>
#include <cmath>
#include <vector>

#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 );

  
  class FilePos{
  public:
    

  private:
    

  };

  vector<FilePos> frameStart;
}

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