src/io/DumpReader.cpp

#define _LARGEFILE_SOURCE64
#define _FILE_OFFSET_BITS 64

#include <sys/types.h>
#include <sys/stat.h>

#include <iostream>
#include <math.h>

#include <stdio.h>
#include <stdlib.h>
#include <string.h>


#include "io/ReadWrite.hpp"
#include "utils/simError.h"

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


DumpReader :: DumpReader(const char *in_name ){

  isScanned = false;

#ifdef IS_MPI
  if (worldRank == 0) {
#endif

  inFile = fopen(in_name, "r");
  if(inFile == NULL){
    sprintf(painCave.errMsg,
            "Cannot open file: %s\n", in_name);
    painCave.isFatal = 1;
    simError();
  }
  
  inFileName = 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
  vector<fpos_t*>::iterator i;

  int error;
  error = fclose( inFile );
  if( error ){
    sprintf( painCave.errMsg,
             "Error closing %s\n", inFileName.c_str());
    simError();
  }

  for(i = framePos.begin(); i != framePos.end(); ++i)
    delete *i;
  framePos.clear();
  
#ifdef IS_MPI
  }
  strcpy( checkPointMsg, "Dump file closed successfully." );
  MPIcheckPoint();
#endif

  return;
}

int DumpReader::getNframes( void ){

  if( !isScanned )
    scanFile();
  return framePos.size();
}

void DumpReader::scanFile( void ){

  int i, j;
  int lineNum = 0;
  char readBuffer[2000];
  fpos_t *currPos;

#ifdef IS_MPI
  if( worldRank == 0 ){
#endif // is_mpi
    
    rewind( inFile );
    
    currPos = new fpos_t;
    fgetpos( inFile, currPos );
    fgets( readBuffer, sizeof( readBuffer ), inFile );
    lineNum++;
    if( feof( inFile ) ){
      sprintf( painCave.errMsg,
               "File \"%s\" ended unexpectedly at line %d\n", 
               inFileName.c_str(),
               lineNum );
      painCave.isFatal = 1;
      simError();
    }

    while( !feof( inFile ) ){
      
      framePos.push_back(currPos);

      i = atoi(readBuffer);
      
      fgets( readBuffer, sizeof( readBuffer ), inFile );
      lineNum++;    
      if( feof( inFile ) ){
        sprintf( painCave.errMsg,
                 "File \"%s\" ended unexpectedly at line %d\n",
                 inFileName.c_str(),
                 lineNum );
        painCave.isFatal = 1;
        simError();
      }
            
      for(j=0; j<i; j++){
        
        fgets( readBuffer, sizeof( readBuffer ), inFile );
        lineNum++;    
        if( feof( inFile ) ){
          sprintf( painCave.errMsg,
                   "File \"%s\" ended unexpectedly at line %d,"
                   " with atom %d\n", 
                   inFileName.c_str(),
                   lineNum, 
                   j );
          painCave.isFatal = 1;
          simError();
        }
        
      }
      
      currPos = new fpos_t;
      fgetpos( inFile, currPos );
      fgets( readBuffer, sizeof( readBuffer ), inFile );
      lineNum++;
    }
    
    delete currPos;
    rewind( inFile );
     
    isScanned = true;

#ifdef IS_MPI
  }
  strcpy( checkPointMsg, "Successfully scanned DumpFile\n" );
  MPIcheckPoint();
#endif // is_mpi
}

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

  simnfo = the_simnfo;

  this->readSet( whichFrame );
}


