1 /* -*- Mode:C++; c-file-style:"gnu"; indent-tabs-mode:nil; -*- */
2 /*
3  * Copyright (c) 2009-12 University of Washington
4  *
5  * This program is free software; you can redistribute it and/or modify
6  * it under the terms of the GNU General Public License version 2 as
7  * published by the Free Software Foundation;
8  *
9  * This program is distributed in the hope that it will be useful,
10  * but WITHOUT ANY WARRANTY; without even the implied warranty of
11  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
12  * GNU General Public License for more details.
13  *
14  * You should have received a copy of the GNU General Public License
15  * along with this program; if not, write to the Free Software
16  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
17  *
18  * This file is based on rng-test-suite.cc.
19  *
20  * Modified by Mitch Watrous <watrous@u.washington.edu>
21  *
22  */
23 
24 
25 #include <gsl/gsl_cdf.h>
26 #include <gsl/gsl_histogram.h>
27 #include <gsl/gsl_sf_zeta.h>
28 #include <ctime>
29 #include <fstream>
30 #include <cmath>
31 
32 #include "ns3/boolean.h"
33 #include "ns3/double.h"
34 #include "ns3/string.h"
35 #include "ns3/integer.h"
36 #include "ns3/test.h"
37 #include "ns3/log.h"
38 #include "ns3/rng-seed-manager.h"
39 #include "ns3/random-variable-stream.h"
40 
41 using namespace ns3;
42 
43 NS_LOG_COMPONENT_DEFINE ("RandomVariableStreamGenerators");
44 
45 namespace ns3 {
46 
47 namespace test {
48 
49 namespace RandomVariable {
50 
51 /**
52  * Base class for RandomVariableStream test suites.
53  */
54 class TestCaseBase : public TestCase
55 {
56 public:
57   /** Number of bins for sampling the distributions. */
58   static const uint32_t N_BINS {50};
59   /** Number of samples to draw when populating the distributions. */
60   static const uint32_t N_MEASUREMENTS {1000000};
61   /** Number of retry attempts to pass a chi-square test. */
62   static const uint32_t N_RUNS {5};
63 
64   /**
65    * Constructor
66    * \param [in] name The test case name.
67    */
TestCaseBase(std::string name)68   TestCaseBase (std::string name)
69     : TestCase (name)
70   {}
71 
72   /**
73    * Configure a GSL histogram with uniform bins, with optional
74    * under/over-flow bins.
75    * \param [in,out] h The GSL histogram to configure.
76    * \param [in] start The minimum value of the lowest bin.
77    * \param [in] end The maximum value of the last bin.
78    * \param [in] underflow If \c true the lowest bin should contain the underflow,
79    * \param [in] overflow If \c ture the highest bin should contain the overflow.
80    * \returns A vector of the bin edges, including the top of the highest bin.
81    * This vector has one more entry than the number of bins in the histogram.
82    */
83   std::vector<double>
UniformHistogramBins(gsl_histogram * h,double start,double end,bool underflow=true,bool overflow=true) const84   UniformHistogramBins (gsl_histogram *h, double start, double end,
85                       bool underflow = true, bool overflow = true) const
86   {
87     NS_LOG_FUNCTION (this << h << start << end);
88     std::size_t nBins = gsl_histogram_bins (h);
89     double increment = (end - start) / (nBins - 1.);
90     double d = start;
91 
92     std::vector<double> range (nBins + 1);
93 
94     for (auto & r : range)
95       {
96         r = d;
97         d += increment;
98       }
99     if (underflow)
100       {
101         range[0] = -std::numeric_limits<double>::max ();
102       }
103     if (overflow)
104       {
105         range[nBins] = std::numeric_limits<double>::max ();
106       }
107 
108     gsl_histogram_set_ranges (h, range.data (), nBins + 1);
109     return range;
110   }
111 
112   /**
113    * Compute the average of a random variable.
114    * \param [in] rng The random variable to sample.
115    * \returns The average of \c N_MEASUREMENTS samples.
116    */
117   double
Average(Ptr<RandomVariableStream> rng) const118   Average (Ptr<RandomVariableStream> rng) const
119   {
120     NS_LOG_FUNCTION (this << rng);
121     double sum = 0.0;
122     for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
123       {
124         double value = rng->GetValue ();
125         sum += value;
126       }
127     double valueMean = sum / N_MEASUREMENTS;
128     return valueMean;
129   }
130 
131   /** A factory base class to create new instances of a random variable. */
132   class RngGeneratorBase
133   {
134   public:
135     /**
136      * Create a new instance of a random variable stream
137      * \returns The new random variable stream instance.
138      */
139     virtual Ptr<RandomVariableStream> Create (void) const = 0;
140   };
141 
142   /**
143    * Factory class to create new instances of a particular random variable stream.
144    *
145    * \tparam RNG The type of random variable generator to create.
146    */
147   template <typename RNG>
148   class RngGenerator : public RngGeneratorBase
149   {
150   public:
151     /**
152      * Constructor.
153      * \param [in] anti Create antithetic streams if \c true.
154      */
RngGenerator(bool anti=false)155     RngGenerator (bool anti = false)
156       : m_anti (anti)
157     {}
158 
159     // Inherited
160     Ptr<RandomVariableStream>
Create(void) const161     Create (void) const
162     {
163       auto rng = CreateObject<RNG> ();
164       rng->SetAttribute ("Antithetic", BooleanValue (m_anti));
165       return rng;
166     }
167 
168   private:
169     /** Whether to create antithetic random variable streams. */
170     bool m_anti;
171   };
172 
173   /**
174    * Compute the chi squared value of a sampled distribution
175    * compared to the expected distribution.
176    *
177    * This function captures the actual computation of the chi square,
178    * given an expected distribution.
179    *
180    * The random variable is sampled \c N_MEASUREMENTS times, filling
181    * a histogram. The chi square value is formed by comparing to the
182    * expected distribution.
183    * \param [in,out] h The histogram, which defines the binning for sampling.
184    * \param [in] expected The expected distribution.
185    * \param [in] rng The random variable to sample.
186    * \returns The chi square value.
187    */
188   double
ChiSquared(gsl_histogram * h,const std::vector<double> & expected,Ptr<RandomVariableStream> rng) const189   ChiSquared (gsl_histogram * h,
190                 const std::vector<double> & expected,
191                 Ptr<RandomVariableStream> rng) const
192   {
193     NS_LOG_FUNCTION (this << h << expected.size () << rng);
194     NS_ASSERT_MSG (gsl_histogram_bins (h) == expected.size (),
195                    "Histogram and expected vector have different sizes.");
196 
197     // Sample the rng into the histogram
198     for (std::size_t i = 0; i < N_MEASUREMENTS; ++i)
199       {
200         double value = rng->GetValue ();
201         gsl_histogram_increment (h, value);
202       }
203 
204     // Compute the chi square value
205     double chiSquared = 0;
206     std::size_t nBins = gsl_histogram_bins (h);
207     for (std::size_t i = 0; i < nBins; ++i)
208       {
209         double hbin = gsl_histogram_get (h, i);
210         double tmp = hbin - expected[i];
211         tmp *= tmp;
212         tmp /= expected[i];
213         chiSquared += tmp;
214       }
215 
216     return chiSquared;
217   }
218 
219   /**
220    * Compute the chi square value from a random variable.
221    *
222    * This function sets up the binning and expected distribution
223    * needed to actually compute the chi squared value, which
224    * should be done by a call to ChiSquared.
225    *
226    * This is the point of customization expected to be implemented
227    * in derived classes with the appropriate histogram binning and
228    * expected distribution.  For example
229    *
230    *    SomeRngTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
231    *    {
232    *      gsl_histogram * h = gsl_histogram_alloc (N_BINS);
233    *      auto range = UniformHistogramBins (h, -4., 4.);
234    *      std::vector<double> expected (N_BINS);
235    *      // Populated expected
236    *      for (std::size_t i = 0; i < N_BINS; ++i)
237    *        {
238    *          expected[i] = ...;
239    *          expected[i] *= N_MEASUREMENTS;
240    *        }
241    *      double chiSquared = ChiSquared (h, expected, rng);
242    *      gsl_histogram_free (h);
243    *      return chiSquared;
244    *    }
245    *
246    * \param [in] rng The random number generator to test.
247    * \returns The chi squared value.
248    */
249   virtual double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const250   ChiSquaredTest (Ptr<RandomVariableStream> rng) const
251   {
252     return 0;
253   }
254 
255   /**
256    * Average the chi squared value over some number of runs,
257    * each run with a new instance of the random number generator.
258    * \param [in] generator The factory to create instances of the
259    *             random number generator.
260    * \param [in] nRuns The number of runs to average over.
261    * \returns The average chi square over the number of runs.
262    */
263   double
ChiSquaredsAverage(const RngGeneratorBase * generator,std::size_t nRuns) const264   ChiSquaredsAverage (const RngGeneratorBase * generator,
265                         std::size_t nRuns) const
266   {
267     NS_LOG_FUNCTION (this << generator << nRuns);
268 
269     double sum = 0.;
270     for (std::size_t i = 0; i < nRuns; ++i)
271       {
272         auto rng = generator->Create ();
273         double result = ChiSquaredTest (rng);
274         sum += result;
275       }
276     sum /= (double)nRuns;
277     return sum;
278   }
279 
280   /**
281    * Set the seed used for this test suite.
282    *
283    * This test suite is designed to be run with both deterministic and
284    * random seed and run number values.  Deterministic values can be used
285    * for basic regression testing; random values can be used to more
286    * exhaustively test the generated streams, with the side effect of
287    * occasional test failures.
288    *
289    * By default, this test suite will use the default values of RngSeed = 1
290    * and RngRun = 1.  Users can configure any other seed and run number
291    * in the usual way, but the special value of RngRun = 0 results in
292    * selecting a RngSeed value that corresponds to the seconds since epoch
293    * (\c time (0) from \c ctime).  Note: this is not a recommended practice for
294    * seeding normal simulations, as described in the ns-3 manual, but
295    * allows the test to be exposed to a wider range of seeds.
296    *
297    * In either case, the values produced will be checked with a chi-squared
298    * test.
299    *
300    * For example, this command will cause this test suite to use the
301    * deterministic value of seed=3 and default run number=1 every time:
302    *   NS_GLOBAL_VALUE="RngSeed=3" ./test.py -s random-variable-stream-generators
303    * or equivalently (to see log output):
304    *   NS_LOG="RandomVariableStreamGenerators" NS_GLOBAL_VALUE="RngSeed=3" ./waf --run "test-runner --suite=random-variable-stream-generators"
305    *
306    * Conversely, this command will cause this test suite to use a seed
307    * based on time-of-day, and run number=0:
308    *   NS_GLOBAL_VALUE="RngRun=0" ./test.py -s random-variable-stream-generators
309    */
310   void
SetTestSuiteSeed(void)311   SetTestSuiteSeed (void)
312   {
313     if (m_seedSet == false)
314       {
315         uint32_t seed;
316         if (RngSeedManager::GetRun () == 0)
317           {
318             seed = static_cast<uint32_t> (time (0));
319             m_seedSet = true;
320             NS_LOG_DEBUG ("Special run number value of zero; seeding with time of day: " << seed);
321 
322           }
323         else
324           {
325             seed = RngSeedManager::GetSeed ();
326             m_seedSet = true;
327             NS_LOG_DEBUG ("Using the values seed: " <<
328                           seed << " and run: " << RngSeedManager::GetRun ());
329           }
330         SeedManager::SetSeed (seed);
331       }
332   }
333 
334 private:
335   /** \c true if we've already set the seed the correctly. */
336   bool m_seedSet = false;
337 
338 }; // class TestCaseBase
339 
340 
341 /**
342  * Test case for uniform distribution random variable stream generator.
343  */
344 class UniformTestCase : public TestCaseBase
345 {
346 public:
347 
348   // Constructor
349   UniformTestCase ();
350 
351   // Inherited
352   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
353 
354 private:
355   // Inherited
356   virtual void DoRun (void);
357 };
358 
UniformTestCase()359 UniformTestCase::UniformTestCase ()
360   : TestCaseBase ("Uniform Random Variable Stream Generator")
361 {}
362 
363 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const364 UniformTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
365 {
366   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
367 
368   // Note that this assumes that the range for u is [0,1], which is
369   // the default range for this distribution.
370   gsl_histogram_set_ranges_uniform (h, 0., 1.);
371 
372   std::vector<double> expected (N_BINS, ((double)N_MEASUREMENTS / (double)N_BINS));
373 
374   double chiSquared = ChiSquared (h, expected, rng);
375   gsl_histogram_free (h);
376   return chiSquared;
377 }
378 
379 void
DoRun(void)380 UniformTestCase::DoRun (void)
381 {
382   NS_LOG_FUNCTION (this);
383   SetTestSuiteSeed ();
384 
385   double confidence = 0.99;
386   double maxStatistic = gsl_cdf_chisq_Pinv (confidence, (N_BINS - 1));
387   NS_LOG_DEBUG ("Chi square required at " << confidence << " confidence for " << N_BINS << " bins is " << maxStatistic);
388 
389   double result = maxStatistic;
390   // If chi-squared test fails, re-try it up to N_RUNS times
391   for (uint32_t i = 0; i < N_RUNS; ++i)
392     {
393       Ptr<UniformRandomVariable> rng = CreateObject<UniformRandomVariable> ();
394       result = ChiSquaredTest (rng);
395       NS_LOG_DEBUG ("Chi square result is " << result);
396       if (result < maxStatistic)
397         {
398           break;
399         }
400     }
401 
402   NS_TEST_ASSERT_MSG_LT (result, maxStatistic, "Chi-squared statistic out of range");
403 
404   double min = 0.0;
405   double max = 10.0;
406   double value;
407 
408   // Create the RNG with the specified range.
409   Ptr<UniformRandomVariable> x = CreateObject<UniformRandomVariable> ();
410 
411   x->SetAttribute ("Min", DoubleValue (min));
412   x->SetAttribute ("Max", DoubleValue (max));
413 
414   // Test that values are always within the range:
415   //
416   //     [min, max)
417   //
418   for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
419     {
420       value = x->GetValue ();
421       NS_TEST_ASSERT_MSG_EQ ((value >= min), true, "Value less than minimum.");
422       NS_TEST_ASSERT_MSG_LT (value, max, "Value greater than or equal to maximum.");
423     }
424 
425   // Boundary checking on GetInteger; should be [min,max]; from bug 1964
426   static const uint32_t UNIFORM_INTEGER_MIN {0};
427   static const uint32_t UNIFORM_INTEGER_MAX {4294967295U};
428   // [0,0] should return 0
429   uint32_t intValue;
430   intValue = x->GetInteger (UNIFORM_INTEGER_MIN, UNIFORM_INTEGER_MIN);
431   NS_TEST_ASSERT_MSG_EQ (intValue, UNIFORM_INTEGER_MIN, "Uniform RV GetInteger boundary testing");
432   // [UNIFORM_INTEGER_MAX, UNIFORM_INTEGER_MAX] should return UNIFORM_INTEGER_MAX
433   intValue = x->GetInteger (UNIFORM_INTEGER_MAX, UNIFORM_INTEGER_MAX);
434   NS_TEST_ASSERT_MSG_EQ (intValue, UNIFORM_INTEGER_MAX, "Uniform RV GetInteger boundary testing");
435   // [0,1] should return mix of 0 or 1
436   intValue = 0;
437   for (int i = 0; i < 20; i++)
438     {
439       intValue += x->GetInteger (UNIFORM_INTEGER_MIN, UNIFORM_INTEGER_MIN + 1);
440     }
441   NS_TEST_ASSERT_MSG_GT (intValue, 0, "Uniform RV GetInteger boundary testing");
442   NS_TEST_ASSERT_MSG_LT (intValue, 20, "Uniform RV GetInteger boundary testing");
443   // [MAX-1,MAX] should return mix of MAX-1 or MAX
444   uint32_t count = 0;
445   for (int i = 0; i < 20; i++)
446     {
447       intValue = x->GetInteger (UNIFORM_INTEGER_MAX - 1, UNIFORM_INTEGER_MAX);
448       if (intValue == UNIFORM_INTEGER_MAX)
449         {
450           count++;
451         }
452     }
453   NS_TEST_ASSERT_MSG_GT (count, 0, "Uniform RV GetInteger boundary testing");
454   NS_TEST_ASSERT_MSG_LT (count, 20, "Uniform RV GetInteger boundary testing");
455   // multiple [0,UNIFORM_INTEGER_MAX] should return non-zero
456   intValue = x->GetInteger (UNIFORM_INTEGER_MIN, UNIFORM_INTEGER_MAX);
457   uint32_t intValue2 = x->GetInteger (UNIFORM_INTEGER_MIN, UNIFORM_INTEGER_MAX);
458   NS_TEST_ASSERT_MSG_GT (intValue + intValue2, 0, "Uniform RV GetInteger boundary testing");
459 
460 }
461 
462 /**
463  * Test case for antithetic uniform distribution random variable stream generator
464  */
465 class UniformAntitheticTestCase : public TestCaseBase
466 {
467 public:
468 
469   // Constructor
470   UniformAntitheticTestCase ();
471 
472   // Inherited
473   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
474 
475 private:
476   // Inherited
477   virtual void DoRun (void);
478 };
479 
UniformAntitheticTestCase()480 UniformAntitheticTestCase::UniformAntitheticTestCase ()
481   : TestCaseBase ("Antithetic Uniform Random Variable Stream Generator")
482 {}
483 
484 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const485 UniformAntitheticTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
486 {
487   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
488 
489   // Note that this assumes that the range for u is [0,1], which is
490   // the default range for this distribution.
491   gsl_histogram_set_ranges_uniform (h, 0., 1.);
492 
493   std::vector<double> expected (N_BINS, ((double)N_MEASUREMENTS / (double)N_BINS));
494 
495   double chiSquared = ChiSquared (h, expected, rng);
496   gsl_histogram_free (h);
497   return chiSquared;
498 }
499 
500 void
DoRun(void)501 UniformAntitheticTestCase::DoRun (void)
502 {
503   NS_LOG_FUNCTION (this);
504   SetTestSuiteSeed ();
505 
506   auto generator = RngGenerator<UniformRandomVariable> (true);
507   double sum = ChiSquaredsAverage (&generator, N_RUNS);
508   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
509   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
510 
511   double min = 0.0;
512   double max = 10.0;
513   double value;
514 
515   // Create the RNG with the specified range.
516   Ptr<UniformRandomVariable> x = CreateObject<UniformRandomVariable> ();
517 
518   // Make this generate antithetic values.
519   x->SetAttribute ("Antithetic", BooleanValue (true));
520 
521   x->SetAttribute ("Min", DoubleValue (min));
522   x->SetAttribute ("Max", DoubleValue (max));
523 
524   // Test that values are always within the range:
525   //
526   //     [min, max)
527   //
528   for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
529     {
530       value = x->GetValue ();
531       NS_TEST_ASSERT_MSG_EQ ((value >= min), true, "Value less than minimum.");
532       NS_TEST_ASSERT_MSG_LT (value, max, "Value greater than or equal to maximum.");
533     }
534 
535 
536 }
537 /**
538  * Test case for constant random variable stream generator
539  */
540 class ConstantTestCase : public TestCaseBase
541 {
542 public:
543   // Constructor
544   ConstantTestCase ();
545 
546 private:
547   // Inherited
548   virtual void DoRun (void);
549 
550   /** Tolerance for testing rng values against expectation. */
551   static constexpr double TOLERANCE {1e-8};
552 };
553 
ConstantTestCase()554 ConstantTestCase::ConstantTestCase ()
555   : TestCaseBase ("Constant Random Variable Stream Generator")
556 {}
557 
558 void
DoRun(void)559 ConstantTestCase::DoRun (void)
560 {
561   NS_LOG_FUNCTION (this);
562   SetTestSuiteSeed ();
563 
564   Ptr<ConstantRandomVariable> c = CreateObject<ConstantRandomVariable> ();
565 
566   double constant;
567 
568   // Test that the constant value can be changed using its attribute.
569   constant = 10.0;
570   c->SetAttribute ("Constant", DoubleValue (constant));
571   NS_TEST_ASSERT_MSG_EQ_TOL (c->GetValue (), constant, TOLERANCE, "Constant value changed");
572   c->SetAttribute ("Constant", DoubleValue (20.0));
573   NS_TEST_ASSERT_MSG_NE (c->GetValue (), constant, "Constant value not changed");
574 
575   // Test that the constant value does not change.
576   constant = c->GetValue ();
577   for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
578     {
579       NS_TEST_ASSERT_MSG_EQ_TOL (c->GetValue (), constant, TOLERANCE, "Constant value changed in loop");
580     }
581 }
582 
583 /**
584  * Test case for sequential random variable stream generator
585  */
586 class SequentialTestCase : public TestCaseBase
587 {
588 public:
589   // Constructor
590   SequentialTestCase ();
591 
592 private:
593   // Inherited
594   virtual void DoRun (void);
595 
596   /** Tolerance for testing rng values against expectation. */
597   static constexpr double TOLERANCE {1e-8};
598 };
599 
SequentialTestCase()600 SequentialTestCase::SequentialTestCase ()
601   : TestCaseBase ("Sequential Random Variable Stream Generator")
602 {}
603 
604 void
DoRun(void)605 SequentialTestCase::DoRun (void)
606 {
607   NS_LOG_FUNCTION (this);
608   SetTestSuiteSeed ();
609 
610   Ptr<SequentialRandomVariable> s = CreateObject<SequentialRandomVariable> ();
611 
612   // The following four attributes should give the sequence
613   //
614   //    4, 4, 7, 7, 10, 10
615   //
616   s->SetAttribute ("Min", DoubleValue (4));
617   s->SetAttribute ("Max", DoubleValue (11));
618   s->SetAttribute ("Increment", StringValue ("ns3::UniformRandomVariable[Min=3.0|Max=3.0]"));
619   s->SetAttribute ("Consecutive", IntegerValue (2));
620 
621   double value;
622 
623   // Test that the sequencet is correct.
624   value = s->GetValue ();
625   NS_TEST_ASSERT_MSG_EQ_TOL (value, 4, TOLERANCE, "Sequence value 1 wrong.");
626   value = s->GetValue ();
627   NS_TEST_ASSERT_MSG_EQ_TOL (value, 4, TOLERANCE, "Sequence value 2 wrong.");
628   value = s->GetValue ();
629   NS_TEST_ASSERT_MSG_EQ_TOL (value, 7, TOLERANCE, "Sequence value 3 wrong.");
630   value = s->GetValue ();
631   NS_TEST_ASSERT_MSG_EQ_TOL (value, 7, TOLERANCE, "Sequence value 4 wrong.");
632   value = s->GetValue ();
633   NS_TEST_ASSERT_MSG_EQ_TOL (value, 10, TOLERANCE, "Sequence value 5 wrong.");
634   value = s->GetValue ();
635   NS_TEST_ASSERT_MSG_EQ_TOL (value, 10, TOLERANCE, "Sequence value 6 wrong.");
636 
637 }
638 
639 /**
640  * Test case for normal distribution random variable stream generator
641  */
642 class NormalTestCase : public TestCaseBase
643 {
644 public:
645   // Constructor
646   NormalTestCase ();
647 
648   // Inherited
649   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
650 
651 private:
652   // Inherited
653   virtual void DoRun (void);
654 
655   /** Tolerance for testing rng values against expectation, in rms. */
656   static constexpr double TOLERANCE {5};
657 };
658 
NormalTestCase()659 NormalTestCase::NormalTestCase ()
660   : TestCaseBase ("Normal Random Variable Stream Generator")
661 {}
662 
663 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const664 NormalTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
665 {
666   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
667   auto range = UniformHistogramBins (h, -4., 4.);
668 
669   std::vector<double> expected (N_BINS);
670 
671   // Note that this assumes that n has mean equal to zero and standard
672   // deviation equal to one, which are their default values for this
673   // distribution.
674   double sigma = 1.;
675 
676   for (std::size_t i = 0; i < N_BINS; ++i)
677     {
678       expected[i] = gsl_cdf_gaussian_P (range[i + 1], sigma) - gsl_cdf_gaussian_P (range[i], sigma);
679       expected[i] *= N_MEASUREMENTS;
680     }
681 
682   double chiSquared = ChiSquared (h, expected, rng);
683   gsl_histogram_free (h);
684   return chiSquared;
685 }
686 
687 void
DoRun(void)688 NormalTestCase::DoRun (void)
689 {
690   NS_LOG_FUNCTION (this);
691   SetTestSuiteSeed ();
692 
693   auto generator = RngGenerator<NormalRandomVariable> ();
694   auto rng = generator.Create ();
695 
696   double sum = ChiSquaredsAverage (&generator, N_RUNS);
697   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
698   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
699 
700   double mean = 5.0;
701   double variance = 2.0;
702 
703   // Create the RNG with the specified range.
704   Ptr<NormalRandomVariable> x = CreateObject<NormalRandomVariable> ();
705   x->SetAttribute ("Mean", DoubleValue (mean));
706   x->SetAttribute ("Variance", DoubleValue (variance));
707 
708   // Calculate the mean of these values.
709   double valueMean = Average (x);
710 
711   // The expected value for the mean of the values returned by a
712   // normally distributed random variable is equal to mean.
713   double expectedMean = mean;
714   double expectedRms  = mean / std::sqrt (variance * N_MEASUREMENTS);
715 
716   // Test that values have approximately the right mean value.
717   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedRms * TOLERANCE, "Wrong mean value.");
718 }
719 
720 /**
721  * Test case for antithetic normal distribution random variable stream generator
722  */
723 class NormalAntitheticTestCase : public TestCaseBase
724 {
725 public:
726   // Constructor
727   NormalAntitheticTestCase ();
728 
729   // Inherited
730   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
731 
732 private:
733   // Inherited
734   virtual void DoRun (void);
735 
736   /** Tolerance for testing rng values against expectation, in rms. */
737   static constexpr double TOLERANCE {5};
738 };
739 
NormalAntitheticTestCase()740 NormalAntitheticTestCase::NormalAntitheticTestCase ()
741   : TestCaseBase ("Antithetic Normal Random Variable Stream Generator")
742 {}
743 
744 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const745 NormalAntitheticTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
746 {
747   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
748   auto range = UniformHistogramBins (h, -4, 4);
749 
750   std::vector<double> expected (N_BINS);
751 
752   // Note that this assumes that n has mean equal to zero and standard
753   // deviation equal to one, which are their default values for this
754   // distribution.
755   double sigma = 1.;
756 
757   for (std::size_t i = 0; i < N_BINS; ++i)
758     {
759       expected[i] = gsl_cdf_gaussian_P (range[i + 1], sigma) - gsl_cdf_gaussian_P (range[i], sigma);
760       expected[i] *= N_MEASUREMENTS;
761     }
762 
763   double chiSquared = ChiSquared (h, expected, rng);
764 
765   gsl_histogram_free (h);
766   return chiSquared;
767 }
768 
769 void
DoRun(void)770 NormalAntitheticTestCase::DoRun (void)
771 {
772   NS_LOG_FUNCTION (this);
773   SetTestSuiteSeed ();
774 
775   auto generator = RngGenerator<NormalRandomVariable> (true);
776   double sum = ChiSquaredsAverage (&generator, N_RUNS);
777   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
778   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
779 
780   double mean = 5.0;
781   double variance = 2.0;
782 
783   // Create the RNG with the specified range.
784   Ptr<NormalRandomVariable> x = CreateObject<NormalRandomVariable> ();
785   x->SetAttribute ("Mean", DoubleValue (mean));
786   x->SetAttribute ("Variance", DoubleValue (variance));
787 
788   // Make this generate antithetic values.
789   x->SetAttribute ("Antithetic", BooleanValue (true));
790 
791   // Calculate the mean of these values.
792   double valueMean = Average (x);
793 
794   // The expected value for the mean of the values returned by a
795   // normally distributed random variable is equal to mean.
796   double expectedMean = mean;
797   double expectedRms  = mean / std::sqrt (variance * N_MEASUREMENTS);
798 
799   // Test that values have approximately the right mean value.
800   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedRms * TOLERANCE, "Wrong mean value.");
801 }
802 
803 /**
804  * Test case for exponential distribution random variable stream generator
805  */
806 class ExponentialTestCase : public TestCaseBase
807 {
808 public:
809   // Constructor
810   ExponentialTestCase ();
811 
812   // Inherited
813   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
814 
815 private:
816   // Inherited
817   virtual void DoRun (void);
818 
819   /** Tolerance for testing rng values against expectation, in rms. */
820   static constexpr double TOLERANCE {5};
821 };
822 
ExponentialTestCase()823 ExponentialTestCase::ExponentialTestCase ()
824   : TestCaseBase ("Exponential Random Variable Stream Generator")
825 {}
826 
827 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const828 ExponentialTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
829 {
830   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
831   auto range = UniformHistogramBins (h, 0, 10, false);
832 
833   std::vector<double> expected (N_BINS);
834 
835   // Note that this assumes that e has mean equal to one, which is the
836   // default value for this distribution.
837   double mu = 1.;
838 
839   for (std::size_t i = 0; i < N_BINS; ++i)
840     {
841       expected[i] = gsl_cdf_exponential_P (range[i + 1], mu) - gsl_cdf_exponential_P (range[i], mu);
842       expected[i] *= N_MEASUREMENTS;
843     }
844 
845   double chiSquared = ChiSquared (h, expected, rng);
846 
847   gsl_histogram_free (h);
848   return chiSquared;
849 }
850 
851 void
DoRun(void)852 ExponentialTestCase::DoRun (void)
853 {
854   NS_LOG_FUNCTION (this);
855   SetTestSuiteSeed ();
856 
857   auto generator = RngGenerator<ExponentialRandomVariable> ();
858   double sum = ChiSquaredsAverage (&generator, N_RUNS);
859   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
860   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
861 
862   double mean = 3.14;
863   double bound = 0.0;
864 
865   // Create the RNG with the specified range.
866   Ptr<ExponentialRandomVariable> x = CreateObject<ExponentialRandomVariable> ();
867   x->SetAttribute ("Mean", DoubleValue (mean));
868   x->SetAttribute ("Bound", DoubleValue (bound));
869 
870   // Calculate the mean of these values.
871   double valueMean = Average (x);
872   double expectedMean = mean;
873   double expectedRms  = std::sqrt (mean / N_MEASUREMENTS);
874 
875   // Test that values have approximately the right mean value.
876   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedRms * TOLERANCE, "Wrong mean value.");
877 }
878 
879 /**
880  * Test case for antithetic exponential distribution random variable stream generator
881  */
882 class ExponentialAntitheticTestCase : public TestCaseBase
883 {
884 public:
885   // Constructor
886   ExponentialAntitheticTestCase ();
887 
888   // Inherited
889   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
890 
891 private:
892   // Inherited
893   virtual void DoRun (void);
894 
895   /** Tolerance for testing rng values against expectation, in rms. */
896   static constexpr double TOLERANCE {5};
897 };
898 
ExponentialAntitheticTestCase()899 ExponentialAntitheticTestCase::ExponentialAntitheticTestCase ()
900   : TestCaseBase ("Antithetic Exponential Random Variable Stream Generator")
901 {}
902 
903 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const904 ExponentialAntitheticTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
905 {
906   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
907   auto range = UniformHistogramBins (h, 0, 10, false);
908 
909   std::vector<double> expected (N_BINS);
910 
911   // Note that this assumes that e has mean equal to one, which is the
912   // default value for this distribution.
913   double mu = 1.;
914 
915   for (std::size_t i = 0; i < N_BINS; ++i)
916     {
917       expected[i] = gsl_cdf_exponential_P (range[i + 1], mu) - gsl_cdf_exponential_P (range[i], mu);
918       expected[i] *= N_MEASUREMENTS;
919     }
920 
921   double chiSquared = ChiSquared (h, expected, rng);
922 
923   gsl_histogram_free (h);
924   return chiSquared;
925 }
926 
927 void
DoRun(void)928 ExponentialAntitheticTestCase::DoRun (void)
929 {
930   NS_LOG_FUNCTION (this);
931   SetTestSuiteSeed ();
932 
933   auto generator = RngGenerator<ExponentialRandomVariable> (true);
934   double sum = ChiSquaredsAverage (&generator, N_RUNS);
935   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
936   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
937 
938   double mean = 3.14;
939   double bound = 0.0;
940 
941   // Create the RNG with the specified range.
942   Ptr<ExponentialRandomVariable> x = CreateObject<ExponentialRandomVariable> ();
943   x->SetAttribute ("Mean", DoubleValue (mean));
944   x->SetAttribute ("Bound", DoubleValue (bound));
945 
946   // Make this generate antithetic values.
947   x->SetAttribute ("Antithetic", BooleanValue (true));
948 
949   // Calculate the mean of these values.
950   double valueMean = Average (x);
951   double expectedMean = mean;
952   double expectedRms  = std::sqrt (mean / N_MEASUREMENTS);
953 
954   // Test that values have approximately the right mean value.
955   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedRms * TOLERANCE, "Wrong mean value.");
956 }
957 
958 /**
959  * Test case for Pareto distribution random variable stream generator
960  */
961 class ParetoTestCase : public TestCaseBase
962 {
963 public:
964   // Constructor
965   ParetoTestCase ();
966 
967   // Inherited
968   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
969 
970 private:
971   // Inherited
972   virtual void DoRun (void);
973 
974   /**
975    * Tolerance for testing rng values against expectation,
976    * as a fraction of mean value.
977    */
978   static constexpr double TOLERANCE {1e-2};
979 };
980 
ParetoTestCase()981 ParetoTestCase::ParetoTestCase ()
982   : TestCaseBase ("Pareto Random Variable Stream Generator")
983 {}
984 
985 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const986 ParetoTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
987 {
988   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
989   auto range = UniformHistogramBins (h, 1, 10, false);
990 
991   std::vector<double> expected (N_BINS);
992 
993   double shape = 2.0;
994   double scale = 1.0;
995 
996   for (std::size_t i = 0; i < N_BINS; ++i)
997     {
998       expected[i] = gsl_cdf_pareto_P (range[i + 1], shape, scale) - gsl_cdf_pareto_P (range[i], shape, scale);
999       expected[i] *= N_MEASUREMENTS;
1000     }
1001 
1002   double chiSquared = ChiSquared (h, expected, rng);
1003 
1004   gsl_histogram_free (h);
1005   return chiSquared;
1006 }
1007 
1008 void
DoRun(void)1009 ParetoTestCase::DoRun (void)
1010 {
1011   NS_LOG_FUNCTION (this);
1012   SetTestSuiteSeed ();
1013 
1014   auto generator = RngGenerator<ParetoRandomVariable> ();
1015   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1016   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1017   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1018 
1019   double shape = 2.0;
1020   double scale = 1.0;
1021 
1022   // Create the RNG with the specified range.
1023   Ptr<ParetoRandomVariable> x = CreateObject<ParetoRandomVariable> ();
1024   x->SetAttribute ("Shape", DoubleValue (shape));
1025   x->SetAttribute ("Scale", DoubleValue (scale));
1026 
1027   // Calculate the mean of these values.
1028   double valueMean = Average (x);
1029 
1030   // The expected value for the mean is given by
1031   //
1032   //                   shape * scale
1033   //     E[value]  =  ---------------  ,
1034   //                     shape - 1
1035   //
1036   // where
1037   //
1038   //     scale  =  mean * (shape - 1.0) / shape .
1039   double expectedMean = (shape * scale) / (shape - 1.0);
1040 
1041   // Test that values have approximately the right mean value.
1042   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1043 }
1044 
1045 /**
1046  * Test case for antithetic Pareto distribution random variable stream generator
1047  */
1048 class ParetoAntitheticTestCase : public TestCaseBase
1049 {
1050 public:
1051   // Constructor
1052   ParetoAntitheticTestCase ();
1053 
1054   // Inherited
1055   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
1056 
1057 private:
1058   // Inherited
1059   virtual void DoRun (void);
1060 
1061   /**
1062    * Tolerance for testing rng values against expectation,
1063    * as a fraction of mean value.
1064    */
1065   static constexpr double TOLERANCE {1e-2};
1066 };
1067 
ParetoAntitheticTestCase()1068 ParetoAntitheticTestCase::ParetoAntitheticTestCase ()
1069   : TestCaseBase ("Antithetic Pareto Random Variable Stream Generator")
1070 {}
1071 
1072 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const1073 ParetoAntitheticTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
1074 {
1075   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
1076   auto range = UniformHistogramBins (h, 1, 10, false);
1077 
1078   std::vector<double> expected (N_BINS);
1079 
1080   double shape = 2.0;
1081   double scale = 1.0;
1082 
1083   for (std::size_t i = 0; i < N_BINS; ++i)
1084     {
1085       expected[i] = gsl_cdf_pareto_P (range[i + 1], shape, scale) - gsl_cdf_pareto_P (range[i], shape, scale);
1086       expected[i] *= N_MEASUREMENTS;
1087     }
1088 
1089   double chiSquared = ChiSquared (h, expected, rng);
1090 
1091   gsl_histogram_free (h);
1092   return chiSquared;
1093 }
1094 
1095 void
DoRun(void)1096 ParetoAntitheticTestCase::DoRun (void)
1097 {
1098   NS_LOG_FUNCTION (this);
1099   SetTestSuiteSeed ();
1100 
1101   auto generator = RngGenerator<ParetoRandomVariable> (true);
1102   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1103   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1104   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1105 
1106   double shape = 2.0;
1107   double scale = 1.0;
1108 
1109   // Create the RNG with the specified range.
1110   Ptr<ParetoRandomVariable> x = CreateObject<ParetoRandomVariable> ();
1111   x->SetAttribute ("Shape", DoubleValue (shape));
1112   x->SetAttribute ("Scale", DoubleValue (scale));
1113 
1114   // Make this generate antithetic values.
1115   x->SetAttribute ("Antithetic", BooleanValue (true));
1116 
1117   // Calculate the mean of these values.
1118   double valueMean = Average (x);
1119 
1120   // The expected value for the mean is given by
1121   //
1122   //                   shape * scale
1123   //     E[value]  =  ---------------  ,
1124   //                     shape - 1
1125   //
1126   // where
1127   //
1128   //     scale  =  mean * (shape - 1.0) / shape .
1129   //
1130   double expectedMean = (shape * scale) / (shape - 1.0);
1131 
1132   // Test that values have approximately the right mean value.
1133   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1134 }
1135 
1136 /**
1137  * Test case for Weibull distribution random variable stream generator
1138  */
1139 class WeibullTestCase : public TestCaseBase
1140 {
1141 public:
1142   // Constructor
1143   WeibullTestCase ();
1144 
1145   // Inherited
1146   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
1147 
1148 private:
1149   // Inherited
1150   virtual void DoRun (void);
1151 
1152   /**
1153    * Tolerance for testing rng values against expectation,
1154    * as a fraction of mean value.
1155    */
1156   static constexpr double TOLERANCE {1e-2};
1157 };
1158 
WeibullTestCase()1159 WeibullTestCase::WeibullTestCase ()
1160   : TestCaseBase ("Weibull Random Variable Stream Generator")
1161 {}
1162 
1163 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const1164 WeibullTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
1165 {
1166   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
1167   auto range = UniformHistogramBins (h, 1, 10, false);
1168 
1169   std::vector<double> expected (N_BINS);
1170 
1171   // Note that this assumes that p has shape equal to one and scale
1172   // equal to one, which are their default values for this
1173   // distribution.
1174   double a = 1.0;
1175   double b = 1.0;
1176 
1177   for (std::size_t i = 0; i < N_BINS; ++i)
1178     {
1179       expected[i] = gsl_cdf_weibull_P (range[i + 1], a, b) - gsl_cdf_weibull_P (range[i], a, b);
1180       expected[i] *= N_MEASUREMENTS;
1181       NS_LOG_INFO ("weibull: " << expected[i]);
1182     }
1183 
1184   double chiSquared = ChiSquared (h, expected, rng);
1185 
1186   gsl_histogram_free (h);
1187   return chiSquared;
1188 }
1189 
1190 void
DoRun(void)1191 WeibullTestCase::DoRun (void)
1192 {
1193   NS_LOG_FUNCTION (this);
1194   SetTestSuiteSeed ();
1195 
1196   auto generator = RngGenerator<WeibullRandomVariable> ();
1197   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1198   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1199   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1200 
1201   double scale = 5.0;
1202   double shape = 1.0;
1203 
1204   // Create the RNG with the specified range.
1205   Ptr<WeibullRandomVariable> x = CreateObject<WeibullRandomVariable> ();
1206   x->SetAttribute ("Scale", DoubleValue (scale));
1207   x->SetAttribute ("Shape", DoubleValue (shape));
1208 
1209   // Calculate the mean of these values.
1210   double valueMean = Average (x);
1211 
1212   // The expected value for the mean of the values returned by a
1213   // Weibull distributed random variable is
1214   //
1215   //     E[value]  =  scale * Gamma(1 + 1 / shape)  ,
1216   //
1217   // where Gamma() is the Gamma function.  Note that
1218   //
1219   //     Gamma(n)  =  (n - 1)!
1220   //
1221   // if n is a positive integer.
1222   //
1223   // For this test,
1224   //
1225   //     Gamma(1 + 1 / shape)  =  Gamma(1 + 1 / 1)
1226   //                           =  Gamma(2)
1227   //                           =  (2 - 1)!
1228   //                           =  1
1229   //
1230   // which means
1231   //
1232   //     E[value]  =  scale  .
1233   //
1234   double expectedMean = scale;
1235 
1236   // Test that values have approximately the right mean value.
1237   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1238 }
1239 
1240 /**
1241  * Test case for antithetic Weibull distribution random variable stream generator
1242  */
1243 class WeibullAntitheticTestCase : public TestCaseBase
1244 {
1245 public:
1246   // Constructor
1247   WeibullAntitheticTestCase ();
1248 
1249   // Inherited
1250   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
1251 
1252 private:
1253   // Inherited
1254   virtual void DoRun (void);
1255 
1256   /**
1257    * Tolerance for testing rng values against expectation,
1258    * as a fraction of mean value.
1259    */
1260   static constexpr double TOLERANCE {1e-2};
1261 };
1262 
WeibullAntitheticTestCase()1263 WeibullAntitheticTestCase::WeibullAntitheticTestCase ()
1264   : TestCaseBase ("Antithetic Weibull Random Variable Stream Generator")
1265 {}
1266 
1267 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const1268 WeibullAntitheticTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
1269 {
1270   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
1271   auto range = UniformHistogramBins (h, 1, 10, false);
1272 
1273   std::vector<double> expected (N_BINS);
1274 
1275   // Note that this assumes that p has shape equal to one and scale
1276   // equal to one, which are their default values for this
1277   // distribution.
1278   double a = 1.0;
1279   double b = 1.0;
1280 
1281   for (std::size_t i = 0; i < N_BINS; ++i)
1282     {
1283       expected[i] = gsl_cdf_weibull_P (range[i + 1], a, b) - gsl_cdf_weibull_P (range[i], a, b);
1284       expected[i] *= N_MEASUREMENTS;
1285     }
1286 
1287   double chiSquared = ChiSquared (h, expected, rng);
1288 
1289   gsl_histogram_free (h);
1290   return chiSquared;
1291 }
1292 
1293 void
DoRun(void)1294 WeibullAntitheticTestCase::DoRun (void)
1295 {
1296   NS_LOG_FUNCTION (this);
1297   SetTestSuiteSeed ();
1298 
1299   auto generator = RngGenerator<WeibullRandomVariable> (true);
1300   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1301   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1302   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1303 
1304   double scale = 5.0;
1305   double shape = 1.0;
1306 
1307   // Create the RNG with the specified range.
1308   Ptr<WeibullRandomVariable> x = CreateObject<WeibullRandomVariable> ();
1309   x->SetAttribute ("Scale", DoubleValue (scale));
1310   x->SetAttribute ("Shape", DoubleValue (shape));
1311 
1312   // Make this generate antithetic values.
1313   x->SetAttribute ("Antithetic", BooleanValue (true));
1314 
1315   // Calculate the mean of these values.
1316   double valueMean = Average (x);
1317 
1318   // The expected value for the mean of the values returned by a
1319   // Weibull distributed random variable is
1320   //
1321   //     E[value]  =  scale * Gamma(1 + 1 / shape)  ,
1322   //
1323   // where Gamma() is the Gamma function.  Note that
1324   //
1325   //     Gamma(n)  =  (n - 1)!
1326   //
1327   // if n is a positive integer.
1328   //
1329   // For this test,
1330   //
1331   //     Gamma(1 + 1 / shape)  =  Gamma(1 + 1 / 1)
1332   //                           =  Gamma(2)
1333   //                           =  (2 - 1)!
1334   //                           =  1
1335   //
1336   // which means
1337   //
1338   //     E[value]  =  scale  .
1339   //
1340   double expectedMean = scale;
1341 
1342   // Test that values have approximately the right mean value.
1343   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1344 }
1345 
1346 /**
1347  * Test case for log-normal distribution random variable stream generator
1348  */
1349 class LogNormalTestCase : public TestCaseBase
1350 {
1351 public:
1352   // Constructor
1353   LogNormalTestCase ();
1354 
1355   // Inherited
1356   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
1357 
1358 private:
1359   // Inherited
1360   virtual void DoRun (void);
1361 
1362   /**
1363    * Tolerance for testing rng values against expectation,
1364    * as a fraction of mean value.
1365    */
1366   static constexpr double TOLERANCE {3e-2};
1367 };
1368 
LogNormalTestCase()1369 LogNormalTestCase::LogNormalTestCase ()
1370   : TestCaseBase ("Log-Normal Random Variable Stream Generator")
1371 {}
1372 
1373 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const1374 LogNormalTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
1375 {
1376   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
1377   auto range = UniformHistogramBins (h, 0, 10, false);
1378 
1379   std::vector<double> expected (N_BINS);
1380 
1381   // Note that this assumes that n has mu equal to zero and sigma
1382   // equal to one, which are their default values for this
1383   // distribution.
1384   double mu = 0.0;
1385   double sigma = 1.0;
1386 
1387   for (std::size_t i = 0; i < N_BINS; ++i)
1388     {
1389       expected[i] = gsl_cdf_lognormal_P (range[i + 1], mu, sigma) - gsl_cdf_lognormal_P (range[i], mu, sigma);
1390       expected[i] *= N_MEASUREMENTS;
1391     }
1392 
1393   double chiSquared = ChiSquared (h, expected, rng);
1394 
1395   gsl_histogram_free (h);
1396   return chiSquared;
1397 }
1398 
1399 void
DoRun(void)1400 LogNormalTestCase::DoRun (void)
1401 {
1402   NS_LOG_FUNCTION (this);
1403   SetTestSuiteSeed ();
1404 
1405   auto generator = RngGenerator<LogNormalRandomVariable> ();
1406   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1407   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1408 
1409   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1410 
1411   double mu = 5.0;
1412   double sigma = 2.0;
1413 
1414   // Create the RNG with the specified range.
1415   Ptr<LogNormalRandomVariable> x = CreateObject<LogNormalRandomVariable> ();
1416   x->SetAttribute ("Mu", DoubleValue (mu));
1417   x->SetAttribute ("Sigma", DoubleValue (sigma));
1418 
1419   // Calculate the mean of these values.
1420   double valueMean = Average (x);
1421 
1422   // The expected value for the mean of the values returned by a
1423   // log-normally distributed random variable is equal to
1424   //
1425   //                             2
1426   //                   mu + sigma  / 2
1427   //     E[value]  =  e                 .
1428   //
1429   double expectedMean = std::exp (mu + sigma * sigma / 2.0);
1430 
1431   // Test that values have approximately the right mean value.
1432   //
1433   /**
1434    * \todo This test fails sometimes if the required tolerance is less
1435    * than 3%, which may be because there is a bug in the
1436    * implementation or that the mean of this distribution is more
1437    * sensitive to its parameters than the others are.
1438    */
1439   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1440 }
1441 
1442 /**
1443  * Test case for antithetic log-normal distribution random variable stream generator
1444  */
1445 class LogNormalAntitheticTestCase : public TestCaseBase
1446 {
1447 public:
1448   // Constructor
1449   LogNormalAntitheticTestCase ();
1450 
1451   // Inherited
1452   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
1453 
1454 private:
1455   // Inherited
1456   virtual void DoRun (void);
1457 
1458   /**
1459    * Tolerance for testing rng values against expectation,
1460    * as a fraction of mean value.
1461    */
1462   static constexpr double TOLERANCE {3e-2};
1463 };
1464 
LogNormalAntitheticTestCase()1465 LogNormalAntitheticTestCase::LogNormalAntitheticTestCase ()
1466   : TestCaseBase ("Antithetic Log-Normal Random Variable Stream Generator")
1467 {}
1468 
1469 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const1470 LogNormalAntitheticTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
1471 {
1472   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
1473   auto range = UniformHistogramBins (h, 0, 10, false);
1474 
1475   std::vector<double> expected (N_BINS);
1476 
1477   // Note that this assumes that n has mu equal to zero and sigma
1478   // equal to one, which are their default values for this
1479   // distribution.
1480   double mu = 0.0;
1481   double sigma = 1.0;
1482 
1483   for (std::size_t i = 0; i < N_BINS; ++i)
1484     {
1485       expected[i] = gsl_cdf_lognormal_P (range[i + 1], mu, sigma) - gsl_cdf_lognormal_P (range[i], mu, sigma);
1486       expected[i] *= N_MEASUREMENTS;
1487     }
1488 
1489   double chiSquared = ChiSquared (h, expected, rng);
1490 
1491   gsl_histogram_free (h);
1492   return chiSquared;
1493 }
1494 
1495 void
DoRun(void)1496 LogNormalAntitheticTestCase::DoRun (void)
1497 {
1498   NS_LOG_FUNCTION (this);
1499   SetTestSuiteSeed ();
1500 
1501   auto generator = RngGenerator<LogNormalRandomVariable> (true);
1502   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1503   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1504   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1505 
1506   double mu = 5.0;
1507   double sigma = 2.0;
1508 
1509   // Create the RNG with the specified range.
1510   Ptr<LogNormalRandomVariable> x = CreateObject<LogNormalRandomVariable> ();
1511   x->SetAttribute ("Mu", DoubleValue (mu));
1512   x->SetAttribute ("Sigma", DoubleValue (sigma));
1513 
1514   // Make this generate antithetic values.
1515   x->SetAttribute ("Antithetic", BooleanValue (true));
1516 
1517   // Calculate the mean of these values.
1518   double valueMean = Average (x);
1519 
1520   // The expected value for the mean of the values returned by a
1521   // log-normally distributed random variable is equal to
1522   //
1523   //                             2
1524   //                   mu + sigma  / 2
1525   //     E[value]  =  e                 .
1526   //
1527   double expectedMean = std::exp (mu + sigma * sigma / 2.0);
1528 
1529   // Test that values have approximately the right mean value.
1530   //
1531   /**
1532    * \todo This test fails sometimes if the required tolerance is less
1533    * than 3%, which may be because there is a bug in the
1534    * implementation or that the mean of this distribution is more
1535    * sensitive to its parameters than the others are.
1536    */
1537   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1538 }
1539 
1540 /**
1541  * Test case for gamma distribution random variable stream generator
1542  */
1543 class GammaTestCase : public TestCaseBase
1544 {
1545 public:
1546   // Constructor
1547   GammaTestCase ();
1548 
1549   // Inherited
1550   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
1551 
1552 private:
1553   // Inherited
1554   virtual void DoRun (void);
1555 
1556   /**
1557    * Tolerance for testing rng values against expectation,
1558    * as a fraction of mean value.
1559    */
1560   static constexpr double TOLERANCE {1e-2};
1561 };
1562 
GammaTestCase()1563 GammaTestCase::GammaTestCase ()
1564   : TestCaseBase ("Gamma Random Variable Stream Generator")
1565 {}
1566 
1567 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const1568 GammaTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
1569 {
1570   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
1571   auto range = UniformHistogramBins (h, 0, 10, false);
1572 
1573   std::vector<double> expected (N_BINS);
1574 
1575   // Note that this assumes that n has alpha equal to one and beta
1576   // equal to one, which are their default values for this
1577   // distribution.
1578   double alpha = 1.0;
1579   double beta = 1.0;
1580 
1581   for (std::size_t i = 0; i < N_BINS; ++i)
1582     {
1583       expected[i] = gsl_cdf_gamma_P (range[i + 1], alpha, beta) - gsl_cdf_gamma_P (range[i], alpha, beta);
1584       expected[i] *= N_MEASUREMENTS;
1585     }
1586 
1587   double chiSquared = ChiSquared (h, expected, rng);
1588 
1589   gsl_histogram_free (h);
1590   return chiSquared;
1591 }
1592 
1593 void
DoRun(void)1594 GammaTestCase::DoRun (void)
1595 {
1596   NS_LOG_FUNCTION (this);
1597   SetTestSuiteSeed ();
1598 
1599   auto generator = RngGenerator<GammaRandomVariable> ();
1600   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1601   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1602   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1603 
1604   double alpha = 5.0;
1605   double beta = 2.0;
1606 
1607   // Create the RNG with the specified range.
1608   Ptr<GammaRandomVariable> x = CreateObject<GammaRandomVariable> ();
1609   x->SetAttribute ("Alpha", DoubleValue (alpha));
1610   x->SetAttribute ("Beta", DoubleValue (beta));
1611 
1612   // Calculate the mean of these values.
1613   double valueMean = Average (x);
1614 
1615   // The expected value for the mean of the values returned by a
1616   // gammaly distributed random variable is equal to
1617   //
1618   //     E[value]  =  alpha * beta  .
1619   //
1620   double expectedMean = alpha * beta;
1621 
1622   // Test that values have approximately the right mean value.
1623   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1624 }
1625 
1626 /**
1627  * Test case for antithetic gamma distribution random variable stream generator
1628  */
1629 class GammaAntitheticTestCase : public TestCaseBase
1630 {
1631 public:
1632   // Constructor
1633   GammaAntitheticTestCase ();
1634 
1635   // Inherited
1636   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
1637 
1638 private:
1639   // Inherited
1640   virtual void DoRun (void);
1641 
1642   /**
1643    * Tolerance for testing rng values against expectation,
1644    * as a fraction of mean value.
1645    */
1646   static constexpr double TOLERANCE {1e-2};
1647 };
1648 
GammaAntitheticTestCase()1649 GammaAntitheticTestCase::GammaAntitheticTestCase ()
1650   : TestCaseBase ("Antithetic Gamma Random Variable Stream Generator")
1651 {}
1652 
1653 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const1654 GammaAntitheticTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
1655 {
1656   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
1657   auto range = UniformHistogramBins (h, 0, 10, false);
1658 
1659   std::vector<double> expected (N_BINS);
1660 
1661   // Note that this assumes that n has alpha equal to one and beta
1662   // equal to one, which are their default values for this
1663   // distribution.
1664   double alpha = 1.0;
1665   double beta = 1.0;
1666 
1667   for (std::size_t i = 0; i < N_BINS; ++i)
1668     {
1669       expected[i] = gsl_cdf_gamma_P (range[i + 1], alpha, beta) - gsl_cdf_gamma_P (range[i], alpha, beta);
1670       expected[i] *= N_MEASUREMENTS;
1671     }
1672 
1673   double chiSquared = ChiSquared (h, expected, rng);
1674 
1675   gsl_histogram_free (h);
1676   return chiSquared;
1677 }
1678 
1679 void
DoRun(void)1680 GammaAntitheticTestCase::DoRun (void)
1681 {
1682   NS_LOG_FUNCTION (this);
1683   SetTestSuiteSeed ();
1684 
1685   auto generator = RngGenerator<GammaRandomVariable> (true);
1686   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1687   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1688   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1689 
1690   double alpha = 5.0;
1691   double beta = 2.0;
1692 
1693   // Create the RNG with the specified range.
1694   Ptr<GammaRandomVariable> x = CreateObject<GammaRandomVariable> ();
1695 
1696   // Make this generate antithetic values.
1697   x->SetAttribute ("Antithetic", BooleanValue (true));
1698 
1699   x->SetAttribute ("Alpha", DoubleValue (alpha));
1700   x->SetAttribute ("Beta", DoubleValue (beta));
1701 
1702   // Calculate the mean of these values.
1703   double valueMean = Average (x);
1704 
1705   // The expected value for the mean of the values returned by a
1706   // gammaly distributed random variable is equal to
1707   //
1708   //     E[value]  =  alpha * beta  .
1709   //
1710   double expectedMean = alpha * beta;
1711 
1712   // Test that values have approximately the right mean value.
1713   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1714 }
1715 
1716 /**
1717  * Test case for Erlang distribution random variable stream generator
1718  */
1719 class ErlangTestCase : public TestCaseBase
1720 {
1721 public:
1722   // Constructor
1723   ErlangTestCase ();
1724 
1725   // Inherited
1726   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
1727 
1728 private:
1729   // Inherited
1730   virtual void DoRun (void);
1731 
1732   /**
1733    * Tolerance for testing rng values against expectation,
1734    * as a fraction of mean value.
1735    */
1736   static constexpr double TOLERANCE {1e-2};
1737 };
1738 
ErlangTestCase()1739 ErlangTestCase::ErlangTestCase ()
1740   : TestCaseBase ("Erlang Random Variable Stream Generator")
1741 {}
1742 
1743 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const1744 ErlangTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
1745 {
1746   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
1747   auto range = UniformHistogramBins (h, 0, 10, false);
1748 
1749   std::vector<double> expected (N_BINS);
1750 
1751   // Note that this assumes that n has k equal to one and lambda
1752   // equal to one, which are their default values for this
1753   // distribution.
1754   uint32_t k = 1;
1755   double lambda = 1.0;
1756 
1757   // Note that Erlang distribution is equal to the gamma distribution
1758   // when k is an iteger, which is why the gamma distribution's cdf
1759   // function can be used here.
1760   for (std::size_t i = 0; i < N_BINS; ++i)
1761     {
1762       expected[i] = gsl_cdf_gamma_P (range[i + 1], k, lambda) - gsl_cdf_gamma_P (range[i], k, lambda);
1763       expected[i] *= N_MEASUREMENTS;
1764     }
1765 
1766   double chiSquared = ChiSquared (h, expected, rng);
1767 
1768   gsl_histogram_free (h);
1769   return chiSquared;
1770 }
1771 
1772 void
DoRun(void)1773 ErlangTestCase::DoRun (void)
1774 {
1775   NS_LOG_FUNCTION (this);
1776   SetTestSuiteSeed ();
1777 
1778   auto generator = RngGenerator<ErlangRandomVariable> ();
1779   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1780   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1781   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1782 
1783   uint32_t k = 5;
1784   double lambda = 2.0;
1785 
1786   // Create the RNG with the specified range.
1787   Ptr<ErlangRandomVariable> x = CreateObject<ErlangRandomVariable> ();
1788   x->SetAttribute ("K", IntegerValue (k));
1789   x->SetAttribute ("Lambda", DoubleValue (lambda));
1790 
1791   // Calculate the mean of these values.
1792   double valueMean = Average (x);
1793 
1794   // The expected value for the mean of the values returned by a
1795   // Erlangly distributed random variable is equal to
1796   //
1797   //     E[value]  =  k * lambda  .
1798   //
1799   double expectedMean = k * lambda;
1800 
1801   // Test that values have approximately the right mean value.
1802   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1803 }
1804 
1805 /**
1806  * Test case for antithetic Erlang distribution random variable stream generator
1807  */
1808 class ErlangAntitheticTestCase : public TestCaseBase
1809 {
1810 public:
1811   // Constructor
1812   ErlangAntitheticTestCase ();
1813 
1814   // Inherited
1815   double ChiSquaredTest (Ptr<RandomVariableStream> rng) const;
1816 
1817 private:
1818   // Inherited
1819   virtual void DoRun (void);
1820 
1821   /**
1822    * Tolerance for testing rng values against expectation,
1823    * as a fraction of mean value.
1824    */
1825   static constexpr double TOLERANCE {1e-2};
1826 };
1827 
ErlangAntitheticTestCase()1828 ErlangAntitheticTestCase::ErlangAntitheticTestCase ()
1829   : TestCaseBase ("Antithetic Erlang Random Variable Stream Generator")
1830 {}
1831 
1832 double
ChiSquaredTest(Ptr<RandomVariableStream> rng) const1833 ErlangAntitheticTestCase::ChiSquaredTest (Ptr<RandomVariableStream> rng) const
1834 {
1835   gsl_histogram * h = gsl_histogram_alloc (N_BINS);
1836   auto range = UniformHistogramBins (h, 0, 10, false);
1837 
1838   std::vector<double> expected (N_BINS);
1839 
1840   // Note that this assumes that n has k equal to one and lambda
1841   // equal to one, which are their default values for this
1842   // distribution.
1843   uint32_t k = 1;
1844   double lambda = 1.0;
1845 
1846   // Note that Erlang distribution is equal to the gamma distribution
1847   // when k is an iteger, which is why the gamma distribution's cdf
1848   // function can be used here.
1849   for (std::size_t i = 0; i < N_BINS; ++i)
1850     {
1851       expected[i] = gsl_cdf_gamma_P (range[i + 1], k, lambda) - gsl_cdf_gamma_P (range[i], k, lambda);
1852       expected[i] *= N_MEASUREMENTS;
1853     }
1854 
1855   double chiSquared = ChiSquared (h, expected, rng);
1856 
1857   gsl_histogram_free (h);
1858   return chiSquared;
1859 }
1860 
1861 void
DoRun(void)1862 ErlangAntitheticTestCase::DoRun (void)
1863 {
1864   NS_LOG_FUNCTION (this);
1865   SetTestSuiteSeed ();
1866 
1867   auto generator = RngGenerator<ErlangRandomVariable> (true);
1868   double sum = ChiSquaredsAverage (&generator, N_RUNS);
1869   double maxStatistic = gsl_cdf_chisq_Qinv (0.05, N_BINS);
1870   NS_TEST_ASSERT_MSG_LT (sum, maxStatistic, "Chi-squared statistic out of range");
1871 
1872   uint32_t k = 5;
1873   double lambda = 2.0;
1874 
1875   // Create the RNG with the specified range.
1876   Ptr<ErlangRandomVariable> x = CreateObject<ErlangRandomVariable> ();
1877 
1878   // Make this generate antithetic values.
1879   x->SetAttribute ("Antithetic", BooleanValue (true));
1880 
1881   x->SetAttribute ("K", IntegerValue (k));
1882   x->SetAttribute ("Lambda", DoubleValue (lambda));
1883 
1884   // Calculate the mean of these values.
1885   double valueMean = Average (x);
1886 
1887   // The expected value for the mean of the values returned by a
1888   // Erlangly distributed random variable is equal to
1889   //
1890   //     E[value]  =  k * lambda  .
1891   //
1892   double expectedMean = k * lambda;
1893 
1894   // Test that values have approximately the right mean value.
1895   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1896 }
1897 
1898 /**
1899  * Test case for Zipf distribution random variable stream generator
1900  */
1901 class ZipfTestCase : public TestCaseBase
1902 {
1903 public:
1904   // Constructor
1905   ZipfTestCase ();
1906 
1907 private:
1908   // Inherited
1909   virtual void DoRun (void);
1910 
1911   /**
1912    * Tolerance for testing rng values against expectation,
1913    * as a fraction of mean value.
1914    */
1915   static constexpr double TOLERANCE {1e-2};
1916 };
1917 
ZipfTestCase()1918 ZipfTestCase::ZipfTestCase ()
1919   : TestCaseBase ("Zipf Random Variable Stream Generator")
1920 {}
1921 
1922 void
DoRun(void)1923 ZipfTestCase::DoRun (void)
1924 {
1925   NS_LOG_FUNCTION (this);
1926   SetTestSuiteSeed ();
1927 
1928   uint32_t n = 1;
1929   double alpha = 2.0;
1930 
1931   // Create the RNG with the specified range.
1932   Ptr<ZipfRandomVariable> x = CreateObject<ZipfRandomVariable> ();
1933   x->SetAttribute ("N", IntegerValue (n));
1934   x->SetAttribute ("Alpha", DoubleValue (alpha));
1935 
1936   // Calculate the mean of these values.
1937   double valueMean = Average (x);
1938 
1939   // The expected value for the mean of the values returned by a
1940   // Zipfly distributed random variable is equal to
1941   //
1942   //                   H
1943   //                    N, alpha - 1
1944   //     E[value]  =  ---------------
1945   //                     H
1946   //                      N, alpha
1947   //
1948   // where
1949   //
1950   //                    N
1951   //                   ---
1952   //                   \     -alpha
1953   //     H          =  /    m        .
1954   //      N, alpha     ---
1955   //                   m=1
1956   //
1957   // For this test,
1958   //
1959   //                      -(alpha - 1)
1960   //                     1
1961   //     E[value]  =  ---------------
1962   //                      -alpha
1963   //                     1
1964   //
1965   //               =  1  .
1966   //
1967   double expectedMean = 1.0;
1968 
1969   // Test that values have approximately the right mean value.
1970   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
1971 }
1972 
1973 /**
1974  * Test case for antithetic Zipf distribution random variable stream generator
1975  */
1976 class ZipfAntitheticTestCase : public TestCaseBase
1977 {
1978 public:
1979   // Constructor
1980   ZipfAntitheticTestCase ();
1981 
1982 private:
1983   // Inherited
1984   virtual void DoRun (void);
1985 
1986   /**
1987    * Tolerance for testing rng values against expectation,
1988    * as a fraction of mean value.
1989    */
1990   static constexpr double TOLERANCE {1e-2};
1991 };
1992 
ZipfAntitheticTestCase()1993 ZipfAntitheticTestCase::ZipfAntitheticTestCase ()
1994   : TestCaseBase ("Antithetic Zipf Random Variable Stream Generator")
1995 {}
1996 
1997 void
DoRun(void)1998 ZipfAntitheticTestCase::DoRun (void)
1999 {
2000   NS_LOG_FUNCTION (this);
2001   SetTestSuiteSeed ();
2002 
2003   uint32_t n = 1;
2004   double alpha = 2.0;
2005 
2006   // Create the RNG with the specified range.
2007   Ptr<ZipfRandomVariable> x = CreateObject<ZipfRandomVariable> ();
2008   x->SetAttribute ("N", IntegerValue (n));
2009   x->SetAttribute ("Alpha", DoubleValue (alpha));
2010 
2011   // Make this generate antithetic values.
2012   x->SetAttribute ("Antithetic", BooleanValue (true));
2013 
2014   // Calculate the mean of these values.
2015   double valueMean = Average (x);
2016 
2017   // The expected value for the mean of the values returned by a
2018   // Zipfly distributed random variable is equal to
2019   //
2020   //                   H
2021   //                    N, alpha - 1
2022   //     E[value]  =  ---------------
2023   //                     H
2024   //                      N, alpha
2025   //
2026   // where
2027   //
2028   //                    N
2029   //                   ---
2030   //                   \     -alpha
2031   //     H          =  /    m        .
2032   //      N, alpha     ---
2033   //                   m=1
2034   //
2035   // For this test,
2036   //
2037   //                      -(alpha - 1)
2038   //                     1
2039   //     E[value]  =  ---------------
2040   //                      -alpha
2041   //                     1
2042   //
2043   //               =  1  .
2044   //
2045   double expectedMean = 1.0;
2046 
2047   // Test that values have approximately the right mean value.
2048   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
2049 }
2050 
2051 /**
2052  * Test case for Zeta distribution random variable stream generator
2053  */
2054 class ZetaTestCase : public TestCaseBase
2055 {
2056 public:
2057   // Constructor
2058   ZetaTestCase ();
2059 
2060 private:
2061   // Inherited
2062   virtual void DoRun (void);
2063 
2064   /**
2065    * Tolerance for testing rng values against expectation,
2066    * as a fraction of mean value.
2067    */
2068   static constexpr double TOLERANCE {1e-2};
2069 };
2070 
ZetaTestCase()2071 ZetaTestCase::ZetaTestCase ()
2072   : TestCaseBase ("Zeta Random Variable Stream Generator")
2073 {}
2074 
2075 void
DoRun(void)2076 ZetaTestCase::DoRun (void)
2077 {
2078   NS_LOG_FUNCTION (this);
2079   SetTestSuiteSeed ();
2080 
2081   double alpha = 5.0;
2082 
2083   // Create the RNG with the specified range.
2084   Ptr<ZetaRandomVariable> x = CreateObject<ZetaRandomVariable> ();
2085   x->SetAttribute ("Alpha", DoubleValue (alpha));
2086 
2087   // Calculate the mean of these values.
2088   double valueMean = Average (x);
2089 
2090   // The expected value for the mean of the values returned by a
2091   // zetaly distributed random variable is equal to
2092   //
2093   //                   zeta(alpha - 1)
2094   //     E[value]  =  ---------------   for alpha > 2 ,
2095   //                     zeta(alpha)
2096   //
2097   // where zeta(alpha) is the Riemann zeta function.
2098   //
2099   // There are no simple analytic forms for the Riemann zeta function,
2100   // which is why the gsl library is used in this test to calculate
2101   // the known mean of the values.
2102   double expectedMean =
2103     gsl_sf_zeta_int (static_cast<int> (alpha - 1)) /
2104     gsl_sf_zeta_int (static_cast<int> (alpha)       );
2105 
2106   // Test that values have approximately the right mean value.
2107   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
2108 }
2109 
2110 /**
2111  * Test case for antithetic Zeta distribution random variable stream generator
2112  */
2113 class ZetaAntitheticTestCase : public TestCaseBase
2114 {
2115 public:
2116   // Constructor
2117   ZetaAntitheticTestCase ();
2118 
2119 private:
2120   // Inherited
2121   virtual void DoRun (void);
2122 
2123   /**
2124    * Tolerance for testing rng values against expectation,
2125    * as a fraction of mean value.
2126    */
2127   static constexpr double TOLERANCE {1e-2};
2128 };
2129 
ZetaAntitheticTestCase()2130 ZetaAntitheticTestCase::ZetaAntitheticTestCase ()
2131   : TestCaseBase ("Antithetic Zeta Random Variable Stream Generator")
2132 {}
2133 
2134 void
DoRun(void)2135 ZetaAntitheticTestCase::DoRun (void)
2136 {
2137   NS_LOG_FUNCTION (this);
2138   SetTestSuiteSeed ();
2139 
2140   double alpha = 5.0;
2141 
2142   // Create the RNG with the specified range.
2143   Ptr<ZetaRandomVariable> x = CreateObject<ZetaRandomVariable> ();
2144   x->SetAttribute ("Alpha", DoubleValue (alpha));
2145 
2146   // Make this generate antithetic values.
2147   x->SetAttribute ("Antithetic", BooleanValue (true));
2148 
2149   // Calculate the mean of these values.
2150   double valueMean = Average (x);
2151 
2152   // The expected value for the mean of the values returned by a
2153   // zetaly distributed random variable is equal to
2154   //
2155   //                   zeta(alpha - 1)
2156   //     E[value]  =  ---------------   for alpha > 2 ,
2157   //                     zeta(alpha)
2158   //
2159   // where zeta(alpha) is the Riemann zeta function.
2160   //
2161   // There are no simple analytic forms for the Riemann zeta function,
2162   // which is why the gsl library is used in this test to calculate
2163   // the known mean of the values.
2164   double expectedMean =
2165     gsl_sf_zeta_int (static_cast<int> (alpha) - 1) /
2166     gsl_sf_zeta_int (static_cast<int> (alpha)       );
2167 
2168   // Test that values have approximately the right mean value.
2169   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
2170 }
2171 
2172 /**
2173  * Test case for deterministic random variable stream generator
2174  */
2175 class DeterministicTestCase : public TestCaseBase
2176 {
2177 public:
2178   // Constructor
2179   DeterministicTestCase ();
2180 
2181 private:
2182   // Inherited
2183   virtual void DoRun (void);
2184 
2185   /** Tolerance for testing rng values against expectation. */
2186   static constexpr double TOLERANCE {1e-8};
2187 };
2188 
DeterministicTestCase()2189 DeterministicTestCase::DeterministicTestCase ()
2190   : TestCaseBase ("Deterministic Random Variable Stream Generator")
2191 {}
2192 
2193 void
DoRun(void)2194 DeterministicTestCase::DoRun (void)
2195 {
2196   NS_LOG_FUNCTION (this);
2197   SetTestSuiteSeed ();
2198 
2199   Ptr<DeterministicRandomVariable> s = CreateObject<DeterministicRandomVariable> ();
2200 
2201   // The following array should give the sequence
2202   //
2203   //    4, 4, 7, 7, 10, 10 .
2204   //
2205   double array1 [] = { 4, 4, 7, 7, 10, 10};
2206   std::size_t count1 = 6;
2207   s->SetValueArray (array1, count1);
2208 
2209   double value;
2210 
2211   // Test that the first sequence is correct.
2212   value = s->GetValue ();
2213   NS_TEST_ASSERT_MSG_EQ_TOL (value, 4, TOLERANCE, "Sequence 1 value 1 wrong.");
2214   value = s->GetValue ();
2215   NS_TEST_ASSERT_MSG_EQ_TOL (value, 4, TOLERANCE, "Sequence 1 value 2 wrong.");
2216   value = s->GetValue ();
2217   NS_TEST_ASSERT_MSG_EQ_TOL (value, 7, TOLERANCE, "Sequence 1 value 3 wrong.");
2218   value = s->GetValue ();
2219   NS_TEST_ASSERT_MSG_EQ_TOL (value, 7, TOLERANCE, "Sequence 1 value 4 wrong.");
2220   value = s->GetValue ();
2221   NS_TEST_ASSERT_MSG_EQ_TOL (value, 10, TOLERANCE, "Sequence 1 value 5 wrong.");
2222   value = s->GetValue ();
2223   NS_TEST_ASSERT_MSG_EQ_TOL (value, 10, TOLERANCE, "Sequence 1 value 6 wrong.");
2224 
2225   // The following array should give the sequence
2226   //
2227   //    1000, 2000, 7, 7 .
2228   //
2229   double array2 [] = { 1000, 2000, 3000, 4000};
2230   std::size_t count2 = 4;
2231   s->SetValueArray (array2, count2);
2232 
2233   // Test that the second sequence is correct.
2234   value = s->GetValue ();
2235   NS_TEST_ASSERT_MSG_EQ_TOL (value, 1000, TOLERANCE, "Sequence 2 value 1 wrong.");
2236   value = s->GetValue ();
2237   NS_TEST_ASSERT_MSG_EQ_TOL (value, 2000, TOLERANCE, "Sequence 2 value 2 wrong.");
2238   value = s->GetValue ();
2239   NS_TEST_ASSERT_MSG_EQ_TOL (value, 3000, TOLERANCE, "Sequence 2 value 3 wrong.");
2240   value = s->GetValue ();
2241   NS_TEST_ASSERT_MSG_EQ_TOL (value, 4000, TOLERANCE, "Sequence 2 value 4 wrong.");
2242   value = s->GetValue ();
2243 }
2244 
2245 /**
2246  * Test case for empirical distribution random variable stream generator
2247  */
2248 class EmpiricalTestCase : public TestCaseBase
2249 {
2250 public:
2251   // Constructor
2252   EmpiricalTestCase ();
2253 
2254 private:
2255   // Inherited
2256   virtual void DoRun (void);
2257 
2258   /**
2259    * Tolerance for testing rng values against expectation,
2260    * as a fraction of mean value.
2261    */
2262   static constexpr double TOLERANCE {1e-2};
2263 };
2264 
EmpiricalTestCase()2265 EmpiricalTestCase::EmpiricalTestCase ()
2266   : TestCaseBase ("Empirical Random Variable Stream Generator")
2267 {}
2268 
2269 void
DoRun(void)2270 EmpiricalTestCase::DoRun (void)
2271 {
2272   NS_LOG_FUNCTION (this);
2273   SetTestSuiteSeed ();
2274 
2275   // Create the RNG with a uniform distribution between 0 and 10.
2276   Ptr<EmpiricalRandomVariable> x = CreateObject<EmpiricalRandomVariable> ();
2277   x->SetInterpolate (false);
2278   x->CDF ( 0.0,  0.0);
2279   x->CDF ( 5.0,  0.25);
2280   x->CDF (10.0,  1.0);
2281 
2282   // Check that only the correct values are returned
2283   for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
2284     {
2285       double value = x->GetValue ();
2286       NS_TEST_EXPECT_MSG_EQ ( (value == 5) || (value == 10), true,
2287                               "Incorrect value returned, expected only 5 or 10.");
2288     }
2289 
2290   // Calculate the mean of the sampled values.
2291   double valueMean = Average (x);
2292 
2293   // The expected distribution with sampled values is
2294   //     Value     Probability
2295   //      5        25%
2296   //     10        75%
2297   //
2298   // The expected mean is
2299   //
2300   //     E[value]  =  5 * 25%  +  10 * 75%  =  8.75
2301   //
2302   // Test that values have approximately the right mean value.
2303   double expectedMean = 8.75;
2304   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
2305 
2306 
2307   // Calculate the mean of the interpolated values.
2308   x->SetInterpolate (true);
2309   valueMean = Average (x);
2310 
2311   // The expected distribution (with interpolation) is
2312   //     Bin     Probability
2313   //     [0, 5)     25%
2314   //     [5, 10)    75%
2315   //
2316   // Each bin is uniformly sampled, so the average of the samples in the
2317   // bin is the center of the bin.
2318   //
2319   // The expected mean is
2320   //
2321   //     E[value]  =  2.5 * 25% + 7.5 * 75% = 6.25
2322   //
2323   expectedMean = 6.25;
2324 
2325   // Test that values have approximately the right mean value.
2326   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
2327 
2328   // Bug 2082: Create the RNG with a uniform distribution between -1 and 1.
2329   Ptr<EmpiricalRandomVariable> y = CreateObject<EmpiricalRandomVariable> ();
2330   y->SetInterpolate (false);
2331   y->CDF (-1.0,  0.0);
2332   y->CDF (0.0,  0.5);
2333   y->CDF (1.0,  1.0);
2334   NS_TEST_ASSERT_MSG_LT (y->GetValue (), 2, "Empirical variable with negative domain");
2335 }
2336 
2337 /**
2338  * Test case for antithetic empirical distribution random variable stream generator
2339  */
2340 class EmpiricalAntitheticTestCase : public TestCaseBase
2341 {
2342 public:
2343   // Constructor
2344   EmpiricalAntitheticTestCase ();
2345 
2346 private:
2347   // Inherited
2348   virtual void DoRun (void);
2349 
2350   /**
2351    * Tolerance for testing rng values against expectation,
2352    * as a fraction of mean value.
2353    */
2354   static constexpr double TOLERANCE {1e-2};
2355 };
2356 
EmpiricalAntitheticTestCase()2357 EmpiricalAntitheticTestCase::EmpiricalAntitheticTestCase ()
2358   : TestCaseBase ("EmpiricalAntithetic Random Variable Stream Generator")
2359 {}
2360 
2361 void
DoRun(void)2362 EmpiricalAntitheticTestCase::DoRun (void)
2363 {
2364   NS_LOG_FUNCTION (this);
2365   SetTestSuiteSeed ();
2366 
2367   // Create the RNG with a uniform distribution between 0 and 10.
2368   Ptr<EmpiricalRandomVariable> x = CreateObject<EmpiricalRandomVariable> ();
2369   x->SetInterpolate (false);
2370   x->CDF ( 0.0,  0.0);
2371   x->CDF ( 5.0,  0.25);
2372   x->CDF (10.0,  1.0);
2373 
2374   // Make this generate antithetic values.
2375   x->SetAttribute ("Antithetic", BooleanValue (true));
2376 
2377   // Check that only the correct values are returned
2378   for (uint32_t i = 0; i < N_MEASUREMENTS; ++i)
2379     {
2380       double value = x->GetValue ();
2381       NS_TEST_EXPECT_MSG_EQ ( (value == 5) || (value == 10), true,
2382                               "Incorrect value returned, expected only 5 or 10.");
2383     }
2384 
2385   // Calculate the mean of these values.
2386   double valueMean = Average (x);
2387   // Expected
2388   //    E[value] = 5 * 25%  + 10 * 75%  = 8.75
2389   double expectedMean = 8.75;
2390   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
2391 
2392   // Check interpolated sampling
2393   x->SetInterpolate (true);
2394   valueMean = Average (x);
2395 
2396   // The expected value for the mean of the values returned by this
2397   // empirical distribution with interpolation is
2398   //
2399   //     E[value]  =  2.5 * 25% + 7.5 * 75% = 6.25
2400   //
2401   expectedMean = 6.25;
2402 
2403   // Test that values have approximately the right mean value.
2404   NS_TEST_ASSERT_MSG_EQ_TOL (valueMean, expectedMean, expectedMean * TOLERANCE, "Wrong mean value.");
2405 }
2406 
2407 /**
2408  * Test case for caching of Normal RV parameters (see issue #302)
2409  */
2410 class NormalCachingTestCase : public TestCaseBase
2411 {
2412 public:
2413   // Constructor
2414   NormalCachingTestCase ();
2415 
2416 private:
2417   // Inherited
2418   virtual void DoRun (void);
2419 };
2420 
NormalCachingTestCase()2421 NormalCachingTestCase::NormalCachingTestCase ()
2422   : TestCaseBase ("NormalRandomVariable caching of parameters")
2423 {}
2424 
2425 void
DoRun(void)2426 NormalCachingTestCase::DoRun (void)
2427 {
2428   NS_LOG_FUNCTION (this);
2429   SetTestSuiteSeed ();
2430 
2431   Ptr<NormalRandomVariable> n = CreateObject<NormalRandomVariable> ();
2432   double v1 = n->GetValue (-10, 1, 10);  // Mean -10, variance 1, bounded to [-20,0]
2433   double v2 = n->GetValue (10, 1, 10);   // Mean 10, variance 1, bounded to [0,20]
2434 
2435   NS_TEST_ASSERT_MSG_LT (v1, 0, "Incorrect value returned, expected < 0");
2436   NS_TEST_ASSERT_MSG_GT (v2, 0, "Incorrect value returned, expected > 0");
2437 }
2438 
2439 /**
2440  * RandomVariableStream test suite, covering all random number variable
2441  * stream generator types.
2442  */
2443 class RandomVariableSuite : public TestSuite
2444 {
2445 public:
2446   // Constructor
2447   RandomVariableSuite ();
2448 };
2449 
RandomVariableSuite()2450 RandomVariableSuite::RandomVariableSuite ()
2451   : TestSuite ("random-variable-stream-generators", UNIT)
2452 {
2453   AddTestCase (new UniformTestCase);
2454   AddTestCase (new UniformAntitheticTestCase);
2455   AddTestCase (new ConstantTestCase);
2456   AddTestCase (new SequentialTestCase);
2457   AddTestCase (new NormalTestCase);
2458   AddTestCase (new NormalAntitheticTestCase);
2459   AddTestCase (new ExponentialTestCase);
2460   AddTestCase (new ExponentialAntitheticTestCase);
2461   AddTestCase (new ParetoTestCase);
2462   AddTestCase (new ParetoAntitheticTestCase);
2463   AddTestCase (new WeibullTestCase);
2464   AddTestCase (new WeibullAntitheticTestCase);
2465   AddTestCase (new LogNormalTestCase);
2466   /// \todo This test is currently disabled because it fails sometimes.
2467   /// A possible reason for the failure is that the antithetic code is
2468   /// not implemented properly for this log-normal case.
2469   /*
2470   AddTestCase (new LogNormalAntitheticTestCase);
2471   */
2472   AddTestCase (new GammaTestCase);
2473   /// \todo This test is currently disabled because it fails sometimes.
2474   /// A possible reason for the failure is that the antithetic code is
2475   /// not implemented properly for this gamma case.
2476   /*
2477   AddTestCase (new GammaAntitheticTestCase);
2478   */
2479   AddTestCase (new ErlangTestCase);
2480   AddTestCase (new ErlangAntitheticTestCase);
2481   AddTestCase (new ZipfTestCase);
2482   AddTestCase (new ZipfAntitheticTestCase);
2483   AddTestCase (new ZetaTestCase);
2484   AddTestCase (new ZetaAntitheticTestCase);
2485   AddTestCase (new DeterministicTestCase);
2486   AddTestCase (new EmpiricalTestCase);
2487   AddTestCase (new EmpiricalAntitheticTestCase);
2488   /// Issue #302:  NormalRandomVariable produces stale values
2489   AddTestCase (new NormalCachingTestCase);
2490 }
2491 
2492 static RandomVariableSuite randomVariableSuite;
2493 
2494 }  // namespace RandomVariable
2495 
2496 }  // namespace test
2497 
2498 }  // namespace ns3
2499 
2500