Go back to household_stochastic_simulation

This file defines the main program.

<project-file type=“source”/>

<content> #include “networkstructure.h” #include “netprocess.h” #include “parameters.h” #include “eventloggers.h” #include “genrand2.h”

#include <iostream> #include <fstream> #include <time.h>

int main(int argc, const char *argv[]) {

  Parameters params;

StringArgument ArgNet(“network”, “network data file name (output of the netgen program)”, Argument::REQUIRED);

  StringArgument ArgI0("initial-infections", "the initial infections data file name (a list of node IDs)", Argument::NOT_REQUIRED);
  LongArgument ArgSeed("seed", "the seed to initialize the random number generator", Argument::NOT_REQUIRED);	
  StringArgument ArgLog("log", "log file name", Argument::REQUIRED);
  FloatArgument ArgBetaHousehold("household-contact-rate", "the transmission rate", Argument::REQUIRED);
  FloatArgument ArgBetaExternal("external-contact-rate", "the transmission rate", Argument::REQUIRED);
  FloatArgument ArgBeta("contact-rate", "the transmission rate", Argument::REQUIRED);
  DistributionArgument ArgTinf("infectious-period", "the infectious period", Argument::REQUIRED);
  DistributionArgument ArgTlat("latent-period", "the latent period", Argument::REQUIRED);
  FloatArgument ArgQuit("quit-time", "the simulation end time", Argument::REQUIRED);
  BoolArgument ArgSIR("has-immunity", "SIR/SIS model", Argument::REQUIRED);
  IntArgument ArgRuns("runs", "the number of simulation runs", Argument::NOT_REQUIRED, 1);
  BoolArgument ArgDumpNet("dump-network", "whether to dump the contact network in the log file", Argument::NOT_REQUIRED, false);
  BoolArgument ArgCases("count-cases", "if true, count the number of new cases in each time unit, otherwise, count the number of infected individual in each time unit", Argument::NOT_REQUIRED, false);

  params.addArgument(ArgNet);
  params.addArgument(ArgI0);
  params.addArgument(ArgSeed);
  params.addArgument(ArgLog);
  params.addArgument(ArgBetaHousehold);
  params.addArgument(ArgBetaExternal);
  params.addArgument(ArgTinf);
  params.addArgument(ArgTlat);
  params.addArgument(ArgQuit);
  params.addArgument(ArgSIR);
  params.addArgument(ArgRuns);
  params.addArgument(ArgDumpNet);
  params.addArgument(ArgCases);
  
  if (!params.parse(argc, argv)) {
  	params.usage(argv[0]);
  	return 1;
  }
  
  // initialize the random number generator
  time_t seed;
  if (ArgSeed.set())
  	seed = ArgSeed.value();
  else
  	time(&seed);
  sgenrand2(seed);

  NetworkStructure net(ArgNet.value());

  // set initial infections
  std::vector<int> I0;
  int n = net.nodes();
  if (ArgI0.set()) {
      FILE *fI0 = fopen(ArgI0.value(), "r");
      if (fI0 == NULL) {
          std::cerr<<"Error: initial infections data file not readable!"<<std::endl;
          return 1;
      }
      int a;
      while (fscanf(fI0, "%d", &a) > 0)
          I0.push_back(a);
      fclose(fI0);
  }
  else
      I0.push_back(genrand2i() % n + 1);
  std::ofstream log(ArgLog.value());

  net.infect(I0);
  if (ArgDumpNet.value()) net.dumpNet(log);
  net.dumpHouseholdCompartments(log);
  net.dumpEffectiveDegrees(log);

  // SimStreamLogger::registerEventLogger(new SimDieLogger());
  // SimProcess::registerEventLogger(new SimContactLogger());
  // SimProcess::registerEventLogger(new SimInfectionStateLogger());
  Counter *counter = NULL;
  if (ArgCases.value()) 
  	counter = new CaseCounter();
  else
  	counter = new Counter();
  SimProcess::registerEventLogger(counter);
  // SimStreamLogger::registerEventLogger(new SimChangeOwnerLogger());
  NetProcess proc(log, ArgBetaHousehold.value(), ArgBetaExternal.value(),
                  ArgTlat.value(), ArgTinf.value(), ArgSIR.value(), net);

  proc.initialInfections(I0);

  // std::cout<<"running ..."<<std::endl;
  for (int i = 1; i <= ArgRuns.value(); i++) {
      log << "Run : " << i << std::endl;
      proc.schedule(new SimQuitEvent(ArgQuit.value()));
      proc.run();
      proc.countAffectedHouseholds();
      if (i < ArgRuns.value()) {
          counter->reset();
          proc.reset();
      }
  }

  log.close();
  return 0;

} </content>