diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DNBeta.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DNBeta.java
index c357b66df5d85cd3cc187e6fde1abc17771e07b1..4410c96c63e0d267354d4ef18a9c17cd699e8d21 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DNBeta.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DNBeta.java
@@ -5,13 +5,13 @@
*
* Copyright (C) 1998 Ross Ihaka
* Copyright (c) 2000-12, The R Core Team
- * Copyright (c) 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
package com.oracle.truffle.r.library.stats;
-import static com.oracle.truffle.r.library.stats.GammaFunctions.dpoisRaw;
+import static com.oracle.truffle.r.library.stats.DPois.dpoisRaw;
import com.oracle.truffle.r.library.stats.StatsFunctions.Function4_1;
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DNChisq.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DNChisq.java
index 8fe5845b3441e24e8921e601459456e7347eb6ac..656bb8daa0157f7256dc75c23c5bd30664f9e152 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DNChisq.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DNChisq.java
@@ -6,13 +6,13 @@
* Copyright (C) 1998 Ross Ihaka
* Copyright (c) 2000-15, The R Core Team
* Copyright (c) 2004-15, The R Foundation
- * Copyright (c) 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
package com.oracle.truffle.r.library.stats;
-import static com.oracle.truffle.r.library.stats.GammaFunctions.dpoisRaw;
+import static com.oracle.truffle.r.library.stats.DPois.dpoisRaw;
import com.oracle.truffle.r.library.stats.Chisq.DChisq;
import com.oracle.truffle.r.library.stats.StatsFunctions.Function3_1;
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DPois.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DPois.java
index bc14734a766afbf46ea03e5df58082a6f4701d47..82dc7bb0e0c751e8c4724700d5567b260890ceec 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DPois.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DPois.java
@@ -4,7 +4,7 @@
* http://www.gnu.org/licenses/gpl-2.0.html
*
* Copyright (c) 2000--2014, The R Core Team
- * Copyright (c) 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
@@ -12,11 +12,17 @@
// Author: Catherine Loader, catherine@research.bell-labs.com, October 23, 2000.
package com.oracle.truffle.r.library.stats;
-import static com.oracle.truffle.r.library.stats.GammaFunctions.dpoisRaw;
+import static com.oracle.truffle.r.library.stats.DPQ.rd0;
+import static com.oracle.truffle.r.library.stats.DPQ.rd1;
+import static com.oracle.truffle.r.library.stats.DPQ.rdexp;
+import static com.oracle.truffle.r.library.stats.DPQ.rdfexp;
+import static com.oracle.truffle.r.library.stats.MathConstants.DBL_MIN;
+import static com.oracle.truffle.r.library.stats.MathConstants.M_2PI;
import static com.oracle.truffle.r.library.stats.RMath.forceint;
import com.oracle.truffle.r.library.stats.DPQ.EarlyReturn;
import com.oracle.truffle.r.library.stats.StatsFunctions.Function2_1;
+import com.oracle.truffle.r.runtime.RRuntime;
public final class DPois implements Function2_1 {
@Override
@@ -39,4 +45,26 @@ public final class DPois implements Function2_1 {
return dpoisRaw(forceint(x), lambda, giveLog);
}
+
+ public static double dpoisRaw(double x, double lambda, boolean giveLog) {
+ /*
+ * x >= 0 ; integer for dpois(), but not e.g. for pgamma()! lambda >= 0
+ */
+ if (lambda == 0) {
+ return (x == 0) ? rd1(giveLog) : rd0(giveLog);
+ }
+ if (!RRuntime.isFinite(lambda)) {
+ return rd0(giveLog);
+ }
+ if (x < 0) {
+ return rd0(giveLog);
+ }
+ if (x <= lambda * DBL_MIN) {
+ return (rdexp(-lambda, giveLog));
+ }
+ if (lambda < x * DBL_MIN) {
+ return (rdexp(-lambda + x * Math.log(lambda) - GammaFunctions.lgammafn(x + 1), giveLog));
+ }
+ return rdfexp(M_2PI * x, -RMath.stirlerr(x) - RMath.bd0(x, lambda), giveLog);
+ }
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dbinom.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dbinom.java
index 0486d311a7cce855fd617c28224dae4a39454e27..a44f26278278ef964622f2c3a0cbc547e0d1018f 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dbinom.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dbinom.java
@@ -5,7 +5,7 @@
*
* Copyright (c) 2000-2014, The R Core Team
* Copyright (c) 2008, The R Foundation
- * Copyright (c) 2016, 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
@@ -38,11 +38,11 @@ public final class Dbinom implements StatsFunctions.Function3_1 {
if (n == 0) {
return DPQ.rd1(giveLog);
}
- double lc = (p < 0.1) ? -GammaFunctions.bd0(n, n * q) - n * p : n * Math.log(q);
+ double lc = (p < 0.1) ? -RMath.bd0(n, n * q) - n * p : n * Math.log(q);
return DPQ.rdexp(lc, giveLog);
}
if (x == n) {
- double lc = (q < 0.1) ? -GammaFunctions.bd0(n, n * p) - n * q : n * Math.log(p);
+ double lc = (q < 0.1) ? -RMath.bd0(n, n * p) - n * q : n * Math.log(p);
return DPQ.rdexp(lc, giveLog);
}
if (x < 0 || x > n) {
@@ -53,7 +53,7 @@ public final class Dbinom implements StatsFunctions.Function3_1 {
* n*p or n*q can underflow to zero if n and p or q are small. This used to occur in dbeta,
* and gives NaN as from R 2.3.0.
*/
- double lc = GammaFunctions.stirlerr(n) - GammaFunctions.stirlerr(x) - GammaFunctions.stirlerr(n - x) - GammaFunctions.bd0(x, n * p) - GammaFunctions.bd0(n - x, n * q);
+ double lc = RMath.stirlerr(n) - RMath.stirlerr(x) - RMath.stirlerr(n - x) - RMath.bd0(x, n * p) - RMath.bd0(n - x, n * q);
/* f = (M_2PI*x*(n-x))/n; could overflow or underflow */
/*
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dt.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dt.java
index 45eb01fe6b8bd4ead373dccfced2c307d3e50321..86cc49430547bea99ff6f3fcef4f56a7d84b9503 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dt.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dt.java
@@ -4,7 +4,7 @@
* http://www.gnu.org/licenses/gpl-2.0.html
*
* Copyright (c) 2000--2014, The R Core Team
- * Copyright (c) 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
@@ -12,11 +12,11 @@
// Author: Catherine Loader, catherine@research.bell-labs.com, October 23, 2000.
package com.oracle.truffle.r.library.stats;
-import static com.oracle.truffle.r.library.stats.GammaFunctions.bd0;
-import static com.oracle.truffle.r.library.stats.GammaFunctions.stirlerr;
import static com.oracle.truffle.r.library.stats.MathConstants.DBL_EPSILON;
import static com.oracle.truffle.r.library.stats.MathConstants.M_1_SQRT_2PI;
import static com.oracle.truffle.r.library.stats.MathConstants.M_LN_SQRT_2PI;
+import static com.oracle.truffle.r.library.stats.RMath.bd0;
+import static com.oracle.truffle.r.library.stats.RMath.stirlerr;
import com.oracle.truffle.r.library.stats.StatsFunctions.Function2_1;
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 8b221ce9340d4d5f9701d510a3ac14de224879ac..a1c57d865d60e600618a1acbc8e1e3f6203e760b 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
@@ -8,7 +8,7 @@
* Copyright (c) 1998--2014, The R Core Team
* Copyright (c) 2002--2010, The R Foundation
* Copyright (C) 2005--2006, Morten Welinder
- * Copyright (c) 2014, 2016, Oracle and/or its affiliates
+ * Copyright (c) 2014, 2017, Oracle and/or its affiliates
*
* based on AS 91 (C) 1979 Royal Statistical Society
* and on AS 111 (C) 1977 Royal Statistical Society
@@ -21,7 +21,6 @@ package com.oracle.truffle.r.library.stats;
import static com.oracle.truffle.r.library.stats.DPQ.rd0;
import static com.oracle.truffle.r.library.stats.DPQ.rd1;
import static com.oracle.truffle.r.library.stats.DPQ.rdexp;
-import static com.oracle.truffle.r.library.stats.DPQ.rdfexp;
import static com.oracle.truffle.r.library.stats.DPQ.rdt0;
import static com.oracle.truffle.r.library.stats.DPQ.rdt1;
import static com.oracle.truffle.r.library.stats.DPQ.rdtclog;
@@ -30,7 +29,6 @@ import static com.oracle.truffle.r.library.stats.DPQ.rlog1exp;
import static com.oracle.truffle.r.library.stats.MathConstants.DBL_EPSILON;
import static com.oracle.truffle.r.library.stats.MathConstants.DBL_MIN;
import static com.oracle.truffle.r.library.stats.MathConstants.M_1_SQRT_2PI;
-import static com.oracle.truffle.r.library.stats.MathConstants.M_2PI;
import static com.oracle.truffle.r.library.stats.MathConstants.M_LN2;
import static com.oracle.truffle.r.library.stats.MathConstants.M_SQRT_32;
import static com.oracle.truffle.r.library.stats.RMath.fmax2;
@@ -76,87 +74,6 @@ public abstract class GammaFunctions {
// TODO Many of the functions below that are not directly supporting qgamma should eventually be
// factored out to support those R functions they represent.
- // TODO for all occurrences of ML_ERR_return_NAN;
- // this expands to:
- // ML_ERROR(ME_DOMAIN, ""); return ML_NAN;
- // i.e., raise "argument out of domain" error and return NaN
-
- //
- // stirlerr
- //
-
- private static final double S0 = 0.083333333333333333333; /* 1/12 */
- private static final double S1 = 0.00277777777777777777778; /* 1/360 */
- private static final double S2 = 0.00079365079365079365079365; /* 1/1260 */
- private static final double S3 = 0.000595238095238095238095238; /* 1/1680 */
- private static final double S4 = 0.0008417508417508417508417508; /* 1/1188 */
-
- /*
- * error for 0, 0.5, 1.0, 1.5, ..., 14.5, 15.0.
- */
- @CompilationFinal private static final double[] sferr_halves = new double[]{0.0, /*
- * n=0 - wrong,
- * place holder
- * only
- */
- 0.1534264097200273452913848, /* 0.5 */
- 0.0810614667953272582196702, /* 1.0 */
- 0.0548141210519176538961390, /* 1.5 */
- 0.0413406959554092940938221, /* 2.0 */
- 0.03316287351993628748511048, /* 2.5 */
- 0.02767792568499833914878929, /* 3.0 */
- 0.02374616365629749597132920, /* 3.5 */
- 0.02079067210376509311152277, /* 4.0 */
- 0.01848845053267318523077934, /* 4.5 */
- 0.01664469118982119216319487, /* 5.0 */
- 0.01513497322191737887351255, /* 5.5 */
- 0.01387612882307074799874573, /* 6.0 */
- 0.01281046524292022692424986, /* 6.5 */
- 0.01189670994589177009505572, /* 7.0 */
- 0.01110455975820691732662991, /* 7.5 */
- 0.010411265261972096497478567, /* 8.0 */
- 0.009799416126158803298389475, /* 8.5 */
- 0.009255462182712732917728637, /* 9.0 */
- 0.008768700134139385462952823, /* 9.5 */
- 0.008330563433362871256469318, /* 10.0 */
- 0.007934114564314020547248100, /* 10.5 */
- 0.007573675487951840794972024, /* 11.0 */
- 0.007244554301320383179543912, /* 11.5 */
- 0.006942840107209529865664152, /* 12.0 */
- 0.006665247032707682442354394, /* 12.5 */
- 0.006408994188004207068439631, /* 13.0 */
- 0.006171712263039457647532867, /* 13.5 */
- 0.005951370112758847735624416, /* 14.0 */
- 0.005746216513010115682023589, /* 14.5 */
- 0.005554733551962801371038690 /* 15.0 */
- };
-
- static double stirlerr(double n) {
-
- double nn;
-
- if (n <= 15.0) {
- nn = n + n;
- if (nn == (int) nn) {
- return (sferr_halves[(int) nn]);
- }
- return (lgammafn(n + 1.) - (n + 0.5) * Math.log(n) + n - M_LN_SQRT_2PI);
- }
-
- nn = n * n;
- if (n > 500) {
- return ((S0 - S1 / nn) / n);
- }
- if (n > 80) {
- return ((S0 - (S1 - S2 / nn) / nn) / n);
- }
- if (n > 35) {
- return ((S0 - (S1 - (S2 - S3 / nn) / nn) / nn) / n);
- }
- /* 15 < n <= 35 : */
- return ((S0 - (S1 - (S2 - (S3 - S4 / nn) / nn) / nn) / nn) / n);
- }
-
//
// lgammacor
//
@@ -178,8 +95,7 @@ public abstract class GammaFunctions {
double tmp;
if (x < 10) {
- // TODO ML_ERR_return_NAN
- return Double.NaN;
+ return RMathError.defaultError();
} else if (x >= lgc_xmax) {
/* allow to underflow below */
RMathError.error(MLError.UNDERFLOW, "lgammacor");
@@ -319,7 +235,7 @@ public abstract class GammaFunctions {
value *= i;
}
} else { /* normal case */
- value = Math.exp((y - 0.5) * Math.log(y) - y + M_LN_SQRT_2PI + ((2 * y == (int) (2 * y)) ? stirlerr(y) : lgammacor(y)));
+ value = Math.exp((y - 0.5) * Math.log(y) - y + M_LN_SQRT_2PI + ((2 * y == (int) (2 * y)) ? RMath.stirlerr(y) : lgammacor(y)));
}
if (x > 0) {
return value;
@@ -845,12 +761,12 @@ public abstract class GammaFunctions {
return rd0(giveLog);
}
if (xplus1 > 1) {
- return dpoisRaw(xplus1 - 1, lambda, giveLog);
+ return DPois.dpoisRaw(xplus1 - 1, lambda, giveLog);
}
if (lambda > Math.abs(xplus1 - 1) * M_cutoff) {
return rdexp(-lambda - lgammafn(xplus1), giveLog);
} else {
- double d = dpoisRaw(xplus1, lambda, giveLog);
+ double d = DPois.dpoisRaw(xplus1, lambda, giveLog);
return giveLog ? d + Math.log(xplus1 / lambda) : d * (xplus1 / lambda);
}
}
@@ -880,7 +796,7 @@ public abstract class GammaFunctions {
double f1 = logp ? RMath.log1p(sum) : 1 + sum;
double f2;
if (alph > 1) {
- f2 = dpoisRaw(alph, x, logp);
+ f2 = DPois.dpoisRaw(alph, x, logp);
f2 = logp ? f2 + x : f2 * Math.exp(x);
} else if (logp) {
f2 = alph * Math.log(x) - lgamma1p(alph);
@@ -1082,7 +998,7 @@ public abstract class GammaFunctions {
return 1 / sum;
} else {
- double d = dnorm(localX, 0., 1., false);
+ double d = new DNorm().evaluate(localX, 0., 1., false);
return d / Math.exp(lp);
}
}
@@ -1158,7 +1074,7 @@ public abstract class GammaFunctions {
double ndOverP = dpnorm(s2pt, !lowerTail, np);
return np + RMath.log1p(f * ndOverP);
} else {
- double nd = dnorm(s2pt, 0., 1., logp);
+ double nd = new DNorm().evaluate(s2pt, 0., 1., logp);
return np + f * nd;
}
} /* ppois_asymp() */
@@ -1241,96 +1157,6 @@ public abstract class GammaFunctions {
return pgammaRaw(localX, alph, lowerTail, logp);
}
- //
- // dpois
- //
-
- public static double dpoisRaw(double x, double lambda, boolean giveLog) {
- /*
- * x >= 0 ; integer for dpois(), but not e.g. for pgamma()! lambda >= 0
- */
- if (lambda == 0) {
- return (x == 0) ? rd1(giveLog) : rd0(giveLog);
- }
- if (!RRuntime.isFinite(lambda)) {
- return rd0(giveLog);
- }
- if (x < 0) {
- return rd0(giveLog);
- }
- if (x <= lambda * DBL_MIN) {
- return (rdexp(-lambda, giveLog));
- }
- if (lambda < x * DBL_MIN) {
- return (rdexp(-lambda + x * Math.log(lambda) - lgammafn(x + 1), giveLog));
- }
- return rdfexp(M_2PI * x, -stirlerr(x) - bd0(x, lambda), giveLog);
- }
-
- //
- // bd0
- //
-
- static double bd0(double x, double np) {
- double ej;
- double s;
- double s1;
- double v;
- int j;
-
- if (!RRuntime.isFinite(x) || !RRuntime.isFinite(np) || np == 0.0) {
- return RMathError.defaultError();
- }
-
- if (Math.abs(x - np) < 0.1 * (x + np)) {
- v = (x - np) / (x + np);
- s = (x - np) * v; /* s using v -- change by MM */
- ej = 2 * x * v;
- v = v * v;
- for (j = 1;; j++) { /* Taylor series */
- ej *= v;
- s1 = s + ej / ((j << 1) + 1);
- if (s1 == s) { /* last term was effectively 0 */
- return s1;
- }
- s = s1;
- }
- }
- /* else: | x - np | is not too small */
- return x * Math.log(x / np) + np - x;
- }
-
- //
- // dnorm
- //
-
- static double dnorm(double x, double mu, double sigma, boolean giveLog) {
- double localX = x;
- if (Double.isNaN(localX) || Double.isNaN(mu) || Double.isNaN(sigma)) {
- return localX + mu + sigma;
- }
- if (!RRuntime.isFinite(sigma)) {
- return rd0(giveLog);
- }
- if (!RRuntime.isFinite(localX) && mu == localX) {
- return Double.NaN; /* x-mu is NaN */
- }
- if (sigma <= 0) {
- if (sigma < 0) {
- return RMathError.defaultError();
- }
- /* sigma == 0 */
- return (localX == mu) ? Double.POSITIVE_INFINITY : rd0(giveLog);
- }
- localX = (localX - mu) / sigma;
-
- if (!RRuntime.isFinite(localX)) {
- return rd0(giveLog);
- }
- return giveLog ? -(M_LN_SQRT_2PI + 0.5 * localX * localX + Math.log(sigma)) : M_1_SQRT_2PI * Math.exp(-0.5 * localX * localX) / sigma;
- /* M_1_SQRT_2PI = 1 / sqrt(2 * pi) */
- }
-
//
// dgamma
//
@@ -1361,11 +1187,11 @@ public abstract class GammaFunctions {
}
if (shape < 1) {
- pr = dpoisRaw(shape, x / scale, giveLog);
+ pr = DPois.dpoisRaw(shape, x / scale, giveLog);
return giveLog ? pr + Math.log(shape / x) : pr * shape / x;
}
/* else shape >= 1 */
- pr = dpoisRaw(shape - 1, x / scale, giveLog);
+ pr = DPois.dpoisRaw(shape - 1, x / scale, giveLog);
return giveLog ? pr - Math.log(scale) : pr / scale;
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pbeta.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pbeta.java
index 4c12e76f64cc9dc11324231c7a825de34b62abb3..c195a4bbc295502935d69141b53f601d80fc244d 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pbeta.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pbeta.java
@@ -4,7 +4,7 @@
* http://www.gnu.org/licenses/gpl-2.0.html
*
* Copyright (c) 2006, The R Core Team
- * Copyright (c) 2016, 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
@@ -14,7 +14,6 @@ import com.oracle.truffle.api.CompilerDirectives.TruffleBoundary;
import com.oracle.truffle.api.profiles.BranchProfile;
import com.oracle.truffle.r.library.stats.StatsFunctions.Function3_2;
import com.oracle.truffle.r.library.stats.TOMS708.Bratio;
-import com.oracle.truffle.r.runtime.RError;
import com.oracle.truffle.r.runtime.RError.Message;
// transcribed from pbeta.c
@@ -68,7 +67,7 @@ public final class Pbeta implements Function3_2 {
@TruffleBoundary
private static void doWarning(double x, double a, double b, int ierr) {
- RError.warning(RError.SHOW_CALLER, Message.GENERIC, String.format("pbeta_raw(%g, a=%g, b=%g, ..) -> bratio() gave error code %d", x, a, b, ierr));
+ RMathError.warning(Message.GENERIC, String.format("pbeta_raw(%g, a=%g, b=%g, ..) -> bratio() gave error code %d", x, a, b, ierr));
}
static double pbeta(double x, double a, double b, boolean lowerTail, boolean logP, BranchProfile nanProfile) {
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pnorm.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pnorm.java
index 2d73851c472337e07e46861bed2700d4ccb9a2d2..bc15917ca04e4378a1e1b51191a6ecdc0453246b 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pnorm.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Pnorm.java
@@ -6,7 +6,7 @@
* Copyright (C) 1998 Ross Ihaka
* Copyright (c) 2000--2013, The R Core Team
* Copyright (c) 2003, The R Foundation
- * Copyright (c) 2016, 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
@@ -223,7 +223,6 @@ public final class Pnorm implements StatsFunctions.Function3_2 {
// if(*ccum < min) *ccum = 0.;
// }
// #endif
- return;
}
}
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QBeta.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QBeta.java
index 8e8c7512edb90b564a5677756a68258f040b6c3a..ce7e709e319af2771b6851be27feaa387ad77a3a 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QBeta.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QBeta.java
@@ -5,7 +5,7 @@
*
* Copyright (c) 1995, 1996, Robert Gentleman and Ross Ihaka
* Copyright (c) 1998-2015, The R Core Team
- * Copyright (c) 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
@@ -117,11 +117,11 @@ public final class QBeta implements Function3_2 {
/**
*
- * @param alpha
- * @param p
- * @param q
- * @param lower_tail
- * @param log_p
+ * @param alpha parameter of beta distribution
+ * @param p probability
+ * @param q quantile
+ * @param lower_tail whether to use lower tail of the distribution function
+ * @param log_p is the probability given as log(p)
* @param swap_01 {true, NA, false}: if NA, algorithm decides swap_tail
* @param log_q_cut if == Inf: return Math.log(qbeta(..)); otherwise, if finite: the bound
* for switching to Math.log(x)-scale; see use_log_x
@@ -490,7 +490,7 @@ public final class QBeta implements Function3_2 {
// (cancellation in (u_n -u) => may differ from adj:
D = RMath.fmin2(Math.abs(adj), Math.abs(u_n - u));
/* R_ifDEBUG_printf(" delta(u)=%g\n", u_n - u); */
- debugPrintf(" it{in}=%d, delta(u)=%9.3g, D/|.|=%.3g\n",
+ debugPrintf(" it{in}=%d, delta(u)=%.3g, D/|.|=%.3g\n",
i_inn, u_n - u, D / Math.abs(u_n + u));
if (D <= 4e-16 * Math.abs(u_n + u)) {
converged(log_p, qb);
@@ -524,7 +524,7 @@ public final class QBeta implements Function3_2 {
: (y - a) * Math.exp(logbeta + r * Math.log(xinbta) + t * Math.log1p(-xinbta));
if (i_pb >= n_N && w * wprev <= 0.)
prev = RMath.fmax2(Math.abs(adj), fpu);
- debugPrintf("N(i=%2d): x0=%#17.15g, pb(x0)=%#17.15g, w=%#17.15g, %s prev=%g,",
+ debugPrintf("N(i=%d): x0=%.15g, pb(x0)=%.15g, w=%.15g, %s prev=%g,",
i_pb, xinbta, y, w, (w * wprev <= 0.) ? "new" : "old", prev);
g = 1;
int i_inn;
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QuantileSearch.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QuantileSearch.java
index 8ccd4a00705284edf04e43af5b520616cfa275b1..05ea19680faafbbdfca665731df700181267d4c7 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QuantileSearch.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/QuantileSearch.java
@@ -6,7 +6,7 @@
* Copyright (C) 1998 Ross Ihaka
* Copyright (c) 2000-2016, The R Core Team
* Copyright (c) 2003-2016, The R Foundation
- * Copyright (c) 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
@@ -15,8 +15,8 @@ package com.oracle.truffle.r.library.stats;
/**
* Searches for a quantile of given random variable using it's distribution function. The search
* takes steps of given size until it reaches the quantile or until it steps over it. This class and
- * its {@code {@link #simpleSearch(double, double, double)}} method correspond to several
- * {@code do_search} functions in GnuR.
+ * its {@link #simpleSearch(double, double, double)} method correspond to several {@code do_search}
+ * functions in GnuR.
*/
public final class QuantileSearch {
/**
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RBeta.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RBeta.java
index cc59560caa6b01c6559d67e1f37e6559e0e4e3ea..0aceb6a7997c63167a5c806f15d410d56a2c7b72 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RBeta.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RBeta.java
@@ -6,7 +6,7 @@
* Copyright (c) 1995, 1996, 1997 Robert Gentleman and Ross Ihaka
* Copyright (c) 1998-2013, The R Core Team
* Copyright (c) 2003-2015, The R Foundation
- * Copyright (c) 2016, 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
@@ -27,7 +27,6 @@ public final class RBeta extends RandFunction2_Double {
// TODO: state variables
private static double beta = 0;
private static double gamma = 1;
- private static double delta;
private static double k1 = 0;
private static double k2 = 0;
private static double olda = -1.0;
@@ -70,7 +69,7 @@ public final class RBeta extends RandFunction2_Double {
/* changed notation, now also a <= b (was reversed) */
if (!qsame) { /* initialize */
beta = 1.0 / a;
- delta = 1.0 + b - a;
+ double delta = 1.0 + b - a;
k1 = delta * (0.0138889 + 0.0416667 * a) / (b * beta - 0.777778);
k2 = 0.25 + (0.5 + 0.25 / delta) * a;
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RMath.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RMath.java
index c3f9b54793a7e349070491f7060797f3e76fe35c..57e33d4af095f5d632f5cddd1d0fbafc668dbe76 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RMath.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RMath.java
@@ -6,16 +6,17 @@
* Copyright (C) 1998 Ross Ihaka
* Copyright (c) 1998--2012, The R Core Team
* Copyright (c) 2004, The R Foundation
- * Copyright (c) 2013, 2016, Oracle and/or its affiliates
+ * Copyright (c) 2013, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
package com.oracle.truffle.r.library.stats;
import static com.oracle.truffle.r.library.stats.LBeta.lbeta;
+import static com.oracle.truffle.r.library.stats.MathConstants.M_LN_SQRT_2PI;
import com.oracle.truffle.api.CompilerDirectives.CompilationFinal;
-import com.oracle.truffle.api.CompilerDirectives.TruffleBoundary;
+import com.oracle.truffle.r.library.stats.RMathError.MLError;
import com.oracle.truffle.r.runtime.RRuntime;
/**
@@ -193,7 +194,7 @@ public final class RMath {
/* else */
if (x < xmin) {
/* answer less than half precision because x too near -1 */
- fail("ERROR: ML_ERROR(ME_PRECISION, \"log1p\")");
+ RMathError.error(MLError.PRECISION, "log1p");
}
return Math.log(1 + x);
}
@@ -210,11 +211,11 @@ public final class RMath {
int i;
if (n < 1 || n > 1000) {
- return Double.NaN; // ML_ERR_return_NAN;
+ return RMathError.defaultError();
}
if (x < -1.1 || x > 1.1) {
- return Double.NaN; // ML_ERR_return_NAN;
+ return RMathError.defaultError();
}
twox = x * 2;
@@ -228,8 +229,111 @@ public final class RMath {
return (b0 - b2) * 0.5;
}
- @TruffleBoundary
- private static void fail(String message) {
- throw new RuntimeException(message);
+ //
+ // stirlerr
+ //
+
+ private static final double S0 = 0.083333333333333333333; /* 1/12 */
+ private static final double S1 = 0.00277777777777777777778; /* 1/360 */
+ private static final double S2 = 0.00079365079365079365079365; /* 1/1260 */
+ private static final double S3 = 0.000595238095238095238095238; /* 1/1680 */
+ private static final double S4 = 0.0008417508417508417508417508; /* 1/1188 */
+
+ /*
+ * error for 0, 0.5, 1.0, 1.5, ..., 14.5, 15.0.
+ */
+ @CompilationFinal private static final double[] sferr_halves = new double[]{
+ 0.0, /* n=0 - wrong, place holder only */
+ 0.1534264097200273452913848, /* 0.5 */
+ 0.0810614667953272582196702, /* 1.0 */
+ 0.0548141210519176538961390, /* 1.5 */
+ 0.0413406959554092940938221, /* 2.0 */
+ 0.03316287351993628748511048, /* 2.5 */
+ 0.02767792568499833914878929, /* 3.0 */
+ 0.02374616365629749597132920, /* 3.5 */
+ 0.02079067210376509311152277, /* 4.0 */
+ 0.01848845053267318523077934, /* 4.5 */
+ 0.01664469118982119216319487, /* 5.0 */
+ 0.01513497322191737887351255, /* 5.5 */
+ 0.01387612882307074799874573, /* 6.0 */
+ 0.01281046524292022692424986, /* 6.5 */
+ 0.01189670994589177009505572, /* 7.0 */
+ 0.01110455975820691732662991, /* 7.5 */
+ 0.010411265261972096497478567, /* 8.0 */
+ 0.009799416126158803298389475, /* 8.5 */
+ 0.009255462182712732917728637, /* 9.0 */
+ 0.008768700134139385462952823, /* 9.5 */
+ 0.008330563433362871256469318, /* 10.0 */
+ 0.007934114564314020547248100, /* 10.5 */
+ 0.007573675487951840794972024, /* 11.0 */
+ 0.007244554301320383179543912, /* 11.5 */
+ 0.006942840107209529865664152, /* 12.0 */
+ 0.006665247032707682442354394, /* 12.5 */
+ 0.006408994188004207068439631, /* 13.0 */
+ 0.006171712263039457647532867, /* 13.5 */
+ 0.005951370112758847735624416, /* 14.0 */
+ 0.005746216513010115682023589, /* 14.5 */
+ 0.005554733551962801371038690 /* 15.0 */
+ };
+
+ static double stirlerr(double n) {
+ double nn;
+
+ if (n <= 15.0) {
+ nn = n + n;
+ if (nn == (int) nn) {
+ return (sferr_halves[(int) nn]);
+ }
+ return (GammaFunctions.lgammafn(n + 1.) - (n + 0.5) * Math.log(n) + n - M_LN_SQRT_2PI);
+ }
+
+ nn = n * n;
+ if (n > 500) {
+ return ((S0 - S1 / nn) / n);
+ }
+ if (n > 80) {
+ return ((S0 - (S1 - S2 / nn) / nn) / n);
+ }
+ if (n > 35) {
+ return ((S0 - (S1 - (S2 - S3 / nn) / nn) / nn) / n);
+ }
+ /* 15 < n <= 35 : */
+ return ((S0 - (S1 - (S2 - (S3 - S4 / nn) / nn) / nn) / nn) / n);
+ }
+
+ //
+ // bd0
+ //
+
+ /**
+ * Evaluates the "deviance part". Transcribed from GnuR.
+ */
+ static double bd0(double x, double np) {
+ double ej;
+ double s;
+ double s1;
+ double v;
+ int j;
+
+ if (!RRuntime.isFinite(x) || !RRuntime.isFinite(np) || np == 0.0) {
+ return RMathError.defaultError();
+ }
+
+ if (Math.abs(x - np) < 0.1 * (x + np)) {
+ v = (x - np) / (x + np);
+ s = (x - np) * v; /* s using v -- change by MM */
+ ej = 2 * x * v;
+ v = v * v;
+ for (j = 1;; j++) { /* Taylor series */
+ ej *= v;
+ s1 = s + ej / ((j << 1) + 1);
+ if (s1 == s) { /* last term was effectively 0 */
+ return s1;
+ }
+ s = s1;
+ }
+ }
+ /* else: | x - np | is not too small */
+ return x * Math.log(x / np) + np - x;
}
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RMathError.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RMathError.java
index 31ab7145dccb5391f74d7fc05f509526ba73d2e7..0644a4c0123d9a6fe33918f804fa0a12771634cd 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RMathError.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RMathError.java
@@ -4,12 +4,13 @@
* http://www.gnu.org/licenses/gpl-2.0.html
*
* Copyright (c) 1998-2016, The R Core Team
- * Copyright (c) 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
package com.oracle.truffle.r.library.stats;
+import com.oracle.truffle.api.CompilerDirectives.TruffleBoundary;
import com.oracle.truffle.r.runtime.RError;
import com.oracle.truffle.r.runtime.RError.Message;
@@ -36,6 +37,7 @@ public final class RMathError {
this.message = message;
}
+ @TruffleBoundary
public void warning(String arg) {
RError.warning(RError.SHOW_CALLER, message, arg);
}
@@ -64,6 +66,7 @@ public final class RMathError {
/**
* Corresponds to macros {@code MATHLIB_WARNINGX} in GnuR.
*/
+ @TruffleBoundary
public static void warning(RError.Message message, Object... args) {
RError.warning(RError.SHOW_CALLER, message, args);
}
diff --git a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/TOMS708.java b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/TOMS708.java
index 9f0ee17d7b44b2a40a399ec10a63b0beccd7c870..a569fe86ca2604564b6d98764950dbca7fba85b8 100644
--- a/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/TOMS708.java
+++ b/com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/TOMS708.java
@@ -5,7 +5,7 @@
*
* Copyright (c) 1995-2012, The R Core Team
* Copyright (c) 2003, The R Foundation
- * Copyright (c) 2016, 2016, Oracle and/or its affiliates
+ * Copyright (c) 2016, 2017, Oracle and/or its affiliates
*
* All rights reserved.
*/
@@ -20,6 +20,7 @@ import com.oracle.truffle.r.runtime.RError;
import com.oracle.truffle.r.runtime.RError.Message;
import com.oracle.truffle.r.runtime.RInternalError;
import com.oracle.truffle.r.runtime.RRuntime;
+import com.oracle.truffle.r.runtime.Utils;
/*
* transcribed from toms708.c - as the original file contains no copyright header, we assume that it is copyright R code and R foundation.
@@ -555,7 +556,7 @@ public final class TOMS708 {
if (b * z == 0.0) { // should not happen, but does, e.g.,
// for pbeta(1e-320, 1e-5, 0.5) i.e., _subnormal_ x,
// Warning ... bgrat(a=20.5, b=1e-05, x=1, y=9.99989e-321): ..
- RError.warning(RError.SHOW_CALLER, Message.GENERIC, String.format(
+ RMathError.warning(Message.GENERIC, Utils.stringFormat(
"bgrat(a=%f, b=%f, x=%f, y=%f): b*z == 0 underflow, hence inaccurate pbeta()",
a, b, x, y));
/* L_Error: THE EXPANSION CANNOT BE COMPUTED */
@@ -627,7 +628,7 @@ public final class TOMS708 {
break;
} else if (n == n_terms_bgrat) {
// never? ; please notify R-core if seen:
- RError.warning(RError.SHOW_CALLER, Message.GENERIC, String.format("bgrat(a=%f, b=%f, x=%f,..): did *not* converge; dj=%f, rel.err=%f\n", a, b, x, dj, Math.abs(dj) / (sum + l)));
+ RMathError.warning(Message.GENERIC, Utils.stringFormat("bgrat(a=%f, b=%f, x=%f,..): did *not* converge; dj=%f, rel.err=%f\n", a, b, x, dj, Math.abs(dj) / (sum + l)));
}
}
diff --git a/mx.fastr/copyrights/overrides b/mx.fastr/copyrights/overrides
index 602b12154d1ac2624bea8b5ecd71c5a008d05cad..53a014c74828085892acc561a70efeafe9399759 100644
--- a/mx.fastr/copyrights/overrides
+++ b/mx.fastr/copyrights/overrides
@@ -228,7 +228,6 @@ com.oracle.truffle.r.nodes/src/com/oracle/truffle/r/nodes/access/UpdateSlotNode.
com.oracle.truffle.r.nodes/src/com/oracle/truffle/r/nodes/binary/CastTypeNode.java,purdue.copyright
com.oracle.truffle.r.nodes/src/com/oracle/truffle/r/nodes/function/CallMatcherNode.java,purdue.copyright
com.oracle.truffle.r.nodes/src/com/oracle/truffle/r/nodes/function/S3FunctionLookupNode.java,purdue.copyright
-com.oracle.truffle.r.nodes/src/com/oracle/truffle/r/nodes/function/UseMethodInternalNode.java,purdue.copyright
com.oracle.truffle.r.nodes/src/com/oracle/truffle/r/nodes/objects/AsS4.java,gnu_r_gentleman_ihaka.copyright
com.oracle.truffle.r.nodes/src/com/oracle/truffle/r/nodes/objects/DispatchGeneric.java,gnu_r_gentleman_ihaka.copyright
com.oracle.truffle.r.nodes/src/com/oracle/truffle/r/nodes/objects/ExecuteMethod.java,gnu_r_gentleman_ihaka.copyright