void DumpReader :: readSet( int whichFrame ){

  int i;
  unsigned int j;

#ifdef IS_MPI
  int done, which_node, which_atom; // loop counter
#endif //is_mpi

  const int BUFFERSIZE = 2000; // size of the read buffer
  int nTotObjs; // the number of atoms
  char read_buffer[BUFFERSIZE]; //the line buffer for reading

  char *eof_test; // ptr to see when we reach the end of the file
  char *parseErr;

  vector<StuntDouble*> integrableObjects;


#ifndef IS_MPI

  fsetpos(inFile, framePos[whichFrame]);
  eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
  if( eof_test == NULL ){
    sprintf( painCave.errMsg,
             "DumpReader error: error reading 1st line of \"%s\"\n",
             inFileName.c_str() );
    painCave.isFatal = 1;
    simError();
  }

  nTotObjs = atoi( read_buffer );

  if( nTotObjs != simnfo->getTotIntegrableObjects() ){
    sprintf( painCave.errMsg,
             "DumpReader error. %s nIntegrable, %d, "
             "does not match the meta-data file's nIntegrable, %d.\n",
             inFileName.c_str(), nTotObjs, simnfo->getTotIntegrableObjects());
    painCave.isFatal = 1;
    simError();
  }

  //read the box mat from the comment line

  eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
  if(eof_test == NULL){
    sprintf( painCave.errMsg,
             "error in reading commment in %s\n", inFileName.c_str());
    painCave.isFatal = 1;
    simError();
  }

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

  //parse dump lines

  for( i=0; i < simnfo->n_mol; i++){

    integrableObjects = (simnfo->molecules[i]).getIntegrableObjects();

    for(j = 0; j < integrableObjects.size(); j++){

      eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
      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",
              inFileName.c_str(), nTotObjs, i );
        painCave.isFatal = 1;
        simError();
      }
      
      parseErr = parseDumpLine( read_buffer, integrableObjects[j]);
      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 *MolToProcMap = mpiSim->getMolToProcMap();
  int localIndex;
  int nCurObj;
  int nitems;

  nTotObjs = simnfo->getTotIntegrableObjects();
  haveError = 0;
  if (worldRank == 0) {
     fsetpos(inFile,  framePos[whichFrame]);

    eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
    if( eof_test == NULL ){
      sprintf( painCave.errMsg,
               "Error reading 1st line of %s \n ",inFileName.c_str());
      haveError = 1;
      simError();
    }

    nitems = atoi( read_buffer );

    // Check to see that the number of integrable objects in the
    // intial configuration file is the same as derived from the
    // meta-data file.

    if( nTotObjs != nitems){
      sprintf( painCave.errMsg,
               "DumpReader Error. %s nIntegrable, %d, "
               "does not match the meta-data file's nIntegrable, %d.\n",
               inFileName.c_str(), nTotObjs, simnfo->getTotIntegrableObjects());
      haveError= 1;
      simError();
    }

    //read the boxMat from the comment line

    eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
    if(eof_test == NULL){
      sprintf( painCave.errMsg,
               "error in reading commment in %s\n", inFileName.c_str());
      haveError = 1;
      simError();
    }

    //Every single processor will parse the comment line by itself
    //By using this way, we might lose some efficiency, but if we want to add
    //more parameters into comment line, we only need to modify function
    //parseCommentLine

    MPI_Bcast(read_buffer, BUFFERSIZE, MPI_CHAR, 0, MPI_COMM_WORLD);

    parseErr = parseCommentLine( read_buffer, simnfo);

    if( parseErr != NULL ){
      strcpy( painCave.errMsg, parseErr );
      haveError = 1;
      simError();
    }

    for (i=0 ; i < mpiSim->getNMolGlobal(); i++) {
      which_node = MolToProcMap[i];
      if(which_node == 0){
       //molecules belong to master node

      localIndex = mpiSim->getGlobalToLocalMol(i);

      if(localIndex == -1) {
        strcpy(painCave.errMsg, "Molecule not found on node 0!");
        haveError = 1;
        simError();
      }

       integrableObjects = (simnfo->molecules[localIndex]).getIntegrableObjects();
       for(j=0; j < integrableObjects.size(); j++){
        
          eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
          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",
                    inFileName.c_str(), nTotObjs, i );
                haveError= 1;
                simError();
          }
          
          if(haveError) nodeZeroError();

          parseDumpLine(read_buffer, integrableObjects[j]);
           
       }


      }
      else{
      //molecule belongs to slave nodes

        MPI_Recv(&nCurObj, 1, MPI_INT, which_node,
               TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus);
      
       for(j=0; j < nCurObj; j++){
        
          eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
          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",
                    inFileName.c_str(), nTotObjs, i );
                haveError= 1;
                simError();
          }
          
          if(haveError) nodeZeroError();

            MPI_Send(read_buffer, BUFFERSIZE, MPI_CHAR, which_node,
                      TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD);
           
       }

      }
      
    }
    
  }
  else{
  //actions taken at slave nodes
    MPI_Bcast(read_buffer, BUFFERSIZE, MPI_CHAR, 0, MPI_COMM_WORLD);

    parseErr = parseCommentLine( read_buffer, simnfo);

    if( parseErr != NULL ){
      strcpy( painCave.errMsg, parseErr );
      haveError = 1;
      simError();
    }
  
    for (i=0 ; i < mpiSim->getNMolGlobal(); i++) {
      which_node = MolToProcMap[i];
      
      if(which_node == worldRank){
      //molecule with global index i belongs to this processor
      
        localIndex = mpiSim->getGlobalToLocalMol(i);

        if(localIndex == -1) {
          sprintf(painCave.errMsg, "Molecule not found on node %d\n", worldRank);
          haveError = 1;
          simError();
        }

        integrableObjects = (simnfo->molecules[localIndex]).getIntegrableObjects();        

        nCurObj = integrableObjects.size();
        
        MPI_Send(&nCurObj, 1, MPI_INT, 0,
                        TAKE_THIS_TAG_INT, MPI_COMM_WORLD);

        for(j = 0; j < integrableObjects.size(); j++){

          MPI_Recv(read_buffer, BUFFERSIZE, MPI_CHAR, 0,
                              TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD, &istatus);

          parseErr = parseDumpLine(read_buffer, integrableObjects[j]);

          if( parseErr != NULL ){
                strcpy( painCave.errMsg, parseErr );
                simError();
          }

        }
          
      }
      
    }

  }

#endif
}

