Take the 2-minute tour ×
Mathematics Stack Exchange is a question and answer site for people studying math at any level and professionals in related fields. It's 100% free, no registration required.

I'm facing a practical problem where I've calculated a formula that, with the help of some programming, can bring me to my final answer. However, the numbers involved are so big that it takes ages to compute, and I think it might not be necessary. I have the following formula,

$\sum\limits^{k}_{i=m}(N-i)^{k-i}(\frac{1}{N})^k\frac{k!}{(k-i)!i!} \leq a$

from which I need to calculate N. The rest of the values are constants, but in the range of the 100,000s. The factorial there is giving me a headache, since the values involved are too large; what simplifications could I make that will loosen the bounds slightly and thereby simplify the calculation? Are there any standard tricks? Or perhaps a way to calculate this in matlab / octave?

share|improve this question

4 Answers 4

up vote 8 down vote accepted

You need Stirling's approximation. It is very accurate for large factorials.

share|improve this answer

You don't need to compute the individual factorials in order to compute $k!/(k-i)!i!$, since that's the binomial coefficient $\binom{k}{i}$. A simple algorithm for computing binomial coefficients can be found on Wikipedia. A more sophisticated algorithm is due to Goetgheluck (JSTOR); implementations can be found here and here.

Of course, with numbers of the size that you have, this might still not be feasible, and in this case I also recommend Stirling's formula.

share|improve this answer

Here is a good writeup about implementing Stirling's approximation, along with a reference implementation:

http://threebrothers.org/brendan/blog/stirlings-approximation-formula-clojure/

share|improve this answer
 import java.util.HashMap;

/**
 * Convenience methods for deterministic combinations.
 * 
 * Usage:
 * <pre> 
 * Factorial f = new Factorial();
 * long ten_over_two = f.over(10,2);     //takes a while, because no cache yet
 * System.out.println( f.factorial( 8)); //instantaneous
 * </pre>
 *  
 * Notes:
 * <li>instantiate and use, it will reuse its own cache, so after a while it is supposed to run in O(1).
 *     Code is super fast (0.0ms)
 * <li>For very big numbers, use Stirling's approximation, @see http://en.wikipedia.org/wiki/Stirling%27s_approximation
 * <li>Overflow is a big problem. 
 *     the max value is 9223372036854775807 < 21!, 
 *     1.7976931348623157E308 <171!
 *     Either it is super fast or it breaks, the combinations grow too fast, and for 2k~n it is biggest; 
 *     it cannot even handle (29 14) but for k/n small it can function barely.
 * 
 */
public class Factorial {

    static private HashMap<Integer, Double> cache = new HashMap<Integer, Double> (); 

    public Double factorial(int n)
    {
    Integer N = Integer.valueOf( n );
    if(cache.containsKey( N) ) return cache.get( N );

    if(n<1) {
        cache.put( N, 0D );
        return 0D;
    }
    if(n==1) {
        cache.put( N, 1D );
        return 1D;
    }

    Double temp = 1D;
    for (int i = 1; i < n; i++) {
        temp*=i;
        cache.put( i, temp );
    }

    Double prev = factorial(n-1);
    if(prev>Double.MAX_VALUE/n) throw new IllegalArgumentException("factorial("+n+") overflow");

    return n*prev;
    }

    // 52!/47! = 52*51*50*49*48 so extra(52,47) = 52*(51,47) ... and (47,47)=1
    private Double extra(int n, int minus)
    {
        if(n<=minus) return 1D;
        Double result = extra(n-1,minus);
        if(result>Double.MAX_VALUE/n) throw new IllegalArgumentException(  "overflow");
        return n* result;
    }

    /**
     * 
     * This method computes (n k) = n!/(k! * (n-k)!) 
     * <p/>
     * Example usage:
     * <li>over(52,2) = 2598960
     * <li>over(5,2) = 10
     * <li>over(2,2) = 1 //ignoring swaps
     *
     * @param n the big set
     * @param k how much to draw out of n
     * @return all combinations of drawing k out of n, ignoring swaps. Simply double the amount for regarding swaps.
     *
     */
    public Double over(int n, int k){
        if(k==n) return 1D;
        if(k>n) return 0D;
        if(k<=0) return 0D;
        //k>=1 and n>k
        if(k<n-k) return over(n,n-k); 
        Double over = extra(n,k);
        Double divisor = factorial(k);      
        return over/divisor;            
    }

    /**
     * basic test and usage
     */
    public static void main(String[] args) 
    {
        Double fuc= new Factorial().factorial( 170);     //takes a while, because no cache yet
        System.out.println(fuc);


        {
            int n = 100;
            int k = 50;
            double start = System.currentTimeMillis(); 
            //Factorial f = new Factorial();
            Double ten_over_two = new Factorial().over(n,k);     //takes a while, because no cache yet
            double delta = System.currentTimeMillis()-start;
            System.out.println( "ten_over_two=" +ten_over_two + ": in " + delta + "ms");
        }

      {
            int n = 8;
            double start = System.currentTimeMillis(); 
        Factorial f = new Factorial();
        System.out.println( f.factorial( n)); //instantaneous
            double delta = System.currentTimeMillis()-start;

            System.out.println( "factorial("+n+") in " + delta + "ms");
      }

    }

}//~class}
share|improve this answer
    
Why would that be faster ? –  The Game Oct 9 at 13:40
    
because it reuses already computed factors. –  user182100 Oct 9 at 13:52

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Not the answer you're looking for? Browse other questions tagged or ask your own question.