Gamma { public static double trigamma(double x) { if (x > 0 && x <= S_LIMIT) { return 1 / (x * x); } if (x >= C_LIMIT) { double inv = 1 / (x * x); return 1 / x + inv / 2 + inv / x * (1.0 / 6 - inv * (1.0 / 30 + inv / 42)); } return trigamma(x + 1) + 1 / (x * x); } private  Gamma(); static double logGamma(double x); static double regularizedGammaP(double a, double x); static double regularizedGammaP(double a,
                                           double x,
                                           double epsilon,
                                           int maxIterations); static double regularizedGammaQ(double a, double x); static double regularizedGammaQ(final double a,
                                           double x,
                                           double epsilon,
                                           int maxIterations); static double digamma(double x); static double trigamma(double x); static double lanczos(final double x); static double invGamma1pm1(final double x); static double logGamma1p(final double x); static double gamma(final double x); static final double GAMMA; static final double LANCZOS_G; }
@Test public void testTrigamma() { double eps = 1e-8; double[] data = { 1e-4, 1.0000000164469368793e8, 1e-3, 1.0000016425331958690e6, 1e-2, 10001.621213528313220, 1e-1, 101.43329915079275882, 1, 1.6449340668482264365, 2, 0.64493406684822643647, 3, 0.39493406684822643647, 4, 0.28382295573711532536, 5, 0.22132295573711532536, 10, 0.10516633568168574612, 20, 0.051270822935203119832, 50, 0.020201333226697125806, 100, 0.010050166663333571395 }; for (int i = data.length - 2; i >= 0; i -= 2) { Assert.assertEquals(String.format("trigamma %.0f", data[i]), data[i + 1], Gamma.trigamma(data[i]), eps); } }