You are viewing a plain text version of this content. The canonical link for it is here.
Posted to commits@mahout.apache.org by ra...@apache.org on 2018/06/29 16:10:51 UTC
[07/18] mahout git commit: MAHOUT-2033 Fixed Map-Reduce Refactor
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/math/Polynomial.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/math/Polynomial.java b/core/src/main/java/org/apache/mahout/math/jet/math/Polynomial.java
deleted file mode 100644
index 723e7d0..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/math/Polynomial.java
+++ /dev/null
@@ -1,98 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-/*
-Copyright 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.math;
-
-/**
- * Polynomial functions.
- */
-public final class Polynomial {
-
- private Polynomial() {
- }
-
- /**
- * Evaluates the given polynomial of degree <tt>N</tt> at <tt>x</tt>, assuming coefficient of N is 1.0. Otherwise same
- * as <tt>polevl()</tt>.
- * <pre>
- * 2 N
- * y = C + C x + C x +...+ C x
- * 0 1 2 N
- *
- * where C = 1 and hence is omitted from the array.
- * N
- *
- * Coefficients are stored in reverse order:
- *
- * coef[0] = C , ..., coef[N-1] = C .
- * N-1 0
- *
- * Calling arguments are otherwise the same as polevl().
- * </pre>
- * In the interest of speed, there are no checks for out of bounds arithmetic.
- *
- * @param x argument to the polynomial.
- * @param coef the coefficients of the polynomial.
- * @param N the degree of the polynomial.
- */
- public static double p1evl(double x, double[] coef, int N) {
-
- double ans = x + coef[0];
-
- for (int i = 1; i < N; i++) {
- ans = ans * x + coef[i];
- }
-
- return ans;
- }
-
- /**
- * Evaluates the given polynomial of degree <tt>N</tt> at <tt>x</tt>.
- * <pre>
- * 2 N
- * y = C + C x + C x +...+ C x
- * 0 1 2 N
- *
- * Coefficients are stored in reverse order:
- *
- * coef[0] = C , ..., coef[N] = C .
- * N 0
- * </pre>
- * In the interest of speed, there are no checks for out of bounds arithmetic.
- *
- * @param x argument to the polynomial.
- * @param coef the coefficients of the polynomial.
- * @param N the degree of the polynomial.
- */
- public static double polevl(double x, double[] coef, int N) {
- double ans = coef[0];
-
- for (int i = 1; i <= N; i++) {
- ans = ans * x + coef[i];
- }
-
- return ans;
- }
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/math/package-info.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/math/package-info.java b/core/src/main/java/org/apache/mahout/math/jet/math/package-info.java
deleted file mode 100644
index 3cda850..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/math/package-info.java
+++ /dev/null
@@ -1,5 +0,0 @@
-/**
- * Tools for basic and advanced mathematics: Arithmetics and Algebra, Polynomials and Chebyshev series, Bessel and Airy
- * functions, Function Objects for generic function evaluation, etc.
- */
-package org.apache.mahout.math.jet.math;
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/AbstractContinousDistribution.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/AbstractContinousDistribution.java b/core/src/main/java/org/apache/mahout/math/jet/random/AbstractContinousDistribution.java
deleted file mode 100644
index 8ca03d0..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/AbstractContinousDistribution.java
+++ /dev/null
@@ -1,51 +0,0 @@
-/**
- * Licensed to the Apache Software Foundation (ASF) under one
- * or more contributor license agreements. See the NOTICE file
- * distributed with this work for additional information
- * regarding copyright ownership. The ASF licenses this file
- * to you under the Apache License, Version 2.0 (the
- * "License"); you may not use this file except in compliance
- * with the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing,
- * software distributed under the License is distributed on an
- * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
- * KIND, either express or implied. See the License for the
- * specific language governing permissions and limitations
- * under the License.
- */
-/*
-Copyright 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random;
-
-/**
- * Abstract base class for all continuous distributions. Continuous distributions have
- * probability density and a cumulative distribution functions.
- *
- */
-public abstract class AbstractContinousDistribution extends AbstractDistribution {
- public double cdf(double x) {
- throw new UnsupportedOperationException("Can't compute pdf for " + this.getClass().getName());
- }
-
- public double pdf(double x) {
- throw new UnsupportedOperationException("Can't compute pdf for " + this.getClass().getName());
- }
-
- /**
- * @return A random number from the distribution; returns <tt>(int) Math.round(nextDouble())</tt>.
- * Override this method if necessary.
- */
- @Override
- public int nextInt() {
- return (int) Math.round(nextDouble());
- }
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/AbstractDiscreteDistribution.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/AbstractDiscreteDistribution.java b/core/src/main/java/org/apache/mahout/math/jet/random/AbstractDiscreteDistribution.java
deleted file mode 100644
index d93d76c..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/AbstractDiscreteDistribution.java
+++ /dev/null
@@ -1,27 +0,0 @@
-/*
-Copyright 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random;
-
-/**
- * Abstract base class for all discrete distributions.
- *
- */
-public abstract class AbstractDiscreteDistribution extends AbstractDistribution {
-
- /** Makes this class non instantiable, but still let's others inherit from it. */
- protected AbstractDiscreteDistribution() {
- }
-
- /** Returns a random number from the distribution; returns <tt>(double) nextInt()</tt>. */
- @Override
- public double nextDouble() {
- return nextInt();
- }
-
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/AbstractDistribution.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/AbstractDistribution.java b/core/src/main/java/org/apache/mahout/math/jet/random/AbstractDistribution.java
deleted file mode 100644
index 8e9cb0e..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/AbstractDistribution.java
+++ /dev/null
@@ -1,87 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-/*
-Copyright � 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random;
-
-import java.util.Random;
-
-import org.apache.mahout.math.function.DoubleFunction;
-import org.apache.mahout.math.function.IntFunction;
-
-public abstract class AbstractDistribution extends DoubleFunction implements IntFunction {
-
- private Random randomGenerator;
-
- /** Makes this class non instantiable, but still let's others inherit from it. */
- protected AbstractDistribution() {
- }
-
- protected Random getRandomGenerator() {
- return randomGenerator;
- }
-
- protected double randomDouble() {
- return randomGenerator.nextDouble();
- }
-
- /**
- * Equivalent to <tt>nextDouble()</tt>. This has the effect that distributions can now be used as function objects,
- * returning a random number upon function evaluation.
- */
- @Override
- public double apply(double dummy) {
- return nextDouble();
- }
-
- /**
- * Equivalent to <tt>nextInt()</tt>. This has the effect that distributions can now be used as function objects,
- * returning a random number upon function evaluation.
- */
- @Override
- public int apply(int dummy) {
- return nextInt();
- }
-
- /**
- * Returns a random number from the distribution.
- * @return A new sample from this distribution.
- */
- public abstract double nextDouble();
-
- /**
- * @return
- * A random number from the distribution; returns <tt>(int) Math.round(nextDouble())</tt>. Override this
- * method if necessary.
- */
- public abstract int nextInt();
-
- /**
- * Sets the uniform random generator internally used.
- * @param randomGenerator the new PRNG
- */
- public void setRandomGenerator(Random randomGenerator) {
- this.randomGenerator = randomGenerator;
- }
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/Exponential.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/Exponential.java b/core/src/main/java/org/apache/mahout/math/jet/random/Exponential.java
deleted file mode 100644
index 06472c2..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/Exponential.java
+++ /dev/null
@@ -1,81 +0,0 @@
-/*
-Copyright � 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random;
-
-import java.util.Locale;
-import java.util.Random;
-
-public class Exponential extends AbstractContinousDistribution {
- // rate parameter for the distribution. Mean is 1/lambda.
- private double lambda;
-
- /**
- * Provides a negative exponential distribution given a rate parameter lambda and an underlying
- * random number generator. The mean of this distribution will be equal to 1/lambda.
- *
- * @param lambda The rate parameter of the distribution.
- * @param randomGenerator The PRNG that is used to generate values.
- */
- public Exponential(double lambda, Random randomGenerator) {
- setRandomGenerator(randomGenerator);
- this.lambda = lambda;
- }
-
- /**
- * Returns the cumulative distribution function.
- * @param x The point at which the cumulative distribution function is to be evaluated.
- * @return Returns the integral from -infinity to x of the PDF, also known as the cumulative distribution
- * function.
- */
- @Override
- public double cdf(double x) {
- if (x <= 0.0) {
- return 0.0;
- }
- return 1.0 - Math.exp(-x * lambda);
- }
-
- /**
- * Returns a random number from the distribution.
- */
- @Override
- public double nextDouble() {
- return -Math.log1p(-randomDouble()) / lambda;
- }
-
- /**
- * Returns the value of the probability density function at a particular point.
- * @param x The point at which the probability density function is to be evaluated.
- * @return The value of the probability density function at the specified point.
- */
- @Override
- public double pdf(double x) {
- if (x < 0.0) {
- return 0.0;
- }
- return lambda * Math.exp(-x * lambda);
- }
-
- /**
- * Sets the rate parameter.
- * @param lambda The new value of the rate parameter.
- */
- public void setState(double lambda) {
- this.lambda = lambda;
- }
-
- /**
- * Returns a String representation of the receiver.
- */
- @Override
- public String toString() {
- return String.format(Locale.ENGLISH, "%s(%.4f)", this.getClass().getName(), lambda);
- }
-
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/Gamma.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/Gamma.java b/core/src/main/java/org/apache/mahout/math/jet/random/Gamma.java
deleted file mode 100644
index 914157b..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/Gamma.java
+++ /dev/null
@@ -1,302 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-/*
-Copyright � 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random;
-
-import org.apache.mahout.math.jet.stat.Probability;
-
-import java.util.Random;
-
-public class Gamma extends AbstractContinousDistribution {
- // shape
- private final double alpha;
-
- // rate
- private final double rate;
-
- /**
- * Constructs a Gamma distribution with a given shape (alpha) and rate (beta).
- *
- * @param alpha The shape parameter.
- * @param rate The rate parameter.
- * @param randomGenerator The random number generator that generates bits for us.
- * @throws IllegalArgumentException if <tt>alpha <= 0.0 || alpha <= 0.0</tt>.
- */
- public Gamma(double alpha, double rate, Random randomGenerator) {
- this.alpha = alpha;
- this.rate = rate;
- setRandomGenerator(randomGenerator);
- }
-
- /**
- * Returns the cumulative distribution function.
- * @param x The end-point where the cumulation should end.
- */
- @Override
- public double cdf(double x) {
- return Probability.gamma(alpha, rate, x);
- }
-
- /** Returns a random number from the distribution. */
- @Override
- public double nextDouble() {
- return nextDouble(alpha, rate);
- }
-
- /** Returns a random number from the distribution; bypasses the internal state.
- * *
- * Gamma Distribution - Acceptance Rejection combined with *
- * Acceptance Complement *
- * *
- ******************************************************************
- * *
- * FUNCTION: - gds samples a random number from the standard *
- * gamma distribution with parameter a > 0. *
- * Acceptance Rejection gs for a < 1 , *
- * Acceptance Complement gd for a >= 1 . *
- * REFERENCES: - J.H. Ahrens, U. Dieter (1974): Computer methods *
- * for sampling from gamma, beta, Poisson and *
- * binomial distributions, Computing 12, 223-246. *
- * - J.H. Ahrens, U. Dieter (1982): Generating gamma *
- * variates by a modified rejection technique, *
- * Communications of the ACM 25, 47-54. *
- * SUBPROGRAMS: - drand(seed) ... (0,1)-Uniform generator with *
- * unsigned long integer *seed *
- * - NORMAL(seed) ... Normal generator N(0,1). *
- * *
- * @param alpha Shape parameter.
- * @param rate Rate parameter (=1/scale).
- * @return A gamma distributed sample.
- */
- public double nextDouble(double alpha, double rate) {
- if (alpha <= 0.0) {
- throw new IllegalArgumentException();
- }
- if (rate <= 0.0) {
- throw new IllegalArgumentException();
- }
-
- double gds;
- double b = 0.0;
- if (alpha < 1.0) { // CASE A: Acceptance rejection algorithm gs
- b = 1.0 + 0.36788794412 * alpha; // Step 1
- while (true) {
- double p = b * randomDouble();
- if (p <= 1.0) { // Step 2. Case gds <= 1
- gds = Math.exp(Math.log(p) / alpha);
- if (Math.log(randomDouble()) <= -gds) {
- return gds / rate;
- }
- } else { // Step 3. Case gds > 1
- gds = -Math.log((b - p) / alpha);
- if (Math.log(randomDouble()) <= (alpha - 1.0) * Math.log(gds)) {
- return gds / rate;
- }
- }
- }
- } else { // CASE B: Acceptance complement algorithm gd (gaussian distribution, box muller transformation)
- double ss = 0.0;
- double s = 0.0;
- double d = 0.0;
- if (alpha != -1.0) { // Step 1. Preparations
- ss = alpha - 0.5;
- s = Math.sqrt(ss);
- d = 5.656854249 - 12.0 * s;
- }
- // Step 2. Normal deviate
- double v12;
- double v1;
- do {
- v1 = 2.0 * randomDouble() - 1.0;
- double v2 = 2.0 * randomDouble() - 1.0;
- v12 = v1 * v1 + v2 * v2;
- } while (v12 > 1.0);
- double t = v1 * Math.sqrt(-2.0 * Math.log(v12) / v12);
- double x = s + 0.5 * t;
- gds = x * x;
- if (t >= 0.0) {
- return gds / rate;
- } // Immediate acceptance
-
- double u = randomDouble();
- if (d * u <= t * t * t) {
- return gds / rate;
- } // Squeeze acceptance
-
- double q0 = 0.0;
- double si = 0.0;
- double c = 0.0;
- if (alpha != -1.0) { // Step 4. Set-up for hat case
- double r = 1.0 / alpha;
- double q9 = 0.0001710320;
- double q8 = -0.0004701849;
- double q7 = 0.0006053049;
- double q6 = 0.0003340332;
- double q5 = -0.0003349403;
- double q4 = 0.0015746717;
- double q3 = 0.0079849875;
- double q2 = 0.0208333723;
- double q1 = 0.0416666664;
- q0 = ((((((((q9 * r + q8) * r + q7) * r + q6) * r + q5) * r + q4) * r + q3) * r + q2) * r + q1) * r;
- if (alpha > 3.686) {
- if (alpha > 13.022) {
- b = 1.77;
- si = 0.75;
- c = 0.1515 / s;
- } else {
- b = 1.654 + 0.0076 * ss;
- si = 1.68 / s + 0.275;
- c = 0.062 / s + 0.024;
- }
- } else {
- b = 0.463 + s - 0.178 * ss;
- si = 1.235;
- c = 0.195 / s - 0.079 + 0.016 * s;
- }
- }
- double v;
- double q;
- double a9 = 0.104089866;
- double a8 = -0.112750886;
- double a7 = 0.110368310;
- double a6 = -0.124385581;
- double a5 = 0.142873973;
- double a4 = -0.166677482;
- double a3 = 0.199999867;
- double a2 = -0.249999949;
- double a1 = 0.333333333;
- if (x > 0.0) { // Step 5. Calculation of q
- v = t / (s + s); // Step 6.
- if (Math.abs(v) > 0.25) {
- q = q0 - s * t + 0.25 * t * t + (ss + ss) * Math.log1p(v);
- } else {
- q = q0 + 0.5 * t * t * ((((((((a9 * v + a8) * v + a7) * v + a6)
- * v + a5) * v + a4) * v + a3) * v + a2) * v + a1) * v;
- } // Step 7. Quotient acceptance
- if (Math.log1p(-u) <= q) {
- return gds / rate;
- }
- }
-
- double e7 = 0.000247453;
- double e6 = 0.001353826;
- double e5 = 0.008345522;
- double e4 = 0.041664508;
- double e3 = 0.166666848;
- double e2 = 0.499999994;
- double e1 = 1.000000000;
- while (true) { // Step 8. Double exponential deviate t
- double sign_u;
- double e;
- do {
- e = -Math.log(randomDouble());
- u = randomDouble();
- u = u + u - 1.0;
- sign_u = u > 0 ? 1.0 : -1.0;
- t = b + e * si * sign_u;
- } while (t <= -0.71874483771719); // Step 9. Rejection of t
- v = t / (s + s); // Step 10. New q(t)
- if (Math.abs(v) > 0.25) {
- q = q0 - s * t + 0.25 * t * t + (ss + ss) * Math.log1p(v);
- } else {
- q = q0 + 0.5 * t * t * ((((((((a9 * v + a8) * v + a7) * v + a6)
- * v + a5) * v + a4) * v + a3) * v + a2) * v + a1) * v;
- }
- if (q <= 0.0) {
- continue;
- } // Step 11.
- double w;
- if (q > 0.5) {
- w = Math.exp(q) - 1.0;
- } else {
- w = ((((((e7 * q + e6) * q + e5) * q + e4) * q + e3) * q + e2) * q + e1) * q;
- } // Step 12. Hat acceptance
- if (c * u * sign_u <= w * Math.exp(e - 0.5 * t * t)) {
- x = s + 0.5 * t;
- return x * x / rate;
- }
- }
- }
- }
-
- /** Returns the probability distribution function.
- * @param x Where to compute the density function.
- *
- * @return The value of the gamma density at x.
- */
- @Override
- public double pdf(double x) {
- if (x < 0) {
- throw new IllegalArgumentException();
- }
- if (x == 0) {
- if (alpha == 1.0) {
- return rate;
- } else if (alpha < 1) {
- return Double.POSITIVE_INFINITY;
- } else {
- return 0;
- }
- }
- if (alpha == 1.0) {
- return rate * Math.exp(-x * rate);
- }
- return rate * Math.exp((alpha - 1.0) * Math.log(x * rate) - x * rate - logGamma(alpha));
- }
-
- @Override
- public String toString() {
- return this.getClass().getName() + '(' + rate + ',' + alpha + ')';
- }
-
- /** Returns a quick approximation of <tt>log(gamma(x))</tt>. */
- public static double logGamma(double x) {
-
- if (x <= 0.0 /* || x > 1.3e19 */) {
- return -999;
- }
-
- double z;
- for (z = 1.0; x < 11.0; x++) {
- z *= x;
- }
-
- double r = 1.0 / (x * x);
- double c6 = -1.9175269175269175e-03;
- double c5 = 8.4175084175084175e-04;
- double c4 = -5.9523809523809524e-04;
- double c3 = 7.9365079365079365e-04;
- double c2 = -2.7777777777777777e-03;
- double c1 = 8.3333333333333333e-02;
- double g = c1 + r * (c2 + r * (c3 + r * (c4 + r * (c5 + r + c6))));
- double c0 = 9.1893853320467274e-01;
- g = (x - 0.5) * Math.log(x) - x + c0 + g / x;
- if (z == 1.0) {
- return g;
- }
- return g - Math.log(z);
- }
-
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/NegativeBinomial.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/NegativeBinomial.java b/core/src/main/java/org/apache/mahout/math/jet/random/NegativeBinomial.java
deleted file mode 100644
index 1e577eb..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/NegativeBinomial.java
+++ /dev/null
@@ -1,106 +0,0 @@
-/*
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-/*
-Copyright � 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random;
-
-import org.apache.mahout.math.jet.math.Arithmetic;
-import org.apache.mahout.math.jet.stat.Probability;
-
-import java.util.Random;
-
-/** Mostly deprecated until unit tests are in place. Until this time, this class/interface is unsupported. */
-public final class NegativeBinomial extends AbstractDiscreteDistribution {
-
- private final int r;
- private final double p;
-
- private final Gamma gamma;
- private final Poisson poisson;
-
- /**
- * Constructs a Negative Binomial distribution which describes the probability of getting
- * a particular number of negative trials (k) before getting a fixed number of positive
- * trials (r) where each positive trial has probability (p) of being successful.
- *
- * @param r the required number of positive trials.
- * @param p the probability of success.
- * @param randomGenerator a uniform random number generator.
- */
- public NegativeBinomial(int r, double p, Random randomGenerator) {
- setRandomGenerator(randomGenerator);
- this.r = r;
- this.p = p;
- this.gamma = new Gamma(r, 1, randomGenerator);
- this.poisson = new Poisson(0.0, randomGenerator);
- }
-
- /**
- * Returns the cumulative distribution function.
- */
- public double cdf(int k) {
- return Probability.negativeBinomial(k, r, p);
- }
-
- /**
- * Returns the probability distribution function.
- */
- public double pdf(int k) {
- return Arithmetic.binomial(k + r - 1, r - 1) * Math.pow(p, r) * Math.pow(1.0 - p, k);
- }
-
- @Override
- public int nextInt() {
- return nextInt(r, p);
- }
-
- /**
- * Returns a sample from this distribution. The value returned will
- * be the number of negative samples required before achieving r
- * positive samples. Each successive sample is taken independently
- * from a Bernouli process with probability p of success.
- *
- * The algorithm used is taken from J.H. Ahrens, U. Dieter (1974):
- * Computer methods for sampling from gamma, beta, Poisson and
- * binomial distributions, Computing 12, 223--246.
- *
- * This algorithm is essentially the same as described at
- * http://en.wikipedia.org/wiki/Negative_binomial_distribution#Gamma.E2.80.93Poisson_mixture
- * except that the notion of positive and negative outcomes is uniformly
- * inverted. Because the inversion is complete and consistent, this
- * definition is effectively identical to that defined on wikipedia.
- */
- public int nextInt(int r, double p) {
- return this.poisson.nextInt(gamma.nextDouble(r, p / (1.0 - p)));
- }
-
- /**
- * Returns a String representation of the receiver.
- */
- @Override
- public String toString() {
- return this.getClass().getName() + '(' + r + ',' + p + ')';
- }
-
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/Normal.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/Normal.java b/core/src/main/java/org/apache/mahout/math/jet/random/Normal.java
deleted file mode 100644
index 7ceac22..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/Normal.java
+++ /dev/null
@@ -1,110 +0,0 @@
-/*
-Copyright � 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random;
-
-import org.apache.mahout.math.jet.stat.Probability;
-
-import java.util.Locale;
-import java.util.Random;
-
-/**
- * Implements a normal distribution specified mean and standard deviation.
- */
-public class Normal extends AbstractContinousDistribution {
-
- private double mean;
- private double variance;
- private double standardDeviation;
-
- private double cache; // cache for Box-Mueller algorithm
- private boolean cacheFilled; // Box-Mueller
-
- private double normalizer; // performance cache
-
- /**
- * @param mean The mean of the resulting distribution.
- * @param standardDeviation The standard deviation of the distribution.
- * @param randomGenerator The random number generator to use. This can be null if you don't
- * need to generate any numbers.
- */
- public Normal(double mean, double standardDeviation, Random randomGenerator) {
- setRandomGenerator(randomGenerator);
- setState(mean, standardDeviation);
- }
-
- /**
- * Returns the cumulative distribution function.
- */
- @Override
- public double cdf(double x) {
- return Probability.normal(mean, variance, x);
- }
-
- /** Returns the probability density function. */
- @Override
- public double pdf(double x) {
- double diff = x - mean;
- return normalizer * Math.exp(-(diff * diff) / (2.0 * variance));
- }
-
- /**
- * Returns a random number from the distribution.
- */
- @Override
- public double nextDouble() {
- // Uses polar Box-Muller transformation.
- if (cacheFilled) {
- cacheFilled = false;
- return cache;
- }
-
- double x;
- double y;
- double r;
- do {
- x = 2.0 * randomDouble() - 1.0;
- y = 2.0 * randomDouble() - 1.0;
- r = x * x + y * y;
- } while (r >= 1.0);
-
- double z = Math.sqrt(-2.0 * Math.log(r) / r);
- cache = this.mean + this.standardDeviation * x * z;
- cacheFilled = true;
- return this.mean + this.standardDeviation * y * z;
- }
-
- /** Sets the uniform random generator internally used. */
- @Override
- public final void setRandomGenerator(Random randomGenerator) {
- super.setRandomGenerator(randomGenerator);
- this.cacheFilled = false;
- }
-
- /**
- * Sets the mean and variance.
- * @param mean The new value for the mean.
- * @param standardDeviation The new value for the standard deviation.
- */
- public final void setState(double mean, double standardDeviation) {
- if (mean != this.mean || standardDeviation != this.standardDeviation) {
- this.mean = mean;
- this.standardDeviation = standardDeviation;
- this.variance = standardDeviation * standardDeviation;
- this.cacheFilled = false;
-
- this.normalizer = 1.0 / Math.sqrt(2.0 * Math.PI * variance);
- }
- }
-
- /** Returns a String representation of the receiver. */
- @Override
- public String toString() {
- return String.format(Locale.ENGLISH, "%s(m=%f, sd=%f)", this.getClass().getName(), mean, standardDeviation);
- }
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/Poisson.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/Poisson.java b/core/src/main/java/org/apache/mahout/math/jet/random/Poisson.java
deleted file mode 100644
index 497691e..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/Poisson.java
+++ /dev/null
@@ -1,296 +0,0 @@
-/*
-Copyright � 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random;
-
-import org.apache.mahout.math.jet.math.Arithmetic;
-
-import java.util.Random;
-
-/** Partially deprecated until unit tests are in place. Until this time, this class/interface is unsupported. */
-public final class Poisson extends AbstractDiscreteDistribution {
-
- private final double mean;
-
- // precomputed and cached values (for performance only)
- // cache for < SWITCH_MEAN
- private double myOld = -1.0;
- private double p;
- private double q;
- private double p0;
- private final double[] pp = new double[36];
- private int llll;
-
- // cache for >= SWITCH_MEAN
- private double myLast = -1.0;
- private double ll;
- private int k2;
- private int k4;
- private int k1;
- private int k5;
- private double dl;
- private double dr;
- private double r1;
- private double r2;
- private double r4;
- private double r5;
- private double lr;
- private double lMy;
- private double cPm;
- private double f1;
- private double f2;
- private double f4;
- private double f5;
- private double p1;
- private double p2;
- private double p3;
- private double p4;
- private double p5;
- private double p6;
-
- // cache for both;
-
-
- private static final double MEAN_MAX = Integer.MAX_VALUE;
- // for all means larger than that, we don't try to compute a poisson deviation, but return the mean.
- private static final double SWITCH_MEAN = 10.0; // switch from method A to method B
-
-
- /** Constructs a poisson distribution. Example: mean=1.0. */
- public Poisson(double mean, Random randomGenerator) {
- setRandomGenerator(randomGenerator);
- this.mean = mean;
- }
-
- private static double f(int k, double lNu, double cPm) {
- return Math.exp(k * lNu - Arithmetic.logFactorial(k) - cPm);
- }
-
- @Override
- public int nextInt() {
- return nextInt(mean);
- }
-
- /** Returns a random number from the distribution; bypasses the internal state. */
- public int nextInt(double theMean) {
- /******************************************************************
- * *
- * Poisson Distribution - Patchwork Rejection/Inversion *
- * *
- ******************************************************************
- * *
- * For parameter my < 10 Tabulated Inversion is applied. *
- * For my >= 10 Patchwork Rejection is employed: *
- * The area below the histogram function f(x) is rearranged in *
- * its body by certain point reflections. Within a large center *
- * interval variates are sampled efficiently by rejection from *
- * uniform hats. Rectangular immediate acceptance regions speed *
- * up the generation. The remaining tails are covered by *
- * exponential functions. *
- * *
- *****************************************************************/
- Random gen = getRandomGenerator();
-
- //double t, g, my_k;
-
- //double gx, gy, px, py, e, x, xx, delta, v;
- //int sign;
-
- //static double p,q,p0,pp[36];
- //static long ll,m;
-
- int m;
- if (theMean < SWITCH_MEAN) { // CASE B: Inversion- start new table and calculate p0
- if (theMean != myOld) {
- myOld = theMean;
- llll = 0;
- p = Math.exp(-theMean);
- q = p;
- p0 = p;
- //for (k=pp.length; --k >=0;) pp[k] = 0;
- }
- m = theMean > 1.0 ? (int) theMean : 1;
- while (true) {
- double u = gen.nextDouble();
- int k = 0;
- if (u <= p0) {
- return k;
- }
- if (llll != 0) { // Step T. Table comparison
- int i = u > 0.458 ? Math.min(llll, m) : 1;
- for (k = i; k <= llll; k++) {
- if (u <= pp[k]) {
- return k;
- }
- }
- if (llll == 35) {
- continue;
- }
- }
- for (k = llll + 1; k <= 35; k++) { // Step C. Creation of new prob.
- p *= theMean / k;
- q += p;
- pp[k] = q;
- if (u <= q) {
- llll = k;
- return k;
- }
- }
- llll = 35;
- }
- // end my < SWITCH_MEAN
- } else if (theMean < MEAN_MAX) { // CASE A: acceptance complement
- //static double my_last = -1.0;
- //static long int m, k2, k4, k1, k5;
- //static double dl, dr, r1, r2, r4, r5, ll, lr, l_my, c_pm,
- // f1, f2, f4, f5, p1, p2, p3, p4, p5, p6;
-
- m = (int) theMean;
- if (theMean != myLast) { // set-up
- myLast = theMean;
-
- // approximate deviation of reflection points k2, k4 from my - 1/2
- double Ds = Math.sqrt(theMean + 0.25);
-
- // mode m, reflection points k2 and k4, and points k1 and k5, which
- // delimit the centre region of h(x)
- k2 = (int) Math.ceil(theMean - 0.5 - Ds);
- k4 = (int) (theMean - 0.5 + Ds);
- k1 = k2 + k2 - m + 1;
- k5 = k4 + k4 - m;
-
- // range width of the critical left and right centre region
- dl = k2 - k1;
- dr = k5 - k4;
-
- // recurrence constants r(k) = p(k)/p(k-1) at k = k1, k2, k4+1, k5+1
- r1 = theMean / k1;
- r2 = theMean / k2;
- r4 = theMean / (k4 + 1);
- r5 = theMean / (k5 + 1);
-
- // reciprocal values of the scale parameters of expon. tail envelopes
- ll = Math.log(r1); // expon. tail left
- lr = -Math.log(r5); // expon. tail right
-
- // Poisson constants, necessary for computing function values f(k)
- lMy = Math.log(theMean);
- cPm = m * lMy - Arithmetic.logFactorial(m);
-
- // function values f(k) = p(k)/p(m) at k = k2, k4, k1, k5
- f2 = f(k2, lMy, cPm);
- f4 = f(k4, lMy, cPm);
- f1 = f(k1, lMy, cPm);
- f5 = f(k5, lMy, cPm);
-
- // area of the two centre and the two exponential tail regions
- // area of the two immediate acceptance regions between k2, k4
- p1 = f2 * (dl + 1.0); // immed. left
- p2 = f2 * dl + p1; // centre left
- p3 = f4 * (dr + 1.0) + p2; // immed. right
- p4 = f4 * dr + p3; // centre right
- p5 = f1 / ll + p4; // expon. tail left
- p6 = f5 / lr + p5; // expon. tail right
- } // end set-up
-
- while (true) {
- // generate uniform number U -- U(0, p6)
- // case distinction corresponding to U
- double W;
- double V;
- double U;
- int Y;
- int X;
- int Dk;
- if ((U = gen.nextDouble() * p6) < p2) { // centre left
-
- // immediate acceptance region R2 = [k2, m) *[0, f2), X = k2, ... m -1
- if ((V = U - p1) < 0.0) {
- return k2 + (int) (U / f2);
- }
- // immediate acceptance region R1 = [k1, k2)*[0, f1), X = k1, ... k2-1
- if ((W = V / dl) < f1) {
- return k1 + (int) (V / f1);
- }
-
- // computation of candidate X < k2, and its counterpart Y > k2
- // either squeeze-acceptance of X or acceptance-rejection of Y
- Dk = gen.nextInt((int) dl) + 1;
- if (W <= f2 - Dk * (f2 - f2 / r2)) { // quick accept of
- return k2 - Dk; // X = k2 - Dk
- }
- if ((V = f2 + f2 - W) < 1.0) { // quick reject of Y
- Y = k2 + Dk;
- if (V <= f2 + Dk * (1.0 - f2) / (dl + 1.0)) { // quick accept of
- return Y; // Y = k2 + Dk
- }
- if (V <= f(Y, lMy, cPm)) {
- return Y;
- } // final accept of Y
- }
- X = k2 - Dk;
- } else if (U < p4) { // centre right
- // immediate acceptance region R3 = [m, k4+1)*[0, f4), X = m, ... k4
- if ((V = U - p3) < 0.0) {
- return k4 - (int) ((U - p2) / f4);
- }
- // immediate acceptance region R4 = [k4+1, k5+1)*[0, f5)
- if ((W = V / dr) < f5) {
- return k5 - (int) (V / f5);
- }
-
- // computation of candidate X > k4, and its counterpart Y < k4
- // either squeeze-acceptance of X or acceptance-rejection of Y
- Dk = gen.nextInt((int) dr) + 1;
- if (W <= f4 - Dk * (f4 - f4 * r4)) { // quick accept of
- return k4 + Dk; // X = k4 + Dk
- }
- if ((V = f4 + f4 - W) < 1.0) { // quick reject of Y
- Y = k4 - Dk;
- if (V <= f4 + Dk * (1.0 - f4) / dr) { // quick accept of
- return Y; // Y = k4 - Dk
- }
- if (V <= f(Y, lMy, cPm)) {
- return Y;
- } // final accept of Y
- }
- X = k4 + Dk;
- } else {
- W = gen.nextDouble();
- if (U < p5) { // expon. tail left
- Dk = (int) (1.0 - Math.log(W) / ll);
- if ((X = k1 - Dk) < 0) {
- continue;
- } // 0 <= X <= k1 - 1
- W *= (U - p4) * ll; // W -- U(0, h(x))
- if (W <= f1 - Dk * (f1 - f1 / r1)) {
- return X;
- } // quick accept of X
- } else { // expon. tail right
- Dk = (int) (1.0 - Math.log(W) / lr);
- X = k5 + Dk; // X >= k5 + 1
- W *= (U - p5) * lr; // W -- U(0, h(x))
- if (W <= f5 - Dk * (f5 - f5 * r5)) {
- return X;
- } // quick accept of X
- }
- }
-
- // acceptance-rejection test of candidate X from the original area
- // test, whether W <= f(k), with W = U*h(x) and U -- U(0, 1)
- // log f(X) = (X - m)*log(my) - log X! + log m!
- if (Math.log(W) <= X * lMy - Arithmetic.logFactorial(X) - cPm) {
- return X;
- }
- }
- } else { // mean is too large
- return (int) theMean;
- }
- }
-
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/Uniform.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/Uniform.java b/core/src/main/java/org/apache/mahout/math/jet/random/Uniform.java
deleted file mode 100644
index 32c8b90..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/Uniform.java
+++ /dev/null
@@ -1,164 +0,0 @@
-/*
-Copyright � 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random;
-
-import org.apache.mahout.common.RandomUtils;
-
-import java.util.Random;
-
-public class Uniform extends AbstractContinousDistribution {
-
- private double min;
- private double max;
-
- /**
- * Constructs a uniform distribution with the given minimum and maximum, using a {@link
- * org.apache.mahout.math.jet.random.engine.MersenneTwister} seeded with the given seed.
- */
- public Uniform(double min, double max, int seed) {
- this(min, max, RandomUtils.getRandom(seed));
- }
-
- /** Constructs a uniform distribution with the given minimum and maximum. */
- public Uniform(double min, double max, Random randomGenerator) {
- setRandomGenerator(randomGenerator);
- setState(min, max);
- }
-
- /** Constructs a uniform distribution with <tt>min=0.0</tt> and <tt>max=1.0</tt>. */
- public Uniform(Random randomGenerator) {
- this(0, 1, randomGenerator);
- }
-
- /** Returns the cumulative distribution function (assuming a continous uniform distribution). */
- @Override
- public double cdf(double x) {
- if (x <= min) {
- return 0.0;
- }
- if (x >= max) {
- return 1.0;
- }
- return (x - min) / (max - min);
- }
-
- /** Returns a uniformly distributed random <tt>boolean</tt>. */
- public boolean nextBoolean() {
- return randomDouble() > 0.5;
- }
-
- /**
- * Returns a uniformly distributed random number in the open interval <tt>(min,max)</tt> (excluding <tt>min</tt> and
- * <tt>max</tt>).
- */
- @Override
- public double nextDouble() {
- return min + (max - min) * randomDouble();
- }
-
- /**
- * Returns a uniformly distributed random number in the open interval <tt>(from,to)</tt> (excluding <tt>from</tt> and
- * <tt>to</tt>). Pre conditions: <tt>from <= to</tt>.
- */
- public double nextDoubleFromTo(double from, double to) {
- return from + (to - from) * randomDouble();
- }
-
- /**
- * Returns a uniformly distributed random number in the open interval <tt>(from,to)</tt> (excluding <tt>from</tt> and
- * <tt>to</tt>). Pre conditions: <tt>from <= to</tt>.
- */
- public float nextFloatFromTo(float from, float to) {
- return (float) nextDoubleFromTo(from, to);
- }
-
- /**
- * Returns a uniformly distributed random number in the closed interval
- * <tt>[from,to]</tt> (including <tt>from</tt>
- * and <tt>to</tt>). Pre conditions: <tt>from <= to</tt>.
- */
- public int nextIntFromTo(int from, int to) {
- return (int) (from + (long) ((1L + to - from) * randomDouble()));
- }
-
- /**
- * Returns a uniformly distributed random number in the closed interval <tt>[from,to]</tt> (including <tt>from</tt>
- * and <tt>to</tt>). Pre conditions: <tt>from <= to</tt>.
- */
- public long nextLongFromTo(long from, long to) {
- /* Doing the thing turns out to be more tricky than expected.
- avoids overflows and underflows.
- treats cases like from=-1, to=1 and the like right.
- the following code would NOT solve the problem: return (long) (Doubles.randomFromTo(from,to));
-
- rounding avoids the unsymmetric behaviour of casts from double to long: (long) -0.7 = 0, (long) 0.7 = 0.
- checking for overflows and underflows is also necessary.
- */
-
- // first the most likely and also the fastest case.
- if (from >= 0 && to < Long.MAX_VALUE) {
- return from + (long) nextDoubleFromTo(0.0, to - from + 1);
- }
-
- // would we get a numeric overflow?
- // if not, we can still handle the case rather efficient.
- double diff = (double) to - (double) from + 1.0;
- if (diff <= Long.MAX_VALUE) {
- return from + (long) nextDoubleFromTo(0.0, diff);
- }
-
- // now the pathologic boundary cases.
- // they are handled rather slow.
- long random;
- if (from == Long.MIN_VALUE) {
- if (to == Long.MAX_VALUE) {
- //return Math.round(nextDoubleFromTo(from,to));
- int i1 = nextIntFromTo(Integer.MIN_VALUE, Integer.MAX_VALUE);
- int i2 = nextIntFromTo(Integer.MIN_VALUE, Integer.MAX_VALUE);
- return ((i1 & 0xFFFFFFFFL) << 32) | (i2 & 0xFFFFFFFFL);
- }
- random = Math.round(nextDoubleFromTo(Long.MIN_VALUE, to + 1));
- if (random > to) {
- random = Long.MIN_VALUE;
- }
- } else {
- random = Math.round(nextDoubleFromTo(from - 1, to));
- if (random < from) {
- random = to;
- }
- }
- return random;
- }
-
- /** Returns the probability distribution function (assuming a continous uniform distribution). */
- @Override
- public double pdf(double x) {
- if (x <= min || x >= max) {
- return 0.0;
- }
- return 1.0 / (max - min);
- }
-
- /** Sets the internal state. */
- public void setState(double min, double max) {
- if (max < min) {
- setState(max, min);
- return;
- }
- this.min = min;
- this.max = max;
- }
-
-
- /** Returns a String representation of the receiver. */
- @Override
- public String toString() {
- return this.getClass().getName() + '(' + min + ',' + max + ')';
- }
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/engine/MersenneTwister.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/engine/MersenneTwister.java b/core/src/main/java/org/apache/mahout/math/jet/random/engine/MersenneTwister.java
deleted file mode 100644
index 8bca895..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/engine/MersenneTwister.java
+++ /dev/null
@@ -1,275 +0,0 @@
-/**
- * Licensed to the Apache Software Foundation (ASF) under one or more
- * contributor license agreements. See the NOTICE file distributed with
- * this work for additional information regarding copyright ownership.
- * The ASF licenses this file to You under the Apache License, Version 2.0
- * (the "License"); you may not use this file except in compliance with
- * the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing, software
- * distributed under the License is distributed on an "AS IS" BASIS,
- * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
- * See the License for the specific language governing permissions and
- * limitations under the License.
- */
-
-/*
-Copyright 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-
-package org.apache.mahout.math.jet.random.engine;
-
-import java.util.Date;
-/**
- MersenneTwister (MT19937) is one of the strongest uniform pseudo-random number generators
- known so far; at the same time it is quick.
- Produces uniformly distributed <tt>int</tt>'s and <tt>long</tt>'s in the closed intervals
- <tt>[Integer.MIN_VALUE,Integer.MAX_VALUE]</tt> and <tt>[Long.MIN_VALUE,Long.MAX_VALUE]</tt>, respectively,
- as well as <tt>float</tt>'s and <tt>double</tt>'s in the open unit intervals <tt>(0.0f,1.0f)</tt>
- and <tt>(0.0,1.0)</tt>, respectively.
- The seed can be any 32-bit integer except <tt>0</tt>. Shawn J. Cokus commented that perhaps the
- seed should preferably be odd.
- <p>
- <b>Quality:</b> MersenneTwister is designed to pass the k-distribution test. It has an
- astronomically large period of 2<sup>19937</sup>-1 (=10<sup>6001</sup>) and 623-dimensional
- equidistribution up to 32-bit accuracy.
- It passes many stringent statistical tests, including the
- <A HREF="http://stat.fsu.edu/~geo/diehard.html">diehard</A> test of G. Marsaglia
- and the load test of P. Hellekalek and S. Wegenkittl.
- <p>
- <b>Performance:</b> Its speed is comparable to other modern generators (in particular,
- as fast as <tt>java.util.Random.nextFloat()</tt>).
- 2.5 million calls to <tt>raw()</tt> per second (Pentium Pro 200 Mhz, JDK 1.2, NT).
- Be aware, however, that there is a non-negligible amount of overhead required to initialize the data
- structures used by a MersenneTwister. Code like
- {@code
- double sum = 0.0;
- for (int i=0; i<100000; ++i) {
- RandomElement twister = new MersenneTwister(new Date());
- sum += twister.raw();
- }
- }
- will be wildly inefficient. Consider using
- {@code
- double sum = 0.0;
- RandomElement twister = new MersenneTwister(new Date());
- for (int i=0; i<100000; ++i) {
- sum += twister.raw();
- }
- }
- instead. This allows the cost of constructing the MersenneTwister object
- to be borne only once, rather than once for each iteration in the loop.
- <p>
- <b>Implementation:</b> After M. Matsumoto and T. Nishimura,
- "Mersenne Twister: A 623-Dimensionally Equidistributed Uniform Pseudo-Random Number Generator",
- ACM Transactions on Modeling and Computer Simulation,
- Vol. 8, No. 1, January 1998, pp 3--30.
- <dl>
- <dt>More info on <a HREF="http://www.math.keio.ac.jp/~matumoto/eindex.html"> Masumoto's homepage</a>.</dt>
- <dt>More info on <a HREF="http://www.ncsa.uiuc.edu/Apps/CMP/RNG/www-rng.html"> Pseudo-random number
- generators is on the Web</a>.</dt>
- <dt>Yet <a HREF="http://nhse.npac.syr.edu/random"> some more info</a>.</dt>
- <p>
- The correctness of this implementation has been verified against the published output sequence
- <a href="http://www.math.keio.ac.jp/~nisimura/random/real2/mt19937-2.out">mt19937-2.out</a> of the C-implementation
- <a href="http://www.math.keio.ac.jp/~nisimura/random/real2/mt19937-2.c">mt19937-2.c</a>.
- (Call <tt>test(1000)</tt> to print the sequence).
- <dt>
- Note that this implementation is <b>not synchronized</b>.</dt>
- </dl>
- <p>
- <b>Details:</b> MersenneTwister is designed with consideration of the flaws of various existing generators in mind.
- It is an improved version of TT800, a very successful generator.
- MersenneTwister is based on linear recurrences modulo 2.
- Such generators are very fast, have extremely long periods, and appear quite robust.
- MersenneTwister produces 32-bit numbers, and every <tt>k</tt>-dimensional vector of such
- numbers appears the same number of times as <tt>k</tt> successive values over the
- period length, for each <tt>k <= 623</tt> (except for the zero vector, which appears one time less).
- If one looks at only the first <tt>n <= 16</tt> bits of each number, then the property holds
- for even larger <tt>k</tt>, as shown in the following table (taken from the publication cited above):
- <table width="75%" border="1" cellspacing="0" cellpadding="0" summary="property table">
- <tr>
- <td width="2%" align="center"> <div>n</div> </td>
- <td width="6%" align="center"> <div>1</div> </td>
- <td width="5%" align="center"> <div>2</div> </td>
- <td width="5%" align="center"> <div>3</div> </td>
- <td width="5%" align="center"> <div>4</div> </td>
- <td width="5%" align="center"> <div>5</div> </td>
- <td width="5%" align="center"> <div>6</div> </td>
- <td width="5%" align="center"> <div>7</div> </td>
- <td width="5%" align="center"> <div>8</div> </td>
- <td width="5%" align="center"> <div>9</div> </td>
- <td width="5%" align="center"> <div>10</div> </td>
- <td width="5%" align="center"> <div>11</div> </td>
- <td width="10%" align="center"> <div>12 .. 16</div> </td>
- <td width="10%" align="center"> <div>17 .. 32</div> </td>
- </tr>
- <tr>
- <td width="2%" align="center"> <div>k</div> </td>
- <td width="6%" align="center"> <div>19937</div> </td>
- <td width="5%" align="center"> <div>9968</div> </td>
- <td width="5%" align="center"> <div>6240</div> </td>
- <td width="5%" align="center"> <div>4984</div> </td>
- <td width="5%" align="center"> <div>3738</div> </td>
- <td width="5%" align="center"> <div>3115</div> </td>
- <td width="5%" align="center"> <div>2493</div> </td>
- <td width="5%" align="center"> <div>2492</div> </td>
- <td width="5%" align="center"> <div>1869</div> </td>
- <td width="5%" align="center"> <div>1869</div> </td>
- <td width="5%" align="center"> <div>1248</div> </td>
- <td width="10%" align="center"> <div>1246</div> </td>
- <td width="10%" align="center"> <div>623</div> </td>
- </tr>
- </table>
- <p>
- MersenneTwister generates random numbers in batches of 624 numbers at a time, so
- the caching and pipelining of modern systems is exploited.
- The generator is implemented to generate the output by using the fastest arithmetic
- operations only: 32-bit additions and bit operations (no division, no multiplication, no mod).
- These operations generate sequences of 32 random bits (<tt>int</tt>'s).
- <tt>long</tt>'s are formed by concatenating two 32 bit <tt>int</tt>'s.
- <tt>float</tt>'s are formed by dividing the interval <tt>[0.0,1.0]</tt> into 2<sup>32</sup>
- sub intervals, then randomly choosing one subinterval.
- <tt>double</tt>'s are formed by dividing the interval <tt>[0.0,1.0]</tt> into 2<sup>64</sup>
- sub intervals, then randomly choosing one subinterval.
- <p>
- @author wolfgang.hoschek@cern.ch
- @version 1.0, 09/24/99
- @see java.util.Random
- */
-public final class MersenneTwister extends RandomEngine {
-
- /* Period parameters */
- private static final int N = 624;
- private static final int M = 397;
- private static final int MATRIX_A = 0x9908b0df; /* constant vector a */
- private static final int UPPER_MASK = 0x80000000; /* most significant w-r bits */
- private static final int LOWER_MASK = 0x7fffffff; /* least significant r bits */
-
- /* for tempering */
- private static final int TEMPERING_MASK_B = 0x9d2c5680;
- private static final int TEMPERING_MASK_C = 0xefc60000;
-
- private static final int MAG0 = 0x0;
- private static final int MAG1 = MATRIX_A;
- //private static final int[] mag01=new int[] {0x0, MATRIX_A};
- /* mag01[x] = x * MATRIX_A for x=0,1 */
-
- private static final int DEFAULT_SEED = 4357;
-
- private int mti;
- private final int[] mt = new int[N]; /* set initial seeds: N = 624 words */
-
- /**
- * Constructs and returns a random number generator with a default seed, which is a <b>constant</b>. Thus using this
- * constructor will yield generators that always produce exactly the same sequence. This method is mainly intended to
- * ease testing and debugging.
- */
- public MersenneTwister() {
- this(DEFAULT_SEED);
- }
-
- /** Constructs and returns a random number generator with the given seed.
- * @param seed A number that is used to initialize the internal state of the generator.
- */
- public MersenneTwister(int seed) {
- setSeed(seed);
- }
-
- /**
- * Constructs and returns a random number generator seeded with the given date.
- *
- * @param d typically <tt>new Date()</tt>
- */
- public MersenneTwister(Date d) {
- this((int) d.getTime());
- }
-
- /** Generates N words at one time */
- void nextBlock() {
- int y;
- int kk;
-
- for (kk = 0; kk < N - M; kk++) {
- y = (mt[kk] & UPPER_MASK) | (mt[kk + 1] & LOWER_MASK);
- mt[kk] = mt[kk + M] ^ (y >>> 1) ^ ((y & 0x1) == 0 ? MAG0 : MAG1);
- }
- for (; kk < N - 1; kk++) {
- y = (mt[kk] & UPPER_MASK) | (mt[kk + 1] & LOWER_MASK);
- mt[kk] = mt[kk + (M - N)] ^ (y >>> 1) ^ ((y & 0x1) == 0 ? MAG0 : MAG1);
- }
- y = (mt[N - 1] & UPPER_MASK) | (mt[0] & LOWER_MASK);
- mt[N - 1] = mt[M - 1] ^ (y >>> 1) ^ ((y & 0x1) == 0 ? MAG0 : MAG1);
-
- this.mti = 0;
- }
-
- /**
- * Returns a 32 bit uniformly distributed random number in the closed interval
- * <tt>[Integer.MIN_VALUE,Integer.MAX_VALUE]</tt>
- * (including <tt>Integer.MIN_VALUE</tt> and <tt>Integer.MAX_VALUE</tt>).
- */
- @Override
- public int nextInt() {
- /* Each single bit including the sign bit will be random */
- if (mti == N) {
- nextBlock();
- } // generate N ints at one time
-
- int y = mt[mti++];
- y ^= y >>> 11; // y ^= TEMPERING_SHIFT_U(y );
- y ^= (y << 7) & TEMPERING_MASK_B; // y ^= TEMPERING_SHIFT_S(y) & TEMPERING_MASK_B;
- y ^= (y << 15) & TEMPERING_MASK_C; // y ^= TEMPERING_SHIFT_T(y) & TEMPERING_MASK_C;
- // y &= 0xffffffff; //you may delete this line if word size = 32
- y ^= y >>> 18; // y ^= TEMPERING_SHIFT_L(y);
-
- return y;
- }
-
- /** Sets the receiver's seed. This method resets the receiver's entire internal state.
- * @param seed An integer that is used to reset the internal state of the generator */
- void setSeed(int seed) {
- mt[0] = seed;
- for (int i = 1; i < N; i++) {
- mt[i] = 1812433253 * (mt[i - 1] ^ (mt[i - 1] >> 30)) + i;
- /* See Knuth TAOCP Vol2. 3rd Ed. P.106 for multiplier. */
- /* In the previous versions, MSBs of the seed affect */
- /* only MSBs of the array mt[]. */
- /* 2002/01/09 modified by Makoto Matsumoto */
- //mt[i] &= 0xffffffff;
- /* for >32 bit machines */
- }
- //log.info("init done");
- mti = N;
- }
-
- /**
- * Sets the receiver's seed in a fashion compatible with the
- * reference C implementation. See
- * http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/VERSIONS/C-LANG/980409/mt19937int.c
- *
- * This method isn't as good as the default method due to poor distribution of the
- * resulting states.
- *
- * @param seed An integer that is used to reset the internal state in the same way as
- * done in the 1999 reference implementation. Should only be used for testing, not
- * actual coding.
- */
- void setReferenceSeed(int seed) {
- for (int i = 0; i < N; i++) {
- mt[i] = seed & 0xffff0000;
- seed = 69069 * seed + 1;
- mt[i] |= (seed & 0xffff0000) >>> 16;
- seed = 69069 * seed + 1;
- }
- //log.info("init done");
- mti = N;
- }
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/engine/RandomEngine.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/engine/RandomEngine.java b/core/src/main/java/org/apache/mahout/math/jet/random/engine/RandomEngine.java
deleted file mode 100644
index f832b1d..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/engine/RandomEngine.java
+++ /dev/null
@@ -1,169 +0,0 @@
-/**
- * Licensed to the Apache Software Foundation (ASF) under one
- * or more contributor license agreements. See the NOTICE file
- * distributed with this work for additional information
- * regarding copyright ownership. The ASF licenses this file
- * to you under the Apache License, Version 2.0 (the
- * "License"); you may not use this file except in compliance
- * with the License. You may obtain a copy of the License at
- *
- * http://www.apache.org/licenses/LICENSE-2.0
- *
- * Unless required by applicable law or agreed to in writing,
- * software distributed under the License is distributed on an
- * "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
- * KIND, either express or implied. See the License for the
- * specific language governing permissions and limitations
- * under the License.
- */
-/*
-Copyright � 1999 CERN - European Organization for Nuclear Research.
-Permission to use, copy, modify, distribute and sell this software and its documentation for any purpose
-is hereby granted without fee, provided that the above copyright notice appear in all copies and
-that both that copyright notice and this permission notice appear in supporting documentation.
-CERN makes no representations about the suitability of this software for any purpose.
-It is provided "as is" without expressed or implied warranty.
-*/
-package org.apache.mahout.math.jet.random.engine;
-
-import org.apache.mahout.math.function.DoubleFunction;
-import org.apache.mahout.math.function.IntFunction;
-
-/**
- * Abstract base class for uniform pseudo-random number generating engines.
- * <p>
- * Most probability distributions are obtained by using a <b>uniform</b> pseudo-random number generation engine
- * followed by a transformation to the desired distribution.
- * Thus, subclasses of this class are at the core of computational statistics, simulations, Monte Carlo methods, etc.
- * <p>
- * Subclasses produce uniformly distributed <tt>int</tt>'s and <tt>long</tt>'s in the closed intervals
- * <tt>[Integer.MIN_VALUE,Integer.MAX_VALUE]</tt> and <tt>[Long.MIN_VALUE,Long.MAX_VALUE]</tt>, respectively,
- * as well as <tt>float</tt>'s and <tt>double</tt>'s in the open unit intervals <tt>(0.0f,1.0f)</tt> and
- * <tt>(0.0,1.0)</tt>, respectively.
- * <p>
- * Subclasses need to override one single method only: <tt>nextInt()</tt>.
- * All other methods generating different data types or ranges are usually layered upon <tt>nextInt()</tt>.
- * <tt>long</tt>'s are formed by concatenating two 32 bit <tt>int</tt>'s.
- * <tt>float</tt>'s are formed by dividing the interval <tt>[0.0f,1.0f]</tt> into 2<sup>32</sup> sub intervals,
- * then randomly choosing one subinterval.
- * <tt>double</tt>'s are formed by dividing the interval <tt>[0.0,1.0]</tt> into 2<sup>64</sup> sub intervals,
- * then randomly choosing one subinterval.
- * <p>
- * Note that this implementation is <b>not synchronized</b>.
- *
- * @see MersenneTwister
- * @see java.util.Random
- */
-public abstract class RandomEngine extends DoubleFunction implements IntFunction {
-
- /**
- * Equivalent to <tt>raw()</tt>. This has the effect that random engines can now be used as function objects,
- * returning a random number upon function evaluation.
- */
- @Override
- public double apply(double dummy) {
- return raw();
- }
-
- /**
- * Equivalent to <tt>nextInt()</tt>. This has the effect that random engines can now be used as function objects,
- * returning a random number upon function evaluation.
- */
- @Override
- public int apply(int dummy) {
- return nextInt();
- }
-
- /**
- * @return a 64 bit uniformly distributed random number in the open unit interval {@code (0.0,1.0)} (excluding
- * 0.0 and 1.0).
- */
- public double nextDouble() {
- double nextDouble;
-
- do {
- // -9.223372036854776E18 == (double) Long.MIN_VALUE
- // 5.421010862427522E-20 == 1 / Math.pow(2,64) == 1 / ((double) Long.MAX_VALUE - (double) Long.MIN_VALUE);
- nextDouble = (nextLong() - -9.223372036854776E18) * 5.421010862427522E-20;
- }
- // catch loss of precision of long --> double conversion
- while (!(nextDouble > 0.0 && nextDouble < 1.0));
-
- // --> in (0.0,1.0)
- return nextDouble;
-
- /*
- nextLong == Long.MAX_VALUE --> 1.0
- nextLong == Long.MIN_VALUE --> 0.0
- nextLong == Long.MAX_VALUE-1 --> 1.0
- nextLong == Long.MAX_VALUE-100000L --> 0.9999999999999946
- nextLong == Long.MIN_VALUE+1 --> 0.0
- nextLong == Long.MIN_VALUE-100000L --> 0.9999999999999946
- nextLong == 1L --> 0.5
- nextLong == -1L --> 0.5
- nextLong == 2L --> 0.5
- nextLong == -2L --> 0.5
- nextLong == 2L+100000L --> 0.5000000000000054
- nextLong == -2L-100000L --> 0.49999999999999456
- */
- }
-
- /**
- * @return a 32 bit uniformly distributed random number in the open unit interval {@code (0.0f, 1.0f)} (excluding
- * 0.0f and 1.0f).
- */
- public float nextFloat() {
- // catch loss of precision of double --> float conversion which could result in a value == 1.0F
- float nextFloat;
- do {
- nextFloat = (float) raw();
- }
- while (nextFloat >= 1.0f);
-
- // --> in [0.0f,1.0f)
- return nextFloat;
- }
-
- /**
- * @return a 32 bit uniformly distributed random number in the closed interval
- * <tt>[Integer.MIN_VALUE,Integer.MAX_VALUE]</tt>
- * (including <tt>Integer.MIN_VALUE</tt> and <tt>Integer.MAX_VALUE</tt>);
- */
- public abstract int nextInt();
-
- /**
- * @return a 64 bit uniformly distributed random number in the closed interval
- * <tt>[Long.MIN_VALUE,Long.MAX_VALUE]</tt>
- * (including <tt>Long.MIN_VALUE</tt> and <tt>Long.MAX_VALUE</tt>).
- */
- public long nextLong() {
- // concatenate two 32-bit strings into one 64-bit string
- return ((nextInt() & 0xFFFFFFFFL) << 32) | (nextInt() & 0xFFFFFFFFL);
- }
-
- /**
- * @return a 32 bit uniformly distributed random number in the open unit interval {@code (0.0, 1.0)} (excluding
- * 0.0 and 1.0).
- */
- public double raw() {
- int nextInt;
- do { // accept anything but zero
- nextInt = nextInt(); // in [Integer.MIN_VALUE,Integer.MAX_VALUE]-interval
- } while (nextInt == 0);
-
- // transform to (0.0,1.0)-interval
- // 2.3283064365386963E-10 == 1.0 / Math.pow(2,32)
- return (nextInt & 0xFFFFFFFFL) * 2.3283064365386963E-10;
-
- /*
- nextInt == Integer.MAX_VALUE --> 0.49999999976716936
- nextInt == Integer.MIN_VALUE --> 0.5
- nextInt == Integer.MAX_VALUE-1 --> 0.4999999995343387
- nextInt == Integer.MIN_VALUE+1 --> 0.5000000002328306
- nextInt == 1 --> 2.3283064365386963E-10
- nextInt == -1 --> 0.9999999997671694
- nextInt == 2 --> 4.6566128730773926E-10
- nextInt == -2 --> 0.9999999995343387
- */
- }
-}
http://git-wip-us.apache.org/repos/asf/mahout/blob/49ad8cb4/core/src/main/java/org/apache/mahout/math/jet/random/engine/package-info.java
----------------------------------------------------------------------
diff --git a/core/src/main/java/org/apache/mahout/math/jet/random/engine/package-info.java b/core/src/main/java/org/apache/mahout/math/jet/random/engine/package-info.java
deleted file mode 100644
index e092010..0000000
--- a/core/src/main/java/org/apache/mahout/math/jet/random/engine/package-info.java
+++ /dev/null
@@ -1,7 +0,0 @@
-/**
- * Engines generating strong uniformly distributed pseudo-random numbers;
- * Needed by all JET probability distributions since they rely on uniform random numbers to generate random
- * numbers from their own distribution.
- * Thus, the classes of this package are at the core of computational statistics, simulations, Monte Carlo methods, etc.
- */
-package org.apache.mahout.math.jet.random.engine;