1 package nl.tudelft.simulation.jstats.distributions;
2
3 import org.djutils.exceptions.Throw;
4
5 import nl.tudelft.simulation.jstats.math.ProbMath;
6 import nl.tudelft.simulation.jstats.streams.StreamInterface;
7
8 /**
9 * The Erlang distribution. For more information on this distribution see
10 * <a href="https://mathworld.wolfram.com/ErlangDistribution.html"> http://mathworld.wolfram.com/ErlangDistribution.html
11 * </a><br>
12 * The Erlang distribution is the distribution of a sum of k independent exponential variables with the scale parameter as the
13 * mean. The scale parameter is equal to 1/rate or 1/λ, giving the entire Erlang distribution a mean of k*scale.
14 * <p>
15 * Copyright (c) 2002-2024 Delft University of Technology, Jaffalaan 5, 2628 BX Delft, the Netherlands. All rights reserved. See
16 * for project information <a href="https://simulation.tudelft.nl/" target="_blank"> https://simulation.tudelft.nl</a>. The DSOL
17 * project is distributed under a three-clause BSD-style license, which can be found at
18 * <a href="https://https://simulation.tudelft.nl/dsol/docs/latest/license.html" target="_blank">
19 * https://https://simulation.tudelft.nl/dsol/docs/latest/license.html</a>.
20 * </p>
21 * @author <a href="https://www.linkedin.com/in/peterhmjacobs">Peter Jacobs </a>
22 * @author <a href="https://www.tudelft.nl/averbraeck">Alexander Verbraeck</a>
23 */
24 public class DistErlang extends DistContinuous
25 {
26 /** */
27 private static final long serialVersionUID = 1L;
28
29 /**
30 * k is the shape parameter of the Erlang distribution. The shape k is the number of times a drawing is done from the
31 * exponential distribution, where the Erlang distribution is the sum of these k independent exponential variables.
32 */
33 private final int k;
34
35 /** scale is the mean of a single exponential distribution (1/rate), of which k are summed. */
36 private final double scale;
37
38 /** the rate value of the Erlang distribution (1 / scale). */
39 private final double lambda;
40
41 /** distGamma is the underlying gamma distribution. */
42 private final DistGamma distGamma;
43
44 /** GAMMATHRESHOLD is the threshold above which we use a gamma function and below repeated drawing. */
45 private static final short GAMMATHRESHOLD = 10;
46
47 /**
48 * Construct a new Erlang distribution with k and a mean (so not k and a rate) as parameters. It is the distribution of a
49 * sum of k independent exponential variables with the scale parameter as the mean. The scale parameter is equal to 1/rate
50 * or 1/λ, giving the entire Erlang distribution a mean of k*scale.
51 * @param stream StreamInterface; the random number stream
52 * @param scale double; the mean of a single sample from the exponential distribution, of which k are summed. Equal to
53 * 1/rate or 1/λ.
54 * @param k int; the shape parameter of the Erlang distribution. The shape k is the number of times a drawing is done from
55 * the exponential distribution, where the Erlang distribution is the sum of these k independent exponential
56 * variables.
57 * @throws IllegalArgumentException when k <= 0 or scale <= 0
58 */
59 public DistErlang(final StreamInterface stream, final double scale, final int k)
60 {
61 super(stream);
62 Throw.when(k <= 0 || scale <= 0.0, IllegalArgumentException.class, "Error Erlang - k <= 0 or scale <= 0");
63 this.k = k;
64 this.scale = scale;
65 this.lambda = 1.0 / scale;
66 this.distGamma = this.k <= DistErlang.GAMMATHRESHOLD ? null : new DistGamma(stream, this.k, this.scale);
67 }
68
69 /** {@inheritDoc} */
70 @Override
71 public double draw()
72 {
73 if (this.k <= DistErlang.GAMMATHRESHOLD)
74 {
75 // according to Law and Kelton, Simulation Modeling and Analysis
76 // repeated drawing and composition is usually faster for k<=10
77 double product = 1.0;
78 for (int i = 1; i <= this.k; i++)
79 {
80 product = product * this.stream.nextDouble();
81 }
82 return -this.scale * Math.log(product);
83 }
84 // and using the gamma distribution is faster for k>10
85 return this.distGamma.draw();
86 }
87
88 /** {@inheritDoc} */
89 @Override
90 public double getProbabilityDensity(final double x)
91 {
92 if (x < 0)
93 {
94 return 0;
95 }
96 return this.lambda * Math.exp(-this.lambda * x) * Math.pow(this.lambda * x, this.k - 1)
97 / ProbMath.factorial(this.k - 1);
98 }
99
100 /**
101 * @return k
102 */
103 public int getK()
104 {
105 return this.k;
106 }
107
108 /**
109 * @return scale parameter
110 */
111 public double getScale()
112 {
113 return this.scale;
114 }
115
116 /** {@inheritDoc} */
117 @Override
118 public void setStream(final StreamInterface stream)
119 {
120 super.setStream(stream);
121 if (this.distGamma != null)
122 {
123 this.distGamma.setStream(stream);
124 }
125 }
126
127 /** {@inheritDoc} */
128 @Override
129 public String toString()
130 {
131 return "Erlang(" + this.scale + "," + this.k + ")";
132 }
133 }