1 /* -*- Mode:C++; c-file-style:"gnu"; indent-tabs-mode:nil; -*- */
2 /*
3 * Copyright (c) 2006 Georgia Tech Research Corporation
4 * Copyright (c) 2011 Mathieu Lacage
5 *
6 * This program is free software; you can redistribute it and/or modify
7 * it under the terms of the GNU General Public License version 2 as
8 * published by the Free Software Foundation;
9 *
10 * This program is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
13 * GNU General Public License for more details.
14 *
15 * You should have received a copy of the GNU General Public License
16 * along with this program; if not, write to the Free Software
17 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
18 *
19 * Authors: Rajib Bhattacharjea<raj.b@gatech.edu>
20 * Hadi Arbabi<marbabi@cs.odu.edu>
21 * Mathieu Lacage <mathieu.lacage@gmail.com>
22 *
23 * Modified by Mitch Watrous <watrous@u.washington.edu>
24 *
25 */
26 #include "random-variable-stream.h"
27 #include "assert.h"
28 #include "boolean.h"
29 #include "double.h"
30 #include "integer.h"
31 #include "string.h"
32 #include "pointer.h"
33 #include "log.h"
34 #include "rng-stream.h"
35 #include "rng-seed-manager.h"
36 #include "unused.h"
37 #include <cmath>
38 #include <iostream>
39 #include <algorithm> // upper_bound
40
41 /**
42 * \file
43 * \ingroup randomvariable
44 * ns3::RandomVariableStream and related implementations
45 */
46
47 namespace ns3 {
48
49 NS_LOG_COMPONENT_DEFINE ("RandomVariableStream");
50
51 NS_OBJECT_ENSURE_REGISTERED (RandomVariableStream);
52
53 TypeId
GetTypeId(void)54 RandomVariableStream::GetTypeId (void)
55 {
56 static TypeId tid = TypeId ("ns3::RandomVariableStream")
57 .SetParent<Object> ()
58 .SetGroupName ("Core")
59 .AddAttribute ("Stream",
60 "The stream number for this RNG stream. -1 means \"allocate a stream automatically\". "
61 "Note that if -1 is set, Get will return -1 so that it is not possible to know which "
62 "value was automatically allocated.",
63 IntegerValue (-1),
64 MakeIntegerAccessor (&RandomVariableStream::SetStream,
65 &RandomVariableStream::GetStream),
66 MakeIntegerChecker<int64_t>())
67 .AddAttribute ("Antithetic", "Set this RNG stream to generate antithetic values",
68 BooleanValue (false),
69 MakeBooleanAccessor (&RandomVariableStream::SetAntithetic,
70 &RandomVariableStream::IsAntithetic),
71 MakeBooleanChecker ())
72 ;
73 return tid;
74 }
75
RandomVariableStream()76 RandomVariableStream::RandomVariableStream ()
77 : m_rng (0)
78 {
79 NS_LOG_FUNCTION (this);
80 }
~RandomVariableStream()81 RandomVariableStream::~RandomVariableStream ()
82 {
83 NS_LOG_FUNCTION (this);
84 delete m_rng;
85 }
86
87 void
SetAntithetic(bool isAntithetic)88 RandomVariableStream::SetAntithetic (bool isAntithetic)
89 {
90 NS_LOG_FUNCTION (this << isAntithetic);
91 m_isAntithetic = isAntithetic;
92 }
93 bool
IsAntithetic(void) const94 RandomVariableStream::IsAntithetic (void) const
95 {
96 NS_LOG_FUNCTION (this);
97 return m_isAntithetic;
98 }
99 void
SetStream(int64_t stream)100 RandomVariableStream::SetStream (int64_t stream)
101 {
102 NS_LOG_FUNCTION (this << stream);
103 // negative values are not legal.
104 NS_ASSERT (stream >= -1);
105 delete m_rng;
106 if (stream == -1)
107 {
108 // The first 2^63 streams are reserved for automatic stream
109 // number assignment.
110 uint64_t nextStream = RngSeedManager::GetNextStreamIndex ();
111 NS_ASSERT (nextStream <= ((1ULL) << 63));
112 m_rng = new RngStream (RngSeedManager::GetSeed (),
113 nextStream,
114 RngSeedManager::GetRun ());
115 }
116 else
117 {
118 // The last 2^63 streams are reserved for deterministic stream
119 // number assignment.
120 uint64_t base = ((1ULL) << 63);
121 uint64_t target = base + stream;
122 m_rng = new RngStream (RngSeedManager::GetSeed (),
123 target,
124 RngSeedManager::GetRun ());
125 }
126 m_stream = stream;
127 }
128 int64_t
GetStream(void) const129 RandomVariableStream::GetStream (void) const
130 {
131 NS_LOG_FUNCTION (this);
132 return m_stream;
133 }
134
135 RngStream *
Peek(void) const136 RandomVariableStream::Peek (void) const
137 {
138 NS_LOG_FUNCTION (this);
139 return m_rng;
140 }
141
142 NS_OBJECT_ENSURE_REGISTERED (UniformRandomVariable);
143
144 TypeId
GetTypeId(void)145 UniformRandomVariable::GetTypeId (void)
146 {
147 static TypeId tid = TypeId ("ns3::UniformRandomVariable")
148 .SetParent<RandomVariableStream>()
149 .SetGroupName ("Core")
150 .AddConstructor<UniformRandomVariable> ()
151 .AddAttribute ("Min", "The lower bound on the values returned by this RNG stream.",
152 DoubleValue (0),
153 MakeDoubleAccessor (&UniformRandomVariable::m_min),
154 MakeDoubleChecker<double>())
155 .AddAttribute ("Max", "The upper bound on the values returned by this RNG stream.",
156 DoubleValue (1.0),
157 MakeDoubleAccessor (&UniformRandomVariable::m_max),
158 MakeDoubleChecker<double>())
159 ;
160 return tid;
161 }
UniformRandomVariable()162 UniformRandomVariable::UniformRandomVariable ()
163 {
164 // m_min and m_max are initialized after constructor by attributes
165 NS_LOG_FUNCTION (this);
166 }
167
168 double
GetMin(void) const169 UniformRandomVariable::GetMin (void) const
170 {
171 NS_LOG_FUNCTION (this);
172 return m_min;
173 }
174 double
GetMax(void) const175 UniformRandomVariable::GetMax (void) const
176 {
177 NS_LOG_FUNCTION (this);
178 return m_max;
179 }
180
181 double
GetValue(double min,double max)182 UniformRandomVariable::GetValue (double min, double max)
183 {
184 NS_LOG_FUNCTION (this << min << max);
185 double v = min + Peek ()->RandU01 () * (max - min);
186 if (IsAntithetic ())
187 {
188 v = min + (max - v);
189 }
190 return v;
191 }
192 uint32_t
GetInteger(uint32_t min,uint32_t max)193 UniformRandomVariable::GetInteger (uint32_t min, uint32_t max)
194 {
195 NS_LOG_FUNCTION (this << min << max);
196 NS_ASSERT (min <= max);
197 return static_cast<uint32_t> ( GetValue ((double) (min), (double) (max) + 1.0) );
198 }
199
200 double
GetValue(void)201 UniformRandomVariable::GetValue (void)
202 {
203 NS_LOG_FUNCTION (this);
204 return GetValue (m_min, m_max);
205 }
206 uint32_t
GetInteger(void)207 UniformRandomVariable::GetInteger (void)
208 {
209 NS_LOG_FUNCTION (this);
210 return (uint32_t)GetValue (m_min, m_max + 1);
211 }
212
213 NS_OBJECT_ENSURE_REGISTERED (ConstantRandomVariable);
214
215 TypeId
GetTypeId(void)216 ConstantRandomVariable::GetTypeId (void)
217 {
218 static TypeId tid = TypeId ("ns3::ConstantRandomVariable")
219 .SetParent<RandomVariableStream>()
220 .SetGroupName ("Core")
221 .AddConstructor<ConstantRandomVariable> ()
222 .AddAttribute ("Constant", "The constant value returned by this RNG stream.",
223 DoubleValue (0),
224 MakeDoubleAccessor (&ConstantRandomVariable::m_constant),
225 MakeDoubleChecker<double>())
226 ;
227 return tid;
228 }
ConstantRandomVariable()229 ConstantRandomVariable::ConstantRandomVariable ()
230 {
231 // m_constant is initialized after constructor by attributes
232 NS_LOG_FUNCTION (this);
233 }
234
235 double
GetConstant(void) const236 ConstantRandomVariable::GetConstant (void) const
237 {
238 NS_LOG_FUNCTION (this);
239 return m_constant;
240 }
241
242 double
GetValue(double constant)243 ConstantRandomVariable::GetValue (double constant)
244 {
245 NS_LOG_FUNCTION (this << constant);
246 return constant;
247 }
248 uint32_t
GetInteger(uint32_t constant)249 ConstantRandomVariable::GetInteger (uint32_t constant)
250 {
251 NS_LOG_FUNCTION (this << constant);
252 return constant;
253 }
254
255 double
GetValue(void)256 ConstantRandomVariable::GetValue (void)
257 {
258 NS_LOG_FUNCTION (this);
259 return GetValue (m_constant);
260 }
261 uint32_t
GetInteger(void)262 ConstantRandomVariable::GetInteger (void)
263 {
264 NS_LOG_FUNCTION (this);
265 return (uint32_t)GetValue (m_constant);
266 }
267
268 NS_OBJECT_ENSURE_REGISTERED (SequentialRandomVariable);
269
270 TypeId
GetTypeId(void)271 SequentialRandomVariable::GetTypeId (void)
272 {
273 static TypeId tid = TypeId ("ns3::SequentialRandomVariable")
274 .SetParent<RandomVariableStream>()
275 .SetGroupName ("Core")
276 .AddConstructor<SequentialRandomVariable> ()
277 .AddAttribute ("Min", "The first value of the sequence.",
278 DoubleValue (0),
279 MakeDoubleAccessor (&SequentialRandomVariable::m_min),
280 MakeDoubleChecker<double>())
281 .AddAttribute ("Max", "One more than the last value of the sequence.",
282 DoubleValue (0),
283 MakeDoubleAccessor (&SequentialRandomVariable::m_max),
284 MakeDoubleChecker<double>())
285 .AddAttribute ("Increment", "The sequence random variable increment.",
286 StringValue ("ns3::ConstantRandomVariable[Constant=1]"),
287 MakePointerAccessor (&SequentialRandomVariable::m_increment),
288 MakePointerChecker<RandomVariableStream> ())
289 .AddAttribute ("Consecutive", "The number of times each member of the sequence is repeated.",
290 IntegerValue (1),
291 MakeIntegerAccessor (&SequentialRandomVariable::m_consecutive),
292 MakeIntegerChecker<uint32_t>());
293 return tid;
294 }
SequentialRandomVariable()295 SequentialRandomVariable::SequentialRandomVariable ()
296 :
297 m_current (0),
298 m_currentConsecutive (0),
299 m_isCurrentSet (false)
300 {
301 // m_min, m_max, m_increment, and m_consecutive are initialized
302 // after constructor by attributes.
303 NS_LOG_FUNCTION (this);
304 }
305
306 double
GetMin(void) const307 SequentialRandomVariable::GetMin (void) const
308 {
309 NS_LOG_FUNCTION (this);
310 return m_min;
311 }
312
313 double
GetMax(void) const314 SequentialRandomVariable::GetMax (void) const
315 {
316 NS_LOG_FUNCTION (this);
317 return m_max;
318 }
319
320 Ptr<RandomVariableStream>
GetIncrement(void) const321 SequentialRandomVariable::GetIncrement (void) const
322 {
323 NS_LOG_FUNCTION (this);
324 return m_increment;
325 }
326
327 uint32_t
GetConsecutive(void) const328 SequentialRandomVariable::GetConsecutive (void) const
329 {
330 NS_LOG_FUNCTION (this);
331 return m_consecutive;
332 }
333
334 double
GetValue(void)335 SequentialRandomVariable::GetValue (void)
336 {
337 // Set the current sequence value if it hasn't been set.
338 NS_LOG_FUNCTION (this);
339 if (!m_isCurrentSet)
340 {
341 // Start the sequence at its minimium value.
342 m_current = m_min;
343 m_isCurrentSet = true;
344 }
345
346 // Return a sequential series of values
347 double r = m_current;
348 if (++m_currentConsecutive == m_consecutive)
349 { // Time to advance to next
350 m_currentConsecutive = 0;
351 m_current += m_increment->GetValue ();
352 if (m_current >= m_max)
353 {
354 m_current = m_min + (m_current - m_max);
355 }
356 }
357 return r;
358 }
359
360 uint32_t
GetInteger(void)361 SequentialRandomVariable::GetInteger (void)
362 {
363 NS_LOG_FUNCTION (this);
364 return (uint32_t)GetValue ();
365 }
366
367 NS_OBJECT_ENSURE_REGISTERED (ExponentialRandomVariable);
368
369 TypeId
GetTypeId(void)370 ExponentialRandomVariable::GetTypeId (void)
371 {
372 static TypeId tid = TypeId ("ns3::ExponentialRandomVariable")
373 .SetParent<RandomVariableStream>()
374 .SetGroupName ("Core")
375 .AddConstructor<ExponentialRandomVariable> ()
376 .AddAttribute ("Mean", "The mean of the values returned by this RNG stream.",
377 DoubleValue (1.0),
378 MakeDoubleAccessor (&ExponentialRandomVariable::m_mean),
379 MakeDoubleChecker<double>())
380 .AddAttribute ("Bound", "The upper bound on the values returned by this RNG stream.",
381 DoubleValue (0.0),
382 MakeDoubleAccessor (&ExponentialRandomVariable::m_bound),
383 MakeDoubleChecker<double>())
384 ;
385 return tid;
386 }
ExponentialRandomVariable()387 ExponentialRandomVariable::ExponentialRandomVariable ()
388 {
389 // m_mean and m_bound are initialized after constructor by attributes
390 NS_LOG_FUNCTION (this);
391 }
392
393 double
GetMean(void) const394 ExponentialRandomVariable::GetMean (void) const
395 {
396 NS_LOG_FUNCTION (this);
397 return m_mean;
398 }
399 double
GetBound(void) const400 ExponentialRandomVariable::GetBound (void) const
401 {
402 NS_LOG_FUNCTION (this);
403 return m_bound;
404 }
405
406 double
GetValue(double mean,double bound)407 ExponentialRandomVariable::GetValue (double mean, double bound)
408 {
409 NS_LOG_FUNCTION (this << mean << bound);
410 while (1)
411 {
412 // Get a uniform random variable in [0,1].
413 double v = Peek ()->RandU01 ();
414 if (IsAntithetic ())
415 {
416 v = (1 - v);
417 }
418
419 // Calculate the exponential random variable.
420 double r = -mean*std::log (v);
421
422 // Use this value if it's acceptable.
423 if (bound == 0 || r <= bound)
424 {
425 return r;
426 }
427 }
428 }
429 uint32_t
GetInteger(uint32_t mean,uint32_t bound)430 ExponentialRandomVariable::GetInteger (uint32_t mean, uint32_t bound)
431 {
432 NS_LOG_FUNCTION (this << mean << bound);
433 return static_cast<uint32_t> ( GetValue (mean, bound) );
434 }
435
436 double
GetValue(void)437 ExponentialRandomVariable::GetValue (void)
438 {
439 NS_LOG_FUNCTION (this);
440 return GetValue (m_mean, m_bound);
441 }
442 uint32_t
GetInteger(void)443 ExponentialRandomVariable::GetInteger (void)
444 {
445 NS_LOG_FUNCTION (this);
446 return (uint32_t)GetValue (m_mean, m_bound);
447 }
448
449 NS_OBJECT_ENSURE_REGISTERED (ParetoRandomVariable);
450
451 TypeId
GetTypeId(void)452 ParetoRandomVariable::GetTypeId (void)
453 {
454 static TypeId tid = TypeId ("ns3::ParetoRandomVariable")
455 .SetParent<RandomVariableStream>()
456 .SetGroupName ("Core")
457 .AddConstructor<ParetoRandomVariable> ()
458 .AddAttribute ("Scale", "The scale parameter for the Pareto distribution returned by this RNG stream.",
459 DoubleValue (1.0),
460 MakeDoubleAccessor (&ParetoRandomVariable::m_scale),
461 MakeDoubleChecker<double>())
462 .AddAttribute ("Shape", "The shape parameter for the Pareto distribution returned by this RNG stream.",
463 DoubleValue (2.0),
464 MakeDoubleAccessor (&ParetoRandomVariable::m_shape),
465 MakeDoubleChecker<double>())
466 .AddAttribute ("Bound", "The upper bound on the values returned by this RNG stream (if non-zero).",
467 DoubleValue (0.0),
468 MakeDoubleAccessor (&ParetoRandomVariable::m_bound),
469 MakeDoubleChecker<double>())
470 ;
471 return tid;
472 }
ParetoRandomVariable()473 ParetoRandomVariable::ParetoRandomVariable ()
474 {
475 // m_shape, m_shape, and m_bound are initialized after constructor
476 // by attributes
477 NS_LOG_FUNCTION (this);
478 }
479
480 double
GetScale(void) const481 ParetoRandomVariable::GetScale (void) const
482 {
483 NS_LOG_FUNCTION (this);
484 return m_scale;
485 }
486
487 double
GetShape(void) const488 ParetoRandomVariable::GetShape (void) const
489 {
490 NS_LOG_FUNCTION (this);
491 return m_shape;
492 }
493
494 double
GetBound(void) const495 ParetoRandomVariable::GetBound (void) const
496 {
497 NS_LOG_FUNCTION (this);
498 return m_bound;
499 }
500
501 double
GetValue(double scale,double shape,double bound)502 ParetoRandomVariable::GetValue (double scale, double shape, double bound)
503 {
504 // Calculate the scale parameter.
505 NS_LOG_FUNCTION (this << scale << shape << bound);
506
507 while (1)
508 {
509 // Get a uniform random variable in [0,1].
510 double v = Peek ()->RandU01 ();
511 if (IsAntithetic ())
512 {
513 v = (1 - v);
514 }
515
516 // Calculate the Pareto random variable.
517 double r = (scale * ( 1.0 / std::pow (v, 1.0 / shape)));
518
519 // Use this value if it's acceptable.
520 if (bound == 0 || r <= bound)
521 {
522 return r;
523 }
524 }
525 }
526 uint32_t
GetInteger(uint32_t scale,uint32_t shape,uint32_t bound)527 ParetoRandomVariable::GetInteger (uint32_t scale, uint32_t shape, uint32_t bound)
528 {
529 NS_LOG_FUNCTION (this << scale << shape << bound);
530 return static_cast<uint32_t> ( GetValue (scale, shape, bound) );
531 }
532
533 double
GetValue(void)534 ParetoRandomVariable::GetValue (void)
535 {
536 NS_LOG_FUNCTION (this);
537 return GetValue (m_scale, m_shape, m_bound);
538 }
539 uint32_t
GetInteger(void)540 ParetoRandomVariable::GetInteger (void)
541 {
542 NS_LOG_FUNCTION (this);
543 return (uint32_t)GetValue (m_scale, m_shape, m_bound);
544 }
545
546 NS_OBJECT_ENSURE_REGISTERED (WeibullRandomVariable);
547
548 TypeId
GetTypeId(void)549 WeibullRandomVariable::GetTypeId (void)
550 {
551 static TypeId tid = TypeId ("ns3::WeibullRandomVariable")
552 .SetParent<RandomVariableStream>()
553 .SetGroupName ("Core")
554 .AddConstructor<WeibullRandomVariable> ()
555 .AddAttribute ("Scale", "The scale parameter for the Weibull distribution returned by this RNG stream.",
556 DoubleValue (1.0),
557 MakeDoubleAccessor (&WeibullRandomVariable::m_scale),
558 MakeDoubleChecker<double>())
559 .AddAttribute ("Shape", "The shape parameter for the Weibull distribution returned by this RNG stream.",
560 DoubleValue (1),
561 MakeDoubleAccessor (&WeibullRandomVariable::m_shape),
562 MakeDoubleChecker<double>())
563 .AddAttribute ("Bound", "The upper bound on the values returned by this RNG stream.",
564 DoubleValue (0.0),
565 MakeDoubleAccessor (&WeibullRandomVariable::m_bound),
566 MakeDoubleChecker<double>())
567 ;
568 return tid;
569 }
WeibullRandomVariable()570 WeibullRandomVariable::WeibullRandomVariable ()
571 {
572 // m_scale, m_shape, and m_bound are initialized after constructor
573 // by attributes
574 NS_LOG_FUNCTION (this);
575 }
576
577 double
GetScale(void) const578 WeibullRandomVariable::GetScale (void) const
579 {
580 NS_LOG_FUNCTION (this);
581 return m_scale;
582 }
583 double
GetShape(void) const584 WeibullRandomVariable::GetShape (void) const
585 {
586 NS_LOG_FUNCTION (this);
587 return m_shape;
588 }
589 double
GetBound(void) const590 WeibullRandomVariable::GetBound (void) const
591 {
592 NS_LOG_FUNCTION (this);
593 return m_bound;
594 }
595
596 double
GetValue(double scale,double shape,double bound)597 WeibullRandomVariable::GetValue (double scale, double shape, double bound)
598 {
599 NS_LOG_FUNCTION (this << scale << shape << bound);
600 double exponent = 1.0 / shape;
601 while (1)
602 {
603 // Get a uniform random variable in [0,1].
604 double v = Peek ()->RandU01 ();
605 if (IsAntithetic ())
606 {
607 v = (1 - v);
608 }
609
610 // Calculate the Weibull random variable.
611 double r = scale * std::pow ( -std::log (v), exponent);
612
613 // Use this value if it's acceptable.
614 if (bound == 0 || r <= bound)
615 {
616 return r;
617 }
618 }
619 }
620 uint32_t
GetInteger(uint32_t scale,uint32_t shape,uint32_t bound)621 WeibullRandomVariable::GetInteger (uint32_t scale, uint32_t shape, uint32_t bound)
622 {
623 NS_LOG_FUNCTION (this << scale << shape << bound);
624 return static_cast<uint32_t> ( GetValue (scale, shape, bound) );
625 }
626
627 double
GetValue(void)628 WeibullRandomVariable::GetValue (void)
629 {
630 NS_LOG_FUNCTION (this);
631 return GetValue (m_scale, m_shape, m_bound);
632 }
633 uint32_t
GetInteger(void)634 WeibullRandomVariable::GetInteger (void)
635 {
636 NS_LOG_FUNCTION (this);
637 return (uint32_t)GetValue (m_scale, m_shape, m_bound);
638 }
639
640 NS_OBJECT_ENSURE_REGISTERED (NormalRandomVariable);
641
642 const double NormalRandomVariable::INFINITE_VALUE = 1e307;
643
644 TypeId
GetTypeId(void)645 NormalRandomVariable::GetTypeId (void)
646 {
647 static TypeId tid = TypeId ("ns3::NormalRandomVariable")
648 .SetParent<RandomVariableStream>()
649 .SetGroupName ("Core")
650 .AddConstructor<NormalRandomVariable> ()
651 .AddAttribute ("Mean", "The mean value for the normal distribution returned by this RNG stream.",
652 DoubleValue (0.0),
653 MakeDoubleAccessor (&NormalRandomVariable::m_mean),
654 MakeDoubleChecker<double>())
655 .AddAttribute ("Variance", "The variance value for the normal distribution returned by this RNG stream.",
656 DoubleValue (1.0),
657 MakeDoubleAccessor (&NormalRandomVariable::m_variance),
658 MakeDoubleChecker<double>())
659 .AddAttribute ("Bound", "The bound on the values returned by this RNG stream.",
660 DoubleValue (INFINITE_VALUE),
661 MakeDoubleAccessor (&NormalRandomVariable::m_bound),
662 MakeDoubleChecker<double>())
663 ;
664 return tid;
665 }
NormalRandomVariable()666 NormalRandomVariable::NormalRandomVariable ()
667 :
668 m_nextValid (false)
669 {
670 // m_mean, m_variance, and m_bound are initialized after constructor
671 // by attributes
672 NS_LOG_FUNCTION (this);
673 }
674
675 double
GetMean(void) const676 NormalRandomVariable::GetMean (void) const
677 {
678 NS_LOG_FUNCTION (this);
679 return m_mean;
680 }
681 double
GetVariance(void) const682 NormalRandomVariable::GetVariance (void) const
683 {
684 NS_LOG_FUNCTION (this);
685 return m_variance;
686 }
687 double
GetBound(void) const688 NormalRandomVariable::GetBound (void) const
689 {
690 NS_LOG_FUNCTION (this);
691 return m_bound;
692 }
693
694 double
GetValue(double mean,double variance,double bound)695 NormalRandomVariable::GetValue (double mean, double variance, double bound)
696 {
697 NS_LOG_FUNCTION (this << mean << variance << bound);
698 if (m_nextValid)
699 { // use previously generated
700 m_nextValid = false;
701 double x2 = mean + m_v2 * m_y * std::sqrt (variance);
702 if (std::fabs (x2 - mean) <= bound)
703 {
704 return x2;
705 }
706 }
707 while (1)
708 { // See Simulation Modeling and Analysis p. 466 (Averill Law)
709 // for algorithm; basically a Box-Muller transform:
710 // http://en.wikipedia.org/wiki/Box-Muller_transform
711 double u1 = Peek ()->RandU01 ();
712 double u2 = Peek ()->RandU01 ();
713 if (IsAntithetic ())
714 {
715 u1 = (1 - u1);
716 u2 = (1 - u2);
717 }
718 double v1 = 2 * u1 - 1;
719 double v2 = 2 * u2 - 1;
720 double w = v1 * v1 + v2 * v2;
721 if (w <= 1.0)
722 { // Got good pair
723 double y = std::sqrt ((-2 * std::log (w)) / w);
724 double x1 = mean + v1 * y * std::sqrt (variance);
725 // if x1 is in bounds, return it, cache v2 and y
726 if (std::fabs (x1 - mean) <= bound)
727 {
728 m_nextValid = true;
729 m_y = y;
730 m_v2 = v2;
731 return x1;
732 }
733 // otherwise try and return the other if it is valid
734 double x2 = mean + v2 * y * std::sqrt (variance);
735 if (std::fabs (x2 - mean) <= bound)
736 {
737 m_nextValid = false;
738 return x2;
739 }
740 // otherwise, just run this loop again
741 }
742 }
743 }
744
745 uint32_t
GetInteger(uint32_t mean,uint32_t variance,uint32_t bound)746 NormalRandomVariable::GetInteger (uint32_t mean, uint32_t variance, uint32_t bound)
747 {
748 NS_LOG_FUNCTION (this << mean << variance << bound);
749 return static_cast<uint32_t> ( GetValue (mean, variance, bound) );
750 }
751
752 double
GetValue(void)753 NormalRandomVariable::GetValue (void)
754 {
755 NS_LOG_FUNCTION (this);
756 return GetValue (m_mean, m_variance, m_bound);
757 }
758 uint32_t
GetInteger(void)759 NormalRandomVariable::GetInteger (void)
760 {
761 NS_LOG_FUNCTION (this);
762 return (uint32_t)GetValue (m_mean, m_variance, m_bound);
763 }
764
765 NS_OBJECT_ENSURE_REGISTERED (LogNormalRandomVariable);
766
767 TypeId
GetTypeId(void)768 LogNormalRandomVariable::GetTypeId (void)
769 {
770 static TypeId tid = TypeId ("ns3::LogNormalRandomVariable")
771 .SetParent<RandomVariableStream>()
772 .SetGroupName ("Core")
773 .AddConstructor<LogNormalRandomVariable> ()
774 .AddAttribute ("Mu", "The mu value for the log-normal distribution returned by this RNG stream.",
775 DoubleValue (0.0),
776 MakeDoubleAccessor (&LogNormalRandomVariable::m_mu),
777 MakeDoubleChecker<double>())
778 .AddAttribute ("Sigma", "The sigma value for the log-normal distribution returned by this RNG stream.",
779 DoubleValue (1.0),
780 MakeDoubleAccessor (&LogNormalRandomVariable::m_sigma),
781 MakeDoubleChecker<double>())
782 ;
783 return tid;
784 }
LogNormalRandomVariable()785 LogNormalRandomVariable::LogNormalRandomVariable ()
786 {
787 // m_mu and m_sigma are initialized after constructor by
788 // attributes
789 NS_LOG_FUNCTION (this);
790 }
791
792 double
GetMu(void) const793 LogNormalRandomVariable::GetMu (void) const
794 {
795 NS_LOG_FUNCTION (this);
796 return m_mu;
797 }
798 double
GetSigma(void) const799 LogNormalRandomVariable::GetSigma (void) const
800 {
801 NS_LOG_FUNCTION (this);
802 return m_sigma;
803 }
804
805 // The code from this function was adapted from the GNU Scientific
806 // Library 1.8:
807 /* randist/lognormal.c
808 *
809 * Copyright (C) 1996, 1997, 1998, 1999, 2000 James Theiler, Brian Gough
810 *
811 * This program is free software; you can redistribute it and/or modify
812 * it under the terms of the GNU General Public License as published by
813 * the Free Software Foundation; either version 2 of the License, or (at
814 * your option) any later version.
815 *
816 * This program is distributed in the hope that it will be useful, but
817 * WITHOUT ANY WARRANTY; without even the implied warranty of
818 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
819 * General Public License for more details.
820 *
821 * You should have received a copy of the GNU General Public License
822 * along with this program; if not, write to the Free Software
823 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
824 */
825 /* The lognormal distribution has the form
826
827 p(x) dx = 1/(x * sqrt(2 pi sigma^2)) exp(-(ln(x) - zeta)^2/2 sigma^2) dx
828
829 for x > 0. Lognormal random numbers are the exponentials of
830 gaussian random numbers */
831 double
GetValue(double mu,double sigma)832 LogNormalRandomVariable::GetValue (double mu, double sigma)
833 {
834 double v1, v2, r2, normal, x;
835
836 NS_LOG_FUNCTION (this << mu << sigma);
837
838 do
839 {
840 /* choose x,y in uniform square (-1,-1) to (+1,+1) */
841
842 double u1 = Peek ()->RandU01 ();
843 double u2 = Peek ()->RandU01 ();
844 if (IsAntithetic ())
845 {
846 u1 = (1 - u1);
847 u2 = (1 - u2);
848 }
849
850 v1 = -1 + 2 * u1;
851 v2 = -1 + 2 * u2;
852
853 /* see if it is in the unit circle */
854 r2 = v1 * v1 + v2 * v2;
855 }
856 while (r2 > 1.0 || r2 == 0);
857
858 normal = v1 * std::sqrt (-2.0 * std::log (r2) / r2);
859
860 x = std::exp (sigma * normal + mu);
861
862 return x;
863 }
864
865 uint32_t
GetInteger(uint32_t mu,uint32_t sigma)866 LogNormalRandomVariable::GetInteger (uint32_t mu, uint32_t sigma)
867 {
868 NS_LOG_FUNCTION (this << mu << sigma);
869 return static_cast<uint32_t> ( GetValue (mu, sigma));
870 }
871
872 double
GetValue(void)873 LogNormalRandomVariable::GetValue (void)
874 {
875 NS_LOG_FUNCTION (this);
876 return GetValue (m_mu, m_sigma);
877 }
878 uint32_t
GetInteger(void)879 LogNormalRandomVariable::GetInteger (void)
880 {
881 NS_LOG_FUNCTION (this);
882 return (uint32_t)GetValue (m_mu, m_sigma);
883 }
884
885 NS_OBJECT_ENSURE_REGISTERED (GammaRandomVariable);
886
887 TypeId
GetTypeId(void)888 GammaRandomVariable::GetTypeId (void)
889 {
890 static TypeId tid = TypeId ("ns3::GammaRandomVariable")
891 .SetParent<RandomVariableStream>()
892 .SetGroupName ("Core")
893 .AddConstructor<GammaRandomVariable> ()
894 .AddAttribute ("Alpha", "The alpha value for the gamma distribution returned by this RNG stream.",
895 DoubleValue (1.0),
896 MakeDoubleAccessor (&GammaRandomVariable::m_alpha),
897 MakeDoubleChecker<double>())
898 .AddAttribute ("Beta", "The beta value for the gamma distribution returned by this RNG stream.",
899 DoubleValue (1.0),
900 MakeDoubleAccessor (&GammaRandomVariable::m_beta),
901 MakeDoubleChecker<double>())
902 ;
903 return tid;
904 }
GammaRandomVariable()905 GammaRandomVariable::GammaRandomVariable ()
906 :
907 m_nextValid (false)
908 {
909 // m_alpha and m_beta are initialized after constructor by
910 // attributes
911 NS_LOG_FUNCTION (this);
912 }
913
914 double
GetAlpha(void) const915 GammaRandomVariable::GetAlpha (void) const
916 {
917 NS_LOG_FUNCTION (this);
918 return m_alpha;
919 }
920 double
GetBeta(void) const921 GammaRandomVariable::GetBeta (void) const
922 {
923 NS_LOG_FUNCTION (this);
924 return m_beta;
925 }
926
927 /*
928 The code for the following generator functions was adapted from ns-2
929 tools/ranvar.cc
930
931 Originally the algorithm was devised by Marsaglia in 2000:
932 G. Marsaglia, W. W. Tsang: A simple method for generating Gamma variables
933 ACM Transactions on mathematical software, Vol. 26, No. 3, Sept. 2000
934
935 The Gamma distribution density function has the form
936
937 x^(alpha-1) * exp(-x/beta)
938 p(x; alpha, beta) = ----------------------------
939 beta^alpha * Gamma(alpha)
940
941 for x > 0.
942 */
943 double
GetValue(double alpha,double beta)944 GammaRandomVariable::GetValue (double alpha, double beta)
945 {
946 NS_LOG_FUNCTION (this << alpha << beta);
947 if (alpha < 1)
948 {
949 double u = Peek ()->RandU01 ();
950 if (IsAntithetic ())
951 {
952 u = (1 - u);
953 }
954 return GetValue (1.0 + alpha, beta) * std::pow (u, 1.0 / alpha);
955 }
956
957 double x, v, u;
958 double d = alpha - 1.0 / 3.0;
959 double c = (1.0 / 3.0) / std::sqrt (d);
960
961 while (1)
962 {
963 do
964 {
965 // Get a value from a normal distribution that has mean
966 // zero, variance 1, and no bound.
967 double mean = 0.0;
968 double variance = 1.0;
969 double bound = NormalRandomVariable::INFINITE_VALUE;
970 x = GetNormalValue (mean, variance, bound);
971
972 v = 1.0 + c * x;
973 }
974 while (v <= 0);
975
976 v = v * v * v;
977 u = Peek ()->RandU01 ();
978 if (IsAntithetic ())
979 {
980 u = (1 - u);
981 }
982 if (u < 1 - 0.0331 * x * x * x * x)
983 {
984 break;
985 }
986 if (std::log (u) < 0.5 * x * x + d * (1 - v + std::log (v)))
987 {
988 break;
989 }
990 }
991
992 return beta * d * v;
993 }
994
995 uint32_t
GetInteger(uint32_t alpha,uint32_t beta)996 GammaRandomVariable::GetInteger (uint32_t alpha, uint32_t beta)
997 {
998 NS_LOG_FUNCTION (this << alpha << beta);
999 return static_cast<uint32_t> ( GetValue (alpha, beta));
1000 }
1001
1002 double
GetValue(void)1003 GammaRandomVariable::GetValue (void)
1004 {
1005 NS_LOG_FUNCTION (this);
1006 return GetValue (m_alpha, m_beta);
1007 }
1008 uint32_t
GetInteger(void)1009 GammaRandomVariable::GetInteger (void)
1010 {
1011 NS_LOG_FUNCTION (this);
1012 return (uint32_t)GetValue (m_alpha, m_beta);
1013 }
1014
1015 double
GetNormalValue(double mean,double variance,double bound)1016 GammaRandomVariable::GetNormalValue (double mean, double variance, double bound)
1017 {
1018 NS_LOG_FUNCTION (this << mean << variance << bound);
1019 if (m_nextValid)
1020 { // use previously generated
1021 m_nextValid = false;
1022 double x2 = mean + m_v2 * m_y * std::sqrt (variance);
1023 if (std::fabs (x2 - mean) <= bound)
1024 {
1025 return x2;
1026 }
1027 }
1028 while (1)
1029 { // See Simulation Modeling and Analysis p. 466 (Averill Law)
1030 // for algorithm; basically a Box-Muller transform:
1031 // http://en.wikipedia.org/wiki/Box-Muller_transform
1032 double u1 = Peek ()->RandU01 ();
1033 double u2 = Peek ()->RandU01 ();
1034 if (IsAntithetic ())
1035 {
1036 u1 = (1 - u1);
1037 u2 = (1 - u2);
1038 }
1039 double v1 = 2 * u1 - 1;
1040 double v2 = 2 * u2 - 1;
1041 double w = v1 * v1 + v2 * v2;
1042 if (w <= 1.0)
1043 { // Got good pair
1044 double y = std::sqrt ((-2 * std::log (w)) / w);
1045 double x1 = mean + v1 * y * std::sqrt (variance);
1046 // if x1 is in bounds, return it, cache v2 an y
1047 if (std::fabs (x1 - mean) <= bound)
1048 {
1049 m_nextValid = true;
1050 m_y = y;
1051 m_v2 = v2;
1052 return x1;
1053 }
1054 // otherwise try and return the other if it is valid
1055 double x2 = mean + v2 * y * std::sqrt (variance);
1056 if (std::fabs (x2 - mean) <= bound)
1057 {
1058 m_nextValid = false;
1059 return x2;
1060 }
1061 // otherwise, just run this loop again
1062 }
1063 }
1064 }
1065
1066 NS_OBJECT_ENSURE_REGISTERED (ErlangRandomVariable);
1067
1068 TypeId
GetTypeId(void)1069 ErlangRandomVariable::GetTypeId (void)
1070 {
1071 static TypeId tid = TypeId ("ns3::ErlangRandomVariable")
1072 .SetParent<RandomVariableStream>()
1073 .SetGroupName ("Core")
1074 .AddConstructor<ErlangRandomVariable> ()
1075 .AddAttribute ("K", "The k value for the Erlang distribution returned by this RNG stream.",
1076 IntegerValue (1),
1077 MakeIntegerAccessor (&ErlangRandomVariable::m_k),
1078 MakeIntegerChecker<uint32_t>())
1079 .AddAttribute ("Lambda", "The lambda value for the Erlang distribution returned by this RNG stream.",
1080 DoubleValue (1.0),
1081 MakeDoubleAccessor (&ErlangRandomVariable::m_lambda),
1082 MakeDoubleChecker<double>())
1083 ;
1084 return tid;
1085 }
ErlangRandomVariable()1086 ErlangRandomVariable::ErlangRandomVariable ()
1087 {
1088 // m_k and m_lambda are initialized after constructor by attributes
1089 NS_LOG_FUNCTION (this);
1090 }
1091
1092 uint32_t
GetK(void) const1093 ErlangRandomVariable::GetK (void) const
1094 {
1095 NS_LOG_FUNCTION (this);
1096 return m_k;
1097 }
1098 double
GetLambda(void) const1099 ErlangRandomVariable::GetLambda (void) const
1100 {
1101 NS_LOG_FUNCTION (this);
1102 return m_lambda;
1103 }
1104
1105 /*
1106 The code for the following generator functions was adapted from ns-2
1107 tools/ranvar.cc
1108
1109 The Erlang distribution density function has the form
1110
1111 x^(k-1) * exp(-x/lambda)
1112 p(x; k, lambda) = ---------------------------
1113 lambda^k * (k-1)!
1114
1115 for x > 0.
1116 */
1117 double
GetValue(uint32_t k,double lambda)1118 ErlangRandomVariable::GetValue (uint32_t k, double lambda)
1119 {
1120 NS_LOG_FUNCTION (this << k << lambda);
1121 double mean = lambda;
1122 double bound = 0.0;
1123
1124 double result = 0;
1125 for (unsigned int i = 0; i < k; ++i)
1126 {
1127 result += GetExponentialValue (mean, bound);
1128
1129 }
1130
1131 return result;
1132 }
1133
1134 uint32_t
GetInteger(uint32_t k,uint32_t lambda)1135 ErlangRandomVariable::GetInteger (uint32_t k, uint32_t lambda)
1136 {
1137 NS_LOG_FUNCTION (this << k << lambda);
1138 return static_cast<uint32_t> ( GetValue (k, lambda));
1139 }
1140
1141 double
GetValue(void)1142 ErlangRandomVariable::GetValue (void)
1143 {
1144 NS_LOG_FUNCTION (this);
1145 return GetValue (m_k, m_lambda);
1146 }
1147 uint32_t
GetInteger(void)1148 ErlangRandomVariable::GetInteger (void)
1149 {
1150 NS_LOG_FUNCTION (this);
1151 return (uint32_t)GetValue (m_k, m_lambda);
1152 }
1153
1154 double
GetExponentialValue(double mean,double bound)1155 ErlangRandomVariable::GetExponentialValue (double mean, double bound)
1156 {
1157 NS_LOG_FUNCTION (this << mean << bound);
1158 while (1)
1159 {
1160 // Get a uniform random variable in [0,1].
1161 double v = Peek ()->RandU01 ();
1162 if (IsAntithetic ())
1163 {
1164 v = (1 - v);
1165 }
1166
1167 // Calculate the exponential random variable.
1168 double r = -mean*std::log (v);
1169
1170 // Use this value if it's acceptable.
1171 if (bound == 0 || r <= bound)
1172 {
1173 return r;
1174 }
1175 }
1176 }
1177
1178 NS_OBJECT_ENSURE_REGISTERED (TriangularRandomVariable);
1179
1180 TypeId
GetTypeId(void)1181 TriangularRandomVariable::GetTypeId (void)
1182 {
1183 static TypeId tid = TypeId ("ns3::TriangularRandomVariable")
1184 .SetParent<RandomVariableStream>()
1185 .SetGroupName ("Core")
1186 .AddConstructor<TriangularRandomVariable> ()
1187 .AddAttribute ("Mean", "The mean value for the triangular distribution returned by this RNG stream.",
1188 DoubleValue (0.5),
1189 MakeDoubleAccessor (&TriangularRandomVariable::m_mean),
1190 MakeDoubleChecker<double>())
1191 .AddAttribute ("Min", "The lower bound on the values returned by this RNG stream.",
1192 DoubleValue (0.0),
1193 MakeDoubleAccessor (&TriangularRandomVariable::m_min),
1194 MakeDoubleChecker<double>())
1195 .AddAttribute ("Max", "The upper bound on the values returned by this RNG stream.",
1196 DoubleValue (1.0),
1197 MakeDoubleAccessor (&TriangularRandomVariable::m_max),
1198 MakeDoubleChecker<double>())
1199 ;
1200 return tid;
1201 }
TriangularRandomVariable()1202 TriangularRandomVariable::TriangularRandomVariable ()
1203 {
1204 // m_mean, m_min, and m_max are initialized after constructor by
1205 // attributes
1206 NS_LOG_FUNCTION (this);
1207 }
1208
1209 double
GetMean(void) const1210 TriangularRandomVariable::GetMean (void) const
1211 {
1212 NS_LOG_FUNCTION (this);
1213 return m_mean;
1214 }
1215 double
GetMin(void) const1216 TriangularRandomVariable::GetMin (void) const
1217 {
1218 NS_LOG_FUNCTION (this);
1219 return m_min;
1220 }
1221 double
GetMax(void) const1222 TriangularRandomVariable::GetMax (void) const
1223 {
1224 NS_LOG_FUNCTION (this);
1225 return m_max;
1226 }
1227
1228 double
GetValue(double mean,double min,double max)1229 TriangularRandomVariable::GetValue (double mean, double min, double max)
1230 {
1231 // Calculate the mode.
1232 NS_LOG_FUNCTION (this << mean << min << max);
1233 double mode = 3.0 * mean - min - max;
1234
1235 // Get a uniform random variable in [0,1].
1236 double u = Peek ()->RandU01 ();
1237 if (IsAntithetic ())
1238 {
1239 u = (1 - u);
1240 }
1241
1242 // Calculate the triangular random variable.
1243 if (u <= (mode - min) / (max - min) )
1244 {
1245 return min + std::sqrt (u * (max - min) * (mode - min) );
1246 }
1247 else
1248 {
1249 return max - std::sqrt ( (1 - u) * (max - min) * (max - mode) );
1250 }
1251 }
1252
1253 uint32_t
GetInteger(uint32_t mean,uint32_t min,uint32_t max)1254 TriangularRandomVariable::GetInteger (uint32_t mean, uint32_t min, uint32_t max)
1255 {
1256 NS_LOG_FUNCTION (this << mean << min << max);
1257 return static_cast<uint32_t> ( GetValue (mean, min, max) );
1258 }
1259
1260 double
GetValue(void)1261 TriangularRandomVariable::GetValue (void)
1262 {
1263 NS_LOG_FUNCTION (this);
1264 return GetValue (m_mean, m_min, m_max);
1265 }
1266 uint32_t
GetInteger(void)1267 TriangularRandomVariable::GetInteger (void)
1268 {
1269 NS_LOG_FUNCTION (this);
1270 return (uint32_t)GetValue (m_mean, m_min, m_max);
1271 }
1272
1273 NS_OBJECT_ENSURE_REGISTERED (ZipfRandomVariable);
1274
1275 TypeId
GetTypeId(void)1276 ZipfRandomVariable::GetTypeId (void)
1277 {
1278 static TypeId tid = TypeId ("ns3::ZipfRandomVariable")
1279 .SetParent<RandomVariableStream>()
1280 .SetGroupName ("Core")
1281 .AddConstructor<ZipfRandomVariable> ()
1282 .AddAttribute ("N", "The n value for the Zipf distribution returned by this RNG stream.",
1283 IntegerValue (1),
1284 MakeIntegerAccessor (&ZipfRandomVariable::m_n),
1285 MakeIntegerChecker<uint32_t>())
1286 .AddAttribute ("Alpha", "The alpha value for the Zipf distribution returned by this RNG stream.",
1287 DoubleValue (0.0),
1288 MakeDoubleAccessor (&ZipfRandomVariable::m_alpha),
1289 MakeDoubleChecker<double>())
1290 ;
1291 return tid;
1292 }
ZipfRandomVariable()1293 ZipfRandomVariable::ZipfRandomVariable ()
1294 {
1295 // m_n and m_alpha are initialized after constructor by attributes
1296 NS_LOG_FUNCTION (this);
1297 }
1298
1299 uint32_t
GetN(void) const1300 ZipfRandomVariable::GetN (void) const
1301 {
1302 NS_LOG_FUNCTION (this);
1303 return m_n;
1304 }
1305 double
GetAlpha(void) const1306 ZipfRandomVariable::GetAlpha (void) const
1307 {
1308 NS_LOG_FUNCTION (this);
1309 return m_alpha;
1310 }
1311
1312 double
GetValue(uint32_t n,double alpha)1313 ZipfRandomVariable::GetValue (uint32_t n, double alpha)
1314 {
1315 NS_LOG_FUNCTION (this << n << alpha);
1316 // Calculate the normalization constant c.
1317 m_c = 0.0;
1318 for (uint32_t i = 1; i <= n; i++)
1319 {
1320 m_c += (1.0 / std::pow ((double)i,alpha));
1321 }
1322 m_c = 1.0 / m_c;
1323
1324 // Get a uniform random variable in [0,1].
1325 double u = Peek ()->RandU01 ();
1326 if (IsAntithetic ())
1327 {
1328 u = (1 - u);
1329 }
1330
1331 double sum_prob = 0,zipf_value = 0;
1332 for (uint32_t i = 1; i <= m_n; i++)
1333 {
1334 sum_prob += m_c / std::pow ((double)i,m_alpha);
1335 if (sum_prob > u)
1336 {
1337 zipf_value = i;
1338 break;
1339 }
1340 }
1341 return zipf_value;
1342 }
1343
1344 uint32_t
GetInteger(uint32_t n,uint32_t alpha)1345 ZipfRandomVariable::GetInteger (uint32_t n, uint32_t alpha)
1346 {
1347 NS_LOG_FUNCTION (this << n << alpha);
1348 return static_cast<uint32_t> ( GetValue (n, alpha));
1349 }
1350
1351 double
GetValue(void)1352 ZipfRandomVariable::GetValue (void)
1353 {
1354 NS_LOG_FUNCTION (this);
1355 return GetValue (m_n, m_alpha);
1356 }
1357 uint32_t
GetInteger(void)1358 ZipfRandomVariable::GetInteger (void)
1359 {
1360 NS_LOG_FUNCTION (this);
1361 return (uint32_t)GetValue (m_n, m_alpha);
1362 }
1363
1364 NS_OBJECT_ENSURE_REGISTERED (ZetaRandomVariable);
1365
1366 TypeId
GetTypeId(void)1367 ZetaRandomVariable::GetTypeId (void)
1368 {
1369 static TypeId tid = TypeId ("ns3::ZetaRandomVariable")
1370 .SetParent<RandomVariableStream>()
1371 .SetGroupName ("Core")
1372 .AddConstructor<ZetaRandomVariable> ()
1373 .AddAttribute ("Alpha", "The alpha value for the zeta distribution returned by this RNG stream.",
1374 DoubleValue (3.14),
1375 MakeDoubleAccessor (&ZetaRandomVariable::m_alpha),
1376 MakeDoubleChecker<double>())
1377 ;
1378 return tid;
1379 }
ZetaRandomVariable()1380 ZetaRandomVariable::ZetaRandomVariable ()
1381 {
1382 // m_alpha is initialized after constructor by attributes
1383 NS_LOG_FUNCTION (this);
1384 }
1385
1386 double
GetAlpha(void) const1387 ZetaRandomVariable::GetAlpha (void) const
1388 {
1389 NS_LOG_FUNCTION (this);
1390 return m_alpha;
1391 }
1392
1393 double
GetValue(double alpha)1394 ZetaRandomVariable::GetValue (double alpha)
1395 {
1396 NS_LOG_FUNCTION (this << alpha);
1397 m_b = std::pow (2.0, alpha - 1.0);
1398
1399 double u, v;
1400 double X, T;
1401 double test;
1402
1403 do
1404 {
1405 // Get a uniform random variable in [0,1].
1406 u = Peek ()->RandU01 ();
1407 if (IsAntithetic ())
1408 {
1409 u = (1 - u);
1410 }
1411
1412 // Get a uniform random variable in [0,1].
1413 v = Peek ()->RandU01 ();
1414 if (IsAntithetic ())
1415 {
1416 v = (1 - v);
1417 }
1418
1419 X = std::floor (std::pow (u, -1.0 / (m_alpha - 1.0)));
1420 T = std::pow (1.0 + 1.0 / X, m_alpha - 1.0);
1421 test = v * X * (T - 1.0) / (m_b - 1.0);
1422 }
1423 while ( test > (T / m_b) );
1424
1425 return X;
1426 }
1427
1428 uint32_t
GetInteger(uint32_t alpha)1429 ZetaRandomVariable::GetInteger (uint32_t alpha)
1430 {
1431 NS_LOG_FUNCTION (this << alpha);
1432 return static_cast<uint32_t> ( GetValue (alpha));
1433 }
1434
1435 double
GetValue(void)1436 ZetaRandomVariable::GetValue (void)
1437 {
1438 NS_LOG_FUNCTION (this);
1439 return GetValue (m_alpha);
1440 }
1441 uint32_t
GetInteger(void)1442 ZetaRandomVariable::GetInteger (void)
1443 {
1444 NS_LOG_FUNCTION (this);
1445 return (uint32_t)GetValue (m_alpha);
1446 }
1447
1448 NS_OBJECT_ENSURE_REGISTERED (DeterministicRandomVariable);
1449
1450 TypeId
GetTypeId(void)1451 DeterministicRandomVariable::GetTypeId (void)
1452 {
1453 static TypeId tid = TypeId ("ns3::DeterministicRandomVariable")
1454 .SetParent<RandomVariableStream>()
1455 .SetGroupName ("Core")
1456 .AddConstructor<DeterministicRandomVariable> ()
1457 ;
1458 return tid;
1459 }
DeterministicRandomVariable()1460 DeterministicRandomVariable::DeterministicRandomVariable ()
1461 :
1462 m_count (0),
1463 m_next (0),
1464 m_data (0)
1465 {
1466 NS_LOG_FUNCTION (this);
1467 }
~DeterministicRandomVariable()1468 DeterministicRandomVariable::~DeterministicRandomVariable ()
1469 {
1470 // Delete any values currently set.
1471 NS_LOG_FUNCTION (this);
1472 if (m_data != 0)
1473 {
1474 delete[] m_data;
1475 }
1476 }
1477
1478 void
SetValueArray(double * values,std::size_t length)1479 DeterministicRandomVariable::SetValueArray (double* values, std::size_t length)
1480 {
1481 NS_LOG_FUNCTION (this << values << length);
1482 // Delete any values currently set.
1483 if (m_data != 0)
1484 {
1485 delete[] m_data;
1486 }
1487
1488 // Make room for the values being set.
1489 m_data = new double[length];
1490 m_count = length;
1491 m_next = length;
1492
1493 // Copy the values.
1494 for (std::size_t i = 0; i < m_count; i++)
1495 {
1496 m_data[i] = values[i];
1497 }
1498 }
1499
1500 double
GetValue(void)1501 DeterministicRandomVariable::GetValue (void)
1502 {
1503 NS_LOG_FUNCTION (this);
1504 // Make sure the array has been set.
1505 NS_ASSERT (m_count > 0);
1506
1507 if (m_next == m_count)
1508 {
1509 m_next = 0;
1510 }
1511 return m_data[m_next++];
1512 }
1513
1514 uint32_t
GetInteger(void)1515 DeterministicRandomVariable::GetInteger (void)
1516 {
1517 NS_LOG_FUNCTION (this);
1518 return (uint32_t)GetValue ();
1519 }
1520
1521 NS_OBJECT_ENSURE_REGISTERED (EmpiricalRandomVariable);
1522
1523 // ValueCDF methods
ValueCDF(void)1524 EmpiricalRandomVariable::ValueCDF::ValueCDF (void)
1525 : value (0.0),
1526 cdf (0.0)
1527 {
1528 NS_LOG_FUNCTION (this);
1529 }
1530
ValueCDF(double v,double c)1531 EmpiricalRandomVariable::ValueCDF::ValueCDF (double v, double c)
1532 : value (v),
1533 cdf (c)
1534 {
1535 NS_LOG_FUNCTION (this << v << c);
1536 NS_ASSERT (c >= 0.0 && c <= 1.0);
1537 }
1538
1539 bool
operator <(EmpiricalRandomVariable::ValueCDF a,EmpiricalRandomVariable::ValueCDF b)1540 operator < (EmpiricalRandomVariable::ValueCDF a,
1541 EmpiricalRandomVariable::ValueCDF b)
1542 {
1543 return a.cdf < b.cdf;
1544 }
1545
1546 TypeId
GetTypeId(void)1547 EmpiricalRandomVariable::GetTypeId (void)
1548 {
1549 static TypeId tid = TypeId ("ns3::EmpiricalRandomVariable")
1550 .SetParent<RandomVariableStream>()
1551 .SetGroupName ("Core")
1552 .AddConstructor<EmpiricalRandomVariable> ()
1553 .AddAttribute ("Interpolate",
1554 "Treat the CDF as a smooth distribution and interpolate, "
1555 "default is to treat the CDF as a histogram and sample.",
1556 BooleanValue (false),
1557 MakeBooleanAccessor (&EmpiricalRandomVariable::m_interpolate),
1558 MakeBooleanChecker ())
1559 ;
1560 return tid;
1561 }
EmpiricalRandomVariable(void)1562 EmpiricalRandomVariable::EmpiricalRandomVariable (void)
1563 : m_validated (false)
1564 {
1565 NS_LOG_FUNCTION (this);
1566 }
1567
1568 bool
SetInterpolate(bool interpolate)1569 EmpiricalRandomVariable::SetInterpolate (bool interpolate)
1570 {
1571 NS_LOG_FUNCTION (this << interpolate);
1572 bool prev = m_interpolate;
1573 m_interpolate = interpolate;
1574 return prev;
1575 }
1576
1577 uint32_t
GetInteger(void)1578 EmpiricalRandomVariable::GetInteger (void)
1579 {
1580 NS_LOG_FUNCTION (this);
1581 return static_cast<uint32_t> (GetValue ());
1582 }
1583
1584 bool
PreSample(double & value)1585 EmpiricalRandomVariable::PreSample (double & value)
1586 {
1587 NS_LOG_FUNCTION (this);
1588
1589 if (!m_validated)
1590 {
1591 Validate ();
1592 }
1593
1594 // Get a uniform random variable in [0, 1].
1595 double r = Peek ()->RandU01 ();
1596 if (IsAntithetic ())
1597 {
1598 r = (1 - r);
1599 }
1600
1601 value = r;
1602 bool valid = false;
1603 // check extrema
1604 if (r <= m_emp.front ().cdf)
1605 {
1606 value = m_emp.front ().value; // Less than first
1607 valid = true;
1608 }
1609 else if (r >= m_emp.back ().cdf)
1610 {
1611 value = m_emp.back ().value; // Greater than last
1612 valid = true;
1613 }
1614 return valid;
1615 }
1616
1617 double
GetValue(void)1618 EmpiricalRandomVariable::GetValue (void)
1619 {
1620 NS_LOG_FUNCTION (this);
1621
1622 double value;
1623 if (PreSample (value))
1624 {
1625 return value;
1626 }
1627
1628 // value now has the (unused) URNG selector
1629 if (m_interpolate)
1630 {
1631 value = DoInterpolate (value);
1632 }
1633 else
1634 {
1635 value = DoSampleCDF (value);
1636 }
1637 return value;
1638 }
1639
1640 double
DoSampleCDF(double r)1641 EmpiricalRandomVariable::DoSampleCDF (double r)
1642 {
1643 NS_LOG_FUNCTION (this << r);
1644
1645 ValueCDF selector (0, r);
1646 auto bound = std::upper_bound (m_emp.begin (), m_emp.end (), selector);
1647
1648 return bound->value;
1649 }
1650
1651 double
Interpolate(void)1652 EmpiricalRandomVariable::Interpolate (void)
1653 {
1654 NS_LOG_FUNCTION (this);
1655
1656 double value;
1657 if (PreSample (value))
1658 {
1659 return value;
1660 }
1661
1662 // value now has the (unused) URNG selector
1663 value = DoInterpolate (value);
1664 return value;
1665 }
1666
1667 double
DoInterpolate(double r)1668 EmpiricalRandomVariable::DoInterpolate (double r)
1669 {
1670 NS_LOG_FUNCTION (this << r);
1671
1672 // Return a value from the empirical distribution
1673 // This code based (loosely) on code by Bruce Mah (Thanks Bruce!)
1674
1675 // search
1676 ValueCDF selector (0, r);
1677 auto upper = std::upper_bound (m_emp.begin (), m_emp.end (), selector);
1678 auto lower = std::prev (upper, 1);
1679 if (upper == m_emp.begin ())
1680 {
1681 lower = upper;
1682 }
1683
1684 // Interpolate random value in range [v1..v2) based on [c1 .. r .. c2)
1685 double c1 = lower->cdf;
1686 double c2 = upper->cdf;
1687 double v1 = lower->value;
1688 double v2 = upper->value;
1689
1690 double value = (v1 + ((v2 - v1) / (c2 - c1)) * (r - c1));
1691 return value;
1692 }
1693
1694 void
CDF(double v,double c)1695 EmpiricalRandomVariable::CDF (double v, double c)
1696 {
1697 // Add a new empirical datapoint to the empirical cdf
1698 // NOTE. These MUST be inserted in non-decreasing order
1699 NS_LOG_FUNCTION (this << v << c);
1700 m_emp.push_back (ValueCDF (v, c));
1701 }
1702
1703 void
Validate(void)1704 EmpiricalRandomVariable::Validate (void)
1705 {
1706 NS_LOG_FUNCTION (this);
1707 if (m_emp.empty ())
1708 {
1709 NS_FATAL_ERROR ("CDF is not initialized");
1710 }
1711 ValueCDF prior = m_emp[0];
1712 for (auto current : m_emp)
1713 {
1714 if (current.value < prior.value || current.cdf < prior.cdf)
1715 { // Error
1716 std::cerr << "Empirical Dist error,"
1717 << " current value " << current.value
1718 << " prior value " << prior.value
1719 << " current cdf " << current.cdf
1720 << " prior cdf " << prior.cdf << std::endl;
1721 NS_FATAL_ERROR ("Empirical Dist error");
1722 }
1723 prior = current;
1724 }
1725 if (prior.cdf != 1.0)
1726 {
1727 NS_FATAL_ERROR ("CDF does not cover the whole distribution");
1728 }
1729 m_validated = true;
1730 }
1731
1732 } // namespace ns3
1733