char* DumpReader::parseDumpLine(char* readLine, StuntDouble* sd){

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

  double pos[3]; // position place holders
  double vel[3]; // velocity placeholders
  double q[4]; // the quaternions
  double ji[3]; // angular velocity placeholders;
  double qSqr, qLength; // needed to normalize the quaternion vector.


  // set the string tokenizer

  foo = strtok(readLine, " ,;\t");

  // check the atom name to the current atom

  if( strcmp( foo, sd->getType() ) ){
    sprintf( painCave.errMsg,
             "DumpReader error.  Does not"
             " match the meta-data atom %s.\n",
             sd->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",
             inFileName.c_str());
    return strdup( painCave.errMsg );
  }
  pos[0] = atof( foo );

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading postition y from %s\n",
             inFileName.c_str());
    return strdup( painCave.errMsg );
  }
  pos[1] = atof( foo );

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading postition z from %s\n",
             inFileName.c_str());
    return strdup( painCave.errMsg );
  }
  pos[2] = atof( foo );


  // get the velocities

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading velocity x from %s\n",
             inFileName.c_str() );
    return strdup( painCave.errMsg );
  }
  vel[0] = atof( foo );

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading velocity x from %s\n",
             inFileName.c_str() );
    return strdup( painCave.errMsg );
  }
  vel[1] = atof( foo );

  foo = strtok(NULL, " ,;\t");
  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading velocity x from %s\n",
             inFileName.c_str() );
    return strdup( painCave.errMsg );
  }
  vel[2] = atof( foo );


  // add the positions and velocities to the atom

  sd->setPos( pos );
  sd->setVel( vel );

  if (!sd->isDirectional())
    return NULL;

  // get the quaternions

  if( sd->isDirectional() ){

    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
                     "error in reading velocity x from %s\n",
                      inFileName.c_str() );
      return strdup( painCave.errMsg );
    }
    q[0] = atof( foo );

    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
                     "error in reading velocity x from %s\n",
                      inFileName.c_str() );
      return strdup( painCave.errMsg );
    }
    q[1] = atof( foo );

    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
                     "error in reading velocity x from %s\n",
                      inFileName.c_str() );
      return strdup( painCave.errMsg );
    }
    q[2] = atof( foo );

    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
                     "error in reading velocity x from %s\n",
                      inFileName.c_str() );
      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 velocity x from %s\n",
                      inFileName.c_str() );
      return strdup( painCave.errMsg );
    }
    ji[0] = atof( foo );

    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
                     "error in reading velocity x from %s\n",
                      inFileName.c_str() );
      return strdup( painCave.errMsg );
    }
    ji[1] = atof(foo );

    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
                     "error in reading velocity x from %s\n",
                      inFileName.c_str() );
      return strdup( painCave.errMsg );
    }
    ji[2] = atof( foo );


    // check that the quaternion vector is normalized

    qSqr = (q[0] * q[0]) + (q[1] * q[1]) + (q[2] * q[2]) + (q[3] * q[3]);

    if (fabs(qSqr) < 1e-6) {
      sprintf(painCave.errMsg,
          "initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2 ~ 0).\n");
       return strdup(painCave.errMsg);
    }

    qLength = sqrt( qSqr );
    q[0] = q[0] / qLength;
    q[1] = q[1] / qLength;
    q[2] = q[2] / qLength;
    q[3] = q[3] / qLength;

    // add quaternion and angular velocities

    sd->setQ( q );
    sd->setJ( ji );
  }


  return NULL;
}


char* DumpReader::parseCommentLine(char* readLine, SimInfo* entry_plug){

  double currTime;
  double boxMat[9];
  double theBoxMat3[3][3];
  double chi;
  double integralOfChidt;
  double eta[9];

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

  // set the string tokenizer

  foo = strtok(readLine, " ,;\t");
  // set the timeToken.

  if(foo == NULL){
    sprintf( painCave.errMsg,
             "error in reading Time from %s\n",
             inFileName.c_str() );
    return strdup( painCave.errMsg );
  }

  currTime = atof( foo );
  entry_plug->setTime( currTime );

  //get H-Matrix

  for(int i = 0 ; i < 9; i++){
    foo = strtok(NULL, " ,;\t");
    if(foo == NULL){
      sprintf( painCave.errMsg,
               "error in reading H[%d] from %s\n", i, inFileName.c_str() );
      return strdup( painCave.errMsg );
    }
    boxMat[i] = atof( foo );
  }

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

  //set H-Matrix
  entry_plug->setBoxM( theBoxMat3 );

  //get chi and integralOfChidt, they should appear by pair

  if( entry_plug->useInitXSstate ){
    foo = strtok(NULL, " ,;\t\n");
    if(foo != NULL){
      chi = atof(foo);
      
      foo = strtok(NULL, " ,;\t\n");
      if(foo == NULL){
        sprintf( painCave.errMsg,
                 "chi and integralOfChidt should appear by pair in %s\n", inFileName.c_str() );
        return strdup( painCave.errMsg );
      }
      integralOfChidt = atof( foo );
      
      //push chi and integralOfChidt into SimInfo::properties which can be
      //retrieved by integrator later
      DoubleData* chiValue = new DoubleData();
      chiValue->setID(CHIVALUE_ID);
      chiValue->setData(chi);
      entry_plug->addProperty(chiValue);
      
      DoubleData* integralOfChidtValue = new DoubleData();
      integralOfChidtValue->setID(INTEGRALOFCHIDT_ID);
      integralOfChidtValue->setData(integralOfChidt);
      entry_plug->addProperty(integralOfChidtValue);
      
    }
    else
      return NULL;
    
    //get eta
    foo = strtok(NULL, " ,;\t\n");
    if(foo != NULL ){
  
      for(int i = 0 ; i < 9; i++){
        
        if(foo == NULL){
          sprintf( painCave.errMsg,
                   "error in reading eta[%d] from %s\n", i, inFileName.c_str() );
          return strdup( painCave.errMsg );
        }
        eta[i] = atof( foo );
        foo = strtok(NULL, " ,;\t\n");
      }
    }
    else
      return NULL;
    
    //push eta into SimInfo::properties which can be
    //retrieved by integrator later
    //entry_plug->setBoxM( theBoxMat3 );
    DoubleArrayData* etaValue = new DoubleArrayData();
    etaValue->setID(ETAVALUE_ID);
    etaValue->setData(eta, 9);
    entry_plug->addProperty(etaValue);
  }

  return NULL;
}

