verification.hpp

Go to the documentation of this file.

#ifndef _VERIFICATION_HPP
#define _VERIFICATION_HPP

#include <iostream>
#include <assert.h>
#include "geometry.hpp"
#include "random.hpp"
#include "coloring.hpp"

namespace Arak {

  class CriticalStatsEstimator {

  protected:

    unsigned long int n;

    unsigned long int sumIntVertices;

    unsigned long int sumIntEdges;

    unsigned long int sumBdVertices[4];

    std::vector<Geometry::Segment> testSegments;

    std::vector<unsigned long int> sumCrossings;

  public:

    CriticalStatsEstimator() 
      : n(0), sumIntVertices(0), sumIntEdges(0) {
      sumBdVertices[0] = 
  sumBdVertices[1] = 
  sumBdVertices[2] = 
  sumBdVertices[3] = 0;
    }

    CriticalStatsEstimator(int k, const Coloring& c,
         Arak::Util::Random& random = 
         Arak::Util::default_random) 
      : n(0), sumIntVertices(0), sumIntEdges(0), sumCrossings(k, 0) {
      sumBdVertices[0] = 
  sumBdVertices[1] = 
  sumBdVertices[2] = 
  sumBdVertices[3] = 0;
      for (int i = 0; i < k; i++) {
  Geometry::Point p = c.randomPoint(random);
  Geometry::Point q = c.randomPoint(random);
  testSegments.push_back(Geometry::Segment(p, q));
      }
    }

    void update(const Coloring& c);

    template<typename charT, typename traits>
    void report(const Geometry::Rectangle& boundary,
    double scale,
    std::basic_ostream<charT,traits>& out) const {
      // TODO: I don't know the formula unless the boundary is square
      assert(boundary.xmax() - boundary.xmin() ==
       boundary.ymax() - boundary.ymin());
      double side = CGAL::to_double(boundary.xmax() - boundary.xmin());
      double expNumIntVertices = 4.0 * M_PI * pow(side * scale, 2);
      double avgNumIntVertices = double(sumIntVertices) / double(n);
      double expNumIntEdges = 4.0 * side * scale + expNumIntVertices;
      double avgNumIntEdges = double(sumIntEdges) / double(n);
      out << std::setw(28) << "statistic" 
    << std::setw(12) << "average" 
    << std::setw(12) << "expectation" 
    << std::endl;
      out << std::setw(28) << "----------------------------" 
    << std::setw(12) << "-----------" 
    << std::setw(12) << "-----------" 
    << std::endl;
      out << std::setw(28) << "# interior vertices"
    << std::setw(12) << avgNumIntVertices
    << std::setw(12) << expNumIntVertices
    << std::endl;
      out << std::setw(28) << "# interior edges"
    << std::setw(12) << avgNumIntEdges
    << std::setw(12) << expNumIntEdges
    << std::endl;
      for (int i = 0; i < 4; i++) {
  double expNumBdVertices = 2.0 * scale * side;
  double avgNumBdVertices = double(sumBdVertices[i]) / double(n);
  out << std::setw(26) << "# boundary vertices (" 
      << std::setw(0) << i << ")" 
      << std::setw(12) << avgNumBdVertices
      << std::setw(12) << expNumBdVertices
      << std::endl;
      }
      for (unsigned int i = 0; i < testSegments.size(); i++) {
  const Geometry::Segment& s = testSegments[i];
  double length = sqrt(CGAL::to_double(s.squared_length()));
  double expNumCrossings = 2.0 * scale * length;
  double avgNumCrossings = double(sumCrossings[i]) / double(n);
  out << std::setw(26) << "# crossings of test edge (" 
      << std::setw(0) << i << ")" 
      << std::setw(12) << avgNumCrossings
      << std::setw(12) << expNumCrossings
      << std::endl;
      }
    }

  };

}

#endif

---