OOPSE/libmdtools/InitializeFromFile.cpp

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

#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 initFile{
  void nodeZeroError( void );
  void anonymousNodeDie( void );
}

using namespace initFile;

#endif // is_mpi

InitializeFromFile::InitializeFromFile( 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, "Infile opened for reading successfully." );
  MPIcheckPoint();
#endif
  return;  
}

InitializeFromFile::~InitializeFromFile( ){
#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, "Infile closed successfully." );
  MPIcheckPoint();
#endif

  return;
}


void InitializeFromFile :: readInit( SimInfo* the_simnfo ){

  int i, 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 n_atoms; // 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;

  double currTime;
  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,
             "InitializeFromFile error: error reading 1st line of \"%s\"\n",
             c_in_name );
    painCave.isFatal = 1;
    simError();
  }

  n_atoms = atoi( read_buffer );

  if( n_atoms != simnfo->n_atoms ){
    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 );
    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 = parseBoxLine( read_buffer, boxMat, currTime );
  if( parseErr != NULL ){
    strcpy( painCave.errMsg, parseErr );
    painCave.isFatal = 1;
    simError();
  }

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

  simnfo->setBoxM( theBoxMat3 );
  simnfo->setTime( currTime );


  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 %s \n ",c_in_name);
      haveError = 1;
      simError();
    }
    
    n_atoms = atoi( read_buffer );
    
    // 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 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 = parseBoxLine( read_buffer, boxMat, currTime );
    if( parseErr != NULL ){
      strcpy( painCave.errMsg, parseErr );
      haveError = 1;
      simError();
    }

    MPI_Bcast(boxMat, 9, MPI_DOUBLE, 0, MPI_COMM_WORLD );
    MPI_Bcast(&currTime, 1, 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(boxMat, 9, MPI_DOUBLE, 0, MPI_COMM_WORLD);
    MPI_Bcast(&currTime, 1, 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;

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

  simnfo->setBoxM( theBoxMat3 );
  simnfo->setTime( currTime );

   
#endif
}

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

  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 jx, jy, jz; // angular velocity placeholders;
  double qSqr, qLength; // needed to normalize the quaternion vector.
  
  Atom **atoms = simnfo->atoms;
  DirectionalAtom* dAtom;
  
  int n_atoms, atomIndex;

#ifdef IS_MPI
  int j;

  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 );
  }
  pos[0] = 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 );
  }
  pos[1] = 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 );
  }
  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"
             "natoms  = %d, index = %d\n",
             c_in_name, n_atoms, atomIndex );
    return strdup( painCave.errMsg );
  }
  vel[0] = 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 );
  }
  vel[1] = 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 );
  }
  vel[2] = 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]->setPos( pos );
  atoms[atomIndex]->setVel( vel );

  return NULL;
}


char* InitializeFromFile::parseBoxLine(char* readLine, double boxMat[9],
                                       double &time ){

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