diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/GammaFunctions.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/GammaFunctions.java
index 7c28a860b7c2c9158abd41f07e4360f1753f3d26..bdad479b98efd9546b467c40cfb17ba753dc7e19 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/GammaFunctions.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/GammaFunctions.java
@@ -167,7 +167,7 @@ public abstract class GammaFunctions {
private static final double xbig = 94906265.62425156;
private static final double lgc_xmax = 3.745194030963158e306;
- private static double lgammacor(double x) {
+ static double lgammacor(double x) {
double tmp;
if (x < 10) {
@@ -183,6 +183,10 @@ public abstract class GammaFunctions {
return 1 / (x * 12);
}
+ static double lgamma(double x) {
+ throw RError.nyi(RError.SHOW_CALLER, "lgamma from libc");
+ }
+
//
// gammafn
//
@@ -211,7 +215,7 @@ public abstract class GammaFunctions {
private static final double M_LN_SQRT_2PI = 0.918938533204672741780329736406;
- private static double gammafn(double x) {
+ static double gammafn(double x) {
int i;
int n;
double y;
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/LBeta.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/LBeta.java
new file mode 100644
index 0000000000000000000000000000000000000000..6cc0de5a33d977f8dec7d56bf66456d975a15ee3
--- /dev/null
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/LBeta.java
@@ -0,0 +1,66 @@
+/*
+ * This material is distributed under the GNU General Public License
+ * Version 2. You may review the terms of this license at
+ * http://www.gnu.org/licenses/gpl-2.0.html
+ *
+ * Copyright (C) 1998 Ross Ihaka
+ * Copyright (c) 2000--2012, The R Core Team
+ * Copyright (c) 2003, The R Foundation
+ * Copyright (c) 2016, Oracle and/or its affiliates
+ *
+ * All rights reserved.
+ */
+package com.oracle.truffle.r.library.stats;
+
+import static com.oracle.truffle.r.library.stats.GammaFunctions.gammafn;
+import static com.oracle.truffle.r.library.stats.GammaFunctions.lgamma;
+import static com.oracle.truffle.r.library.stats.GammaFunctions.lgammacor;
+import static com.oracle.truffle.r.library.stats.GammaFunctions.lgammafn;
+import static com.oracle.truffle.r.library.stats.MathConstants.M_LN_SQRT_2PI;
+
+public final class LBeta {
+ public static double lbeta(double a, double b) {
+ double corr;
+ double p;
+ double q;
+
+ if (Double.isNaN(a) || Double.isNaN(b)) {
+ return a + b;
+ }
+
+ p = q = a;
+ if (b < p) {
+ p = b;
+ /* := min(a,b) */ }
+ if (b > q) {
+ q = b;
+ /* := max(a,b) */ }
+
+ /* both arguments must be >= 0 */
+ if (p < 0) {
+ return StatsUtil.mlError();
+ } else if (p == 0) {
+ return Double.POSITIVE_INFINITY;
+ } else if (!Double.isFinite(q)) { /* q == +Inf */
+ return Double.NEGATIVE_INFINITY;
+ }
+
+ if (p >= 10) {
+ /* p and q are big. */
+ corr = lgammacor(p) + lgammacor(q) - lgammacor(p + q);
+ return Math.log(q) * -0.5 + M_LN_SQRT_2PI + corr + (p - 0.5) * Math.log(p / (p + q)) + q * Math.log1p(-p / (p + q));
+ } else if (q >= 10) {
+ /* p is small, but q is big. */
+ corr = lgammacor(q) - lgammacor(p + q);
+ return lgammafn(p) + corr + p - p * Math.log(p + q) + (q - 0.5) * Math.log1p(-p / (p + q));
+ } else {
+ /* p and q are small: p <= q < 10. */
+ /* R change for very small args */
+ if (p < 1e-306) {
+ return lgamma(p) + (lgamma(q) - lgamma(p + q));
+ } else {
+ return Math.log(gammafn(p) * (gammafn(q) / gammafn(p + q)));
+ }
+ }
+ }
+}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/MathConstants.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/MathConstants.java
index 8f0f1a07f4c6afeaa617d46b65a59e9ff9911cc0..1fbefb2a995f10460e01aef039b2ad94135b163c 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/MathConstants.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/MathConstants.java
@@ -77,6 +77,8 @@ public final class MathConstants {
// log(sqrt(pi/2)) == log(pi/2)/2
public static final double M_LN_SQRT_PId2 = 0.225791352644727432363097614947;
+ public static final double DBL_EPSILON = Math.ulp(1.0);
+
/**
* Compute the log of a sum from logs of terms, i.e.,
*
@@ -90,7 +92,9 @@ public final class MathConstants {
}
// R_forceint
- public static long forceint(double x) {
- return Math.round(x);
+ public static double forceint(double x) {
+ // Note: in GnuR this is alias for nearbyint, which may not behave exactly like Math.round,
+ // especially Math.round(-0.5) == 0.0, instead of -0.0, does it matter a lot?
+ return Double.isFinite(x) ? Math.round(x) : x;
}
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QHyper.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QHyper.java
new file mode 100644
index 0000000000000000000000000000000000000000..d3e23b0e73276ac48c384b17a182ec639278afcb
--- /dev/null
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QHyper.java
@@ -0,0 +1,98 @@
+/*
+ * This material is distributed under the GNU General Public License
+ * Version 2. You may review the terms of this license at
+ * http://www.gnu.org/licenses/gpl-2.0.html
+ *
+ * Copyright (C) 1998 Ross Ihaka
+ * Copyright (c) 2000--2014, The R Core Team
+ * Copyright (c) 2016, Oracle and/or its affiliates
+ *
+ * All rights reserved.
+ */
+package com.oracle.truffle.r.library.stats;
+
+import static com.oracle.truffle.r.library.stats.MathConstants.DBL_EPSILON;
+import static com.oracle.truffle.r.library.stats.MathConstants.forceint;
+import static com.oracle.truffle.r.library.stats.StatsUtil.fmax2;
+import static com.oracle.truffle.r.library.stats.StatsUtil.fmin2;
+import static com.oracle.truffle.r.library.stats.StatsUtil.lfastchoose;
+
+import com.oracle.truffle.r.library.stats.DPQ.EarlyReturn;
+
+public final class QHyper {
+ public static double qhyper(double p, double nr, double nb, double n, boolean lowerTail, boolean logP) {
+ /* This is basically the same code as ./phyper.c *used* to be --> FIXME! */
+ double capN;
+ double xstart;
+ double xend;
+ double xr;
+ double xb;
+ double sum;
+ double term;
+ boolean smallN;
+ if (Double.isNaN(p) || Double.isNaN(nr) || Double.isNaN(nb) || Double.isNaN(n)) {
+ return p + nr + nb + n;
+ }
+ if (!Double.isFinite(p) || !Double.isFinite(nr) || !Double.isFinite(nb) || !Double.isFinite(n)) {
+ return StatsUtil.mlError();
+ }
+
+ nr = forceint(nr);
+ nb = forceint(nb);
+ capN = nr + nb;
+ n = forceint(n);
+ if (nr < 0 || nb < 0 || n < 0 || n > capN) {
+ return StatsUtil.mlError();
+ }
+
+ /*
+ * Goal: Find xr (= #{red balls in sample}) such that phyper(xr, NR,NB, n) >= p > phyper(xr
+ * - 1, NR,NB, n)
+ */
+
+ xstart = fmax2(0, n - nb);
+ xend = fmin2(n, nr);
+
+ try {
+ DPQ.qP01Boundaries(p, xstart, xend, lowerTail, logP);
+ } catch (EarlyReturn ex) {
+ return ex.result;
+ }
+
+ xr = xstart;
+ xb = n - xr; /* always ( = #{black balls in sample} ) */
+
+ smallN = (capN < 1000); /* won't have underflow in product below */
+ /*
+ * if N is small, term := product.ratio( bin.coef ); otherwise work with its Math.logarithm
+ * to protect against underflow
+ */
+ term = lfastchoose(nr, xr) + lfastchoose(nb, xb) - lfastchoose(capN, n);
+ if (smallN) {
+ term = Math.exp(term);
+ }
+ nr -= xr;
+ nb -= xb;
+
+ if (!lowerTail || logP) {
+ p = DPQ.dtQIv(p, lowerTail, logP);
+ }
+ p *= 1 - 1000 * DBL_EPSILON; /* was 64, but failed on FreeBSD sometimes */
+ sum = smallN ? term : Math.exp(term);
+
+ while (sum < p && xr < xend) {
+ xr++;
+ nb++;
+ if (smallN) {
+ term *= (nr / xr) * (xb / nb);
+ } else {
+ term += Math.log((nr / xr) * (xb / nb));
+ }
+ sum += smallN ? term : Math.exp(term);
+ xb--;
+ nr--;
+ }
+ return xr;
+ }
+
+}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RChisq.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RChisq.java
index b9658bb23bffefbba09c6b0763049198326fd90b..3b4c21ac102e287b7f0c21c0fdd01a6dd7eff473 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RChisq.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RChisq.java
@@ -11,13 +11,19 @@
*/
package com.oracle.truffle.r.library.stats;
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction1_Double;
import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandomNumberProvider;
-public final class RChisq {
+public final class RChisq implements RandFunction1_Double {
public static double rchisq(double df, RandomNumberProvider rand) {
if (!Double.isFinite(df) || df < 0.0) {
return StatsUtil.mlError();
}
return new RGamma().evaluate(df / 2.0, 2.0, rand);
}
+
+ @Override
+ public double evaluate(double a, RandomNumberProvider rand) {
+ return rchisq(a, rand);
+ }
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RExp.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RExp.java
new file mode 100644
index 0000000000000000000000000000000000000000..078f113dc20a407979ce7f75e2eb634ac663ea82
--- /dev/null
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RExp.java
@@ -0,0 +1,25 @@
+/*
+ * This material is distributed under the GNU General Public License
+ * Version 2. You may review the terms of this license at
+ * http://www.gnu.org/licenses/gpl-2.0.html
+ *
+ * Copyright (C) 1998 Ross Ihaka
+ * Copyright (c) 1998--2008, The R Core Team
+ * Copyright (c) 2016, 2016, Oracle and/or its affiliates
+ *
+ * All rights reserved.
+ */
+package com.oracle.truffle.r.library.stats;
+
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction1_Double;
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandomNumberProvider;
+
+public final class RExp implements RandFunction1_Double {
+ @Override
+ public double evaluate(double scale, RandomNumberProvider rand) {
+ if (!Double.isFinite(scale) || scale <= 0.0) {
+ return scale == 0. ? 0. : StatsUtil.mlError();
+ }
+ return scale * rand.expRand();
+ }
+}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RGeom.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RGeom.java
new file mode 100644
index 0000000000000000000000000000000000000000..555d703a529b0e2021bf70e6ac116b93cdcd9467
--- /dev/null
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RGeom.java
@@ -0,0 +1,25 @@
+/*
+ * This material is distributed under the GNU General Public License
+ * Version 2. You may review the terms of this license at
+ * http://www.gnu.org/licenses/gpl-2.0.html
+ *
+ * Copyright (C) 1998 Ross Ihaka
+ * Copyright (c) 1998--2008, The R Core Team
+ * Copyright (c) 2016, 2016, Oracle and/or its affiliates
+ *
+ * All rights reserved.
+ */
+package com.oracle.truffle.r.library.stats;
+
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction1_Double;
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandomNumberProvider;
+
+public final class RGeom implements RandFunction1_Double {
+ @Override
+ public double evaluate(double p, RandomNumberProvider rand) {
+ if (!Double.isFinite(p) || p <= 0 || p > 1) {
+ return StatsUtil.mlError();
+ }
+ return RPois.rpois(rand.expRand() * ((1 - p) / p), rand);
+ }
+}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RHyper.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RHyper.java
new file mode 100644
index 0000000000000000000000000000000000000000..32d9c6afee1cdf69754660e0d751057aefd13389
--- /dev/null
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RHyper.java
@@ -0,0 +1,364 @@
+/*
+ * This material is distributed under the GNU General Public License
+ * Version 2. You may review the terms of this license at
+ * http://www.gnu.org/licenses/gpl-2.0.html
+ *
+ * Copyright (C) 1998 Ross Ihaka
+ * Copyright (c) 2000--2009, The R Core Team
+ * Copyright (c) 2003--2009, The R Foundation
+ * Copyright (c) 2016, 2016, Oracle and/or its affiliates
+ *
+ * All rights reserved.
+ */
+package com.oracle.truffle.r.library.stats;
+
+import static com.oracle.truffle.r.library.stats.MathConstants.M_LN_SQRT_2PI;
+import static com.oracle.truffle.r.library.stats.MathConstants.forceint;
+
+import com.oracle.truffle.api.CompilerDirectives.TruffleBoundary;
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction3_Double;
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandomNumberProvider;
+import com.oracle.truffle.r.runtime.RError;
+import com.oracle.truffle.r.runtime.RError.Message;
+
+public final class RHyper implements RandFunction3_Double {
+ private static final double[] al = {
+ 0.0, /* ln(0!)=ln(1) */
+ 0.0, /* ln(1!)=ln(1) */
+ 0.69314718055994530941723212145817, /* ln(2) */
+ 1.79175946922805500081247735838070, /* ln(6) */
+ 3.17805383034794561964694160129705, /* ln(24) */
+ 4.78749174278204599424770093452324,
+ 6.57925121201010099506017829290394,
+ 8.52516136106541430016553103634712
+ /*
+ * 10.60460290274525022841722740072165, approx. value below = 10.6046028788027;
+ * rel.error = 2.26 10^{-9}
+ *
+ * FIXME: Use constants and if(n > ..) decisions from ./stirlerr.c ----- will be
+ * even *faster* for n > 500 (or so)
+ */
+ };
+
+ // afc(i) := ln( i! ) [Math.logarithm of the factorial i]
+ private static double afc(int i) {
+ // If (i > 7), use Stirling's approximation, otherwise use table lookup.
+ if (i < 0) {
+ RError.warning(RError.SHOW_CALLER, Message.GENERIC, String.format("RHyper.java: afc(i), i=%d < 0 -- SHOULD NOT HAPPEN!\n", i));
+ return -1;
+ }
+ if (i <= 7) {
+ return al[i];
+ }
+ // else i >= 8 :
+ double di = i;
+ double i2 = di * di;
+ return (di + 0.5) * Math.log(di) - di + M_LN_SQRT_2PI +
+ (0.0833333333333333 - 0.00277777777777778 / i2) / di;
+ }
+
+ /* These should become 'thread_local globals' : */
+ private static int ks = -1;
+ private static int n1s = -1;
+ private static int n2s = -1;
+ private static int m;
+ private static int minjx;
+ private static int maxjx;
+ private static int k;
+ private static int n1;
+ private static int n2; // <- not allowing larger integer par
+ private static double tn;
+
+ // II :
+ private static double w;
+ // III:
+ // Checkstyle: stop
+ private static double a, d, s, xl, xr, kl, kr, lamdl, lamdr, p1, p2, p3;
+ // // Checkstyle: resume
+
+ private static final double scale = 1e25; // scaling factor against (early) underflow
+ private static final double con = 57.5646273248511421;
+
+ private static final double deltal = 0.0078;
+ private static final double deltau = 0.0034;
+
+ private final Rbinom rbinom = new Rbinom();
+
+ // rhyper(NR, NB, n) -- NR 'red', NB 'blue', n drawn, how many are 'red'
+ @Override
+ @TruffleBoundary
+ public double evaluate(double nn1in, double nn2in, double kkin, RandomNumberProvider rand) {
+ /* extern double afc(int); */
+
+ int nn1;
+ int nn2;
+ int kk;
+ int ix; // return value (coerced to double at the very end)
+ boolean setup1;
+ boolean setup2;
+
+ /* check parameter validity */
+
+ if (!Double.isFinite(nn1in) || !Double.isFinite(nn2in) || !Double.isFinite(kkin)) {
+ return StatsUtil.mlError();
+ }
+
+ nn1in = forceint(nn1in);
+ nn2in = forceint(nn2in);
+ kkin = forceint(kkin);
+
+ if (nn1in < 0 || nn2in < 0 || kkin < 0 || kkin > nn1in + nn2in) {
+ return StatsUtil.mlError();
+ }
+ if (nn1in >= Integer.MAX_VALUE || nn2in >= Integer.MAX_VALUE || kkin >= Integer.MAX_VALUE) {
+ /*
+ * large n -- evade integer overflow (and inappropriate algorithms) --------
+ */
+ // FIXME: Much faster to give rbinom() approx when appropriate; -> see Kuensch(1989)
+ // Johnson, Kotz,.. p.258 (top) mention the *four* different binomial approximations
+ if (kkin == 1.) { // Bernoulli
+ return rbinom.evaluate(kkin, nn1in / (nn1in + nn2in), rand);
+ }
+ // Slow, but safe: return F^{-1}(U) where F(.) = phyper(.) and U ~ U[0,1]
+ return QHyper.qhyper(rand.unifRand(), nn1in, nn2in, kkin, false, false);
+ }
+ nn1 = (int) nn1in;
+ nn2 = (int) nn2in;
+ kk = (int) kkin;
+
+ /* if new parameter values, initialize */
+ if (nn1 != n1s || nn2 != n2s) {
+ setup1 = true;
+ setup2 = true;
+ } else if (kk != ks) {
+ setup1 = false;
+ setup2 = true;
+ } else {
+ setup1 = false;
+ setup2 = false;
+ }
+ if (setup1) {
+ n1s = nn1;
+ n2s = nn2;
+ tn = nn1 + nn2;
+ if (nn1 <= nn2) {
+ n1 = nn1;
+ n2 = nn2;
+ } else {
+ n1 = nn2;
+ n2 = nn1;
+ }
+ }
+ if (setup2) {
+ ks = kk;
+ if (kk + kk >= tn) {
+ k = (int) (tn - kk);
+ } else {
+ k = kk;
+ }
+ }
+ if (setup1 || setup2) {
+ m = (int) ((k + 1.) * (n1 + 1.) / (tn + 2.));
+ minjx = Math.max(0, k - n2);
+ maxjx = Math.min(n1, k);
+ }
+ /* generate random variate --- Three basic cases */
+
+ if (minjx == maxjx) { /*
+ * I: degenerate distribution ----------------
+ */
+ ix = maxjx;
+
+ } else if (m - minjx < 10) { // II: (Scaled) algorithm HIN
+ // (inverse transformation)
+ // ----
+ // 25*Math.log(10) = Math.log(scale) { <==> Math.exp(con) == scale }
+ if (setup1 || setup2) {
+ double lw; // Math.log(w); w = Math.exp(lw) * scale = Math.exp(lw +
+ // Math.log(scale)) = Math.exp(lw + con)
+ if (k < n2) {
+ lw = afc(n2) + afc(n1 + n2 - k) - afc(n2 - k) - afc(n1 + n2);
+ } else {
+ lw = afc(n1) + afc(k) - afc(k - n2) - afc(n1 + n2);
+ }
+ w = Math.exp(lw + con);
+ }
+ L10: while (true) {
+ double p = w;
+ ix = minjx;
+ double u = rand.unifRand() * scale;
+ while (u > p) {
+ u -= p;
+ p *= ((double) n1 - ix) * (k - ix);
+ ix++;
+ p = p / ix / (n2 - k + ix);
+ if (ix > maxjx) {
+ continue L10;
+ }
+ // FIXME if(p == 0.) we also "have lost" => goto L10
+ }
+ break L10;
+ }
+ } else { /* III : H2PE Algorithm --------------------------------------- */
+
+ double u;
+ double v;
+
+ if (setup1 || setup2) {
+ s = Math.sqrt((tn - k) * k * n1 * n2 / (tn - 1) / tn / tn);
+
+ /* remark: d is defined in reference without int. */
+ /* the truncation centers the cell boundaries at 0.5 */
+
+ d = (int) (1.5 * s) + .5;
+ xl = m - d + .5;
+ xr = m + d + .5;
+ a = afc(m) + afc(n1 - m) + afc(k - m) + afc(n2 - k + m);
+ kl = Math.exp(a - afc((int) (xl)) - afc((int) (n1 - xl)) - afc((int) (k - xl)) - afc((int) (n2 - k + xl)));
+ kr = Math.exp(a - afc((int) (xr - 1)) - afc((int) (n1 - xr + 1)) - afc((int) (k - xr + 1)) - afc((int) (n2 - k + xr - 1)));
+ lamdl = -Math.log(xl * (n2 - k + xl) / (n1 - xl + 1) / (k - xl + 1));
+ lamdr = -Math.log((n1 - xr + 1) * (k - xr + 1) / xr / (n2 - k + xr));
+ p1 = d + d;
+ p2 = p1 + kl / lamdl;
+ p3 = p2 + kr / lamdr;
+ }
+ int nUv = 0;
+ L30: while (true) {
+
+ u = rand.unifRand() * p3;
+ v = rand.unifRand();
+ nUv++;
+ if (nUv >= 10000) {
+ RError.warning(RError.SHOW_CALLER, Message.GENERIC, String.format("rhyper() branch III: giving up after %d rejections", nUv));
+ return StatsUtil.mlError();
+ }
+
+ if (u < p1) { /* rectangular region */
+ ix = (int) (xl + u);
+ } else if (u <= p2) { /* left tail */
+ ix = (int) (xl + Math.log(v) / lamdl);
+ if (ix < minjx) {
+ continue L30;
+ }
+ v = v * (u - p1) * lamdl;
+ } else { /* right tail */
+ ix = (int) (xr - Math.log(v) / lamdr);
+ if (ix > maxjx) {
+ continue L30;
+ }
+ v = v * (u - p2) * lamdr;
+ }
+
+ /* acceptance/rejection test */
+ boolean reject = true;
+
+ if (m < 100 || ix <= 50) {
+ /* Math.explicit evaluation */
+ /*
+ * The original algorithm (and TOMS 668) have f = f * i * (n2 - k + i) / (n1 -
+ * i) / (k - i); in the (m > ix) case, but the definition of the recurrence
+ * relation on p134 shows that the +1 is needed.
+ */
+ int i;
+ double f = 1.0;
+ if (m < ix) {
+ for (i = m + 1; i <= ix; i++) {
+ f = f * (n1 - i + 1) * (k - i + 1) / (n2 - k + i) / i;
+ }
+ } else if (m > ix) {
+ for (i = ix + 1; i <= m; i++) {
+ f = f * i * (n2 - k + i) / (n1 - i + 1) / (k - i + 1);
+ }
+ }
+ if (v <= f) {
+ reject = false;
+ }
+ } else {
+ // Checkstyle: stop
+ double e, g, r, t, y;
+ double de, dg, dr, ds, dt, gl, gu, nk, nm, ub;
+ double xk, xm, xn, y1, ym, yn, yk, alv;
+ // Checkstyle: resume
+
+ /* squeeze using upper and lower bounds */
+ y = ix;
+ y1 = y + 1.0;
+ ym = y - m;
+ yn = n1 - y + 1.0;
+ yk = k - y + 1.0;
+ nk = n2 - k + y1;
+ r = -ym / y1;
+ s = ym / yn;
+ t = ym / yk;
+ e = -ym / nk;
+ g = yn * yk / (y1 * nk) - 1.0;
+ dg = 1.0;
+ if (g < 0.0) {
+ dg = 1.0 + g;
+ }
+ gu = g * (1.0 + g * (-0.5 + g / 3.0));
+ gl = gu - .25 * (g * g * g * g) / dg;
+ xm = m + 0.5;
+ xn = n1 - m + 0.5;
+ xk = k - m + 0.5;
+ nm = n2 - k + xm;
+ ub = y * gu - m * gl + deltau + xm * r * (1. + r * (-0.5 + r / 3.0)) + xn * s * (1. + s * (-0.5 + s / 3.0)) + xk * t * (1. + t * (-0.5 + t / 3.0)) +
+ nm * e * (1. + e * (-0.5 + e / 3.0));
+ /* test against upper bound */
+ alv = Math.log(v);
+ if (alv > ub) {
+ reject = true;
+ } else {
+ /* test against lower bound */
+ dr = xm * (r * r * r * r);
+ if (r < 0.0) {
+ dr /= (1.0 + r);
+ }
+ ds = xn * (s * s * s * s);
+ if (s < 0.0) {
+ ds /= (1.0 + s);
+ }
+ dt = xk * (t * t * t * t);
+ if (t < 0.0) {
+ dt /= (1.0 + t);
+ }
+ de = nm * (e * e * e * e);
+ if (e < 0.0) {
+ de /= (1.0 + e);
+ }
+ if (alv < ub - 0.25 * (dr + ds + dt + de) + (y + m) * (gl - gu) - deltal) {
+ reject = false;
+ } else {
+ /*
+ * * Stirling's formula to machine accuracy
+ */
+ if (alv <= (a - afc(ix) - afc(n1 - ix) - afc(k - ix) - afc(n2 - k + ix))) {
+ reject = false;
+ } else {
+ reject = true;
+ }
+ }
+ }
+ } // else
+ if (!reject) {
+ break L30; // e.g. if (reject) goto L30;
+ }
+ }
+ }
+
+ /* return appropriate variate */
+ double ix1 = ix;
+ if ((double) kk + (double) kk >= tn) {
+ if ((double) nn1 > (double) nn2) {
+ ix1 = (double) kk - (double) nn2 + ix1;
+ } else {
+ ix1 = (double) nn1 - ix1;
+ }
+ } else {
+ if ((double) nn1 > (double) nn2) {
+ ix1 = (double) kk - ix1;
+ }
+ }
+ return ix1;
+ }
+}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RPois.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RPois.java
index dc2ba3e7b631fe16104485f57e17f8320d8fc164..ffcc999eace1ce6645b72423f3a7ee81f7cbb191 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RPois.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RPois.java
@@ -13,9 +13,10 @@ package com.oracle.truffle.r.library.stats;
import static com.oracle.truffle.r.library.stats.MathConstants.M_1_SQRT_2PI;
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction1_Double;
import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandomNumberProvider;
-public final class RPois {
+public final class RPois implements RandFunction1_Double {
private static final double a0 = -0.5;
private static final double a1 = 0.3333333;
@@ -272,4 +273,9 @@ public final class RPois {
}
return pois;
}
+
+ @Override
+ public double evaluate(double mu, RandomNumberProvider rand) {
+ return rpois(mu, rand);
+ }
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RandGenerationFunctions.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RandGenerationFunctions.java
index 6ad6a060ba1b9524bf728265da6f0f6dfbee6a14..32d44ac1c923010593429876142efcc0008d77ab 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RandGenerationFunctions.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RandGenerationFunctions.java
@@ -19,6 +19,7 @@ import static com.oracle.truffle.r.runtime.RError.SHOW_CALLER;
import java.util.Arrays;
+import com.oracle.truffle.api.CompilerDirectives.CompilationFinal;
import com.oracle.truffle.api.dsl.Cached;
import com.oracle.truffle.api.dsl.Specialization;
import com.oracle.truffle.api.nodes.Node;
@@ -43,7 +44,7 @@ import com.oracle.truffle.r.runtime.rng.RRNG.NormKind;
import com.oracle.truffle.r.runtime.rng.RandomNumberGenerator;
public final class RandGenerationFunctions {
- private static final RDouble DUMMY_VECTOR = RDouble.valueOf(1);
+ @CompilationFinal private static final RDouble DUMMY_VECTOR = RDouble.valueOf(1);
private RandGenerationFunctions() {
// static class
@@ -73,21 +74,26 @@ public final class RandGenerationFunctions {
// inspired by the DEFRAND{X}_REAL and DEFRAND{X}_INT macros in GnuR
- public interface RandFunction3_Int {
- int evaluate(double a, double b, double c, RandomNumberProvider rand);
+ public interface RandFunction3_Double {
+ double evaluate(double a, double b, double c, RandomNumberProvider rand);
}
- public interface RandFunction2_Int extends RandFunction3_Int {
+ public interface RandFunction2_Double extends RandFunction3_Double {
+ double evaluate(double a, double b, RandomNumberProvider rand);
+
@Override
- default int evaluate(double a, double b, double c, RandomNumberProvider rand) {
+ default double evaluate(double a, double b, double c, RandomNumberProvider rand) {
return evaluate(a, b, rand);
}
-
- int evaluate(double a, double b, RandomNumberProvider rand);
}
- public interface RandFunction2_Double {
- double evaluate(double a, double b, RandomNumberProvider rand);
+ public interface RandFunction1_Double extends RandFunction2_Double {
+ double evaluate(double a, RandomNumberProvider rand);
+
+ @Override
+ default double evaluate(double a, double b, RandomNumberProvider rand) {
+ return evaluate(a, rand);
+ }
}
static final class RandGenerationProfiles {
@@ -107,7 +113,7 @@ public final class RandGenerationFunctions {
}
}
- private static RAbstractIntVector evaluate3Int(Node node, RandFunction3_Int function, int lengthIn, RAbstractDoubleVector a, RAbstractDoubleVector b, RAbstractDoubleVector c,
+ private static RAbstractIntVector evaluate3Int(Node node, RandFunction3_Double function, int lengthIn, RAbstractDoubleVector a, RAbstractDoubleVector b, RAbstractDoubleVector c,
RandGenerationProfiles profiles) {
int length = lengthIn;
int aLength = a.getLength();
@@ -115,7 +121,7 @@ public final class RandGenerationFunctions {
int cLength = c.getLength();
if (aLength == 0 || bLength == 0 || cLength == 0) {
profiles.nanResult.enter();
- RError.warning(SHOW_CALLER, RError.Message.NAN_PRODUCED);
+ RError.warning(SHOW_CALLER, RError.Message.NA_PRODUCED);
int[] nansResult = new int[length];
Arrays.fill(nansResult, RRuntime.INT_NA);
return RDataFactory.createIntVector(nansResult, false);
@@ -131,25 +137,29 @@ public final class RandGenerationFunctions {
double aValue = a.getDataAt(i % aLength);
double bValue = b.getDataAt(i % bLength);
double cValue = c.getDataAt(i % cLength);
- int value = function.evaluate(aValue, bValue, cValue, rand);
- if (Double.isNaN(value)) {
+ double value = function.evaluate(aValue, bValue, cValue, rand);
+ if (Double.isNaN(value) || value < Integer.MIN_VALUE || value > Integer.MAX_VALUE) {
profiles.nan.enter();
nans = true;
+ result[i] = RRuntime.INT_NA;
+ } else {
+ result[i] = (int) value;
}
- result[i] = value;
}
RRNG.putRNGState();
if (nans) {
- RError.warning(SHOW_CALLER, RError.Message.NAN_PRODUCED);
+ RError.warning(SHOW_CALLER, RError.Message.NA_PRODUCED);
}
return RDataFactory.createIntVector(result, !nans);
}
- private static RAbstractDoubleVector evaluate2Double(Node node, RandFunction2_Double function, int lengthIn, RAbstractDoubleVector a, RAbstractDoubleVector b, RandGenerationProfiles profiles) {
+ private static RAbstractDoubleVector evaluate3Double(Node node, RandFunction3_Double function, int lengthIn, RAbstractDoubleVector a, RAbstractDoubleVector b, RAbstractDoubleVector c,
+ RandGenerationProfiles profiles) {
int length = lengthIn;
int aLength = a.getLength();
int bLength = b.getLength();
- if (aLength == 0 || bLength == 0) {
+ int cLength = c.getLength();
+ if (aLength == 0 || bLength == 0 || cLength == 0) {
profiles.nanResult.enter();
RError.warning(SHOW_CALLER, RError.Message.NA_PRODUCED);
return createVectorOf(length, RRuntime.DOUBLE_NA);
@@ -165,8 +175,9 @@ public final class RandGenerationFunctions {
for (int i = 0; profiles.loopConditionProfile.inject(i < length); i++) {
double aValue = a.getDataAt(i % aLength);
double bValue = b.getDataAt(i % bLength);
- double value = function.evaluate(aValue, bValue, rand);
- if (Double.isNaN(value) || Double.isNaN(value)) {
+ double cValue = c.getDataAt(i % cLength);
+ double value = function.evaluate(aValue, bValue, cValue, rand);
+ if (Double.isNaN(value) || RRuntime.isNA(value)) {
profiles.nan.enter();
nans = true;
}
@@ -216,10 +227,10 @@ public final class RandGenerationFunctions {
}
public abstract static class Function3_IntNode extends RExternalBuiltinNode.Arg4 {
- private final RandFunction3_Int function;
+ private final RandFunction3_Double function;
@Child private ConvertToLength convertToLength = ConvertToLengthNodeGen.create();
- protected Function3_IntNode(RandFunction3_Int function) {
+ protected Function3_IntNode(RandFunction3_Double function) {
this.function = function;
}
@@ -239,10 +250,10 @@ public final class RandGenerationFunctions {
}
public abstract static class Function2_IntNode extends RExternalBuiltinNode.Arg3 {
- private final RandFunction2_Int function;
+ private final RandFunction2_Double function;
@Child private ConvertToLength convertToLength = ConvertToLengthNodeGen.create();
- protected Function2_IntNode(RandFunction2_Int function) {
+ protected Function2_IntNode(RandFunction2_Double function) {
this.function = function;
}
@@ -260,6 +271,48 @@ public final class RandGenerationFunctions {
}
}
+ public abstract static class Function1_IntNode extends RExternalBuiltinNode.Arg2 {
+ private final RandFunction1_Double function;
+ @Child private ConvertToLength convertToLength = ConvertToLengthNodeGen.create();
+
+ protected Function1_IntNode(RandFunction1_Double function) {
+ this.function = function;
+ }
+
+ @Override
+ protected void createCasts(CastBuilder casts) {
+ ConvertToLength.addLengthCast(casts);
+ casts.arg(1).asDoubleVector();
+ }
+
+ @Specialization
+ protected RAbstractIntVector evaluate(RAbstractVector length, RAbstractDoubleVector a,
+ @Cached("create()") RandGenerationProfiles profiles) {
+ return evaluate3Int(this, function, convertToLength.execute(length), a, DUMMY_VECTOR, DUMMY_VECTOR, profiles);
+ }
+ }
+
+ public abstract static class Function1_DoubleNode extends RExternalBuiltinNode.Arg2 {
+ private final RandFunction1_Double function;
+ @Child private ConvertToLength convertToLength = ConvertToLengthNodeGen.create();
+
+ protected Function1_DoubleNode(RandFunction1_Double function) {
+ this.function = function;
+ }
+
+ @Override
+ protected void createCasts(CastBuilder casts) {
+ ConvertToLength.addLengthCast(casts);
+ casts.arg(1).asDoubleVector();
+ }
+
+ @Specialization
+ protected RAbstractDoubleVector evaluate(RAbstractVector length, RAbstractDoubleVector a,
+ @Cached("create()") RandGenerationProfiles profiles) {
+ return evaluate3Double(this, function, convertToLength.execute(length), a, DUMMY_VECTOR, DUMMY_VECTOR, profiles);
+ }
+ }
+
public abstract static class Function2_DoubleNode extends RExternalBuiltinNode.Arg3 {
private final RandFunction2_Double function;
@Child private ConvertToLength convertToLength = ConvertToLengthNodeGen.create();
@@ -278,7 +331,7 @@ public final class RandGenerationFunctions {
@Specialization
protected RAbstractDoubleVector evaluate(RAbstractVector length, RAbstractDoubleVector a, RAbstractDoubleVector b,
@Cached("create()") RandGenerationProfiles profiles) {
- return evaluate2Double(this, function, convertToLength.execute(length), a, b, profiles);
+ return evaluate3Double(this, function, convertToLength.execute(length), a, b, DUMMY_VECTOR, profiles);
}
}
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Rbinom.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Rbinom.java
index c63e7e44f3249f3bdd66e7595efa23d543fa2ff8..db6545e1e193ce0332d1ef64f35802b55ab55aa2 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Rbinom.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Rbinom.java
@@ -12,18 +12,18 @@
*/
package com.oracle.truffle.r.library.stats;
-import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction2_Int;
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction2_Double;
import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandomNumberProvider;
import com.oracle.truffle.r.runtime.RRuntime;
// transcribed from rbinom.c
-public final class Rbinom implements RandFunction2_Int {
+public final class Rbinom implements RandFunction2_Double {
private final Qbinom qbinom = new Qbinom();
@Override
- public int evaluate(double nin, double pp, RandomNumberProvider rand) {
+ public double evaluate(double nin, double pp, RandomNumberProvider rand) {
double psave = -1.0;
int nsave = -1;
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Rt.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Rt.java
new file mode 100644
index 0000000000000000000000000000000000000000..2ba6a6271371c136fe6a29a544797de96be0513f
--- /dev/null
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Rt.java
@@ -0,0 +1,32 @@
+/*
+ * This material is distributed under the GNU General Public License
+ * Version 2. You may review the terms of this license at
+ * http://www.gnu.org/licenses/gpl-2.0.html
+ *
+ * Copyright (C) 1998 Ross Ihaka
+ * Copyright (c) 1998--2008, The R Core Team
+ * Copyright (c) 2016, 2016, Oracle and/or its affiliates
+ *
+ * All rights reserved.
+ */
+package com.oracle.truffle.r.library.stats;
+
+import static com.oracle.truffle.r.library.stats.RChisq.rchisq;
+
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction1_Double;
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandomNumberProvider;
+
+public final class Rt implements RandFunction1_Double {
+ @Override
+ public double evaluate(double df, RandomNumberProvider rand) {
+ if (Double.isNaN(df) || df <= 0.0) {
+ return StatsUtil.mlError();
+ }
+
+ if (!Double.isFinite(df)) {
+ return rand.normRand();
+ } else {
+ return rand.normRand() / Math.sqrt(rchisq(df, rand) / df);
+ }
+ }
+}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Signrank.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Signrank.java
new file mode 100644
index 0000000000000000000000000000000000000000..422d76cb6f7c006a3e94ad3779e76fc6d9cfb0b2
--- /dev/null
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Signrank.java
@@ -0,0 +1,53 @@
+/*
+ * This material is distributed under the GNU General Public License
+ * Version 2. You may review the terms of this license at
+ * http://www.gnu.org/licenses/gpl-2.0.html
+ *
+ * Copyright (c) 1999--2014, The R Core Team
+ * Copyright (c) 2016, Oracle and/or its affiliates
+ *
+ * All rights reserved.
+ */
+
+package com.oracle.truffle.r.library.stats;
+
+import static com.oracle.truffle.r.library.stats.MathConstants.forceint;
+
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction1_Double;
+import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandomNumberProvider;
+
+public final class Signrank {
+
+ public static final class RSignrank implements RandFunction1_Double {
+ @Override
+ public double evaluate(double n, RandomNumberProvider rand) {
+ int i;
+ int k;
+ double r;
+
+ if (Double.isNaN(n)) {
+ return n;
+ }
+ if (Double.isInfinite(n)) {
+ // In GnuR these "results" seem to be generated due to the behaviour of R_forceint,
+ // and the "(int) n" cast, which ends up casting +/-infinity to integer...
+ return n < 0 ? StatsUtil.mlError() : 0;
+ }
+
+ n = forceint(n);
+ if (n < 0) {
+ return StatsUtil.mlError();
+ }
+
+ if (n == 0) {
+ return 0;
+ }
+ r = 0.0;
+ k = (int) n;
+ for (i = 0; i < k; i++) {
+ r += (i + 1) * Math.floor(rand.unifRand() + 0.5);
+ }
+ return r;
+ }
+ }
+}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/StatsUtil.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/StatsUtil.java
index 7d509e45bc4d1a5df54d4c882dba89aa283f117d..7a5348110e6dcddb2b05578952f40a334ac2f354 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/StatsUtil.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/StatsUtil.java
@@ -12,6 +12,8 @@
*/
package com.oracle.truffle.r.library.stats;
+import static com.oracle.truffle.r.library.stats.LBeta.lbeta;
+
import com.oracle.truffle.api.CompilerDirectives.CompilationFinal;
import com.oracle.truffle.api.CompilerDirectives.TruffleBoundary;
@@ -119,6 +121,10 @@ public class StatsUtil {
return giveLog ? -0.5 * Math.log(f) + x : Math.exp(x) / Math.sqrt(f);
}
+ public static double lfastchoose(double n, double k) {
+ return -Math.log(n + 1.) - lbeta(n - k + 1., k + 1.);
+ }
+
public static double fsign(double x, double y) {
if (Double.isNaN(x) || Double.isNaN(y)) {
return x + y;
diff --git a/com.oracle.truffle.r.nodes.builtin/src/com/oracle/truffle/r/nodes/builtin/base/foreign/ForeignFunctions.java b/com.oracle.truffle.r.nodes.builtin/src/com/oracle/truffle/r/nodes/builtin/base/foreign/ForeignFunctions.java
index eb498d9e37e77a4be870105a079948f5334602fd..f865b988a9092be62f557cc71c8664c663af3750 100644
--- a/com.oracle.truffle.r.nodes.builtin/src/com/oracle/truffle/r/nodes/builtin/base/foreign/ForeignFunctions.java
+++ b/com.oracle.truffle.r.nodes.builtin/src/com/oracle/truffle/r/nodes/builtin/base/foreign/ForeignFunctions.java
@@ -56,16 +56,23 @@ import com.oracle.truffle.r.library.stats.Qbinom;
import com.oracle.truffle.r.library.stats.Qnorm;
import com.oracle.truffle.r.library.stats.RBeta;
import com.oracle.truffle.r.library.stats.RCauchy;
+import com.oracle.truffle.r.library.stats.RChisq;
+import com.oracle.truffle.r.library.stats.RExp;
import com.oracle.truffle.r.library.stats.RGamma;
+import com.oracle.truffle.r.library.stats.RGeom;
+import com.oracle.truffle.r.library.stats.RHyper;
import com.oracle.truffle.r.library.stats.RLogis;
import com.oracle.truffle.r.library.stats.RNbinomMu;
import com.oracle.truffle.r.library.stats.RNchisq;
+import com.oracle.truffle.r.library.stats.RPois;
import com.oracle.truffle.r.library.stats.RWeibull;
import com.oracle.truffle.r.library.stats.RandGenerationFunctionsFactory;
import com.oracle.truffle.r.library.stats.Rbinom;
import com.oracle.truffle.r.library.stats.Rf;
import com.oracle.truffle.r.library.stats.Rnorm;
+import com.oracle.truffle.r.library.stats.Rt;
import com.oracle.truffle.r.library.stats.Runif;
+import com.oracle.truffle.r.library.stats.Signrank.RSignrank;
import com.oracle.truffle.r.library.stats.SplineFunctionsFactory.SplineCoefNodeGen;
import com.oracle.truffle.r.library.stats.SplineFunctionsFactory.SplineEvalNodeGen;
import com.oracle.truffle.r.library.stats.StatsFunctionsFactory;
@@ -392,7 +399,21 @@ public class ForeignFunctions {
case "rnbinom_mu":
return RandGenerationFunctionsFactory.Function2_DoubleNodeGen.create(new RNbinomMu());
case "rwilcox":
- return RandGenerationFunctionsFactory.Function2_DoubleNodeGen.create(new RWilcox());
+ return RandGenerationFunctionsFactory.Function2_IntNodeGen.create(new RWilcox());
+ case "rchisq":
+ return RandGenerationFunctionsFactory.Function1_DoubleNodeGen.create(new RChisq());
+ case "rexp":
+ return RandGenerationFunctionsFactory.Function1_DoubleNodeGen.create(new RExp());
+ case "rgeom":
+ return RandGenerationFunctionsFactory.Function1_IntNodeGen.create(new RGeom());
+ case "rpois":
+ return RandGenerationFunctionsFactory.Function1_IntNodeGen.create(new RPois());
+ case "rt":
+ return RandGenerationFunctionsFactory.Function1_DoubleNodeGen.create(new Rt());
+ case "rsignrank":
+ return RandGenerationFunctionsFactory.Function1_IntNodeGen.create(new RSignrank());
+ case "rhyper":
+ return RandGenerationFunctionsFactory.Function3_IntNodeGen.create(new RHyper());
case "qgamma":
return StatsFunctionsFactory.Function3_2NodeGen.create(new QgammaFunc());
case "dbinom":
diff --git a/com.oracle.truffle.r.test/src/com/oracle/truffle/r/test/ExpectedTestOutput.test b/com.oracle.truffle.r.test/src/com/oracle/truffle/r/test/ExpectedTestOutput.test
index 95e5030c9f694b7236f4659956ff11cd084181cf..44ab3379493cd615579b7342e03d5f1f9037021a 100644
--- a/com.oracle.truffle.r.test/src/com/oracle/truffle/r/test/ExpectedTestOutput.test
+++ b/com.oracle.truffle.r.test/src/com/oracle/truffle/r/test/ExpectedTestOutput.test
@@ -113472,6 +113472,55 @@ $variables
list(y, z)
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions1#Output.IgnoreWhitespace#
+#set.seed(10); rsignrank(12, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.6, 0.0653, 0.000123, 32e-80, 10))
+ [1] NA NA 0 NA NA 1 0 NA 0 0 0 21
+Warning message:
+In rsignrank(12, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.6, 0.0653, :
+ NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions1#Output.IgnoreWhitespace#
+#set.seed(2); rchisq(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, 32e-80, 8833, 79e70))
+ [1] NaN NaN NaN NaN NaN
+ [6] 9.582518e-02 7.012801e-02 NaN 2.294588e-10 0.000000e+00
+[11] 0.000000e+00 8.813564e+03 7.900000e+71
+Warning message:
+In rchisq(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, :
+ NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions1#Output.IgnoreWhitespace#
+#set.seed(2); rexp(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, 32e-80, 8833, 79e70))
+ [1] NaN NaN 0.000000e+00 0.000000e+00 NaN
+ [6] 1.865352e+00 1.349160e+00 NaN 2.245830e+00 1.407081e+04
+[11] 2.797693e+77 7.550069e-05 1.359958e-72
+Warning message:
+In rexp(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, :
+ NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions1#Output.IgnoreWhitespace#
+#set.seed(2); rgeom(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, 32e-80, 8833, 79e70))
+ [1] NA NA NA NA NA 0 1 NA 32 11643 NA NA
+[13] NA
+Warning message:
+In rgeom(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, :
+ NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions1#Output.IgnoreWhitespace#
+#set.seed(2); rpois(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, 32e-80, 8833, 79e70))
+ [1] NA NA NA NA NA 0 0 NA 0 0 0 8981 NA
+Warning message:
+In rpois(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, :
+ NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions1#Output.IgnoreWhitespace#
+#set.seed(2); rt(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, 32e-80, 8833, 79e70))
+ [1] NaN NaN -0.89691455 NaN NaN 0.19975785
+ [7] -0.06927937 NaN 0.96073767 -Inf Inf 0.87841766
+[13] 1.01282869
+Warning message:
+In rt(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, :
+ NAs produced
+
##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions2#Output.IgnoreWhitespace#
#set.seed(1); rbeta(10, 10, 10)
[1] 0.4202441 0.5231868 0.3929161 0.7104015 0.5416782 0.3949132 0.5618393
@@ -113817,6 +113866,54 @@ In runif(2, 2, numeric()) : NAs produced
Warning message:
In runif(2, numeric(), 2) : NAs produced
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions3#
+#rhyper(1, -5, 5, 20)
+[1] NA
+Warning message:
+In rhyper(1, -5, 5, 20) : NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions3#
+#rhyper(1, 5, 5, -1/0)
+[1] NA
+Warning message:
+In rhyper(1, 5, 5, -1/0) : NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions3#
+#rhyper(1, 5, 5, 1/0)
+[1] NA
+Warning message:
+In rhyper(1, 5, 5, 1/0) : NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions3#
+#rhyper(1, 5, 5, 20)
+[1] NA
+Warning message:
+In rhyper(1, 5, 5, 20) : NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions3#
+#rhyper(1, 5, NaN, 20)
+[1] NA
+Warning message:
+In rhyper(1, 5, NaN, 20) : NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions3#
+#rhyper(1, NA, 5, 20)
+[1] NA
+Warning message:
+In rhyper(1, NA, 5, 20) : NAs produced
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions3#
+#set.seed(3); rhyper(2, 1000, 1000, 5)
+[1] 1 3
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions3#
+#set.seed(3); rhyper(3, 10, 10, 5)
+[1] 2 3 2
+
+##com.oracle.truffle.r.test.library.stats.TestRandGenerationFunctions.testFunctions3#
+#set.seed(3); rhyper(3, 10, 79e70, 2)
+[1] 0 0 0
+
##com.oracle.truffle.r.test.library.stats.TestStats.testCor#
#{ as.integer(cor(c(1,2,3),c(1,2,5))*10000000) }
[1] 9607689
diff --git a/com.oracle.truffle.r.test/src/com/oracle/truffle/r/test/library/stats/TestRandGenerationFunctions.java b/com.oracle.truffle.r.test/src/com/oracle/truffle/r/test/library/stats/TestRandGenerationFunctions.java
index 9a6667e769174c30dd3cc647498b478ce0c0886b..e43514d3faaa56a850ee03305afc137ebca3a397 100644
--- a/com.oracle.truffle.r.test/src/com/oracle/truffle/r/test/library/stats/TestRandGenerationFunctions.java
+++ b/com.oracle.truffle.r.test/src/com/oracle/truffle/r/test/library/stats/TestRandGenerationFunctions.java
@@ -58,4 +58,30 @@ public class TestRandGenerationFunctions extends TestBase {
// wrong parameters
assertEval("runif(-1, 1, 2)");
}
+
+ private static final String[] FUNCTION1_NAMES = {"rchisq", "rexp", "rgeom", "rpois", "rt"};
+
+ @Test
+ public void testFunctions1() {
+ assertEval(Output.IgnoreWhitespace, template("set.seed(2); %0(13, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.5, 0.0653, 0.000123, 32e-80, 8833, 79e70))", FUNCTION1_NAMES));
+ // Note: signrank has loop with 'n' iterations: we have to leave out the large numbers
+ assertEval(Output.IgnoreWhitespace, "set.seed(10); rsignrank(12, c(NA, NaN, 1/0, -1/0, -1, 1, 0.3, -0.6, 0.0653, 0.000123, 32e-80, 10))");
+ }
+
+ @Test
+ public void testFunctions3() {
+ // error: drawn (20) mare than red + blue (5+5)
+ assertEval("rhyper(1, 5, 5, 20)");
+ // error: negative number of balls
+ assertEval("rhyper(1, -5, 5, 20)");
+ // common errors with NA, NaN, Inf
+ assertEval("rhyper(1, NA, 5, 20)");
+ assertEval("rhyper(1, 5, NaN, 20)");
+ assertEval("rhyper(1, 5, 5, 1/0)");
+ assertEval("rhyper(1, 5, 5, -1/0)");
+ // few simple tests (note: rhyper seems to be quite slow even in GnuR).
+ assertEval("set.seed(3); rhyper(3, 10, 10, 5)");
+ assertEval("set.seed(3); rhyper(2, 1000, 1000, 5)");
+ assertEval("set.seed(3); rhyper(3, 10, 79e70, 2)");
+ }
}
diff --git a/mx.fastr/copyrights/overrides b/mx.fastr/copyrights/overrides
index 13c0f05a46af23b63714e3454733f0e4fb35fc86..4d7743770e589cb4ac224e40a635382c54f7265f 100644
--- a/mx.fastr/copyrights/overrides
+++ b/mx.fastr/copyrights/overrides
@@ -35,13 +35,18 @@ com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Cdist.java,g
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/CompleteCases.java,gnu_r_gentleman_ihaka2.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Covcor.java,gnu_r.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Cutree.java,gnu_r.core.copyright
+com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Signrank.java,gnu_r.core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dbinom.java,gnu_r.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DoubleCentre.java,gnu_r.core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DPQ.java,gnu_r.core.copyright
+com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Rt.java,gnu_r_ihaka_core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/GammaFunctions.java,gnu_r_qgamma.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/MathConstants.java,gnu_r.copyright
+com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QHyper.java,gnu_r_ihaka_core.copyright
+com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RHyper.java,gnu_r_ihaka.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/MathInit.java,gnu_r.core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pbeta.java,gnu_r.core.copyright
+com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/LBeta.java,gnu_r_ihaka.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pbinom.java,gnu_r_ihaka.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pf.java,gnu_r_ihaka_core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pnorm.java,gnu_r_ihaka.copyright
@@ -63,6 +68,8 @@ com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RNbinomMu.ja
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RLogis.java,gnu_r_ihaka_core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Rf.java,gnu_r_ihaka_core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RChisq.java,gnu_r_ihaka_core.copyright
+com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RExp.java,gnu_r_ihaka_core.copyright
+com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RGeom.java,gnu_r_ihaka_core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RNchisq.java,gnu_r.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Wilcox.java,gnu_r.core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RWeibull.java,gnu_r_ihaka_core.copyright