OOPSE/libmdtools/DumpReader.cpp

#define _FILE_OFFSET_BITS 64

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