xmdssimulation.cc

XMDS is a code generator that integrates equations. You write them down in human readable form in a

    if(verbose()) {
      printf("simulation stochastic defaulting to 'no'\n");
    }

    myParameters.stochastic=0;  
  }

  // ************************************
  // find 'use_mpi' assignment

  getAssignmentBools(yourElement,"use_mpi",0,1,myBoolList);

  if(myBoolList.size()==1) {

    myParameters.usempi=*myBoolList.begin();

    if(myParameters.usempi&!myParameters.mpiAvailable) {
      throw xmdsException(yourElement,"MPI was not available when xmds was installed");
    }

    if(verbose()) {
      if(myParameters.usempi) {
        printf("use_mpi is 'yes'\n");
      }
      else {
        printf("use_mpi is 'no'\n");
      }
    }
  }
  else {
    if(verbose()) {
      printf("use_mpi defaulting to 'no'\n");
    }

    myParameters.usempi=0;  
  }
        
  if((myParameters.usempi)&(myParameters.nThreads>1)&(!myParameters.stochastic)) {
    throw xmdsException(yourElement,"Non-determinstic MPI and threaded FFTs are not compatible");
  }
        
        
  // ************************************
  // find 'bing' assignment

  getAssignmentBools(yourElement,"bing",0,1,myBoolList);

  if(myBoolList.size()==1) {

    myParameters.bing=*myBoolList.begin();

    if(verbose()) {
      if(myParameters.bing) {
        printf("bing is 'yes'\n");
      }
      else {
        printf("bing is 'no'\n");
      }
    }
  }
  else {
    if(verbose()) {
      printf("bing defaulting to 'no'\n");
    }

    myParameters.bing=0;  
  }

  // ************************************
  // if stochastic get paths, and seed assignments
  if(myParameters.stochastic) {

    // ************************************
    // get paths

    getAssignmentULongs(yourElement,"paths",1,1,myULongList);
    myParameters.nPaths = *myULongList.begin();

    if(myParameters.nPaths < 1) {
      throw xmdsException(yourElement,"number of paths must be > 0");
    }

    if(verbose()) {
      printf("number of paths = %li\n",myParameters.nPaths);
    }


    // ************************************
    // get seeds

    getAssignmentULongs(yourElement,"seed",1,2,myULongList);

    list<unsigned long>::const_iterator pULong = myULongList.begin();
    myParameters.seed[0] = *pULong;
    pULong++;
    myParameters.seed[1] = *pULong;

    if((myParameters.seed[0] < 0)|(myParameters.seed[1] < 0)) {
      throw xmdsException(yourElement,"seeds must be positive");
    }

    if(myParameters.seed[0]==myParameters.seed[1]) {
      throw xmdsException(yourElement,"seeds are indentical");
    }

    if(verbose()) {
      printf("seeds are %li and %li\n",myParameters.seed[0],myParameters.seed[1]);
    }

    // ************************************
    // get noises

    getAssignmentULongs(yourElement,"noises",1,1,myULongList);

    myParameters.nNoises = 2*((*myULongList.begin()+1)/2);

    if(myParameters.nNoises < 2) {
      throw xmdsException(yourElement,"number of noises must be > 0");
    }

    if(verbose()) {
      printf("number of noises = %li\n",myParameters.nNoises);
    }

    // find out what kind of noise is specified, if any
    getAttributeStrings(yourElement, "noises", "kind", NOT_REQD, myParameters.noiseKind);
    // the kind has several values, current possible values are:
    // (this list may be augmented in time)
    // "gaussian"
    // "gaussFast"
    // "poissonian"
    // "uniform"
    if (debugFlag) {
      printf("After getAttributeStrings(), noiseKind = %s\n",myParameters.noiseKind.c_str()); 
    }
     
    if (myParameters.noiseKind == "") {
      if (verbose()) {
        printf("Setting the default noise kind (i.e. \"gaussian\")\n");
      }
      myParameters.noiseKind = "gaussian";  // this is the default option
    }
    else if (myParameters.noiseKind == "gaussian") {
      if (verbose()) {
        printf("Noise kind set to \"gaussian\"\n");
      }
    }
    else if (myParameters.noiseKind == "gaussFast") {
      if (verbose()) {
        printf("Noise kind set to \"gaussFast\"\n");
      }
    }
    else if (myParameters.noiseKind == "poissonian") {
      if (verbose()) {
        printf("Noise kind set to \"poissonian\"\n");
      }
      // if we have poissonian specified, we must have a (nonzero) mean for the
      // distribution specified as well
      XMLString noiseMeanString;
      getAttributeStrings(yourElement, "noises", "mean", REQD, noiseMeanString);
      if (debugFlag) {
	printf("After getAttributeStrings(), noiseMean = %s\n",noiseMeanString.c_str()); 
      }      
      // now convert the string to its relevant floating point value
      myParameters.noiseMean = atof(noiseMeanString.c_str());
      if (verbose()) {
	printf("Mean for Poissonian noise set to %g\n",myParameters.noiseMean);
      }
    }
    else if (myParameters.noiseKind == "uniform") {
      if (verbose()) {
        printf("Noise kind set to \"uniform\"\n");
      }
    }
    else {
      throw xmdsException(yourElement, 
                          "You must specify a correct noise kind.  Current choices are:\n gaussian\n gaussFast\n poissonian\n uniform");
    }

  }

  // ************************************
  // find and process children
  // ************************************

  const NodeList* candidateElements;

  // ************************************
  // find argv element and process

  candidateElements = yourElement->getElementsByTagName("argv",0);

  if(candidateElements->length()>1) {
    throw xmdsException(yourElement,"Multiple <argv> elements!");
  }
  else if(candidateElements->length()==1) {
    myArgv = createxmdsArgv();
    const Element* yourElement = dynamic_cast<const Element*> (candidateElements->item(0));
    myArgv->processElement(yourElement);
  }

  // ************************************
  // find globals element and process if present

  candidateElements = yourElement->getElementsByTagName("globals",0);

  if(candidateElements->length()>1) {
    throw xmdsException(yourElement,"Multiple <globals> elements!");
  }
  else if(candidateElements->length()==1) {
    xmdsElement* newGlobals = createxmdsGlobals();
    const Element* yourElement = dynamic_cast<const Element*> (candidateElements->item(0));
    newGlobals->processElement(yourElement);
  }

  // ************************************
  // find field element and process

  candidateElements = yourElement->getElementsByTagName("field",0);

  if(candidateElements->length()>1) {
    throw xmdsException(yourElement,"Multiple <field> elements!");
  }
  else if(candidateElements->length()==1) {
    myField = createxmdsField();
    const Element* yourElement = dynamic_cast<const Element*> (candidateElements->item(0));
    myField->processElement(yourElement);
  }
  else
    throw xmdsException(yourElement,
                        "Cannot find <field> element");
                                
  // ************************************
  // find output element and process

  candidateElements = yourElement->getElementsByTagName("output",0);

  if(candidateElements->length()>1) {
    throw xmdsException(yourElement,"Multiple <output> elements!");
  }
  else if(candidateElements->length()==1) {
    // now grab the format attribute if available
    XMLString outputFormat = "ascii"; // set it to the default while declaring it
    const Node *thisElement = candidateElements->item(0);
    getAttributeStrings(yourElement, "output", "format", NOT_REQD, outputFormat);
    if (outputFormat == "") {
      outputFormat = "ascii";  // the default option
    }
    else {
      if (outputFormat == "ascii") {
        myParameters.binaryOutput = false;
        if (verbose()) {
          printf("Output format set to 'ascii' in output element\n");
        }
      }
      else if (outputFormat == "binary") {
        myParameters.binaryOutput = true;
        if (verbose()) {
          printf("Output format set to 'binary' in output element\n");
        }
      }
      else {
        throw xmdsException(thisElement, "Output format should be either 'ascii' or 'binary'");
      }
    }

    // now grab the precision attribute if available
    XMLString outputPrecision = "double";  // set it to the default while declaring it
    getAttributeStrings(yourElement, "output", "precision", NOT_REQD, outputPrecision);
    if (outputPrecision == "") {
      outputPrecision = "double";  // the default option
    }
    else {
      if (outputPrecision == "double") {
        myParameters.useDouble = true;
        if (verbose()) {
          printf("Output precision set to 'double' in output element\n");
        }
      }
      else if (outputPrecision == "single") {
        myParameters.useDouble = false;
        if (verbose()) {
          printf("Output precision set to 'single' in output element\n");
        }
      }
      else {
        throw xmdsException(thisElement, "Output precision should be either 'double' or 'single'");
      }
    }

    // now process the output element
    myOutput = createxmdsOutput();
    const Element* yourElement = dynamic_cast<const Element*> (candidateElements->item(0));
    myOutput->processElement(yourElement);
  }
  else
    throw xmdsException(yourElement,"Cannot find <output> element");

  // ************************************
  // find sequence element and process

  candidateElements = yourElement->getElementsByTagName("sequence",0);

  if(candidateElements->length()>1) {
    throw xmdsException(yourElement,"Multiple <sequence> elements at top level");
  }
  else if(candidateElements->length()==1) {
    mySequence = createxmdsSequence();
    const Element* yourElement = dynamic_cast<const Element*> (candidateElements->item(0));
    mySequence->processElement(yourElement);
  }
  else
    throw xmdsException(yourElement,"Cannot find <sequence> element");

  myField->outputSampleCount();
  
  mySequence->outputSampleCount();

  myOutput->finaliseGeometry();

  // ************************************
  // Allow non-stochastic mpi?

  if(myParameters.usempi) {
      if((myField->geometry()->nDims()<2)&(!myParameters.stochastic)) {
        myParameters.usempi=0;
        printf("MPI disabled.  MPI cannot be used for deterministic problems with transverse dimensions < 2\n");
      }
      if((myField->geometry()->nDims()>1)&(!myParameters.stochastic)){
        printf("MPI used for a deterministic problem. Ain't life great? \n");
      }
    }
  
};

