/*
* Open Source Physics software is free software as described near the bottom of this code file.
*
* For additional information and documentation on Open Source Physics please see:
* <http://www.opensourcephysics.org/>
*/
package org.opensourcephysics.numerics.specialfunctions;
/**
* Computes the factorial of an integer and functions closely related to factorials.
*/
public class Factorials {
static double[] cof = {76.18009172947146, -86.50532032941677, 24.01409824083091, -1.231739572450155, 0.1208650973866179e-2, -0.5395239384953e-5};
static long[] fac;
private Factorials() { // prohibit instantiation because all methods are static
}
/**
* Calculates the logarithm of the Gamma function using the Lanczos approximation.
* @param x double
* @return double
*/
public static double gammaln(final double x) {
double y = x, tmp = x+5.5;
tmp -= (x+0.5)*Math.log(tmp);
double sum = 1.000000000190015;
for(int j = 0; j<=5; j++) {
sum += cof[j]/++y;
}
return -tmp+Math.log(2.5066282746310005*sum/x);
}
/**
* Calculates the factorial.
* @param n int
* @return double
*/
public static double factorial(final int n) {
if(n<0) {
throw new IllegalArgumentException(Messages.getString("Factorials.neg_val")); //$NON-NLS-1$
}
if(n<fac.length) {
return fac[n];
}
return Math.exp(gammaln(n+1.0));
}
static {
long val = 1;
int maxN = 1;
while(val*maxN>=val) {
val = val*maxN;
maxN++;
}
fac = new long[maxN]; // maxN should be ~21
fac[0] = 1;
for(int i = 1; i<maxN; i++) {
fac[i] = fac[i-1]*i;
}
}
/**
* Returns log (n!) = log (n * (n-1) * ... 2 * 1)
*
* @param n
* @return log(n!)
*/
public static double logFactorial(int n) {
if(n<0) {
throw new IllegalArgumentException(Messages.getString("Factorials.log_neg_val")); //$NON-NLS-1$
}
return Factorials.gammaln(n+1.0);
}
/**
* Returns the Poisson distribution (nu^n e^(-nu) / n!)
*
* @param nu
* @param n
* @return poisson_nu(n)
*/
public static double poisson(double nu, int n) {
return Math.exp(n*Math.log(nu)-nu-logFactorial(n));
}
/**
* Returns the logarithm of the binomial coefficient (n, k)
* In other notation: log (n choose k)
* (n choose k) represents the number of ways of picking k unordered outcomes from n possibilities
*
* @param n
* @param k
* @return log (n choose k)
*/
public static double logChoose(int n, int k) {
return logFactorial(n)-logFactorial(k)-logFactorial(n-k);
}
/*
* --------------- test code ---------------
* public static void main(String args[]) {
* System.out.println(Factorials.factorial(3));
* System.out.println(Factorials.factorial(20));
* System.out.println(Factorials.factorial(21));
* System.out.println(Factorials.factorial(22));
* }
*/
}
/*
* Open Source Physics software is free software; you can redistribute
* it and/or modify it under the terms of the GNU General Public License (GPL) as
* published by the Free Software Foundation; either version 2 of the License,
* or(at your option) any later version.
* Code that uses any portion of the code in the org.opensourcephysics package
* or any subpackage (subdirectory) of this package must must also be be released
* under the GNU GPL license.
*
* This software is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this; if not, write to the Free Software
* Foundation, Inc., 59 Temple Place, Suite 330, Boston MA 02111-1307 USA
* or view the license online at http://www.gnu.org/copyleft/gpl.html
*
* Copyright (c) 2007 The Open Source Physics project
* http://www.opensourcephysics.org
*/