MarkovChain.java

package umontreal.ssj.markovchainrqmc;
 
import umontreal.ssj.stat.Tally;
import umontreal.ssj.rng.*;
import umontreal.ssj.util.Chrono;
import umontreal.ssj.util.PrintfFormat;
import umontreal.ssj.hups.*;

public abstract class MarkovChain implements Cloneable {
 
   Chrono timer = Chrono.createForSingleThread();
   public int numSteps; // Number of steps
   protected boolean stopped = false;


   public abstract void initialState();

   public abstract void nextStep(RandomStream stream);

   public abstract double getPerformance();



   public Object clone() throws CloneNotSupportedException {
      MarkovChain o = null;
      try {
         o = (MarkovChain) super.clone();
      } catch (CloneNotSupportedException e) {
         System.err.println("This MarkovChain cannot be cloned");
      }
      return o;
   }

   public boolean hasStopped() {
      return stopped;
   }

   public void simulSteps(int numSteps, RandomStream stream) {
      initialState();
      this.numSteps = numSteps;
      int step = 0;
      while (step < numSteps && !hasStopped()) {
         nextStep(stream);
         ++step;
      }
   }

   public void simulSteps(RandomStream stream) {
      simulSteps(Integer.MAX_VALUE, stream);
   }

   public void simulRuns(int n, int numSteps, RandomStream stream, Tally statRuns) {
      statRuns.init();
      for (int i = 0; i < n; i++) {
         simulSteps(numSteps, stream);
         statRuns.add(getPerformance());
      }
   }

   public void simulRunsWithSubstreams(int n, int numSteps, RandomStream stream, Tally statRuns) {
      statRuns.init();
      stream.resetStartStream();
      for (int i = 0; i < n; i++) {
         simulSteps(numSteps, stream);
         statRuns.add(getPerformance());
         stream.resetNextSubstream();
      }
   }

   public double simulMC(int n, int numSteps) {
      Tally statRuns = new Tally();
      simulRunsWithSubstreams(n, numSteps, new MRG32k3a(), statRuns);
      return statRuns.average();
   }

   public double simulMC(int n) {
      return simulMC(n, Integer.MAX_VALUE);
   }

   public void simulRepMC(int n, int numSteps, int m, Tally t) {
      for (int rep = 0; rep < m; ++rep) {
         t.add(simulMC(n, numSteps));
      }
   }

   public void simulRepMC(int n, int m, Tally t) {
      simulRepMC(n, Integer.MAX_VALUE, m, t);
   }

   public void simulRQMC(PointSet p, int m, int numSteps, PointSetRandomization rand, Tally statReps) {
      statReps.init();
      Tally statRuns = new Tally(); // Used within the runs.
      int n = p.getNumPoints(); // Number of points.
      RandomStream stream = p.iterator();
      for (int rep = 0; rep < m; rep++) {
         p.randomize(rand);
         simulRunsWithSubstreams(n, numSteps, stream, statRuns);
         statReps.add(statRuns.average());
      }
   }

   public String simulRunsFormat(int n, int numSteps, RandomStream stream, Tally statRuns) {
      timer.init();
      simulRuns(n, numSteps, stream, statRuns);
      StringBuffer sb = new StringBuffer("----------------------------------------------" + PrintfFormat.NEWLINE);
      sb.append("MC simulations:" + PrintfFormat.NEWLINE);
      sb.append(" Number of runs n          = " + n + PrintfFormat.NEWLINE);
      sb.append(formatResults(statRuns));
      sb.append(" CPU Time = " + timer.format() + PrintfFormat.NEWLINE);
      return sb.toString();
   }

   public String simulRunsWithSubstreamsFormat(int n, int numSteps, RandomStream stream, Tally statRuns) {
      timer.init();
      simulRunsWithSubstreams(n, numSteps, stream, statRuns);
      StringBuffer sb = new StringBuffer("----------------------------------------------" + PrintfFormat.NEWLINE);
      sb.append("MC simulations with substreams:" + PrintfFormat.NEWLINE);
      sb.append(" Number of runs n          = " + n + PrintfFormat.NEWLINE);
      sb.append(formatResults(statRuns));
      sb.append(" CPU Time = " + timer.format() + PrintfFormat.NEWLINE);
      return sb.toString();
   }

   public String simulRQMCFormat(PointSet p, int m, int numSteps, PointSetRandomization rand, Tally statReps) {
      timer.init();
      simulRQMC(p, m, numSteps, rand, statReps);
      int n = p.getNumPoints();
      StringBuffer sb = new StringBuffer("----------------------------------------------" + PrintfFormat.NEWLINE);
      sb.append("RQMC simulations:" + PrintfFormat.NEWLINE + PrintfFormat.NEWLINE);
      sb.append(p.toString());
      sb.append(PrintfFormat.NEWLINE + " Number of indep. randomization, m = " + m + PrintfFormat.NEWLINE);
      sb.append(" Number of points n        = " + n + PrintfFormat.NEWLINE);
      sb.append(formatResultsRQMC(statReps, n));
      sb.append(" CPU Time = " + timer.format() + PrintfFormat.NEWLINE);
      return sb.toString();
   }

   public String testImprovementRQMCFormat(PointSet p, int m, int numSteps, PointSetRandomization rand, double varMC,
         Tally statReps) {
      // Removed next line because numSteps may be infinite!
      // p.randomize (0, numSteps * dimPerStep, noise);
      StringBuffer sb = new StringBuffer(simulRQMCFormat(p, m, numSteps, rand, statReps));
      double var = p.getNumPoints() * statReps.variance();
      sb.append(" Variance ratio: " + PrintfFormat.format(15, 10, 4, varMC / var) + PrintfFormat.NEWLINE);
      return sb.toString();
   }

   public String formatResults(Tally stat) {
      StringBuffer sb = new StringBuffer(" Average value             = ");
      sb.append(PrintfFormat.format(12, 9, 5, stat.average()) + PrintfFormat.NEWLINE);
      sb.append(" Variance                  = ");
      sb.append(PrintfFormat.format(12, 9, 5, stat.variance()) + PrintfFormat.NEWLINE);
      sb.append(stat.formatCIStudent(0.9, 7));
      return sb.toString();
   }

   public String formatResultsRQMC(Tally stat, int numPoints) {
      StringBuffer sb = new StringBuffer(" Average value             = ");
      sb.append(PrintfFormat.format(12, 9, 5, stat.average()) + PrintfFormat.NEWLINE);
      sb.append(" Variance * numPoints      = ");
      sb.append(PrintfFormat.format(12, 9, 5, numPoints * stat.variance()) + PrintfFormat.NEWLINE);
      sb.append(stat.formatCIStudent(0.9, 7));
      return sb.toString();
   }
 
}