1   /*
2    * Licensed to the Apache Software Foundation (ASF) under one or more
3    * contributor license agreements.  See the NOTICE file distributed with
4    * this work for additional information regarding copyright ownership.
5    * The ASF licenses this file to You under the Apache License, Version 2.0
6    * (the "License"); you may not use this file except in compliance with
7    * the License.  You may obtain a copy of the License at
8    * 
9    *      http://www.apache.org/licenses/LICENSE-2.0
10   * 
11   * Unless required by applicable law or agreed to in writing, software
12   * distributed under the License is distributed on an "AS IS" BASIS,
13   * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
14   * See the License for the specific language governing permissions and
15   * limitations under the License.
16   */
17  package org.apache.commons.math.special;
18  
19  import org.apache.commons.math.MathException;
20  import org.apache.commons.math.TestUtils;
21  
22  import junit.framework.TestCase;
23  
24  /**
25   * @version $Revision: 778416 $ $Date: 2009-05-25 09:20:07 -0400 (Mon, 25 May 2009) $
26   */
27  public class GammaTest extends TestCase {
28       
29      public GammaTest(String name) {
30          super(name);
31      }
32  
33      private void testRegularizedGamma(double expected, double a, double x) {
34          try {
35              double actualP = Gamma.regularizedGammaP(a, x);
36              double actualQ = Gamma.regularizedGammaQ(a, x);
37              TestUtils.assertEquals(expected, actualP, 10e-15);
38              TestUtils.assertEquals(actualP, 1.0 - actualQ, 10e-15);
39          } catch(MathException ex){
40              fail(ex.getMessage());
41          }
42      }
43  
44      private void testLogGamma(double expected, double x) {
45          double actual = Gamma.logGamma(x);
46          TestUtils.assertEquals(expected, actual, 10e-15);
47      }
48  
49      public void testRegularizedGammaNanPositive() {
50          testRegularizedGamma(Double.NaN, Double.NaN, 1.0);
51      }
52  
53      public void testRegularizedGammaPositiveNan() {
54          testRegularizedGamma(Double.NaN, 1.0, Double.NaN);
55      }
56      
57      public void testRegularizedGammaNegativePositive() {
58          testRegularizedGamma(Double.NaN, -1.5, 1.0);
59      }
60      
61      public void testRegularizedGammaPositiveNegative() {
62          testRegularizedGamma(Double.NaN, 1.0, -1.0);
63      }
64      
65      public void testRegularizedGammaZeroPositive() {
66          testRegularizedGamma(Double.NaN, 0.0, 1.0);
67      }
68      
69      public void testRegularizedGammaPositiveZero() {
70          testRegularizedGamma(0.0, 1.0, 0.0);
71      }
72      
73      public void testRegularizedGammaPositivePositive() {
74          testRegularizedGamma(0.632120558828558, 1.0, 1.0);
75      }
76      
77      public void testLogGammaNan() {
78          testLogGamma(Double.NaN, Double.NaN);
79      }
80      
81      public void testLogGammaNegative() {
82          testLogGamma(Double.NaN, -1.0);
83      }
84      
85      public void testLogGammaZero() {
86          testLogGamma(Double.NaN, 0.0);
87      }
88      
89      public void testLogGammaPositive() {
90          testLogGamma(0.6931471805599457, 3.0);
91      }
92  
93      public void testDigammaLargeArgs() {
94          double eps = 1e-8;
95          assertEquals(4.6001618527380874002, Gamma.digamma(100), eps);
96          assertEquals(3.9019896734278921970, Gamma.digamma(50), eps);
97          assertEquals(2.9705239922421490509, Gamma.digamma(20), eps);
98          assertEquals(2.9958363947076465821, Gamma.digamma(20.5), eps);
99          assertEquals(2.2622143570941481605, Gamma.digamma(10.1), eps);
100         assertEquals(2.1168588189004379233, Gamma.digamma(8.8), eps);
101         assertEquals(1.8727843350984671394, Gamma.digamma(7), eps);
102         assertEquals(0.42278433509846713939, Gamma.digamma(2), eps);
103         assertEquals(-100.56088545786867450, Gamma.digamma(0.01), eps);
104         assertEquals(-4.0390398965921882955, Gamma.digamma(-0.8), eps);
105         assertEquals(4.2003210041401844726, Gamma.digamma(-6.3), eps);
106     }
107 
108     public void testDigammaSmallArgs() {
109         // values for negative powers of 10 from 1 to 30 as computed by webMathematica with 20 digits
110         // see functions.wolfram.com
111         double[] expected = {-10.423754940411076795, -100.56088545786867450, -1000.5755719318103005,
112                 -10000.577051183514335, -100000.57719921568107, -1.0000005772140199687e6, -1.0000000577215500408e7,
113                 -1.0000000057721564845e8, -1.0000000005772156633e9, -1.0000000000577215665e10, -1.0000000000057721566e11,
114                 -1.0000000000005772157e12, -1.0000000000000577216e13, -1.0000000000000057722e14, -1.0000000000000005772e15, -1e+16,
115                 -1e+17, -1e+18, -1e+19, -1e+20, -1e+21, -1e+22, -1e+23, -1e+24, -1e+25, -1e+26,
116                 -1e+27, -1e+28, -1e+29, -1e+30};
117         for (double n = 1; n < 30; n++) {
118             checkRelativeError(String.format("Test %.0f: ", n), expected[(int) (n - 1)], Gamma.digamma(Math.pow(10.0, -n)), 1e-8);
119         }
120     }
121 
122     public void testTrigamma() {
123         double eps = 1e-8;
124         // computed using webMathematica.  For example, to compute trigamma($i) = Polygamma(1, $i), use
125         //
126         // http://functions.wolfram.com/webMathematica/Evaluated.jsp?name=PolyGamma2&plottype=0&vars={%221%22,%22$i%22}&digits=20
127         double[] data = {
128                 1e-4, 1.0000000164469368793e8,
129                 1e-3, 1.0000016425331958690e6,
130                 1e-2, 10001.621213528313220,
131                 1e-1, 101.43329915079275882,
132                 1, 1.6449340668482264365,
133                 2, 0.64493406684822643647,
134                 3, 0.39493406684822643647,
135                 4, 0.28382295573711532536,
136                 5, 0.22132295573711532536,
137                 10, 0.10516633568168574612,
138                 20, 0.051270822935203119832,
139                 50, 0.020201333226697125806,
140                 100, 0.010050166663333571395
141         };
142         for (int i = data.length - 2; i >= 0; i -= 2) {
143             assertEquals(String.format("trigamma %.0f", data[i]), data[i + 1], Gamma.trigamma(data[i]), eps);
144         }
145     }
146 
147     private void checkRelativeError(String msg, double expected, double actual, double tolerance) {
148         assertEquals(msg, expected, actual, Math.abs(tolerance * actual));
149     }
150 }