Attachment #1978 for bug #1927

View | Details | Raw Unified | Return to bug 1927
Collapse All | Expand All




    }
}

double
ChiSquare (Ptr<RandomVariableStream> rng, const double hMin, const double hMax)
{
  gsl_histogram * h = gsl_histogram_alloc (N_BINS);

  gsl_histogram_set_ranges_uniform (h, hMin, hMax);

  std::vector<double> expected (N_BINS, (double)N_MEASUREMENTS / (double)N_BINS);

  return ChiSquare (rng, h, expected);

}

double
ChiSquare (Ptr<RandomVariableStream> rng, gsl_histogram * h, const std::vector<double> expected)
{
  for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
    {
      gsl_histogram_increment (h, rng->GetValue ());
    }

  // TODO:  use std::valarray for these computations?
  double tmp[N_BINS];

  for (uint32_t i = 0; i < N_BINS; ++i)
    {
      tmp[i] = gsl_histogram_get (h, i);
      tmp[i] -= expected[i];
      tmp[i] *= tmp[i];
      tmp[i] /= expected[i];
    }

  gsl_histogram_free (h);

  double chiSquared = 0;

  for (uint32_t i = 0; i < N_BINS; ++i)
    {
      chiSquared += tmp[i];
    }

  return chiSquared;

}
  

} // anonymous namespace

// ===========================================================================

  RandomVariableStreamUniformTestCase ();
  virtual ~RandomVariableStreamUniformTestCase ();

  double ChiSquaredTest (Ptr<UniformRandomVariable> u);

private:
  virtual void DoRun (void);
};

{
}

double
RandomVariableStreamUniformTestCase::ChiSquaredTest (Ptr<UniformRandomVariable> u)
{
  gsl_histogram * h = gsl_histogram_alloc (N_BINS);

  // Note that this assumes that the range for u is [0,1], which is
  // the default range for this distribution.
  gsl_histogram_set_ranges_uniform (h, 0., 1.);

  for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
    {
      gsl_histogram_increment (h, u->GetValue ());
    }

  double tmp[N_BINS];

  double expected = ((double)N_MEASUREMENTS / (double)N_BINS);

  for (uint32_t i = 0; i < N_BINS; ++i)
    {
      tmp[i] = gsl_histogram_get (h, i);
      tmp[i] -= expected;
      tmp[i] *= tmp[i];
      tmp[i] /= expected;
    }

  gsl_histogram_free (h);

  double chiSquared = 0;

  for (uint32_t i = 0; i < N_BINS; ++i)
    {
      chiSquared += tmp[i];
    }

  return chiSquared;
}

void
RandomVariableStreamUniformTestCase::DoRun (void)
{

  for (uint32_t i = 0; i < N_RUNS; ++i)
    {
      Ptr<UniformRandomVariable> u = CreateObject<UniformRandomVariable> ();
      
      // Note that this assumes that the range for u is [0,1], which is
      // the default range for this distribution.
      double result = ChiSquare (u, 0., 1.);
      sum += result;
    }


  RandomVariableStreamUniformAntitheticTestCase ();
  virtual ~RandomVariableStreamUniformAntitheticTestCase ();

  double ChiSquaredTest (Ptr<UniformRandomVariable> u);

private:
  virtual void DoRun (void);
};

{
}

double
RandomVariableStreamUniformAntitheticTestCase::ChiSquaredTest (Ptr<UniformRandomVariable> u)
{
  gsl_histogram * h = gsl_histogram_alloc (N_BINS);

  // Note that this assumes that the range for u is [0,1], which is
  // the default range for this distribution.
  gsl_histogram_set_ranges_uniform (h, 0., 1.);

  for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
    {
      gsl_histogram_increment (h, u->GetValue ());
    }

  double tmp[N_BINS];

  double expected = ((double)N_MEASUREMENTS / (double)N_BINS);

  for (uint32_t i = 0; i < N_BINS; ++i)
    {
      tmp[i] = gsl_histogram_get (h, i);
      tmp[i] -= expected;
      tmp[i] *= tmp[i];
      tmp[i] /= expected;
    }

  gsl_histogram_free (h);

  double chiSquared = 0;

  for (uint32_t i = 0; i < N_BINS; ++i)
    {
      chiSquared += tmp[i];
    }

  return chiSquared;
}

void
RandomVariableStreamUniformAntitheticTestCase::DoRun (void)
{

      // Make this generate antithetic values.
      u->SetAttribute ("Antithetic", BooleanValue (true));

      // Note that this assumes that the range for u is [0,1], which is
      // the default range for this distribution.
      double result = ChiSquare (u, 0., 1.);
      sum += result;
    }


  RandomVariableStreamNormalTestCase ();
  virtual ~RandomVariableStreamNormalTestCase ();

  double ChiSquaredTest (Ptr<NormalRandomVariable> n);

private:
  virtual void DoRun (void);
};

{
}

void
RandomVariableStreamNormalTestCase::DoRun (void)
{
  SeedManager::SetSeed (time (0));

  double sum = 0.;
  double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);

  gsl_histogram * h = gsl_histogram_alloc (N_BINS);

  double range[N_BINS + 1];

      expected[i] *= N_MEASUREMENTS;
    }

  for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
    {
      gsl_histogram_increment (h, n->GetValue ());
    }

  double tmp[N_BINS];

  for (uint32_t i = 0; i < N_BINS; ++i)
    {
      tmp[i] = gsl_histogram_get (h, i);
      tmp[i] -= expected[i];
      tmp[i] *= tmp[i];
      tmp[i] /= expected[i];
    }

  gsl_histogram_free (h);

  double chiSquared = 0;

  for (uint32_t i = 0; i < N_BINS; ++i)
    {
      chiSquared += tmp[i];
    }

  return chiSquared;
}

void
RandomVariableStreamNormalTestCase::DoRun (void)
{
  SeedManager::SetSeed (time (0));

  double sum = 0.;
  double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);

  for (uint32_t i = 0; i < N_RUNS; ++i)
    {
      Ptr<NormalRandomVariable> n = CreateObject<NormalRandomVariable> ();
      double result = ChiSquare (n, h, expected);
      sum += result;
    }


  AddTestCase (new RandomVariableStreamWeibullTestCase, TestCase::QUICK);
  AddTestCase (new RandomVariableStreamWeibullAntitheticTestCase, TestCase::QUICK);
  AddTestCase (new RandomVariableStreamLogNormalTestCase, TestCase::QUICK);
  /**
   * \todo This test is currently disabled because it fails sometimes.
   * A possible reason for the failure is that the antithetic code is
   * not implemented properly for this log-normal case.
   *
  AddTestCase (new RandomVariableStreamLogNormalAntitheticTestCase, TestCase::QUICK);
  */
  AddTestCase (new RandomVariableStreamGammaTestCase, TestCase::QUICK);
  /** \todo This test is currently disabled because it fails sometimes.
   *  A possible reason for the failure is that the antithetic code is
   *  not implemented properly for this gamma case.
   *
  AddTestCase (new RandomVariableStreamGammaAntitheticTestCase, TestCase::QUICK);
  */
  AddTestCase (new RandomVariableStreamErlangTestCase, TestCase::QUICK);

Return to bug 1927