/*
Copyright (C) 2000-2004
Code contributed by Greg Collecutt, Joseph Hope and Paul Cochrane
This file is part of xmds.
This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License
as published by the Free Software Foundation; either version 2
of the License, or (at your option) any later version.
This program is distributed in the hope that it will be useful,
but WITHOUT ANY WARRANTY; without even the implied warranty of
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
GNU General Public License for more details.
You should have received a copy of the GNU General Public License
along with this program; if not, write to the Free Software
Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
*/
/*
$Id: xmdssequence.cc,v 1.43 2005/10/27 00:58:39 joehope Exp $
*/
/*! @file xmdssequence.cc
@brief Sequence element parsing classes and methods
More detailed explanation...
*/
#include<xmlbasics.h>
#include<dom3.h>
#include<xmdsutils.h>
#include<xmdsclasses.h>
// ******************************************************************************
// *****************************************************************************
// xmdsSequence public
// *****************************************************************************
// ******************************************************************************
extern bool debugFlag;
long nxmdsSequences=0; //!< Number of xmds sequence objects
// ******************************************************************************
xmdsSequence::xmdsSequence(
const xmdsSimulation *const yourSimulation,
const bool& yourVerboseMode) :
xmdsSegment(yourSimulation,yourVerboseMode) {
if(debugFlag) {
nxmdsSequences++;
printf("xmdsSequence::xmdsSequence\n");
printf("nxmdsSequences=%li\n",nxmdsSequences);
}
nCycles=1;
};
// ******************************************************************************
xmdsSequence::~xmdsSequence() {
if(debugFlag) {
nxmdsSequences--;
printf("xmdsSequence::~xmdsSequence\n");
printf("nxmdsSequences=%li\n",nxmdsSequences);
}
};
// ******************************************************************************
void xmdsSequence::processElement(
const Element *const yourElement) {
if(debugFlag) {
printf("xmdsSequence::processElement\n");
}
unsigned long i;
const NodeList* myElements;
const Element* nextElement;
xmdsSegment* newSegment;
if(verbose()) {
printf("Processing sequence element ...\n");
}
// ************************************
// process elements in order
myElements = yourElement->getElementsByTagName("*",0);
i=0;
nextElement = dynamic_cast<const Element*> (myElements->item(0));
if(*nextElement->nodeName() == "cycles") {
if(segmentNumber==0) {
throw xmdsException(nextElement,"Cannot have cycles in top level sequence");
}
if(!nextElement->textContent(0)->asULong(nCycles)) {
throw xmdsException(nextElement,"invalid positive integer format");
}
i++;
}
else
// default number of cycles
nCycles=1;
while(i<myElements->length()) {
nextElement = dynamic_cast<const Element*> (myElements->item(i));
if(*nextElement->nodeName() == "sequence") {
newSegment = createxmdsSequence();
}
else if(*nextElement->nodeName() == "integrate") {
// ************************************
// find 'algorithm' so that we can create the appropriate
// type of xmdsIntegrate class
list<XMLString> tempXMLStringList;
getAssignmentStrings(nextElement,"algorithm",0,1,tempXMLStringList);
if(tempXMLStringList.size()==0) {
// need default algorithm
if(simulation()->parameters()->stochastic) {
printf("Algorithm defaulting to SIEX.\n");
tempXMLStringList.push_back("SIEX");
}
else {
printf("Algorithm defaulting to RK4EX.\n");
tempXMLStringList.push_back("RK4EX");
}
}
if(*tempXMLStringList.begin()=="RK4EX") {
newSegment = createxmdsIntegrateRK4EX();
}
else if(*tempXMLStringList.begin()=="RK4IP") {
newSegment = createxmdsIntegrateRK4IP();
}
else if(*tempXMLStringList.begin()=="ARK45IP") {
newSegment = createxmdsIntegrateARK45IP();
}
else if(*tempXMLStringList.begin()=="ARK45EX") {
newSegment = createxmdsIntegrateARK45EX();
}
else if(*tempXMLStringList.begin()=="SIEX") {
newSegment = createxmdsIntegrateSIEX();
}
else if(*tempXMLStringList.begin()=="SIIP") {
newSegment = createxmdsIntegrateSIIP();
}
else {
sprintf(errorMessage(),"Integrate algorithm '%s' unknown.\nCurrent algorithms are SIEX, SIIP, ARK45IP, ARK45EX, RK4EX, and RK4IP.",
tempXMLStringList.begin()->c_str());
throw xmdsException(yourElement,errorMessage());
}
}
else if(*nextElement->nodeName() == "filter") {
newSegment = createxmdsFilter();
}
else {
sprintf(errorMessage(),"segment <%s> unknown",nextElement->nodeName()->c_str());
throw xmdsException(yourElement,errorMessage());
}
newSegment->processElement(nextElement);
i++;
}
};
// ******************************************************************************
void xmdsSequence::outputSampleCount() const {
if(debugFlag) {
printf("xmdsSequence::outputSampleCount\n");
}
for(unsigned long i=0;i<nCycles;i++) {
for(list<const xmdsSegment*>:: const_iterator ppxmdsSegment
= mySegmentList.begin(); ppxmdsSegment !=
mySegmentList.end(); ppxmdsSegment++) {
(*ppxmdsSegment)->outputSampleCount();
}
}
};
// ******************************************************************************
// ******************************************************************************
// xmdsSequence private
// ******************************************************************************
// ******************************************************************************
// ******************************************************************************
void xmdsSequence::writeDefines(
FILE *const outfile) const {
if(debugFlag) {
printf("xmdsSequence::writeDefines\n");
}
if(verbose()) {
printf("Writing sequence defines ...\n");
}
if(segmentNumber==0) {
fprintf(outfile,
"// ********************************************************\n"
"// segment defines\n"
"\n");
}
xmdsElement::writeDefines(outfile);
fprintf(outfile,"\n");
};
// ******************************************************************************
void xmdsSequence::writeGlobals(
FILE *const outfile) const {
if(debugFlag) {
printf("xmdsSequence::writeGlobals\n");
}
if(verbose()) {
printf("Writing sequence globals ...\n");
}
if(segmentNumber==0) {
fprintf(outfile,
"// ********************************************************\n"
"// segment globals\n"
"\n");
}
xmdsElement::writeGlobals(outfile);
fprintf(outfile,"\n");
};
// ******************************************************************************
void xmdsSequence::writePrototypes(
FILE *const outfile) const {
if(debugFlag) {
printf("xmdsSequence::writePrototypes\n");
}
if(verbose()) {
printf("Writing sequence prototypes ...\n");
}
if(segmentNumber==0) {
fprintf(outfile,
"// ********************************************************\n"
"// segment prototypes\n"
"\n");
if(!simulation()->parameters()->usempi)
fprintf(outfile,"void _segment0(); // top level sequence\n");
else if(simulation()->parameters()->mpiMethod=="Scheduling")
fprintf(outfile,"void _segment0(const int); /*top level sequence*/\n");
else
fprintf(outfile,"void _segment0(const int& size); // top level sequence\n");
fprintf(outfile,"\n");
}
else
fprintf(outfile,"void _segment%li(unsigned long cycle); // sequence\n",segmentNumber);
xmdsElement::writePrototypes(outfile);
fprintf(outfile,"\n");
};
// ******************************************************************************
void xmdsSequence::writeRoutines(
FILE *const outfile) const {
if(debugFlag) {
printf("xmdsSequence::writeRoutines\n");
}
if(segmentNumber==0) {
// I am the main routine
if(verbose()) {
printf("Writing main routine ...\n");
}
fprintf(outfile,
"// ********************************************************\n"
"// main routine\n"
"// ********************************************************\n"
"\n"
"// *************************\n");
if (verbose()) {
printf("argStruct()->nameList.size() = %i\n",simulation()->argStruct()->nameList.size());
printf("argStruct()->typeList.size() = %i\n",simulation()->argStruct()->typeList.size());
printf("argStruct()->defaultValueList.size() = %i\n",simulation()->argStruct()->defaultValueList.size());
printf("argStruct()->nameList.size() = %i\n",simulation()->argStruct()->shortOptionList.size());
printf("argStruct()->typeConversionList.size() = %i\n",simulation()->argStruct()->typeConversionList.size());
list<string>::const_iterator inameList = simulation()->argStruct()->nameList.begin();
list<string>::const_iterator itypeList = simulation()->argStruct()->typeList.begin();
list<string>::const_iterator idefaultValueList = simulation()->argStruct()->defaultValueList.begin();
list<string>::const_iterator ishortOptionList = simulation()->argStruct()->shortOptionList.begin();
list<string>::const_iterator itypeConversionList = simulation()->argStruct()->typeConversionList.begin();
printf("argvStruct()->nameList = \n");
for (unsigned long int i=0; i<simulation()->argStruct()->nameList.size(); i++) {
printf("\t%s\n",inameList->c_str());
inameList++;
}
printf("argvStruct()->typeList = \n");
for (unsigned long int i=0; i<simulation()->argStruct()->typeList.size(); i++) {
printf("\t%s\n",itypeList->c_str());
itypeList++;
}
printf("argvStruct()->defaultValueList = \n");
for (unsigned long int i=0; i<simulation()->argStruct()->defaultValueList.size(); i++) {
printf("\t%s\n",idefaultValueList->c_str());
idefaultValueList++;
}
printf("argvStruct()->shortOptionList = \n");
for (unsigned long int i=0; i<simulation()->argStruct()->shortOptionList.size(); i++) {
printf("\t%s\n",ishortOptionList->c_str());
ishortOptionList++;
}
printf("argvStruct()->typeConversionList = \n");
for (unsigned long int i=0; i<simulation()->argStruct()->typeConversionList.size(); i++) {
printf("\t%s\n",itypeConversionList->c_str());
itypeConversionList++;
}
}
// if(simulation()->parameters()->usempi) {
if(simulation()->argStruct()->nameList.size() == 0) {
fprintf(outfile,"int main(int argc, char *argv[]) {\n");
}
else if (simulation()->argStruct()->nameList.size() != 0) {
fprintf(outfile,
"int main(int argc, char **argv) {\n"
"int resp;\n"
"while (1) {\n"
"\tstatic struct option long_options[] = \n"
"\t\t{\n");
list<string>::const_iterator inameList = simulation()->argStruct()->nameList.begin();
// list<string>::const_iterator itypeList = simulation()->argStruct()->typeList.begin();
// list<string>::const_iterator idefaultValueList = simulation()->argStruct()->defaultValueList.begin();
list<string>::const_iterator ishortOptionList = simulation()->argStruct()->shortOptionList.begin();
list<string>::const_iterator itypeConversionList = simulation()->argStruct()->typeConversionList.begin();
fprintf(outfile,"\t\t\t{\"help\", no_argument, 0, 'h'},\n");
for (unsigned long int i=0; i<simulation()->argStruct()->nameList.size(); i++) {
fprintf(outfile,"\t\t\t{\"%s\", required_argument, 0,'%s'},\n",inameList->c_str(),ishortOptionList->c_str());
inameList++;
ishortOptionList++;
}
fprintf(outfile,
"\t\t\t{0, 0, 0, 0}\n"
"\t\t};\n"
"\tint option_index = 0;\n"
"\tresp = getopt_xmds_long(argc, argv, \"");
ishortOptionList = simulation()->argStruct()->shortOptionList.begin();
for (unsigned long int i=0; i<simulation()->argStruct()->shortOptionList.size(); i++) {
fprintf(outfile,"%s:",ishortOptionList->c_str());
ishortOptionList++;
}
// also need to output the short option version of "help" so that
// '-h' works
fprintf(outfile, "h");
fprintf(outfile,
"\", long_options, &option_index);\n"
"\tif (resp == -1) {\n"
"\t\tbreak;\n"
"\t}\n"
"\tswitch (resp) {\n"
"\t\tcase 'h':\n"
"\t\t\tprintf(\"Usage: %s",
simulation()->parameters()->simulationName.c_str());
ishortOptionList = simulation()->argStruct()->shortOptionList.begin();
list<string>::const_iterator itypeList = simulation()->argStruct()->typeList.begin();
for (unsigned long int i=0; i<simulation()->argStruct()->nameList.size(); i++) {
fprintf(outfile," -%s < %s>",ishortOptionList->c_str(),itypeList->c_str());
ishortOptionList++;
itypeList++;
}
fprintf(outfile,
"\\n\\n\");\n"
"\t\t\tprintf(\"Details:\\n\");\n"
"\t\t\tprintf(\"Option\\t\\tType\\t\\tDefault value\\n\");\n");
inameList = simulation()->argStruct()->nameList.begin();
ishortOptionList = simulation()->argStruct()->shortOptionList.begin();
itypeList = simulation()->argStruct()->typeList.begin();
list<string>::const_iterator idefaultValueList = simulation()->argStruct()->defaultValueList.begin();
for (unsigned long int i=0; i<simulation()->argStruct()->nameList.size(); i++) {
fprintf(outfile,"\t\t\tprintf(\"-%s, --%s\\t%s\\t\\t%s\\n\");\n",ishortOptionList->c_str(),inameList->c_str(),itypeList->c_str(),idefaultValueList->c_str());
inameList++;
ishortOptionList++;
itypeList++;
idefaultValueList++;
}
fprintf(outfile,"\t\t\texit(0);\n");
inameList = simulation()->argStruct()->nameList.begin();
ishortOptionList = simulation()->argStruct()->shortOptionList.begin();
itypeConversionList = simulation()->argStruct()->typeConversionList.begin();
for (unsigned long int i=0; i<simulation()->argStruct()->nameList.size(); i++) {
fprintf(outfile,"\t\tcase '%s':\n",ishortOptionList->c_str());
if (*itypeConversionList == "atof") {
fprintf(outfile,"\t\t\t%s = atof(optarg);\n",inameList->c_str());
fprintf(outfile,"\t\t\tbreak;\n");
}
else if (*itypeConversionList == "atoi") {
fprintf(outfile,"\t\t\t%s = atoi(optarg);\n",inameList->c_str());
fprintf(outfile,"\t\t\tbreak;\n");
}
else if (*itypeConversionList == "") {
fprintf(outfile,"\t\t\t%s = optarg;\n",inameList->c_str());
fprintf(outfile,"\t\t\tbreak;\n");
}
inameList++;
ishortOptionList++;
itypeConversionList++;
}
fprintf(outfile,
"\t\tdefault:\n"
"\t\t\texit(10);\n"
"\t\t}\n"
"\t}\n"
"\tif (optind_xmds < argc) {\n"
"\t\tprintf(\"Sorry, I didn't understand some (or all) of your arguments\\n\");\n"
"\t\t\tprintf(\"Usage: %s",simulation()->parameters()->simulationName.c_str());
ishortOptionList = simulation()->argStruct()->shortOptionList.begin();
itypeList = simulation()->argStruct()->typeList.begin();
for (unsigned long int i=0; i<simulation()->argStruct()->nameList.size(); i++) {
fprintf(outfile," -%s < %s>",ishortOptionList->c_str(),itypeList->c_str());
ishortOptionList++;
itypeList++;
}
fprintf(outfile,"\\n\\n\");\n"
"\t\t\tprintf(\"Details:\\n\");\n"
"\t\t\tprintf(\"Option\\t\\tType\\t\\tDefault value\\n\");\n");
inameList = simulation()->argStruct()->nameList.begin();
ishortOptionList = simulation()->argStruct()->shortOptionList.begin();
itypeList = simulation()->argStruct()->typeList.begin();
idefaultValueList = simulation()->argStruct()->defaultValueList.begin();
for (unsigned long int i=0; i<simulation()->argStruct()->nameList.size(); i++) {
fprintf(outfile,"\t\t\tprintf(\"-%s, --%s\\t%s\\t\\t%s\\n\");\n",ishortOptionList->c_str(),inameList->c_str(),itypeList->c_str(),idefaultValueList->c_str());
inameList++;
ishortOptionList++;
itypeList++;
idefaultValueList++;
}
fprintf(outfile,
"\t\t\texit(0);\n"
"\t}\n\n");
}
else {
fprintf(outfile,"int main() {\n");
}
fprintf(outfile,"\n");
if (simulation()->parameters()->binaryOutput) {
fprintf(outfile,
"#ifndef CPU_IS_BIG_ENDIAN\n"
" #ifndef CPU_IS_LITTLE_ENDIAN\n"
" cout << \"Byte ordering isn't defined for this system!\\n\";\n"
" cout << \"Please define CPU_IS_BIG_ENDIAN and CPU_IS_LITTLE_ENDIAN in config.h\\n\";\n"
" exit(255);\n"
" #endif\n"
"#endif\n");
}
if(simulation()->parameters()->usempi) {
fprintf(outfile,
"MPI_Init(&argc, &argv);\n"
"MPI_Comm_size(MPI_COMM_WORLD, &size);\n"
"MPI_Comm_rank(MPI_COMM_WORLD, &rank);\n"
"\n");
}
// if(simulation()->field()->needsFFTWPlans()) {
// It's possible to need FFTW plans in output but not the main field through post-propagation
fprintf(outfile,"/* set up fftw plans */\n");
if(simulation()->parameters()->nThreads>1) {
fprintf(outfile,"int fftw_threads_init(void);\n");
}
fprintf(outfile,"int _fftw_lattice[64];\n");
fprintf(outfile,"\n");
simulation()->field()->writePlanCreationCalls(outfile,1,simulation()->parameters()->useWisdom);
simulation()->output()->writePlanCreationCalls(outfile,0,simulation()->parameters()->useWisdom);
// }
if(simulation()->parameters()->benchmark) {
if(simulation()->parameters()->usempi){
fprintf(outfile,"double startTime = MPI_Wtime();\n");
}
else{
fprintf(outfile,"time_t startTime = time(NULL);\n");
}
}
if((simulation()->parameters()->usempi)&!(simulation()->parameters()->stochastic)) {
fprintf(outfile,"fftwnd_mpi_local_sizes(_main_forward_plan, &local_nx, &local_x_start, &local_ny_after_transpose, &local_y_start_after_transpose, &total_local_size);\n\n");
list<XMLString> vectorNamesList;
const xmdsVector* tempVector;
const char *typeName;
simulation()->field()->vectorNames(vectorNamesList);
for(list<XMLString>::const_iterator pXMLString = vectorNamesList.begin(); pXMLString != vectorNamesList.end(); pXMLString++) {
if(!simulation()->field()->getVector(*pXMLString,tempVector)) {
throw xmdsException("Lost one of the vectors in xmdssequence::writeRoutines. Bad! Bad dog!\n\n");
}
if(tempVector->vectorType()==DOUBLE) {
typeName="double";
}
else {
typeName="complex";
}
fprintf(outfile,"_%s_%s = new %s[total_local_size*_%s_%s_ncomponents];\n",
simulation()->field()->name()->c_str(),tempVector->name()->c_str(),typeName,
simulation()->field()->name()->c_str(),tempVector->name()->c_str());
fprintf(outfile,"_%s_%s_work = new %s[total_local_size*_%s_%s_ncomponents];\n",
simulation()->field()->name()->c_str(),tempVector->name()->c_str(),typeName,
simulation()->field()->name()->c_str(),tempVector->name()->c_str());
}
fprintf(outfile,"\n");
for(unsigned long i=0; i<simulation()->output()->nMomentGroups(); i++) {
if(!simulation()->output()->momentGroup(i)->needscomplexRawVector()) {
typeName="double";
}
else {
typeName="complex";
}
fprintf(outfile,"_mg%li_raw = new %s[_mg%li_size*_mg%li_raw_ncomponents];\n",i,typeName,i,i);
fprintf(outfile,"_mg%li_fullstep = new double[_mg%li_size*_mg%li_fullstep_ncomponents];\n",i,i,i);
if(simulation()->parameters()->errorCheck) {
fprintf(outfile,"_mg%li_halfstep = new double[_mg%li_size*_mg%li_halfstep_ncomponents];\n",i,i,i);
}
}
fprintf(outfile,"\n");
}
simulation()->field()->assignActiveVectorPointers(outfile,"");
simulation()->output()->assignActiveVectorPointers(outfile,"");
if(simulation()->parameters()->stochastic) {
for(unsigned long i=0;i<simulation()->output()->nMomentGroups();i++) {
fprintf(outfile,"_mg%li_fullstep_initialise();\n",i);
fprintf(outfile,"_mg%li_fullstep_sd_initialise();\n",i);
if(simulation()->parameters()->errorCheck) {
fprintf(outfile,"_mg%li_halfstep_initialise();\n",i);
fprintf(outfile,"_mg%li_halfstep_sd_initialise();\n",i);
}
fprintf(outfile,"\n");
}
}
if((simulation()->parameters()->mpiMethod=="Scheduling")&&(simulation()->parameters()->usempi)&&(simulation()->parameters()->stochastic)) {
if(simulation()->parameters()->errorCheck) {
if(simulation()->parameters()->nPaths>1) {
fprintf(outfile,"printf(\"Rank[%%i] beginning full step paths\\n\",rank);\n");
}
else {
fprintf(outfile,"printf(\"Beginning full step integration ...\\n\");\n");
}
fprintf(outfile,"\n");
fprintf(outfile,"_half_step=0;\n");
fprintf(outfile,"\n");
}
else {
if(simulation()->parameters()->nPaths>1) {
fprintf(outfile,"printf(\"Rank[%%i] beginning paths\\n\",rank);\n");
}
else {
fprintf(outfile,"printf(\"Beginning integration ...\\n\");\n");
}
fprintf(outfile,"\n");
}
if(simulation()->parameters()->errorCheck) {
fprintf(outfile,"_gen1[0] = %li+rank;\n",simulation()->parameters()->seed[0]);
fprintf(outfile,"_gen2[0] = %li+rank;\n",simulation()->parameters()->seed[1]);
fprintf(outfile,"_gen1[1] = 0;\n");
fprintf(outfile,"_gen2[1] = 0;\n");
fprintf(outfile,"_gen1[2] = 0;\n");
fprintf(outfile,"_gen2[2] = 0;\n");
fprintf(outfile,"erand48(_gen1);\n");
fprintf(outfile,"erand48(_gen2);\n");
}
else {
fprintf(outfile,"_gen[0] = %li+rank;\n",simulation()->parameters()->seed[0]);
fprintf(outfile,"_gen[1] = 0;\n");
fprintf(outfile,"_gen[2] = 0;\n");
fprintf(outfile,"erand48(_gen);\n");
}
if(simulation()->parameters()->nPaths>1) {
fprintf(outfile,"\ntotal_paths = _n_paths;\n\n");
}
else {
fprintf(outfile,"\ntotal_paths = 1;\n\n");
}
fprintf(outfile,"//Rank 0 becomes master, other processes become slaves\n"
"if (size > 1){\n"
" if (rank==0)\n"
" master();\n"
" else {\n"
" slave((void *) NULL);\n"
" }\n"
"}\n"
"//Run the simulation using the following if only one processor used\n"
"else {\n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile," _half_step = 0;\n");
fprintf(outfile,
" _segment0(total_paths);\n\n");
if(simulation()->parameters()->errorCheck) {
fprintf(outfile," _gen1[0] = %li+rank;\n",simulation()->parameters()->seed[0]);
fprintf(outfile," _gen2[0] = %li+rank;\n",simulation()->parameters()->seed[1]);
fprintf(outfile," _gen1[1] = 0;\n");
fprintf(outfile," _gen2[1] = 0;\n");
fprintf(outfile," _gen1[2] = 0;\n");
fprintf(outfile," _gen2[2] = 0;\n");
fprintf(outfile," erand48(_gen1);\n");
fprintf(outfile," erand48(_gen2);\n\n"
" _half_step = 1;\n");
if(simulation()->parameters()->nPaths>1) {
fprintf(outfile," printf(\"Rank[%%i] beginning half step paths\\n\",rank);\n");
}
else {
fprintf(outfile," printf(\"Beginning half step integration ...\\n\");\n");
}
fprintf(outfile," _segment0(total_paths);\n\n");
}
fprintf(outfile,"goto write_output;\n"
"}\n\n");
}
else {
if(simulation()->parameters()->errorCheck) {
if(simulation()->parameters()->nPaths>1) {
if(simulation()->parameters()->usempi) {
fprintf(outfile,"printf(\"Rank[%%i] beginning full step paths\\n\",rank);\n");
}
else {
fprintf(outfile,"printf(\"Beginning full step paths\\n\");\n");
}
}
else {
if(simulation()->parameters()->usempi&!simulation()->parameters()->stochastic) {
fprintf(outfile,"printf(\"Rank[%%i] beginning full step integration ...\\n\",rank);\n");
}
else {
fprintf(outfile,"printf(\"Beginning full step integration ...\\n\");\n");
}
}
fprintf(outfile,"\n");
fprintf(outfile,"_half_step=0;\n");
fprintf(outfile,"\n");
}
else {
if(simulation()->parameters()->nPaths>1) {
if(simulation()->parameters()->usempi) {
fprintf(outfile,"printf(\"Rank[%%i] beginning paths\\n\",rank);\n");
}
else {
fprintf(outfile,"printf(\"Beginning paths\\n\");\n");
}
}
else {
if(simulation()->parameters()->usempi&!simulation()->parameters()->stochastic) {
fprintf(outfile,"printf(\"Rank[%%i] beginning integration ...\\n\",rank);\n");
}
else {
fprintf(outfile,"printf(\"Beginning integration ...\\n\");\n");
}
}
fprintf(outfile,"\n");
}
if(simulation()->parameters()->usempi) {
fprintf(outfile,"_segment0(size);\n");
}
else {
fprintf(outfile,"_segment0();\n");
}
fprintf(outfile,"\n");
if(simulation()->parameters()->errorCheck) {
if(simulation()->parameters()->nPaths>1) {
if(simulation()->parameters()->usempi) {
fprintf(outfile,"printf(\"Rank[%%i] beginning half step paths\\n\",rank);\n");
}
else {
fprintf(outfile,"printf(\"Beginning half step paths\\n\");\n");
}
}
else {
if(simulation()->parameters()->usempi&!simulation()->parameters()->stochastic) {
fprintf(outfile,"printf(\"Rank[%%i] beginning half step integration ...\\n\",rank);\n");
}
else {
fprintf(outfile,"printf(\"Beginning half step integration ...\\n\");\n");
}
}
fprintf(outfile,"\n");
fprintf(outfile,"_half_step=1;\n");
fprintf(outfile,"\n");
if(simulation()->parameters()->usempi) {
fprintf(outfile,"_segment0(size);\n");
}
else {
fprintf(outfile,"_segment0();\n");
}
fprintf(outfile,"\n");
}
}
if(simulation()->parameters()->usempi&simulation()->parameters()->stochastic) {
fprintf(outfile,"double* _temp_vector;\n");
fprintf(outfile,"unsigned long _length;\n");
fprintf(outfile,"\n");
for(unsigned long i=0;i<simulation()->output()->nMomentGroups();i++) {
fprintf(outfile,"_length = _mg%li_size*_mg%li_fullstep_ncomponents;\n",i,i);
fprintf(outfile,"_temp_vector = new double[_length];\n");
fprintf(outfile,"\n");
fprintf(outfile,"MPI_Reduce(_mg%li_fullstep,_temp_vector,_length,MPI_DOUBLE,MPI_SUM,0,MPI_COMM_WORLD);\n",i);
fprintf(outfile,"\n");
fprintf(outfile,"for(unsigned long _i0=0; _i0<_length; _i0++)\n");
fprintf(outfile," _mg%li_fullstep[_i0]=_temp_vector[_i0];\n",i);
fprintf(outfile,"\n");
if(simulation()->parameters()->nPaths>1) {
fprintf(outfile,"MPI_Reduce(_mg%li_fullstep_sd,_temp_vector,_length,MPI_DOUBLE,MPI_SUM,0,MPI_COMM_WORLD);\n",i);
fprintf(outfile,"\n");
fprintf(outfile,"for(unsigned long _i0=0; _i0<_length; _i0++)\n");
fprintf(outfile," _mg%li_fullstep_sd[_i0]=_temp_vector[_i0];\n",i);
fprintf(outfile,"\n");
}
if(simulation()->parameters()->errorCheck) {
fprintf(outfile,"MPI_Reduce(_mg%li_halfstep,_temp_vector,_length,MPI_DOUBLE,MPI_SUM,0,MPI_COMM_WORLD);\n",i);
fprintf(outfile,"\n");
fprintf(outfile,"for(unsigned long _i0=0; _i0<_length; _i0++)\n");
fprintf(outfile," _mg%li_halfstep[_i0]=_temp_vector[_i0];\n",i);
fprintf(outfile,"\n");
if(simulation()->parameters()->nPaths>1) {
fprintf(outfile,"MPI_Reduce(_mg%li_halfstep_sd,_temp_vector,_length,MPI_DOUBLE,MPI_SUM,0,MPI_COMM_WORLD);\n",i);
fprintf(outfile,"\n");
fprintf(outfile,"for(unsigned long _i0=0; _i0<_length; _i0++)\n");
fprintf(outfile," _mg%li_halfstep_sd[_i0]=_temp_vector[_i0];\n",i);
fprintf(outfile,"\n");
}
}
fprintf(outfile,"delete[] _temp_vector;\n");
fprintf(outfile,"\n");
}
if(simulation()->parameters()->mpiMethod=="Scheduling") {
fprintf(outfile," /* Adaptive Scheduler single processor break point */\n"
"write_output:\n\n");
}
fprintf(outfile,"if(rank==0)\n");
fprintf(outfile," _write_output();\n");
fprintf(outfile,"\n");
}
else if(simulation()->parameters()->usempi){
fprintf(outfile,"if(rank==0)\n");
fprintf(outfile," _write_output();\n");
fprintf(outfile,"\n");
}
else {
fprintf(outfile,"_write_output();\n");
fprintf(outfile,"\n");
}
if(simulation()->parameters()->benchmark) {
if(simulation()->parameters()->usempi){
fprintf(outfile,"double endTime = MPI_Wtime();\n");
fprintf(outfile,"if(rank==0)\n ");
}
else{
fprintf(outfile,"time_t endTime = time(NULL);\n");
}
fprintf(outfile,"printf(\"Time elapsed for simulation is: %%5.5f seconds\\n\",endTime-startTime);\n");
}
simulation()->field()->writePlanDeletionCalls(outfile);
simulation()->output()->writePlanDeletionCalls(outfile);
if(simulation()->parameters()->usempi) {
fprintf(outfile,"MPI_Finalize();\n");
fprintf(outfile,"\n");
}
if (simulation()->parameters()->bing) {
// this should be the last statement in the main() routine
fprintf(outfile,"printf(\"\\a\");\n");
}
fprintf(outfile,"return 0 ;\n");
fprintf(outfile,"}\n");
fprintf(outfile,"\n");
if((simulation()->parameters()->mpiMethod=="Scheduling")&&(simulation()->parameters()->usempi)&&(simulation()->parameters()->stochastic)) {
fprintf(outfile,
"int master(){ \n"
" \n"
" int outstanding = 0; /*number of slaves that are still processing tasks*/ \n"
" int schedule[size]; /*Batch size scheduled for each slave [resetted every iteration]*/ \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" int sendtag; /*MPI Tag value, fullstep -> 1, halfstep -> 2*/ \n");
else
fprintf(outfile,
" int sendtag = 1; \n");
fprintf(outfile,
" \n"
" double timing[size]; /*Timing function to determine computation to communication ratio*/ \n"
" int partitions[size]; /*Batch size scheduled for each slave [resetted after completion]*/ \n"
" \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" double *active_throughput; /*Active pointers, used for clarity*/ \n"
" double *active_commave; \n"
" int *active_assigned; \n"
" \n");
fprintf(outfile,
" int i, j, k; /*indexes, and size of MPI_COMM_WORLD*/ \n"
" double bufs[size]; /*MPI Input buffer*/ \n"
" int idle = size-1; /*Number of idle processors not including master (size - outstanding -1)*/ \n"
" \n"
" MPI_Status stats[size]; /*MPI Structures*/ \n"
" MPI_Request reqs[size]; \n"
" int indices[size]; \n"
" int ndone; \n"
" \n"
" double schedtime=0.0; /*time spent deciding and dispatching schedules*/ \n"
" double commtime=0.0; /*index for communication latency*/ \n"
" double totaltime=0.0; /*index for seconds per schedule*/ \n"
" double totalcommtime=0.0; /*total communication latency*/ \n"
" double paratime=0.0; /*total parallel walltime for slaves excluding mslave*/ \n"
" \n"
" /************* Scheduling Parameters **************/ \n"
" double calpha = 0.1; /*weighting for communication average*/ \n"
" double talpha = 0.1; /*weighting for throughput average*/ \n"
" double rcalpha = 1- calpha; \n"
" double rtalpha = 1- talpha; \n"
" \n"
" double epsilon = 0.005; /*maximum tolerated communication overhead*/ \n"
" double lower = 2.0; /*minimum tolerated resolution in seconds*/ \n"
" double upper = 10.0; /*maximum tolerated resolution seconds*/ \n"
" /***************************************************/ \n"
" \n"
" pthread_t helper; /*handle for slave thread*/ \n"
" \n"
" double tp1, tp2; \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" int master_half_step; /*This indicates whether the master is sending*/ \n"
" /*out half step paths or not - redundant to the index k*/ \n"
" \n");
fprintf(outfile,
" //Allocate memories for global data structures \n"
" slave_stat = new int[size]; //(int *)malloc(sizeof(int)*size); \n"
" fschedcount = new int[size]; //(int *)malloc(sizeof(int)*size); \n"
" fthroughput = new double[size]; \n"
" fcommave = new double[size]; \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" hschedcount = new int[size]; \n"
" hthroughput = new double[size]; \n"
" hcommave = new double[size]; \n");
fprintf(outfile,
" \n"
" //Initialise slave status arrays \n"
" for (i=0;i<size;i++){ \n"
" slave_stat[i]=0; \n"
" partitions[i]=0; \n"
" \n"
" fschedcount[i]=0; \n"
" fcommave[i]=0.0; \n"
" fthroughput[i]=0.0; \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" \n"
" hschedcount[i]=0; \n"
" hcommave[i]=0.0; \n"
" hthroughput[i]=0.0; \n");
fprintf(outfile,
" } \n"
" \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" //Initially set to full step paths \n"
" _half_step = 0; \n"
" master_half_step = _half_step; \n"
" sendtag = 1; \n");
fprintf(outfile,
" \n"
" //Pthreads are always busy doing something \n"
" slave_stat[0] = 1; \n"
" \n"
" //Initialise mutual exclusion mechanism \n"
" pthread_mutex_init(&tasklock, NULL); \n"
" pthread_mutex_init(&finlock, NULL); \n"
" pthread_mutex_lock(&finlock); \n"
" \n"
" //Create a thread to act as a slave \n"
" if (pthread_create(&helper, NULL, mslave, NULL)!=0) \n"
" printf(\"Thread creation failed\\n\"); \n"
" \n"
" //Listen for messages from all slaves, including the master for some reason \n"
" for (i=0;i<size;i++) \n"
" MPI_Irecv(&bufs[i], 1, MPI_DOUBLE, i, MPI_ANY_TAG,MPI_COMM_WORLD,&reqs[i]); \n"
" \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" //Loop until all paths are finished, and all results recieved \n"
" // -repeat \n"
" // -Test for messages from slaves if \n"
" // -Determine path schedule for idle slaves \n"
" // -Send schedule to slaves \n"
" \n"
" /* The first iteration of this for loop allocates full paths \n"
" The second iteration allocates half paths - master_half_step is a redundant \n"
" indicator for this too. We could use k interchangeably with it.*/ \n"
" \n"
" for (k=0;k<2;k++){ \n"
" \n"
" if (master_half_step == 0){ \n"
" active_throughput = fthroughput; \n"
" active_commave = fcommave; \n"
" active_assigned = &fullpaths_assigned; \n"
" } \n"
" else if (master_half_step == 1){ \n"
" active_throughput = hthroughput; \n"
" active_commave = hcommave; \n"
" active_assigned = &halfpaths_assigned; \n"
" } \n"
" \n"
" \n"
" /**************LISTEN FOR \"TASKS COMPLETED\" MESSAGES FROM SLAVES*******************/ \n"
" \n"
" /*Loop until all fullpaths or, if during the second iteration, halfpaths are assigned \n"
" the outstanding > 0 condition makes the master poll for \"complete\" messages when \n"
" it has actually allocated all its full paths and half paths already.*/ \n"
" \n"
" while ((master_half_step == 0 && fullpaths_assigned < total_paths)|| \n"
" (master_half_step == 1 && (halfpaths_assigned < total_paths || outstanding > 0))){ \n"
" \n"
" /*Test whether there is an incoming message from slaves \n"
" the idle == 0 condition is for the situation when SOME slaves \n"
" have not finished their full paths whilst the master is ready to \n"
" allocate half paths - this allows the master to allocate half paths \n"
" in the meantime, instead of waiting for completion of all full paths*/ \n"
" \n"
" if (outstanding > 0 && (idle == 0 || (fullpaths_assigned + halfpaths_assigned >= 2*total_paths))){ \n");
else fprintf(outfile,
" //Loop until all tasks are finished, and all results recieved \n"
" // -Check for requests from slaves for more tasks \n"
" // -Cyclically allocate a task to idle slaves to thread \n"
" \n"
" while (fullpaths_assigned < total_paths || outstanding > 0){ \n"
" \n"
" if (outstanding > 0){ \n");
fprintf(outfile,
" //Wait for messages from slaves \n"
" \n"
" Testsome: \n"
" MPI_Testsome(size,reqs, &ndone, indices, stats); \n"
" \n"
" if (ndone == 0){ \n"
" /*If no messages, sleep for x secs before checking again*/ \n"
" usleep(10000); \n"
" goto Testsome; \n"
" } \n"
" \n"
" \n"
" for (i=0; i<ndone; i++){ \n"
" //Deal with incoming messages \n"
" j = indices[i]; \n"
" \n"
" //Dynamically determine bandwidth and throughput \n"
" totaltime = MPI_Wtime() - timing[j]; \n"
" commtime = totaltime - bufs[j]; \n"
" \n"
" //Calculate average communication time and average throughput \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" //MPI_TAG == 1 means fullpaths completed \n"
" //MPI_TAG == 2 means halfpaths completed \n"
" \n"
" if (stats[i].MPI_TAG == 1){ \n");
fprintf(outfile,
" if (fcommave[j] == 0.0){ \n"
" fcommave[j] = commtime; \n"
" fthroughput[j] = partitions[j]/totaltime; \n"
" } \n"
" else { \n"
" fcommave[j] = fcommave[j]*rcalpha + (commtime * calpha); \n"
" fthroughput[j] = fthroughput[j]*rtalpha + (partitions[j]/totaltime)*talpha; \n"
" } \n"
" \n"
"#ifdef DEBUG_MPI_Scheduling \n"
" fschedcount[j]+=partitions[j]; \n"
" printf(\"%%d: Rank %%d, Completed %%d tasks in %%5.5f seconds ::\", stats[i].MPI_TAG,j,partitions[j],commtime+bufs[j]);\n"
" printf(\" T %%5.5f : C %%5.5f : R %%5.5f %%%%\\n\",fthroughput[j],fcommave[j],commtime/(bufs[j]+commtime)*100.0); \n"
"#endif \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" } \n"
" else if (stats[i].MPI_TAG == 2){ \n"
" if (hcommave[j] == 0.0){ \n"
" hcommave[j] = commtime; \n"
" hthroughput[j] = partitions[j]/totaltime; \n"
" } \n"
" else { \n"
" hcommave[j] = hcommave[j]*rcalpha + (commtime * calpha); \n"
" hthroughput[j] = hthroughput[j]*rtalpha + (partitions[j]/totaltime)*talpha; \n"
" } \n"
"#ifdef DEBUG_MPI_Scheduling \n"
" hschedcount[j]+=partitions[j]; \n"
" printf(\"%%d: Rank %%d, Completed %%d tasks in %%5.5f seconds ::\", stats[i].MPI_TAG,j,partitions[j],commtime+bufs[j]);\n"
" printf(\" T %%5.5f : C %%5.5f : R %%5.5f %%%%\\n\",hthroughput[j],hcommave[j],commtime/(bufs[j]+commtime)*100.0); \n"
"#endif \n"
" } \n");
fprintf(outfile,
" \n"
" totalcommtime += commtime; \n"
" paratime += totaltime; \n"
" slave_stat[j] = 0; \n"
" idle++; \n"
" outstanding--; \n"
" MPI_Irecv(&bufs[j], 1, MPI_DOUBLE, j, MPI_ANY_TAG,MPI_COMM_WORLD,&reqs[j]); \n"
" } \n"
" } \n"
" \n"
" //Continue to listen for \"completed\" messages from slaves, if no more tasks need to be assigned \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" if (halfpaths_assigned + fullpaths_assigned >= 2*total_paths) \n");
else
fprintf(outfile,
" if (fullpaths_assigned >= total_paths) \n");
fprintf(outfile,
" continue; \n"
" \n"
" /********************SCHEDULE MORE TASKS FOR IDLE SLAVES***********************/ \n"
" \n"
" for (i=0;i<size;i++) \n"
" schedule[i]=0; \n"
" \n"
" tp1 = MPI_Wtime(); \n"
" \n"
" pthread_mutex_lock(&tasklock); \n"
" \n"
" //allocate tasks to free processors \n"
" //other scheduling algorithms may be inserted here if needed \n"
" //this one is a FIFO algorithm \n"
" //scheduling must be mutually exclusive as the slave thread \n"
" //also modifies the global variables below for self-scheduling \n"
" \n"
" for (i=1;i<size;i++){ \n"
" \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" if (*active_assigned >= total_paths) \n");
else
fprintf(outfile,
" if (fullpaths_assigned >= total_paths) \n");
fprintf(outfile,
" break; \n"
" \n"
" //only allocate more tasks to slaves that are idle \n"
" if (slave_stat[i] == 0){ \n"
" //printf(\"allocating tasks to rank %%d\\n\",i); \n"
" slave_stat[i]=1; \n"
" idle--; \n"
" \n"
" //Determine partition size based on slave throughput and \n"
" //communication overhead. Preferable estimated computing times \n"
" //for each schedule is high enough to reduce comm overhead \n"
" //and between upper and lower \n"
" \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" partitions[i] = (int) (MAX((active_commave[i]*active_throughput[i])/epsilon, active_throughput[i]*lower)); \n"
" partitions[i] = (int) (MIN(partitions[i], active_throughput[i]*upper)); \n"
" partitions[i] = (int) (MAX(partitions[i], 1)); \n"
" \n"
" if (*active_assigned + partitions[i] <= total_paths){ \n"
" schedule[i] = partitions[i]; \n"
" *active_assigned += partitions[i]; \n"
" } \n"
" else { \n"
" schedule[i] = total_paths-*active_assigned; \n"
" partitions[i] = schedule[i]; \n"
" *active_assigned += schedule[i]; \n");
else
fprintf(outfile,
" partitions[i] = (int) (MAX((fcommave[i]*fthroughput[i])/epsilon, fthroughput[i]*lower)); \n"
" partitions[i] = (int) (MIN(partitions[i],fthroughput[i]*upper)); \n"
" partitions[i] = (int) (MAX(partitions[i], 1)); \n"
" \n"
" if (fullpaths_assigned + partitions[i] <= total_paths){ \n"
" schedule[i] = partitions[i]; \n"
" fullpaths_assigned += partitions[i]; \n"
" } \n"
" else { \n"
" schedule[i] = total_paths-fullpaths_assigned; \n"
" partitions[i] = schedule[i]; \n"
" fullpaths_assigned += schedule[i]; \n");
fprintf(outfile,
" } \n"
" } \n"
" } \n"
" \n"
" pthread_mutex_unlock(&tasklock); \n"
" \n"
" /**************************SEND SCHEDULE TO SLAVE(S)********************/ \n"
" //send instructions via MPI to slaves that have been shcheduled \n"
" for (i=1;i<size;i++){ \n"
" if (schedule[i]>0){ \n"
" timing[i] = MPI_Wtime(); \n"
" MPI_Send(&schedule[i],1,MPI_INT,i,sendtag,MPI_COMM_WORLD); \n"
" outstanding++; \n"
" } \n"
" } \n"
" tp2 = MPI_Wtime() - tp1; \n"
" schedtime += tp2; \n"
" } \n"
" \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" //Initialise variables to do half steps \n"
" master_half_step = 1; \n"
" sendtag = 2; \n"
" } \n"
" \n");
fprintf(outfile,
" //Block until the thread slave has completed processing \n"
" pthread_mutex_lock(&finlock); \n"
" \n"
" //Tell slave processes to Reduce then exit \n"
" for (i=0;i<size;i++){ \n"
" MPI_Send(NULL, 0, MPI_INT,i,0,MPI_COMM_WORLD); \n"
" } \n"
" \n"
" //Kill slave thread \n"
" pthread_cancel(helper); \n"
" \n"
"#ifdef DEBUG_MPI_Scheduling \n"
" printf(\"\\n**************************************\"); \n"
" printf(\"\\nFullpaths returned by each process:\\n\"); \n"
" print_array2(fschedcount, size); \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" printf(\"\\nHalfpaths returned by each process:\\n\"); \n"
" print_array2(hschedcount, size); \n");
fprintf(outfile,
" \n"
" printf(\"\\nEstimated fullpath throughput and commtime\\n\"); \n"
" print_array(fthroughput, size);printf(\"\\n\"); \n"
" print_array(fcommave, size); \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" printf(\"\\nEstimated halfpath throughput and commtime\\n\"); \n"
" print_array(hthroughput, size);printf(\"\\n\"); \n"
" print_array(hcommave, size); \n");
fprintf(outfile,
" \n"
" printf(\"\\n\"); \n"
" printf(\"Scheduling time = %%5.5f sec\\n\",schedtime); \n"
" printf(\"Communication time = %%5.5f sec\\n\", totalcommtime); \n"
" printf(\"Total parallel time (excl thread) = %%5.5f sec\\n\",paratime); \n"
" printf(\"**************************************\\n\\n\"); \n"
" fflush(stdout); \n"
"#endif \n"
" \n"
" /*Free memory*/ \n"
" delete[] slave_stat; \n"
" delete[] fschedcount; \n"
" delete[] fthroughput; \n"
" delete[] fcommave; \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" delete[] hschedcount; \n"
" delete[] hthroughput; \n"
" delete[] hcommave; \n"
" \n");
fprintf(outfile,
" return 0; \n"
"} \n"
" \n"
" \n"
" \n"
"void *mslave(void *ptr){ \n"
" //This is the threaded slave, which is spawned by the master to do tasks in parallel \n"
" //it self schedules itself more task batches. It takes in one argument, which is \n"
" //the size of the solution array. \n"
" \n"
" int thr_batchsize; \n"
" double thr_throughput=0.0; \n"
" double thr_time_per_batch=2.0; \n"
" double thr_talpha=0.1; \n"
" \n"
" int i; \n"
" struct timeval *tp1,*tp2; \n"
" double tp3; \n"
" \n"
" \n"
" //Initialise timing structure \n"
" tp1 = new struct timeval; //(struct timeval *)malloc(sizeof(struct timeval)); \n"
" tp2 = new struct timeval; //(struct timeval *)malloc(sizeof(struct timeval)); \n"
" \n"
" thr_batchsize = 1; \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" _half_step = 0; \n"
" \n"
" while(fullpaths_assigned + halfpaths_assigned < 2*total_paths){ \n");
else fprintf(outfile,
" while(fullpaths_assigned < total_paths){ \n");
fprintf(outfile,
" //Self schedule more tasks to process \n"
" gettimeofday(tp1,NULL); \n"
" \n"
" /********************SCHEDULE MORE TASKS*************************/ \n"
" pthread_mutex_lock(&tasklock); \n"
" \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" if (halfpaths_assigned + fullpaths_assigned >= 2*total_paths) \n");
else fprintf(outfile,
" if (fullpaths_assigned >= total_paths) \n");
fprintf(outfile,
" break; \n"
" \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" //Set to half step (once), if all full step paths are assigned \n"
" if (fullpaths_assigned >= total_paths && _half_step == 0){ \n"
" printf(\"Rank[%%i] beginning half step paths\\n\",rank); \n"
" fthroughput[0]=thr_throughput; \n"
" _half_step = 1; \n"
" thr_batchsize = 1; \n"
" thr_throughput = 0.0; \n"
" \n"
" _gen1[0] = 346+rank; \n"
" _gen2[0] = 983+rank; \n"
" _gen1[1] = 0; \n"
" _gen2[1] = 0; \n"
" _gen1[2] = 0; \n"
" _gen2[2] = 0; \n"
" erand48(_gen1); \n"
" erand48(_gen2); \n"
" } \n"
" \n"
" if (_half_step==0){ \n");
fprintf(outfile,
" if (fullpaths_assigned + thr_batchsize <= total_paths) \n"
" fullpaths_assigned += thr_batchsize; \n"
" else { \n"
" thr_batchsize = (total_paths - fullpaths_assigned); \n"
" fullpaths_assigned += thr_batchsize; \n"
" } \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" } \n"
" else if (_half_step == 1){ \n"
" if (halfpaths_assigned+thr_batchsize <= total_paths) \n"
" halfpaths_assigned += thr_batchsize; \n"
" else { \n"
" thr_batchsize = (total_paths - halfpaths_assigned); \n"
" halfpaths_assigned += thr_batchsize; \n"
" } \n"
" } \n");
fprintf(outfile,
" \n"
" pthread_mutex_unlock(&tasklock); \n"
" \n"
" /*****************************************************************/ \n"
" \n"
" _segment0(thr_batchsize); \n"
" \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
"#ifdef DEBUG_MPI_Scheduling \n"
" if (_half_step == 0) \n");
fprintf(outfile,
" fschedcount[0]+=thr_batchsize; \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" else \n"
" hschedcount[0]+=thr_batchsize; \n"
"#endif \n");
fprintf(outfile,
" \n"
" gettimeofday(tp2, NULL); \n"
" \n"
" /********************CALCULATE NEW BATCH SIZE*********************/ \n"
" \n"
" tp3 = ((tp2->tv_sec - tp1->tv_sec)+((tp2->tv_usec-tp1->tv_usec)/1000000.0)); \n"
" \n"
" if (thr_throughput == 0.0) \n"
" thr_throughput = thr_batchsize/tp3; \n"
" else \n"
" thr_throughput = (1-thr_talpha)*thr_throughput + thr_talpha * (thr_batchsize/tp3); \n"
" \n"
"#ifdef DEBUG_MPI_Scheduling \n"
" printf(\"Rank 0, Completed %%d tasks in %%5.5f secs, T %%5.5f\\n\",thr_batchsize, tp3,thr_throughput); \n"
"#endif \n"
" \n"
" thr_batchsize = MAX(1,(int) (thr_throughput * thr_time_per_batch)); \n"
" /*****************************************************************/ \n"
" } \n"
" \n"
" //Unlocking indicates that the thread slave has finished processing \n"
" //and allows the master to gather results it has produced \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" hthroughput[0] = thr_throughput; \n"
" delete tp1; \n"
" delete tp2; \n");
else
fprintf(outfile,
" fthroughput[0] = thr_throughput; \n");
fprintf(outfile,
" pthread_mutex_unlock(&finlock); \n"
" return 0; \n"
"} \n"
" \n"
"void *slave(void *ptr){ \n"
" //This is a slave that runs on another process, it requests more tasks from the master \n"
" //when it finishes its current task. It takes the length of the solution array as \n"
" //argument. \n"
" \n"
" MPI_Status stat; \n"
" int count; /*batch size*/ \n"
" int i; \n"
" double tp1, tp2; \n"
" unsigned long _length; \n"
" \n"
" while(1){ \n"
" //Wait for task from master \n"
" MPI_Recv(&count, 1, MPI_INT,0,MPI_ANY_TAG,MPI_COMM_WORLD,&stat); \n"
" \n"
" tp1 = MPI_Wtime(); \n"
" \n"
" //Check message type \n"
" if (stat.MPI_TAG == 0) { \n"
" break; \n"
" } \n");
if(simulation()->parameters()->errorCheck)
fprintf(outfile,
" else if (stat.MPI_TAG == 1){ \n"
" _half_step = 0; \n"
" } \n"
" else if (stat.MPI_TAG == 2 && _half_step == 0){ \n"
" printf(\"Rank[%%i] beginning half step paths\\n\",rank); \n"
" _half_step = 1; \n"
" _gen1[0] = 346+rank; \n"
" _gen2[0] = 983+rank; \n"
" _gen1[1] = 0; \n"
" _gen2[1] = 0; \n"
" _gen1[2] = 0; \n"
" _gen2[2] = 0; \n"
" erand48(_gen1); \n"
" erand48(_gen2); \n"
" } \n");
fprintf(outfile,
" \n"
" _segment0(count); \n"
" \n"
" tp2 = MPI_Wtime() - tp1; \n"
" \n"
" //Send completion notice to the master \n"
" MPI_Send(&tp2,1,MPI_DOUBLE,0,stat.MPI_TAG,MPI_COMM_WORLD); \n"
" } \n"
" //Tell master that slave is done \n"
" MPI_Send(NULL, 0, MPI_DOUBLE,0,0,MPI_COMM_WORLD); \n"
" return 0; \n"
"} \n"
" \n"
"//Debugging functions \n"
" \n"
"int print_array(double *array, int size){ \n"
" int i, j; \n"
" for (i=0, j=0 ; i<size ; i++){ \n"
" printf(\"%%9.2f \", array[i]); \n"
" j++; \n"
" if (j==10){ \n"
" printf(\"\\n\"); \n"
" j = 0; \n"
" } \n"
" } \n"
" return 0; \n"
"} \n"
" \n"
"int print_array2(int *array, int size){ \n"
" int i, j; \n"
" for (i=0, j=0 ; i<size ; i++){ \n"
" printf(\"%%d: %%d\\n\", i, array[i]); \n"
" j++; \n"
" if (j==10){ \n"
" printf(\"\\n\"); \n"
" j = 0; \n"
" } \n"
" } \n"
" return 0; \n"
"} \n\n");
}
fprintf(outfile,"// *************************\n");
if(!simulation()->parameters()->usempi)
fprintf(outfile,"void _segment0() {\n");
else if(simulation()->parameters()->mpiMethod=="Scheduling")
fprintf(outfile,"void _segment0(int count) {\n");
else
fprintf(outfile,"void _segment0(const int& size) {\n");
fprintf(outfile,"\n");
for(list<const xmdsSegment*>:: const_iterator ppxmdsSegment
= mySegmentList.begin(); ppxmdsSegment !=
mySegmentList.end(); ppxmdsSegment++) {
(*ppxmdsSegment)->writeInitialisationCalls(outfile);
}
if((simulation()->parameters()->stochastic)&&(!(simulation()->parameters()->mpiMethod=="Scheduling")||(!simulation()->parameters()->usempi))) {
if(simulation()->parameters()->errorCheck) {
if(simulation()->parameters()->usempi) {
fprintf(outfile,"_gen1[0] = %li+rank;\n",simulation()->parameters()->seed[0]);
fprintf(outfile,"_gen2[0] = %li+rank;\n",simulation()->parameters()->seed[1]);
}
else {
fprintf(outfile,"_gen1[0] = %li;\n",simulation()->parameters()->seed[0]);
fprintf(outfile,"_gen2[0] = %li;\n",simulation()->parameters()->seed[1]);
}
fprintf(outfile,"_gen1[1] = 0;\n");
fprintf(outfile,"_gen2[1] = 0;\n");
fprintf(outfile,"_gen1[2] = 0;\n");
fprintf(outfile,"_gen2[2] = 0;\n");
fprintf(outfile,"erand48(_gen1);\n");
fprintf(outfile,"erand48(_gen2);\n");
}
else {
if(simulation()->parameters()->usempi) {
fprintf(outfile,"_gen[0] = %li+rank;\n",simulation()->parameters()->seed[0]);
}
else {
fprintf(outfile,"_gen[0] = %li;\n",simulation()->parameters()->seed[0]);
}
fprintf(outfile,"_gen[1] = 0;\n");
fprintf(outfile,"_gen[2] = 0;\n");
fprintf(outfile,"erand48(_gen);\n");
}
}
const char* indent;
if((simulation()->parameters()->mpiMethod=="Scheduling")&&(simulation()->parameters()->stochastic)&&(simulation()->parameters()->usempi)) {
fprintf(outfile," for(long _i0=0;_i0<count;_i0++)");
indent=" ";
fprintf(outfile," {\n");
}
else {
if(simulation()->parameters()->nPaths>1) {
if(simulation()->parameters()->usempi) {
fprintf(outfile," for(long _i0=rank;_i0<_n_paths;_i0+=size)");
}
else {
fprintf(outfile," for(long _i0=0;_i0<_n_paths;_i0++)");
}
indent=" ";
fprintf(outfile," {\n");
}
else
indent="";
fprintf(outfile,"\n");
}
fprintf(outfile,"%s%s=0;\n",indent,simulation()->parameters()->propDimName.c_str());
fprintf(outfile,"\n");
simulation()->field()->writeInitialisationCalls(outfile,indent);
for(unsigned long i=0;i<simulation()->output()->nMomentGroups();i++) {
fprintf(outfile,"%s_mg%li_sample_pointer=0;\n",indent,i);
if(simulation()->output()->momentGroup(i)->requiresIntegrations()|(simulation()->parameters()->usempi&!simulation()->parameters()->stochastic)) {
fprintf(outfile,"%s_mg%li_raw_initialise();\n",indent,i);
}
}
fprintf(outfile,"\n");
simulation()->field()->writeSampleCalls(outfile,indent);
if((simulation()->parameters()->mpiMethod=="Scheduling")&&(simulation()->parameters()->stochastic)&&(simulation()->parameters()->usempi)) {
//Write nothing
}
else if(simulation()->parameters()->nPaths>1) {
fprintf(outfile,"%sprintf(\"Starting path %%li\\n\",_i0+1);\n",indent);
}
for(list<const xmdsSegment*>::const_iterator ppxmdsSegment =
mySegmentList.begin(); ppxmdsSegment !=
mySegmentList.end(); ppxmdsSegment++) {
fprintf(outfile,"%s_segment%li(1);\n",indent,(*ppxmdsSegment)->segmentNumber);
}
fprintf(outfile,"\n");
for(unsigned long i=0;i<simulation()->output()->nMomentGroups();i++) {
fprintf(outfile,"%s_mg%li_process();\n",indent,i);
}
if((simulation()->parameters()->mpiMethod=="Scheduling")&&(simulation()->parameters()->stochastic)&&(simulation()->parameters()->usempi)) {
fprintf(outfile," }\n");
}
else if(simulation()->parameters()->nPaths>1) {
fprintf(outfile," }\n");
}
fprintf(outfile,"}\n");
fprintf(outfile,"\n");
}
else {
// I am a sub-sequence
if(verbose()) {
printf("Writing sequence routine ...\n");
}
fprintf(outfile,"// ********************************************************\n");
fprintf(outfile,"// segment %li (sequence) routines\n",segmentNumber);
fprintf(outfile,"\n");
fprintf(outfile,"// *************************\n");
fprintf(outfile,"void _segment%li(unsigned long cycle) {\n",segmentNumber);
for(list<const xmdsSegment*>:: const_iterator ppxmdsSegment
= mySegmentList.begin(); ppxmdsSegment !=
mySegmentList.end(); ppxmdsSegment++) {
(*ppxmdsSegment)->writeInitialisationCalls(outfile);
}
fprintf(outfile,"for(unsigned long i=0; i<%li; i++) {\n",nCycles);
for(list<const xmdsSegment*>::const_iterator ppxmdsSegment =
mySegmentList.begin(); ppxmdsSegment !=
mySegmentList.end(); ppxmdsSegment++) {
fprintf(outfile," _segment%li(i+1);\n",
(*ppxmdsSegment)->segmentNumber);
}
fprintf(outfile," }\n");
fprintf(outfile,"}\n");
fprintf(outfile,"\n");
}
xmdsElement::writeRoutines(outfile);
fprintf(outfile,"\n");
};
// ******************************************************************************
xmdsSegment* xmdsSequence::createxmdsSequence() {
if(debugFlag) {
printf("xmdsSequence::createxmdsSequence\n");
}
xmdsSegment* newxmdsSegment = new xmdsSequence(simulation(),verbose());
addChild((xmdsElement*) newxmdsSegment);
mySegmentList.push_back(newxmdsSegment);
return newxmdsSegment;
}
// ******************************************************************************
xmdsSegment* xmdsSequence::createxmdsIntegrateRK4EX() {
if(debugFlag) {
printf("xmdsSequence::createxmdsIntegrateRK4EX\n");
}
xmdsSegment* newxmdsSegment = new xmdsIntegrateRK4EX(simulation(),verbose());
addChild((xmdsElement*) newxmdsSegment);
mySegmentList.push_back(newxmdsSegment);
return newxmdsSegment;
}
// ******************************************************************************
xmdsSegment* xmdsSequence::createxmdsIntegrateRK4IP() {
if(debugFlag) {
printf("xmdsSequence::createxmdsIntegrateRK4IP\n");
}
xmdsSegment* newxmdsSegment = new xmdsIntegrateRK4IP(simulation(),verbose());
addChild((xmdsElement*) newxmdsSegment);
mySegmentList.push_back(newxmdsSegment);
return newxmdsSegment;
}
// ******************************************************************************
xmdsSegment* xmdsSequence::createxmdsIntegrateARK45EX() {
if(debugFlag) {
printf("xmdsSequence::createxmdsIntegrateARK45EX\n");
}
xmdsSegment* newxmdsSegment = new xmdsIntegrateARK45EX(simulation(),verbose());
addChild((xmdsElement*) newxmdsSegment);
mySegmentList.push_back(newxmdsSegment);
return newxmdsSegment;
}
// ******************************************************************************
xmdsSegment* xmdsSequence::createxmdsIntegrateARK45IP() {
if(debugFlag) {
printf("xmdsSequence::createxmdsIntegrateARK45IP\n");
}
xmdsSegment* newxmdsSegment = new xmdsIntegrateARK45IP(simulation(),verbose());
addChild((xmdsElement*) newxmdsSegment);
mySegmentList.push_back(newxmdsSegment);
return newxmdsSegment;
}
// ******************************************************************************
xmdsSegment* xmdsSequence::createxmdsIntegrateSIEX() {
if(debugFlag) {
printf("xmdsSequence::createxmdsIntegrateSIEX\n");
}
xmdsSegment* newxmdsSegment = new xmdsIntegrateSIEX(simulation(),verbose());
addChild((xmdsElement*) newxmdsSegment);
mySegmentList.push_back(newxmdsSegment);
return newxmdsSegment;
}
// ******************************************************************************
xmdsSegment* xmdsSequence::createxmdsIntegrateSIIP() {
if(debugFlag) {
printf("xmdsSequence::createxmdsIntegrateSIIP\n");
}
xmdsSegment* newxmdsSegment = new xmdsIntegrateSIIP(simulation(),verbose());
addChild((xmdsElement*) newxmdsSegment);
mySegmentList.push_back(newxmdsSegment);
return newxmdsSegment;
}
// ******************************************************************************
xmdsSegment* xmdsSequence::createxmdsFilter() {
if(debugFlag) {
printf("xmdsSequence::createxmdsIntegrateFilter\n");
}
xmdsSegment* newxmdsSegment = new xmdsFilter(simulation(),verbose());
addChild((xmdsElement*) newxmdsSegment);
mySegmentList.push_back(newxmdsSegment);
return newxmdsSegment;
}
syntax highlighted by Code2HTML, v. 0.9.1