// ******************************************************************************
const xmdsSimulation::simulationParametersStruct*
xmdsSimulation::parameters() const {
  return &myParameters;
}

// ******************************************************************************
xmdsField* xmdsSimulation::field() const {
  return myField;
};

// ******************************************************************************
xmdsSimulation::argvStruct* 
xmdsSimulation::argStruct() const {
  return &myArgStruct;
};

// ******************************************************************************
const xmdsOutput* xmdsSimulation::output() const {
  return myOutput;
};

// ******************************************************************************
const xmdsSequence* xmdsSimulation::sequence() const {
  return mySequence;
};

// ******************************************************************************
unsigned long xmdsSimulation::nextSegmentNumber() const {
  if(debugFlag) {
    printf("xmdsSimulation::nextSegmentNumber\n");
  }

  myCurrentSegmentNumber++;

  return myCurrentSegmentNumber-1;
};

// ******************************************************************************
void xmdsSimulation::makeCode(
                              const unsigned long& inFileSplitPoint) const {
  if(debugFlag) {
    printf("xmdsSimulation::makeCode\n");
  }

  if(myOutput != 0) {
    myOutput->setInFileSplitPoint(inFileSplitPoint);
  }

  XMLString fileName = myParameters.simulationName;
  fileName += ".cc";

  FILE* outfile = fopen(fileName.c_str(),"w");

  if(outfile==0) {
    sprintf(errorMessage(),"Unable to open file '%s' for write access",fileName.c_str());
    throw xmdsException(errorMessage());  
  }
  if(verbose()) {
    printf("Beginning to write code ...\n");
  }

  try {
    writeIncludes(outfile);
    writeDefines(outfile);
    writeGlobals(outfile);
    writePrototypes(outfile);
    writeRoutines(outfile);
  }
  catch (xmdsException xmdsExceptionErr) {
    printf("Error: simulation failed to write output code\n");
    printf("due to the following xmdsException:\n");
    printf("%s",xmdsExceptionErr.getError());
    fclose(outfile);
    throw xmdsException("Internal Error"); 
  }
  fclose(outfile);
};

// ******************************************************************************
// ******************************************************************************
//                              xmdsSimulation private
// ******************************************************************************
// ******************************************************************************

// ******************************************************************************
void xmdsSimulation::writeIncludes(
                                   FILE *const outfile) const {
  if(debugFlag) {
    printf("xmdsSimulation::writeIncludes\n");
  }