OOPSE-1.0/libmdtools/InitializeFromFile.cpp

#define _FILE_OFFSET_BITS 64
#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"
#include "GenericData.hpp"

#ifdef IS_MPI
#include <mpi.h>
#include "mpiSimulation.hpp"
#define TAKE_THIS_TAG_CHAR 3134
#define TAKE_THIS_TAG_INT 3135

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
  }
  else{
          sprintf( c_in_name, "mpiNodeParser_%d", worldRank );
  }
  
  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 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;

  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();
  }

  nTotObjs = atoi( read_buffer );

  if( nTotObjs != simnfo->getTotIntegrableObjects() ){
    sprintf( painCave.errMsg,
             "Initialize from File error. %s nIntegrable, %d, "
             "does not match the meta-data file's nIntegrable, %d.\n",
             c_in_name, nTotObjs, simnfo->getTotIntegrableObjects());
    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, 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), 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, 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) {

    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();
    }

    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,
               "Initialize from File error. %s nIntegrable, %d, "
               "does not match the meta-data file's nIntegrable, %d.\n",
               c_in_name, nTotObjs, simnfo->getTotIntegrableObjects());
      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();
    }

    //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), 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, 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), 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, 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* InitializeFromFile::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,
             "Initialize from file 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",
             c_in_name);
    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",
             c_in_name);
    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",
             c_in_name);
    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",
             c_in_name );
    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",
             c_in_name );
    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",
             c_in_name );
    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",
                      c_in_name );
      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",
                      c_in_name );
      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",
                      c_in_name );
      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",
                      c_in_name );
      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",
                      c_in_name );
      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",
                      c_in_name );
      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",
                      c_in_name );
      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* InitializeFromFile::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",
             c_in_name );
    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, c_in_name );
      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", c_in_name );
        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, c_in_name );
          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
    
    DoubleArrayData* etaValue = new DoubleArrayData();
    etaValue->setID(ETAVALUE_ID);
    etaValue->setData(eta, 9);
    entry_plug->addProperty(etaValue);
  }

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