OOPSE/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 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;

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

  //parse dump lines

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

    //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->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(read_buffer, BUFFERSIZE, MPI_CHAR, 0, MPI_COMM_WORLD);

    parseErr = parseCommentLine( read_buffer, simnfo);

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


    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;


#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::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
  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
  for(int i = 0 ; i < 9; i++){
    foo = strtok(NULL, " ,;\t");
    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 );
  }

  //push eta into SimInfo::properties which can be
  //retrieved by integrator later
  //entry_plug->setBoxM( theBoxMat3 );
  DoubleArrayData* etaValue = new DoubleArrayData();
  etaValue->setID(ETAVALUE_ID);
  etaValue->setData(eta, 9);
  entry_plug->addProperty(etaValue);


  return NULL;
}

#ifdef IS_MPI

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