1 /* 2 * Licensed to the Apache Software Foundation (ASF) under one or more 3 * contributor license agreements. See the NOTICE file distributed with 4 * this work for additional information regarding copyright ownership. 5 * The ASF licenses this file to You under the Apache License, Version 2.0 6 * (the "License"); you may not use this file except in compliance with 7 * the License. You may obtain a copy of the License at 8 * 9 * http://www.apache.org/licenses/LICENSE-2.0 10 * 11 * Unless required by applicable law or agreed to in writing, software 12 * distributed under the License is distributed on an "AS IS" BASIS, 13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 14 * See the License for the specific language governing permissions and 15 * limitations under the License. 16 */ 17 18 package org.apache.commons.math.random; 19 20 import java.io.Serializable; 21 import java.security.MessageDigest; 22 import java.security.NoSuchAlgorithmException; 23 import java.security.NoSuchProviderException; 24 import java.security.SecureRandom; 25 import java.util.Collection; 26 27 import org.apache.commons.math.MathException; 28 import org.apache.commons.math.distribution.BetaDistributionImpl; 29 import org.apache.commons.math.distribution.BinomialDistributionImpl; 30 import org.apache.commons.math.distribution.CauchyDistributionImpl; 31 import org.apache.commons.math.distribution.ChiSquaredDistributionImpl; 32 import org.apache.commons.math.distribution.ContinuousDistribution; 33 import org.apache.commons.math.distribution.FDistributionImpl; 34 import org.apache.commons.math.distribution.GammaDistributionImpl; 35 import org.apache.commons.math.distribution.HypergeometricDistributionImpl; 36 import org.apache.commons.math.distribution.IntegerDistribution; 37 import org.apache.commons.math.distribution.PascalDistributionImpl; 38 import org.apache.commons.math.distribution.TDistributionImpl; 39 import org.apache.commons.math.distribution.WeibullDistributionImpl; 40 import org.apache.commons.math.distribution.ZipfDistributionImpl; 41 import org.apache.commons.math.exception.MathInternalError; 42 import org.apache.commons.math.exception.NotStrictlyPositiveException; 43 import org.apache.commons.math.exception.NumberIsTooLargeException; 44 import org.apache.commons.math.exception.util.LocalizedFormats; 45 import org.apache.commons.math.util.FastMath; 46 import org.apache.commons.math.util.MathUtils; 47 48 /** 49 * Implements the {@link RandomData} interface using a {@link RandomGenerator} 50 * instance to generate non-secure data and a {@link java.security.SecureRandom} 51 * instance to provide data for the <code>nextSecureXxx</code> methods. If no 52 * <code>RandomGenerator</code> is provided in the constructor, the default is 53 * to use a generator based on {@link java.util.Random}. To plug in a different 54 * implementation, either implement <code>RandomGenerator</code> directly or 55 * extend {@link AbstractRandomGenerator}. 56 * <p> 57 * Supports reseeding the underlying pseudo-random number generator (PRNG). The 58 * <code>SecurityProvider</code> and <code>Algorithm</code> used by the 59 * <code>SecureRandom</code> instance can also be reset. 60 * </p> 61 * <p> 62 * For details on the default PRNGs, see {@link java.util.Random} and 63 * {@link java.security.SecureRandom}. 64 * </p> 65 * <p> 66 * <strong>Usage Notes</strong>: 67 * <ul> 68 * <li> 69 * Instance variables are used to maintain <code>RandomGenerator</code> and 70 * <code>SecureRandom</code> instances used in data generation. Therefore, to 71 * generate a random sequence of values or strings, you should use just 72 * <strong>one</strong> <code>RandomDataImpl</code> instance repeatedly.</li> 73 * <li> 74 * The "secure" methods are *much* slower. These should be used only when a 75 * cryptographically secure random sequence is required. A secure random 76 * sequence is a sequence of pseudo-random values which, in addition to being 77 * well-dispersed (so no subsequence of values is an any more likely than other 78 * subsequence of the the same length), also has the additional property that 79 * knowledge of values generated up to any point in the sequence does not make 80 * it any easier to predict subsequent values.</li> 81 * <li> 82 * When a new <code>RandomDataImpl</code> is created, the underlying random 83 * number generators are <strong>not</strong> initialized. If you do not 84 * explicitly seed the default non-secure generator, it is seeded with the 85 * current time in milliseconds on first use. The same holds for the secure 86 * generator. If you provide a <code>RandomGenerator</code> to the constructor, 87 * however, this generator is not reseeded by the constructor nor is it reseeded 88 * on first use.</li> 89 * <li> 90 * The <code>reSeed</code> and <code>reSeedSecure</code> methods delegate to the 91 * corresponding methods on the underlying <code>RandomGenerator</code> and 92 * <code>SecureRandom</code> instances. Therefore, <code>reSeed(long)</code> 93 * fully resets the initial state of the non-secure random number generator (so 94 * that reseeding with a specific value always results in the same subsequent 95 * random sequence); whereas reSeedSecure(long) does <strong>not</strong> 96 * reinitialize the secure random number generator (so secure sequences started 97 * with calls to reseedSecure(long) won't be identical).</li> 98 * <li> 99 * This implementation is not synchronized. 100 * </ul> 101 * </p> 102 * 103 * @version $Revision: 1061496 $ $Date: 2011-01-20 21:32:16 +0100 (jeu. 20 janv. 2011) $ 104 */ 105 public class RandomDataImpl implements RandomData, Serializable { 106 107 /** Serializable version identifier */ 108 private static final long serialVersionUID = -626730818244969716L; 109 110 /** underlying random number generator */ 111 private RandomGenerator rand = null; 112 113 /** underlying secure random number generator */ 114 private SecureRandom secRand = null; 115 116 /** 117 * Construct a RandomDataImpl. 118 */ 119 public RandomDataImpl() { 120 } 121 122 /** 123 * Construct a RandomDataImpl using the supplied {@link RandomGenerator} as 124 * the source of (non-secure) random data. 125 * 126 * @param rand 127 * the source of (non-secure) random data 128 * @since 1.1 129 */ 130 public RandomDataImpl(RandomGenerator rand) { 131 super(); 132 this.rand = rand; 133 } 134 135 /** 136 * {@inheritDoc} 137 * <p> 138 * <strong>Algorithm Description:</strong> hex strings are generated using a 139 * 2-step process. 140 * <ol> 141 * <li> 142 * len/2+1 binary bytes are generated using the underlying Random</li> 143 * <li> 144 * Each binary byte is translated into 2 hex digits</li> 145 * </ol> 146 * </p> 147 * 148 * @param len 149 * the desired string length. 150 * @return the random string. 151 * @throws NotStrictlyPositiveException if {@code len <= 0}. 152 */ 153 public String nextHexString(int len) { 154 if (len <= 0) { 155 throw new NotStrictlyPositiveException(LocalizedFormats.LENGTH, len); 156 } 157 158 // Get a random number generator 159 RandomGenerator ran = getRan(); 160 161 // Initialize output buffer 162 StringBuilder outBuffer = new StringBuilder(); 163 164 // Get int(len/2)+1 random bytes 165 byte[] randomBytes = new byte[(len / 2) + 1]; 166 ran.nextBytes(randomBytes); 167 168 // Convert each byte to 2 hex digits 169 for (int i = 0; i < randomBytes.length; i++) { 170 Integer c = Integer.valueOf(randomBytes[i]); 171 172 /* 173 * Add 128 to byte value to make interval 0-255 before doing hex 174 * conversion. This guarantees <= 2 hex digits from toHexString() 175 * toHexString would otherwise add 2^32 to negative arguments. 176 */ 177 String hex = Integer.toHexString(c.intValue() + 128); 178 179 // Make sure we add 2 hex digits for each byte 180 if (hex.length() == 1) { 181 hex = "0" + hex; 182 } 183 outBuffer.append(hex); 184 } 185 return outBuffer.toString().substring(0, len); 186 } 187 188 /** 189 * Generate a random int value uniformly distributed between 190 * <code>lower</code> and <code>upper</code>, inclusive. 191 * 192 * @param lower 193 * the lower bound. 194 * @param upper 195 * the upper bound. 196 * @return the random integer. 197 * @throws NumberIsTooLargeException if {@code lower >= upper}. 198 */ 199 public int nextInt(int lower, int upper) { 200 if (lower >= upper) { 201 throw new NumberIsTooLargeException(LocalizedFormats.LOWER_BOUND_NOT_BELOW_UPPER_BOUND, 202 lower, upper, false); 203 } 204 double r = getRan().nextDouble(); 205 return (int) ((r * upper) + ((1.0 - r) * lower) + r); 206 } 207 208 /** 209 * Generate a random long value uniformly distributed between 210 * <code>lower</code> and <code>upper</code>, inclusive. 211 * 212 * @param lower 213 * the lower bound. 214 * @param upper 215 * the upper bound. 216 * @return the random integer. 217 * @throws NumberIsTooLargeException if {@code lower >= upper}. 218 */ 219 public long nextLong(long lower, long upper) { 220 if (lower >= upper) { 221 throw new NumberIsTooLargeException(LocalizedFormats.LOWER_BOUND_NOT_BELOW_UPPER_BOUND, 222 lower, upper, false); 223 } 224 double r = getRan().nextDouble(); 225 return (long) ((r * upper) + ((1.0 - r) * lower) + r); 226 } 227 228 /** 229 * {@inheritDoc} 230 * <p> 231 * <strong>Algorithm Description:</strong> hex strings are generated in 232 * 40-byte segments using a 3-step process. 233 * <ol> 234 * <li> 235 * 20 random bytes are generated using the underlying 236 * <code>SecureRandom</code>.</li> 237 * <li> 238 * SHA-1 hash is applied to yield a 20-byte binary digest.</li> 239 * <li> 240 * Each byte of the binary digest is converted to 2 hex digits.</li> 241 * </ol> 242 * </p> 243 * 244 * @param len 245 * the length of the generated string 246 * @return the random string 247 * @throws NotStrictlyPositiveException if {@code len <= 0}. 248 */ 249 public String nextSecureHexString(int len) { 250 if (len <= 0) { 251 throw new NotStrictlyPositiveException(LocalizedFormats.LENGTH, len); 252 } 253 254 // Get SecureRandom and setup Digest provider 255 SecureRandom secRan = getSecRan(); 256 MessageDigest alg = null; 257 try { 258 alg = MessageDigest.getInstance("SHA-1"); 259 } catch (NoSuchAlgorithmException ex) { 260 // this should never happen 261 throw new MathInternalError(ex); 262 } 263 alg.reset(); 264 265 // Compute number of iterations required (40 bytes each) 266 int numIter = (len / 40) + 1; 267 268 StringBuilder outBuffer = new StringBuilder(); 269 for (int iter = 1; iter < numIter + 1; iter++) { 270 byte[] randomBytes = new byte[40]; 271 secRan.nextBytes(randomBytes); 272 alg.update(randomBytes); 273 274 // Compute hash -- will create 20-byte binary hash 275 byte hash[] = alg.digest(); 276 277 // Loop over the hash, converting each byte to 2 hex digits 278 for (int i = 0; i < hash.length; i++) { 279 Integer c = Integer.valueOf(hash[i]); 280 281 /* 282 * Add 128 to byte value to make interval 0-255 This guarantees 283 * <= 2 hex digits from toHexString() toHexString would 284 * otherwise add 2^32 to negative arguments 285 */ 286 String hex = Integer.toHexString(c.intValue() + 128); 287 288 // Keep strings uniform length -- guarantees 40 bytes 289 if (hex.length() == 1) { 290 hex = "0" + hex; 291 } 292 outBuffer.append(hex); 293 } 294 } 295 return outBuffer.toString().substring(0, len); 296 } 297 298 /** 299 * Generate a random int value uniformly distributed between 300 * <code>lower</code> and <code>upper</code>, inclusive. This algorithm uses 301 * a secure random number generator. 302 * 303 * @param lower 304 * the lower bound. 305 * @param upper 306 * the upper bound. 307 * @return the random integer. 308 * @throws NumberIsTooLargeException if {@code lower >= upper}. 309 */ 310 public int nextSecureInt(int lower, int upper) { 311 if (lower >= upper) { 312 throw new NumberIsTooLargeException(LocalizedFormats.LOWER_BOUND_NOT_BELOW_UPPER_BOUND, 313 lower, upper, false); 314 } 315 SecureRandom sec = getSecRan(); 316 return lower + (int) (sec.nextDouble() * (upper - lower + 1)); 317 } 318 319 /** 320 * Generate a random long value uniformly distributed between 321 * <code>lower</code> and <code>upper</code>, inclusive. This algorithm uses 322 * a secure random number generator. 323 * 324 * @param lower 325 * the lower bound. 326 * @param upper 327 * the upper bound. 328 * @return the random integer. 329 * @throws NumberIsTooLargeException if {@code lower >= upper}. 330 */ 331 public long nextSecureLong(long lower, long upper) { 332 if (lower >= upper) { 333 throw new NumberIsTooLargeException(LocalizedFormats.LOWER_BOUND_NOT_BELOW_UPPER_BOUND, 334 lower, upper, false); 335 } 336 SecureRandom sec = getSecRan(); 337 return lower + (long) (sec.nextDouble() * (upper - lower + 1)); 338 } 339 340 /** 341 * {@inheritDoc} 342 * <p> 343 * <strong>Algorithm Description</strong>: 344 * <ul><li> For small means, uses simulation of a Poisson process 345 * using Uniform deviates, as described 346 * <a href="http://irmi.epfl.ch/cmos/Pmmi/interactive/rng7.htm"> here.</a> 347 * The Poisson process (and hence value returned) is bounded by 1000 * mean.</li> 348 * 349 * <li> For large means, uses the rejection algorithm described in <br/> 350 * Devroye, Luc. (1981).<i>The Computer Generation of Poisson Random Variables</i> 351 * <strong>Computing</strong> vol. 26 pp. 197-207.</li></ul></p> 352 * 353 * @param mean mean of the Poisson distribution. 354 * @return the random Poisson value. 355 * @throws NotStrictlyPositiveException if {@code mean <= 0}. 356 */ 357 public long nextPoisson(double mean) { 358 if (mean <= 0) { 359 throw new NotStrictlyPositiveException(LocalizedFormats.MEAN, mean); 360 } 361 362 final RandomGenerator generator = getRan(); 363 364 final double pivot = 40.0d; 365 if (mean < pivot) { 366 double p = FastMath.exp(-mean); 367 long n = 0; 368 double r = 1.0d; 369 double rnd = 1.0d; 370 371 while (n < 1000 * mean) { 372 rnd = generator.nextDouble(); 373 r = r * rnd; 374 if (r >= p) { 375 n++; 376 } else { 377 return n; 378 } 379 } 380 return n; 381 } else { 382 final double lambda = FastMath.floor(mean); 383 final double lambdaFractional = mean - lambda; 384 final double logLambda = FastMath.log(lambda); 385 final double logLambdaFactorial = MathUtils.factorialLog((int) lambda); 386 final long y2 = lambdaFractional < Double.MIN_VALUE ? 0 : nextPoisson(lambdaFractional); 387 final double delta = FastMath.sqrt(lambda * FastMath.log(32 * lambda / FastMath.PI + 1)); 388 final double halfDelta = delta / 2; 389 final double twolpd = 2 * lambda + delta; 390 final double a1 = FastMath.sqrt(FastMath.PI * twolpd) * FastMath.exp(1 / 8 * lambda); 391 final double a2 = (twolpd / delta) * FastMath.exp(-delta * (1 + delta) / twolpd); 392 final double aSum = a1 + a2 + 1; 393 final double p1 = a1 / aSum; 394 final double p2 = a2 / aSum; 395 final double c1 = 1 / (8 * lambda); 396 397 double x = 0; 398 double y = 0; 399 double v = 0; 400 int a = 0; 401 double t = 0; 402 double qr = 0; 403 double qa = 0; 404 for (;;) { 405 final double u = nextUniform(0.0, 1); 406 if (u <= p1) { 407 final double n = nextGaussian(0d, 1d); 408 x = n * FastMath.sqrt(lambda + halfDelta) - 0.5d; 409 if (x > delta || x < -lambda) { 410 continue; 411 } 412 y = x < 0 ? FastMath.floor(x) : FastMath.ceil(x); 413 final double e = nextExponential(1d); 414 v = -e - (n * n / 2) + c1; 415 } else { 416 if (u > p1 + p2) { 417 y = lambda; 418 break; 419 } else { 420 x = delta + (twolpd / delta) * nextExponential(1d); 421 y = FastMath.ceil(x); 422 v = -nextExponential(1d) - delta * (x + 1) / twolpd; 423 } 424 } 425 a = x < 0 ? 1 : 0; 426 t = y * (y + 1) / (2 * lambda); 427 if (v < -t && a == 0) { 428 y = lambda + y; 429 break; 430 } 431 qr = t * ((2 * y + 1) / (6 * lambda) - 1); 432 qa = qr - (t * t) / (3 * (lambda + a * (y + 1))); 433 if (v < qa) { 434 y = lambda + y; 435 break; 436 } 437 if (v > qr) { 438 continue; 439 } 440 if (v < y * logLambda - MathUtils.factorialLog((int) (y + lambda)) + logLambdaFactorial) { 441 y = lambda + y; 442 break; 443 } 444 } 445 return y2 + (long) y; 446 } 447 } 448 449 /** 450 * Generate a random value from a Normal (a.k.a. Gaussian) distribution with 451 * the given mean, <code>mu</code> and the given standard deviation, 452 * <code>sigma</code>. 453 * 454 * @param mu 455 * the mean of the distribution 456 * @param sigma 457 * the standard deviation of the distribution 458 * @return the random Normal value 459 * @throws NotStrictlyPositiveException if {@code sigma <= 0}. 460 */ 461 public double nextGaussian(double mu, double sigma) { 462 if (sigma <= 0) { 463 throw new NotStrictlyPositiveException(LocalizedFormats.STANDARD_DEVIATION, sigma); 464 } 465 return sigma * getRan().nextGaussian() + mu; 466 } 467 468 /** 469 * Returns a random value from an Exponential distribution with the given 470 * mean. 471 * <p> 472 * <strong>Algorithm Description</strong>: Uses the <a 473 * href="http://www.jesus.ox.ac.uk/~clifford/a5/chap1/node5.html"> Inversion 474 * Method</a> to generate exponentially distributed random values from 475 * uniform deviates. 476 * </p> 477 * 478 * @param mean the mean of the distribution 479 * @return the random Exponential value 480 * @throws NotStrictlyPositiveException if {@code mean <= 0}. 481 */ 482 public double nextExponential(double mean) { 483 if (mean <= 0.0) { 484 throw new NotStrictlyPositiveException(LocalizedFormats.MEAN, mean); 485 } 486 final RandomGenerator generator = getRan(); 487 double unif = generator.nextDouble(); 488 while (unif == 0.0d) { 489 unif = generator.nextDouble(); 490 } 491 return -mean * FastMath.log(unif); 492 } 493 494 /** 495 * {@inheritDoc} 496 * <p> 497 * <strong>Algorithm Description</strong>: scales the output of 498 * Random.nextDouble(), but rejects 0 values (i.e., will generate another 499 * random double if Random.nextDouble() returns 0). This is necessary to 500 * provide a symmetric output interval (both endpoints excluded). 501 * </p> 502 * 503 * @param lower 504 * the lower bound. 505 * @param upper 506 * the upper bound. 507 * @return a uniformly distributed random value from the interval (lower, 508 * upper) 509 * @throws NumberIsTooLargeException if {@code lower >= upper}. 510 */ 511 public double nextUniform(double lower, double upper) { 512 if (lower >= upper) { 513 throw new NumberIsTooLargeException(LocalizedFormats.LOWER_BOUND_NOT_BELOW_UPPER_BOUND, 514 lower, upper, false); 515 } 516 final RandomGenerator generator = getRan(); 517 518 // ensure nextDouble() isn't 0.0 519 double u = generator.nextDouble(); 520 while (u <= 0.0) { 521 u = generator.nextDouble(); 522 } 523 524 return lower + u * (upper - lower); 525 } 526 527 /** 528 * Generates a random value from the {@link BetaDistributionImpl Beta Distribution}. 529 * This implementation uses {@link #nextInversionDeviate(ContinuousDistribution) inversion} 530 * to generate random values. 531 * 532 * @param alpha first distribution shape parameter 533 * @param beta second distribution shape parameter 534 * @return random value sampled from the beta(alpha, beta) distribution 535 * @throws MathException if an error occurs generating the random value 536 * @since 2.2 537 */ 538 public double nextBeta(double alpha, double beta) throws MathException { 539 return nextInversionDeviate(new BetaDistributionImpl(alpha, beta)); 540 } 541 542 /** 543 * Generates a random value from the {@link BinomialDistributionImpl Binomial Distribution}. 544 * This implementation uses {@link #nextInversionDeviate(ContinuousDistribution) inversion} 545 * to generate random values. 546 * 547 * @param numberOfTrials number of trials of the Binomial distribution 548 * @param probabilityOfSuccess probability of success of the Binomial distribution 549 * @return random value sampled from the Binomial(numberOfTrials, probabilityOfSuccess) distribution 550 * @throws MathException if an error occurs generating the random value 551 * @since 2.2 552 */ 553 public int nextBinomial(int numberOfTrials, double probabilityOfSuccess) throws MathException { 554 return nextInversionDeviate(new BinomialDistributionImpl(numberOfTrials, probabilityOfSuccess)); 555 } 556 557 /** 558 * Generates a random value from the {@link CauchyDistributionImpl Cauchy Distribution}. 559 * This implementation uses {@link #nextInversionDeviate(ContinuousDistribution) inversion} 560 * to generate random values. 561 * 562 * @param median the median of the Cauchy distribution 563 * @param scale the scale parameter of the Cauchy distribution 564 * @return random value sampled from the Cauchy(median, scale) distribution 565 * @throws MathException if an error occurs generating the random value 566 * @since 2.2 567 */ 568 public double nextCauchy(double median, double scale) throws MathException { 569 return nextInversionDeviate(new CauchyDistributionImpl(median, scale)); 570 } 571 572 /** 573 * Generates a random value from the {@link ChiSquaredDistributionImpl ChiSquare Distribution}. 574 * This implementation uses {@link #nextInversionDeviate(ContinuousDistribution) inversion} 575 * to generate random values. 576 * 577 * @param df the degrees of freedom of the ChiSquare distribution 578 * @return random value sampled from the ChiSquare(df) distribution 579 * @throws MathException if an error occurs generating the random value 580 * @since 2.2 581 */ 582 public double nextChiSquare(double df) throws MathException { 583 return nextInversionDeviate(new ChiSquaredDistributionImpl(df)); 584 } 585 586 /** 587 * Generates a random value from the {@link FDistributionImpl F Distribution}. 588 * This implementation uses {@link #nextInversionDeviate(ContinuousDistribution) inversion} 589 * to generate random values. 590 * 591 * @param numeratorDf the numerator degrees of freedom of the F distribution 592 * @param denominatorDf the denominator degrees of freedom of the F distribution 593 * @return random value sampled from the F(numeratorDf, denominatorDf) distribution 594 * @throws MathException if an error occurs generating the random value 595 * @since 2.2 596 */ 597 public double nextF(double numeratorDf, double denominatorDf) throws MathException { 598 return nextInversionDeviate(new FDistributionImpl(numeratorDf, denominatorDf)); 599 } 600 601 /** 602 * Generates a random value from the {@link GammaDistributionImpl Gamma Distribution}. 603 * This implementation uses {@link #nextInversionDeviate(ContinuousDistribution) inversion} 604 * to generate random values. 605 * 606 * @param shape the median of the Gamma distribution 607 * @param scale the scale parameter of the Gamma distribution 608 * @return random value sampled from the Gamma(shape, scale) distribution 609 * @throws MathException if an error occurs generating the random value 610 * @since 2.2 611 */ 612 public double nextGamma(double shape, double scale) throws MathException { 613 return nextInversionDeviate(new GammaDistributionImpl(shape, scale)); 614 } 615 616 /** 617 * Generates a random value from the {@link HypergeometricDistributionImpl Hypergeometric Distribution}. 618 * This implementation uses {@link #nextInversionDeviate(IntegerDistribution) inversion} 619 * to generate random values. 620 * 621 * @param populationSize the population size of the Hypergeometric distribution 622 * @param numberOfSuccesses number of successes in the population of the Hypergeometric distribution 623 * @param sampleSize the sample size of the Hypergeometric distribution 624 * @return random value sampled from the Hypergeometric(numberOfSuccesses, sampleSize) distribution 625 * @throws MathException if an error occurs generating the random value 626 * @since 2.2 627 */ 628 public int nextHypergeometric(int populationSize, int numberOfSuccesses, int sampleSize) throws MathException { 629 return nextInversionDeviate(new HypergeometricDistributionImpl(populationSize, numberOfSuccesses, sampleSize)); 630 } 631 632 /** 633 * Generates a random value from the {@link PascalDistributionImpl Pascal Distribution}. 634 * This implementation uses {@link #nextInversionDeviate(IntegerDistribution) inversion} 635 * to generate random values. 636 * 637 * @param r the number of successes of the Pascal distribution 638 * @param p the probability of success of the Pascal distribution 639 * @return random value sampled from the Pascal(r, p) distribution 640 * @throws MathException if an error occurs generating the random value 641 * @since 2.2 642 */ 643 public int nextPascal(int r, double p) throws MathException { 644 return nextInversionDeviate(new PascalDistributionImpl(r, p)); 645 } 646 647 /** 648 * Generates a random value from the {@link TDistributionImpl T Distribution}. 649 * This implementation uses {@link #nextInversionDeviate(ContinuousDistribution) inversion} 650 * to generate random values. 651 * 652 * @param df the degrees of freedom of the T distribution 653 * @return random value from the T(df) distribution 654 * @throws MathException if an error occurs generating the random value 655 * @since 2.2 656 */ 657 public double nextT(double df) throws MathException { 658 return nextInversionDeviate(new TDistributionImpl(df)); 659 } 660 661 /** 662 * Generates a random value from the {@link WeibullDistributionImpl Weibull Distribution}. 663 * This implementation uses {@link #nextInversionDeviate(ContinuousDistribution) inversion} 664 * to generate random values. 665 * 666 * @param shape the shape parameter of the Weibull distribution 667 * @param scale the scale parameter of the Weibull distribution 668 * @return random value sampled from the Weibull(shape, size) distribution 669 * @throws MathException if an error occurs generating the random value 670 * @since 2.2 671 */ 672 public double nextWeibull(double shape, double scale) throws MathException { 673 return nextInversionDeviate(new WeibullDistributionImpl(shape, scale)); 674 } 675 676 /** 677 * Generates a random value from the {@link ZipfDistributionImpl Zipf Distribution}. 678 * This implementation uses {@link #nextInversionDeviate(IntegerDistribution) inversion} 679 * to generate random values. 680 * 681 * @param numberOfElements the number of elements of the ZipfDistribution 682 * @param exponent the exponent of the ZipfDistribution 683 * @return random value sampled from the Zipf(numberOfElements, exponent) distribution 684 * @throws MathException if an error occurs generating the random value 685 * @since 2.2 686 */ 687 public int nextZipf(int numberOfElements, double exponent) throws MathException { 688 return nextInversionDeviate(new ZipfDistributionImpl(numberOfElements, exponent)); 689 } 690 691 /** 692 * Returns the RandomGenerator used to generate non-secure random data. 693 * <p> 694 * Creates and initializes a default generator if null. 695 * </p> 696 * 697 * @return the Random used to generate random data 698 * @since 1.1 699 */ 700 private RandomGenerator getRan() { 701 if (rand == null) { 702 rand = new JDKRandomGenerator(); 703 rand.setSeed(System.currentTimeMillis()); 704 } 705 return rand; 706 } 707 708 /** 709 * Returns the SecureRandom used to generate secure random data. 710 * <p> 711 * Creates and initializes if null. 712 * </p> 713 * 714 * @return the SecureRandom used to generate secure random data 715 */ 716 private SecureRandom getSecRan() { 717 if (secRand == null) { 718 secRand = new SecureRandom(); 719 secRand.setSeed(System.currentTimeMillis()); 720 } 721 return secRand; 722 } 723 724 /** 725 * Reseeds the random number generator with the supplied seed. 726 * <p> 727 * Will create and initialize if null. 728 * </p> 729 * 730 * @param seed 731 * the seed value to use 732 */ 733 public void reSeed(long seed) { 734 if (rand == null) { 735 rand = new JDKRandomGenerator(); 736 } 737 rand.setSeed(seed); 738 } 739 740 /** 741 * Reseeds the secure random number generator with the current time in 742 * milliseconds. 743 * <p> 744 * Will create and initialize if null. 745 * </p> 746 */ 747 public void reSeedSecure() { 748 if (secRand == null) { 749 secRand = new SecureRandom(); 750 } 751 secRand.setSeed(System.currentTimeMillis()); 752 } 753 754 /** 755 * Reseeds the secure random number generator with the supplied seed. 756 * <p> 757 * Will create and initialize if null. 758 * </p> 759 * 760 * @param seed 761 * the seed value to use 762 */ 763 public void reSeedSecure(long seed) { 764 if (secRand == null) { 765 secRand = new SecureRandom(); 766 } 767 secRand.setSeed(seed); 768 } 769 770 /** 771 * Reseeds the random number generator with the current time in 772 * milliseconds. 773 */ 774 public void reSeed() { 775 if (rand == null) { 776 rand = new JDKRandomGenerator(); 777 } 778 rand.setSeed(System.currentTimeMillis()); 779 } 780 781 /** 782 * Sets the PRNG algorithm for the underlying SecureRandom instance using 783 * the Security Provider API. The Security Provider API is defined in <a 784 * href = 785 * "http://java.sun.com/j2se/1.3/docs/guide/security/CryptoSpec.html#AppA"> 786 * Java Cryptography Architecture API Specification & Reference.</a> 787 * <p> 788 * <strong>USAGE NOTE:</strong> This method carries <i>significant</i> 789 * overhead and may take several seconds to execute. 790 * </p> 791 * 792 * @param algorithm 793 * the name of the PRNG algorithm 794 * @param provider 795 * the name of the provider 796 * @throws NoSuchAlgorithmException 797 * if the specified algorithm is not available 798 * @throws NoSuchProviderException 799 * if the specified provider is not installed 800 */ 801 public void setSecureAlgorithm(String algorithm, String provider) 802 throws NoSuchAlgorithmException, NoSuchProviderException { 803 secRand = SecureRandom.getInstance(algorithm, provider); 804 } 805 806 /** 807 * Generates an integer array of length <code>k</code> whose entries are 808 * selected randomly, without repetition, from the integers 809 * <code>0 through n-1</code> (inclusive). 810 * <p> 811 * Generated arrays represent permutations of <code>n</code> taken 812 * <code>k</code> at a time. 813 * </p> 814 * <p> 815 * <strong>Preconditions:</strong> 816 * <ul> 817 * <li> <code>k <= n</code></li> 818 * <li> <code>n > 0</code></li> 819 * </ul> 820 * If the preconditions are not met, an IllegalArgumentException is thrown. 821 * </p> 822 * <p> 823 * Uses a 2-cycle permutation shuffle. The shuffling process is described <a 824 * href="http://www.maths.abdn.ac.uk/~igc/tch/mx4002/notes/node83.html"> 825 * here</a>. 826 * </p> 827 * 828 * @param n 829 * domain of the permutation (must be positive) 830 * @param k 831 * size of the permutation (must satisfy 0 < k <= n). 832 * @return the random permutation as an int array 833 * @throws NumberIsTooLargeException if {@code k > n}. 834 * @throws NotStrictlyPositiveException if {@code k <= 0}. 835 */ 836 public int[] nextPermutation(int n, int k) { 837 if (k > n) { 838 throw new NumberIsTooLargeException(LocalizedFormats.PERMUTATION_EXCEEDS_N, 839 k, n, true); 840 } 841 if (k == 0) { 842 throw new NotStrictlyPositiveException(LocalizedFormats.PERMUTATION_SIZE, 843 k); 844 } 845 846 int[] index = getNatural(n); 847 shuffle(index, n - k); 848 int[] result = new int[k]; 849 for (int i = 0; i < k; i++) { 850 result[i] = index[n - i - 1]; 851 } 852 853 return result; 854 } 855 856 /** 857 * Uses a 2-cycle permutation shuffle to generate a random permutation. 858 * <strong>Algorithm Description</strong>: Uses a 2-cycle permutation 859 * shuffle to generate a random permutation of <code>c.size()</code> and 860 * then returns the elements whose indexes correspond to the elements of the 861 * generated permutation. This technique is described, and proven to 862 * generate random samples, <a 863 * href="http://www.maths.abdn.ac.uk/~igc/tch/mx4002/notes/node83.html"> 864 * here</a> 865 * 866 * @param c 867 * Collection to sample from. 868 * @param k 869 * sample size. 870 * @return the random sample. 871 * @throws NumberIsTooLargeException if {@code k > c.size()}. 872 * @throws NotStrictlyPositiveException if {@code k <= 0}. 873 */ 874 public Object[] nextSample(Collection<?> c, int k) { 875 int len = c.size(); 876 if (k > len) { 877 throw new NumberIsTooLargeException(LocalizedFormats.SAMPLE_SIZE_EXCEEDS_COLLECTION_SIZE, 878 k, len, true); 879 } 880 if (k <= 0) { 881 throw new NotStrictlyPositiveException(LocalizedFormats.NUMBER_OF_SAMPLES, k); 882 } 883 884 Object[] objects = c.toArray(); 885 int[] index = nextPermutation(len, k); 886 Object[] result = new Object[k]; 887 for (int i = 0; i < k; i++) { 888 result[i] = objects[index[i]]; 889 } 890 return result; 891 } 892 893 /** 894 * Generate a random deviate from the given distribution using the 895 * <a href="http://en.wikipedia.org/wiki/Inverse_transform_sampling"> inversion method.</a> 896 * 897 * @param distribution Continuous distribution to generate a random value from 898 * @return a random value sampled from the given distribution 899 * @throws MathException if an error occurs computing the inverse cumulative distribution function 900 * @since 2.2 901 */ 902 public double nextInversionDeviate(ContinuousDistribution distribution) throws MathException { 903 return distribution.inverseCumulativeProbability(nextUniform(0, 1)); 904 905 } 906 907 /** 908 * Generate a random deviate from the given distribution using the 909 * <a href="http://en.wikipedia.org/wiki/Inverse_transform_sampling"> inversion method.</a> 910 * 911 * @param distribution Integer distribution to generate a random value from 912 * @return a random value sampled from the given distribution 913 * @throws MathException if an error occurs computing the inverse cumulative distribution function 914 * @since 2.2 915 */ 916 public int nextInversionDeviate(IntegerDistribution distribution) throws MathException { 917 final double target = nextUniform(0, 1); 918 final int glb = distribution.inverseCumulativeProbability(target); 919 if (distribution.cumulativeProbability(glb) == 1.0d) { // No mass above 920 return glb; 921 } else { 922 return glb + 1; 923 } 924 } 925 926 // ------------------------Private methods---------------------------------- 927 928 /** 929 * Uses a 2-cycle permutation shuffle to randomly re-order the last elements 930 * of list. 931 * 932 * @param list 933 * list to be shuffled 934 * @param end 935 * element past which shuffling begins 936 */ 937 private void shuffle(int[] list, int end) { 938 int target = 0; 939 for (int i = list.length - 1; i >= end; i--) { 940 if (i == 0) { 941 target = 0; 942 } else { 943 target = nextInt(0, i); 944 } 945 int temp = list[target]; 946 list[target] = list[i]; 947 list[i] = temp; 948 } 949 } 950 951 /** 952 * Returns an array representing n. 953 * 954 * @param n 955 * the natural number to represent 956 * @return array with entries = elements of n 957 */ 958 private int[] getNatural(int n) { 959 int[] natural = new int[n]; 960 for (int i = 0; i < n; i++) { 961 natural[i] = i; 962 } 963 return natural; 964 } 965 966 } 967