1

I have the following function:

private static BigInteger Factorial(int number)
{
     if(number < 2) return BigInteger.One;
     double sum = 0;
     for(int i = 2; i < number; i++)
         sum += Math.Log(i);
     return new BigInteger(Math.Exp(sum));
}

This works fine for small input numbers like 10, however if I pass a bigger number like 50000 it crashes and throws an OverflowException.

I understand why this is happening, the result of Math.Exp(sum) is just too big for a double. But that's why I'm trying to use BigInteger after all, to avoid these type of exceptions.

The problem is that wrapping the result like new BigInteger(Math.Exp(sum)) is useless because Map.Exp(sum) is trying to return a double anyway.

So I decided to use BigInteger.Pow static function:

return BigInteger.Pow(new BigInteger(Math.E), number);

Note the new BigInteger(Math.E) part. BigInteger.Pow takes a BigInteger as the first parameter so I don't have a choice but to wrap my Math.E to a BigInteger.

However, by doing this I am actually truncating the decimal part of my Math.E which ruins my algorithm because I end up with a totally different result.

I was looking for something like BigDouble or something similar, but the only class that I seem to find is BigInteger.

How can I successfully return the correct result as a BigInteger on this function and at the same time have an exception-safe code when big input is received?

Matias Cicero
  • 25,439
  • 13
  • 82
  • 154
  • Do you need an accurate or an adducent value? – Tamas Hegedus Aug 28 '15 at 23:59
  • You will have to find a representation which can hold numbers like `50000!`. Nor `double` nor `BigInteger` is suitable. `BigInteger` can hold the result, but takes pretty much time (several secs) to calculate. – Tamas Hegedus Aug 29 '15 at 00:01
  • You could use the answer [here](http://stackoverflow.com/a/4524254/916299) and add to it. Since they have already implemented multiplication, you can just multiply repeatedly to get powers – James Webster Aug 29 '15 at 00:01
  • @hege_hegedus I really need the accurate value of `Math.E` for this algorithm to work properly. The result of `Math.Pow` or `Math.Exp` can be truncated though. – Matias Cicero Aug 29 '15 at 00:02
  • @JamesWebster I was really trying to avoid multiplication, but I'll take a look at it! – Matias Cicero Aug 29 '15 at 00:03
  • Powers are just multiplication.. Admittedly the compiler can usually make it a bit more efficient by reusing results – James Webster Aug 29 '15 at 00:04
  • @hege_hegedus then how wolframalpha calculates something like [this](http://m.wolframalpha.com/input/?i=5000000000000000%21&x=0&y=0)? – M.kazem Akhgary Aug 29 '15 at 00:42
  • @hege_hegedus if one needs just occasionally compute 50000! several seconds may be ok (like `Enumerable.Range(1, 50000).Aggregate(new BigInteger(1), (i,f)=> i*f)` ). Not clear what OP is actually looking for so... – Alexei Levenkov Aug 29 '15 at 00:48
  • e is not a finite precision value, so you can't get a correct value this way. use otherways like table lookup and/or dynamic programming – phuclv Aug 29 '15 at 02:45

1 Answers1

1

You could use base 2 instead of base e:

private static BigInteger FactorialEstimate(int number)
{
    if (number < 2) return BigInteger.One;
    double sum = 0;
    for (int i = 2; i <= number; i++)
        sum += Math.Log(i, 2);
    return BigInteger.Pow(2, (int)Math.Round(sum));
}

If you need a more accurate answer:

private static BigInteger NthRoot(BigInteger a, int n)
{
    BigInteger min = BigInteger.Zero, max = a, mid = a;
    while (min < max)
    {
        mid = (min + max) >> 1;
        if (BigInteger.Pow(mid, n) < a)
            min = mid + 1;
        else
            max = mid;
    }
    return max;
}

const int Accuracy1 = 16;
const int Accuracy2 = 8;
private static BigInteger FactorialBetterEstimate(int number)
{
    if (number < 2) return BigInteger.One;
    double sumFractPart = 0;
    int sumIntPart = 0, tmpIntPart = 0;
    for (int i = 2; i <= number; i++)
    {
        sumFractPart += Math.Log(i, 2);
        tmpIntPart = (int)sumFractPart;
        sumFractPart -= tmpIntPart;
        sumIntPart += tmpIntPart;
    }

    int correction = Math.Min(Accuracy2, sumIntPart);
    sumIntPart -= correction;
    sumFractPart += correction;

    return NthRoot(BigInteger.Pow(2, Convert.ToInt32(Accuracy1 * sumFractPart)), Accuracy1) << sumIntPart;
}

But you only gain considerable speedup at numbers like 100000 compared to the exact factorial function. And at that scale all the functions get really slow.

Tamas Hegedus
  • 28,755
  • 12
  • 63
  • 97
  • 1
    I'd consider naming the function `FactorialLowEstimate` - using integer types frequently gives callers expectation of exact result which is clearly not the case in this computation. – Alexei Levenkov Aug 29 '15 at 01:10