#ifdef IS_MPI
void DumpReader::nodeZeroError( void ){
  int j, myStatus;

  myStatus = 0;
  for (j = 0; j < mpiSim->getNProcessors(); j++) {
    MPI_Send( &myStatus, 1, MPI_INT, j,
              TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
  }


  MPI_Finalize();
  exit (0);

}

void DumpReader::anonymousNodeDie( void ){

  MPI_Finalize();
  exit (0);
}
#endif
Revision:	3
Committed:	Fri Sep 24 16:27:58 2004 UTC (20 years, 7 months ago) by tim
File size:	18173 byte(s)
Log Message:	change the #include in source files
#	User	Rev	Content
1	gezelter	2	#define _LARGEFILE_SOURCE64
2			#define _FILE_OFFSET_BITS 64
3
4			#include <sys/types.h>
5			#include <sys/stat.h>
6
7			#include <iostream>
8			#include <math.h>
9
10			#include <stdio.h>
11			#include <stdlib.h>
12			#include <string.h>
13
14
15	tim	3	#include "io/ReadWrite.hpp"
16			#include "utils/simError.h"
17	gezelter	2
18			#ifdef IS_MPI
19			#include <mpi.h>
20	tim	3	#include "brains/mpiSimulation.hpp"
21	gezelter	2	#define TAKE_THIS_TAG_CHAR 0
22			#define TAKE_THIS_TAG_INT 1
23			#endif // is_mpi
24
25
26			DumpReader :: DumpReader(const char *in_name ){
27
28			isScanned = false;
29
30			#ifdef IS_MPI
31			if (worldRank == 0) {
32			#endif
33
34			inFile = fopen(in_name, "r");
35			if(inFile == NULL){
36			sprintf(painCave.errMsg,
37			"Cannot open file: %s\n", in_name);
38			painCave.isFatal = 1;
39			simError();
40			}
41
42			inFileName = 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			vector<fpos_t*>::iterator i;
56
57			int error;
58			error = fclose( inFile );
59			if( error ){
60			sprintf( painCave.errMsg,
61			"Error closing %s\n", inFileName.c_str());
62			simError();
63			}
64
65			for(i = framePos.begin(); i != framePos.end(); ++i)
66			delete *i;
67			framePos.clear();
68
69			#ifdef IS_MPI
70			}
71			strcpy( checkPointMsg, "Dump file closed successfully." );
72			MPIcheckPoint();
73			#endif
74
75			return;
76			}
77
78			int DumpReader::getNframes( void ){
79
80			if( !isScanned )
81			scanFile();
82			return framePos.size();
83			}
84
85			void DumpReader::scanFile( void ){
86
87			int i, j;
88			int lineNum = 0;
89			char readBuffer[2000];
90			fpos_t *currPos;
91
92			#ifdef IS_MPI
93			if( worldRank == 0 ){
94			#endif // is_mpi
95
96			rewind( inFile );
97
98			currPos = new fpos_t;
99			fgetpos( inFile, currPos );
100			fgets( readBuffer, sizeof( readBuffer ), inFile );
101			lineNum++;
102			if( feof( inFile ) ){
103			sprintf( painCave.errMsg,
104			"File \"%s\" ended unexpectedly at line %d\n",
105			inFileName.c_str(),
106			lineNum );
107			painCave.isFatal = 1;
108			simError();
109			}
110
111			while( !feof( inFile ) ){
112
113			framePos.push_back(currPos);
114
115			i = atoi(readBuffer);
116
117			fgets( readBuffer, sizeof( readBuffer ), inFile );
118			lineNum++;
119			if( feof( inFile ) ){
120			sprintf( painCave.errMsg,
121			"File \"%s\" ended unexpectedly at line %d\n",
122			inFileName.c_str(),
123			lineNum );
124			painCave.isFatal = 1;
125			simError();
126			}
127
128			for(j=0; j<i; j++){
129
130			fgets( readBuffer, sizeof( readBuffer ), inFile );
131			lineNum++;
132			if( feof( inFile ) ){
133			sprintf( painCave.errMsg,
134			"File \"%s\" ended unexpectedly at line %d,"
135			" with atom %d\n",
136			inFileName.c_str(),
137			lineNum,
138			j );
139			painCave.isFatal = 1;
140			simError();
141			}
142
143			}
144
145			currPos = new fpos_t;
146			fgetpos( inFile, currPos );
147			fgets( readBuffer, sizeof( readBuffer ), inFile );
148			lineNum++;
149			}
150
151			delete currPos;
152			rewind( inFile );
153
154			isScanned = true;
155
156			#ifdef IS_MPI
157			}
158			strcpy( checkPointMsg, "Successfully scanned DumpFile\n" );
159			MPIcheckPoint();
160			#endif // is_mpi
161			}
162
163			void DumpReader :: readFrame( SimInfo* the_simnfo, int whichFrame){
164
165			simnfo = the_simnfo;
166
167			this->readSet( whichFrame );
168			}
169
170
171
172			void DumpReader :: readSet( int whichFrame ){
173
174			int i;
175			unsigned int j;
176
177			#ifdef IS_MPI
178			int done, which_node, which_atom; // loop counter
179			#endif //is_mpi
180
181			const int BUFFERSIZE = 2000; // size of the read buffer
182			int nTotObjs; // the number of atoms
183			char read_buffer[BUFFERSIZE]; //the line buffer for reading
184
185			char *eof_test; // ptr to see when we reach the end of the file
186			char *parseErr;
187
188			vector<StuntDouble*> integrableObjects;
189
190
191			#ifndef IS_MPI
192
193			fsetpos(inFile, framePos[whichFrame]);
194			eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
195			if( eof_test == NULL ){
196			sprintf( painCave.errMsg,
197			"DumpReader error: error reading 1st line of \"%s\"\n",
198			inFileName.c_str() );
199			painCave.isFatal = 1;
200			simError();
201			}
202
203			nTotObjs = atoi( read_buffer );
204
205			if( nTotObjs != simnfo->getTotIntegrableObjects() ){
206			sprintf( painCave.errMsg,
207			"DumpReader error. %s nIntegrable, %d, "
208			"does not match the meta-data file's nIntegrable, %d.\n",
209			inFileName.c_str(), nTotObjs, simnfo->getTotIntegrableObjects());
210			painCave.isFatal = 1;
211			simError();
212			}
213
214			//read the box mat from the comment line
215
216			eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
217			if(eof_test == NULL){
218			sprintf( painCave.errMsg,
219			"error in reading commment in %s\n", inFileName.c_str());
220			painCave.isFatal = 1;
221			simError();
222			}
223
224			parseErr = parseCommentLine( read_buffer, simnfo);
225			if( parseErr != NULL ){
226			strcpy( painCave.errMsg, parseErr );
227			painCave.isFatal = 1;
228			simError();
229			}
230
231			//parse dump lines
232
233			for( i=0; i < simnfo->n_mol; i++){
234
235			integrableObjects = (simnfo->molecules[i]).getIntegrableObjects();
236
237			for(j = 0; j < integrableObjects.size(); j++){
238
239			eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
240			if(eof_test == NULL){
241			sprintf(painCave.errMsg,
242			"error in reading file %s\n"
243			"natoms = %d; index = %d\n"
244			"error reading the line from the file.\n",
245			inFileName.c_str(), nTotObjs, i );
246			painCave.isFatal = 1;
247			simError();
248			}
249
250			parseErr = parseDumpLine( read_buffer, integrableObjects[j]);
251			if( parseErr != NULL ){
252			strcpy( painCave.errMsg, parseErr );
253			painCave.isFatal = 1;
254			simError();
255			}
256			}
257			}
258
259			// MPI Section of code..........
260			#else //IS_MPI
261
262			// first thing first, suspend fatalities.
263			painCave.isEventLoop = 1;
264
265			int myStatus; // 1 = wakeup & success; 0 = error; -1 = AllDone
266			int haveError;
267
268			MPI_Status istatus;
269			int *MolToProcMap = mpiSim->getMolToProcMap();
270			int localIndex;
271			int nCurObj;
272			int nitems;
273
274			nTotObjs = simnfo->getTotIntegrableObjects();
275			haveError = 0;
276			if (worldRank == 0) {
277			fsetpos(inFile, framePos[whichFrame]);
278
279			eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
280			if( eof_test == NULL ){
281			sprintf( painCave.errMsg,
282			"Error reading 1st line of %s \n ",inFileName.c_str());
283			haveError = 1;
284			simError();
285			}
286
287			nitems = atoi( read_buffer );
288
289			// Check to see that the number of integrable objects in the
290			// intial configuration file is the same as derived from the
291			// meta-data file.
292
293			if( nTotObjs != nitems){
294			sprintf( painCave.errMsg,
295			"DumpReader Error. %s nIntegrable, %d, "
296			"does not match the meta-data file's nIntegrable, %d.\n",
297			inFileName.c_str(), nTotObjs, simnfo->getTotIntegrableObjects());
298			haveError= 1;
299			simError();
300			}
301
302			//read the boxMat from the comment line
303
304			eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
305			if(eof_test == NULL){
306			sprintf( painCave.errMsg,
307			"error in reading commment in %s\n", inFileName.c_str());
308			haveError = 1;
309			simError();
310			}
311
312			//Every single processor will parse the comment line by itself
313			//By using this way, we might lose some efficiency, but if we want to add
314			//more parameters into comment line, we only need to modify function
315			//parseCommentLine
316
317			MPI_Bcast(read_buffer, BUFFERSIZE, MPI_CHAR, 0, MPI_COMM_WORLD);
318
319			parseErr = parseCommentLine( read_buffer, simnfo);
320
321			if( parseErr != NULL ){
322			strcpy( painCave.errMsg, parseErr );
323			haveError = 1;
324			simError();
325			}
326
327			for (i=0 ; i < mpiSim->getNMolGlobal(); i++) {
328			which_node = MolToProcMap[i];
329			if(which_node == 0){
330			//molecules belong to master node
331
332			localIndex = mpiSim->getGlobalToLocalMol(i);
333
334			if(localIndex == -1) {
335			strcpy(painCave.errMsg, "Molecule not found on node 0!");
336			haveError = 1;
337			simError();
338			}
339
340			integrableObjects = (simnfo->molecules[localIndex]).getIntegrableObjects();
341			for(j=0; j < integrableObjects.size(); j++){
342
343			eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
344			if(eof_test == NULL){
345			sprintf(painCave.errMsg,
346			"error in reading file %s\n"
347			"natoms = %d; index = %d\n"
348			"error reading the line from the file.\n",
349			inFileName.c_str(), nTotObjs, i );
350			haveError= 1;
351			simError();
352			}
353
354			if(haveError) nodeZeroError();
355
356			parseDumpLine(read_buffer, integrableObjects[j]);
357
358			}
359
360
361			}
362			else{
363			//molecule belongs to slave nodes
364
365			MPI_Recv(&nCurObj, 1, MPI_INT, which_node,
366			TAKE_THIS_TAG_INT, MPI_COMM_WORLD, &istatus);
367
368			for(j=0; j < nCurObj; j++){
369
370			eof_test = fgets(read_buffer, sizeof(read_buffer), inFile);
371			if(eof_test == NULL){
372			sprintf(painCave.errMsg,
373			"error in reading file %s\n"
374			"natoms = %d; index = %d\n"
375			"error reading the line from the file.\n",
376			inFileName.c_str(), nTotObjs, i );
377			haveError= 1;
378			simError();
379			}
380
381			if(haveError) nodeZeroError();
382
383			MPI_Send(read_buffer, BUFFERSIZE, MPI_CHAR, which_node,
384			TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD);
385
386			}
387
388			}
389
390			}
391
392			}
393			else{
394			//actions taken at slave nodes
395			MPI_Bcast(read_buffer, BUFFERSIZE, MPI_CHAR, 0, MPI_COMM_WORLD);
396
397			parseErr = parseCommentLine( read_buffer, simnfo);
398
399			if( parseErr != NULL ){
400			strcpy( painCave.errMsg, parseErr );
401			haveError = 1;
402			simError();
403			}
404
405			for (i=0 ; i < mpiSim->getNMolGlobal(); i++) {
406			which_node = MolToProcMap[i];
407
408			if(which_node == worldRank){
409			//molecule with global index i belongs to this processor
410
411			localIndex = mpiSim->getGlobalToLocalMol(i);
412
413			if(localIndex == -1) {
414			sprintf(painCave.errMsg, "Molecule not found on node %d\n", worldRank);
415			haveError = 1;
416			simError();
417			}
418
419			integrableObjects = (simnfo->molecules[localIndex]).getIntegrableObjects();
420
421			nCurObj = integrableObjects.size();
422
423			MPI_Send(&nCurObj, 1, MPI_INT, 0,
424			TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
425
426			for(j = 0; j < integrableObjects.size(); j++){
427
428			MPI_Recv(read_buffer, BUFFERSIZE, MPI_CHAR, 0,
429			TAKE_THIS_TAG_CHAR, MPI_COMM_WORLD, &istatus);
430
431			parseErr = parseDumpLine(read_buffer, integrableObjects[j]);
432
433			if( parseErr != NULL ){
434			strcpy( painCave.errMsg, parseErr );
435			simError();
436			}
437
438			}
439
440			}
441
442			}
443
444			}
445
446			#endif
447			}
448
449			char* DumpReader::parseDumpLine(char* readLine, StuntDouble* sd){
450
451			char *foo; // the pointer to the current string token
452
453			double pos[3]; // position place holders
454			double vel[3]; // velocity placeholders
455			double q[4]; // the quaternions
456			double ji[3]; // angular velocity placeholders;
457			double qSqr, qLength; // needed to normalize the quaternion vector.
458
459
460			// set the string tokenizer
461
462			foo = strtok(readLine, " ,;\t");
463
464			// check the atom name to the current atom
465
466			if( strcmp( foo, sd->getType() ) ){
467			sprintf( painCave.errMsg,
468			"DumpReader error. Does not"
469			" match the meta-data atom %s.\n",
470			sd->getType() );
471			return strdup( painCave.errMsg );
472			}
473
474			// get the positions
475
476			foo = strtok(NULL, " ,;\t");
477			if(foo == NULL){
478			sprintf( painCave.errMsg,
479			"error in reading postition x from %s\n",
480			inFileName.c_str());
481			return strdup( painCave.errMsg );
482			}
483			pos[0] = atof( foo );
484
485			foo = strtok(NULL, " ,;\t");
486			if(foo == NULL){
487			sprintf( painCave.errMsg,
488			"error in reading postition y from %s\n",
489			inFileName.c_str());
490			return strdup( painCave.errMsg );
491			}
492			pos[1] = atof( foo );
493
494			foo = strtok(NULL, " ,;\t");
495			if(foo == NULL){
496			sprintf( painCave.errMsg,
497			"error in reading postition z from %s\n",
498			inFileName.c_str());
499			return strdup( painCave.errMsg );
500			}
501			pos[2] = atof( foo );
502
503
504			// get the velocities
505
506			foo = strtok(NULL, " ,;\t");
507			if(foo == NULL){
508			sprintf( painCave.errMsg,
509			"error in reading velocity x from %s\n",
510			inFileName.c_str() );
511			return strdup( painCave.errMsg );
512			}
513			vel[0] = atof( foo );
514
515			foo = strtok(NULL, " ,;\t");
516			if(foo == NULL){
517			sprintf( painCave.errMsg,
518			"error in reading velocity x from %s\n",
519			inFileName.c_str() );
520			return strdup( painCave.errMsg );
521			}
522			vel[1] = atof( foo );
523
524			foo = strtok(NULL, " ,;\t");
525			if(foo == NULL){
526			sprintf( painCave.errMsg,
527			"error in reading velocity x from %s\n",
528			inFileName.c_str() );
529			return strdup( painCave.errMsg );
530			}
531			vel[2] = atof( foo );
532
533
534			// add the positions and velocities to the atom
535
536			sd->setPos( pos );
537			sd->setVel( vel );
538
539			if (!sd->isDirectional())
540			return NULL;
541
542			// get the quaternions
543
544			if( sd->isDirectional() ){
545
546			foo = strtok(NULL, " ,;\t");
547			if(foo == NULL){
548			sprintf( painCave.errMsg,
549			"error in reading velocity x from %s\n",
550			inFileName.c_str() );
551			return strdup( painCave.errMsg );
552			}
553			q[0] = atof( foo );
554
555			foo = strtok(NULL, " ,;\t");
556			if(foo == NULL){
557			sprintf( painCave.errMsg,
558			"error in reading velocity x from %s\n",
559			inFileName.c_str() );
560			return strdup( painCave.errMsg );
561			}
562			q[1] = atof( foo );
563
564			foo = strtok(NULL, " ,;\t");
565			if(foo == NULL){
566			sprintf( painCave.errMsg,
567			"error in reading velocity x from %s\n",
568			inFileName.c_str() );
569			return strdup( painCave.errMsg );
570			}
571			q[2] = atof( foo );
572
573			foo = strtok(NULL, " ,;\t");
574			if(foo == NULL){
575			sprintf( painCave.errMsg,
576			"error in reading velocity x from %s\n",
577			inFileName.c_str() );
578			return strdup( painCave.errMsg );
579			}
580			q[3] = atof( foo );
581
582			// get the angular velocities
583
584			foo = strtok(NULL, " ,;\t");
585			if(foo == NULL){
586			sprintf( painCave.errMsg,
587			"error in reading velocity x from %s\n",
588			inFileName.c_str() );
589			return strdup( painCave.errMsg );
590			}
591			ji[0] = atof( foo );
592
593			foo = strtok(NULL, " ,;\t");
594			if(foo == NULL){
595			sprintf( painCave.errMsg,
596			"error in reading velocity x from %s\n",
597			inFileName.c_str() );
598			return strdup( painCave.errMsg );
599			}
600			ji[1] = atof(foo );
601
602			foo = strtok(NULL, " ,;\t");
603			if(foo == NULL){
604			sprintf( painCave.errMsg,
605			"error in reading velocity x from %s\n",
606			inFileName.c_str() );
607			return strdup( painCave.errMsg );
608			}
609			ji[2] = atof( foo );
610
611
612			// check that the quaternion vector is normalized
613
614			qSqr = (q[0] * q[0]) + (q[1] * q[1]) + (q[2] * q[2]) + (q[3] * q[3]);
615
616			if (fabs(qSqr) < 1e-6) {
617			sprintf(painCave.errMsg,
618			"initial quaternion error (q0^2 + q1^2 + q2^2 + q3^2 ~ 0).\n");
619			return strdup(painCave.errMsg);
620			}
621
622			qLength = sqrt( qSqr );
623			q[0] = q[0] / qLength;
624			q[1] = q[1] / qLength;
625			q[2] = q[2] / qLength;
626			q[3] = q[3] / qLength;
627
628			// add quaternion and angular velocities
629
630			sd->setQ( q );
631			sd->setJ( ji );
632			}
633
634
635
636			return NULL;
637			}
638
639
640			char* DumpReader::parseCommentLine(char* readLine, SimInfo* entry_plug){
641
642			double currTime;
643			double boxMat[9];
644			double theBoxMat3[3][3];
645			double chi;
646			double integralOfChidt;
647			double eta[9];
648
649			char *foo; // the pointer to the current string token
650
651			// set the string tokenizer
652
653			foo = strtok(readLine, " ,;\t");
654			// set the timeToken.
655
656			if(foo == NULL){
657			sprintf( painCave.errMsg,
658			"error in reading Time from %s\n",
659			inFileName.c_str() );
660			return strdup( painCave.errMsg );
661			}
662
663			currTime = atof( foo );
664			entry_plug->setTime( currTime );
665
666			//get H-Matrix
667
668			for(int i = 0 ; i < 9; i++){
669			foo = strtok(NULL, " ,;\t");
670			if(foo == NULL){
671			sprintf( painCave.errMsg,
672			"error in reading H[%d] from %s\n", i, inFileName.c_str() );
673			return strdup( painCave.errMsg );
674			}
675			boxMat[i] = atof( foo );
676			}
677
678			for(int i=0;i<3;i++)
679			for(int j=0;j<3;j++) theBoxMat3[i][j] = boxMat[3*j+i];
680
681			//set H-Matrix
682			entry_plug->setBoxM( theBoxMat3 );
683
684			//get chi and integralOfChidt, they should appear by pair
685
686			if( entry_plug->useInitXSstate ){
687			foo = strtok(NULL, " ,;\t\n");
688			if(foo != NULL){
689			chi = atof(foo);
690
691			foo = strtok(NULL, " ,;\t\n");
692			if(foo == NULL){
693			sprintf( painCave.errMsg,
694			"chi and integralOfChidt should appear by pair in %s\n", inFileName.c_str() );
695			return strdup( painCave.errMsg );
696			}
697			integralOfChidt = atof( foo );
698
699			//push chi and integralOfChidt into SimInfo::properties which can be
700			//retrieved by integrator later
701			DoubleData* chiValue = new DoubleData();
702			chiValue->setID(CHIVALUE_ID);
703			chiValue->setData(chi);
704			entry_plug->addProperty(chiValue);
705
706			DoubleData* integralOfChidtValue = new DoubleData();
707			integralOfChidtValue->setID(INTEGRALOFCHIDT_ID);
708			integralOfChidtValue->setData(integralOfChidt);
709			entry_plug->addProperty(integralOfChidtValue);
710
711			}
712			else
713			return NULL;
714
715			//get eta
716			foo = strtok(NULL, " ,;\t\n");
717			if(foo != NULL ){
718
719			for(int i = 0 ; i < 9; i++){
720
721			if(foo == NULL){
722			sprintf( painCave.errMsg,
723			"error in reading eta[%d] from %s\n", i, inFileName.c_str() );
724			return strdup( painCave.errMsg );
725			}
726			eta[i] = atof( foo );
727			foo = strtok(NULL, " ,;\t\n");
728			}
729			}
730			else
731			return NULL;
732
733			//push eta into SimInfo::properties which can be
734			//retrieved by integrator later
735			//entry_plug->setBoxM( theBoxMat3 );
736			DoubleArrayData* etaValue = new DoubleArrayData();
737			etaValue->setID(ETAVALUE_ID);
738			etaValue->setData(eta, 9);
739			entry_plug->addProperty(etaValue);
740			}
741
742			return NULL;
743			}
744
745			#ifdef IS_MPI
746			void DumpReader::nodeZeroError( void ){
747			int j, myStatus;
748
749			myStatus = 0;
750			for (j = 0; j < mpiSim->getNProcessors(); j++) {
751			MPI_Send( &myStatus, 1, MPI_INT, j,
752			TAKE_THIS_TAG_INT, MPI_COMM_WORLD);
753			}
754
755
756			MPI_Finalize();
757			exit (0);
758
759			}
760
761			void DumpReader::anonymousNodeDie( void ){
762
763			MPI_Finalize();
764			exit (0);
765			}
766			#endif