8

I need to calculate the factorial of numbers up to around 100! in order to determine if a series of coin flip-style data is random, as per this Wikipedia entry on Bayesian probability. As you can see there, the necessary formula involves 3 factorial calculations (but, interestingly, two of those factorial calculations are calculated along the way to the third).

I saw this question here, but I'd think that integer is going to get blown out pretty quickly. I could also make a function that is more intelligent about the factorial calculation (ie, if I have 11!/(7!3!), as per the wiki example, I could go to (11*10*9*8)/3!), but that smacks of premature optimization to me, in the sense that I want it to work, but I don't care about speed (yet).

So what's a good C# library I can call to calculate the factorial in order to get that probability? I'm not interested in all the awesomeness that can go into factorial calculation, I just want the result in a way that I can manipulate it. There does not appear to be a factorial function in the Math namespace, hence the question.

Community
  • 1
  • 1
mmr
  • 14,781
  • 29
  • 95
  • 145

6 Answers6

8

You could try Math.NET - I haven't used that library, but they do list Factorial and Logarithmic Factorial.

TrueWill
  • 25,132
  • 10
  • 101
  • 150
  • 2
    See also the [same](http://numerics.mathdotnet.com/api/MathNet.Numerics/SpecialFunctions.htm#Factorial) in [Math.NET Numerics](http://numerics.mathdotnet.com), the successor of the linked Math.NET Iridium library. – Christoph Rüegg Dec 14 '13 at 08:35
4

There has been a previous question on a similar topic. Someone there linked the Fast Factorial Functions web site, which includes some explanations of efficient algorithms and even C# source code.

Community
  • 1
  • 1
bobbymcr
  • 23,769
  • 3
  • 56
  • 67
  • Checking it out now-- man, I really wish there were some best-approximation-may-not-work-for-everyone implementation as part of the base libraries. – mmr Sep 30 '09 at 02:57
  • 1
    The problem with that code is that the downloads appear to be all in java, and I'm just not interested in the work to convert from java to C#. The author probably has C# code to download, but the Math.Net solution puts the result into a double, which is exactly what I need. – mmr Sep 30 '09 at 18:42
3

Do you want to calculate factorials, or binomial coefficients?

It sounds like you want to calculate binomial coefficients - especially as you mention 11!/(7!3!).

There may be a library that can do this for you, but as a (presumably) programmer visiting stack overflow there's no reason not to write one yourself. It's not too complicated.

To avoid memory overflow, don't evaluate the result until all common factors are removed.

This algorithm still needs to be improved, but you have the basis for a good algorithm here. The denominator values need to be split into their prime factors for the best result. As it stands, this will run for n = 50 quite quickly.

float CalculateBinomial(int n, int k)
{
    var numerator = new List<int>();
    var denominator = new List<int>();
    var denominatorOld = new List<int>();

    // again ignore the k! common terms
    for (int i = k + 1; i <= n; i++)
        numerator.Add(i);

    for (int i = 1; i <= (n - k); i++)
    {
        denominator.AddRange(SplitIntoPrimeFactors(i));
    }

    // remove all common factors
    int remainder;                
    for (int i = 0; i < numerator.Count(); i++)
    {
        for (int j = 0; j < denominator.Count() 
            && numerator[i] >= denominator[j]; j++)
        {
            if (denominator[j] > 1)
            {
                int result = Math.DivRem(numerator[i], denominator[j], out remainder);
                if (remainder == 0)
                {
                    numerator[i] = result;
                    denominator[j] = 1;
                }
            }
        }
    }

    float denominatorResult = 1;
    float numeratorResult = 1;

    denominator.RemoveAll(x => x == 1);
    numerator.RemoveAll(x => x == 1);

    denominator.ForEach(d => denominatorResult = denominatorResult * d);
    numerator.ForEach(num => numeratorResult = numeratorResult * num);

    return numeratorResult / denominatorResult;
}

static List<int> Primes = new List<int>() { 2, 3, 5, 7, 11, 13, 17, 19, 
    23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97 };

List<int> SplitIntoPrimeFactors(int x)
{
    var results = new List<int>();
    int remainder = 0;

    int i = 0;
    while (!Primes.Contains(x) && x != 1)
    {
        int result = Math.DivRem(x, Primes[i], out remainder);
        if (remainder == 0)
        {
            results.Add(Primes[i]);
            x = result;
            i = 0;
        }
        else
        {
            i++;
        }
    }
    results.Add(x);
    return results;
}

I can estimate n = 110, k = 50 (returns 6x10^31) but cannot run n = 120, k = 50.

Kirk Broadhurst
  • 27,836
  • 16
  • 104
  • 169
  • And what happens when n ~= 50 and k ~= 2? Overflow! I need some way to handle nontrivial cases. – mmr Sep 30 '09 at 02:55
  • In that case you calculate n = 50 and k = 48 as I pointed out. – Kirk Broadhurst Sep 30 '09 at 02:58
  • ok. What's 48!? Will that fit into an integer? Because that's what your denominator is doing in that case. – mmr Sep 30 '09 at 03:00
  • Regardless of the algorithm, the result of, say, nCr for n = 100 and r = 50 will overflow a 64 bit number. You need to go back to your original problem because you can't store that result. – Kirk Broadhurst Sep 30 '09 at 04:07
  • 1
    Only if I don't go with an algorithm that uses something like Stirling's Approximation or the like, putting the result into a double. The next step in the algorithm is a rudimentary integration, so 'close enough' will be just fine, and the precision loss on a double is ok. – mmr Sep 30 '09 at 18:40
1

The following can calculate the factorial of 5000 in 1 second.

public class Number
{
    #region Fields
    private static long _valueDivision = 1000000000;
    private static int _valueDivisionDigitCount = 9;
    private static string _formatZeros = "000000000";
    List<long> _value;
    #endregion

    #region Properties
    public int ValueCount { get { return _value.Count; } }
    public long ValueAsLong
    {
        get
        {
            return long.Parse(ToString());
        }
        set { SetValue(value.ToString()); }
    }
    #endregion

    #region Constructors
    public Number()
    {
        _value = new List<long>();
    }
    public Number(long value)
        : this()
    {
        SetValue(value.ToString());
    }
    public Number(string value)
        : this()
    {
        SetValue(value);
    }
    private Number(List<long> list)
    {
        _value = list;
    }
    #endregion

    #region Public Methods
    public void SetValue(string value)
    {
        _value.Clear();
        bool finished = false;
        while (!finished)
        {
            if (value.Length > _valueDivisionDigitCount)
            {
                _value.Add(long.Parse(value.Substring(value.Length - _valueDivisionDigitCount)));
                value = value.Remove(value.Length - _valueDivisionDigitCount, _valueDivisionDigitCount);
            }
            else
            {
                _value.Add(long.Parse(value));
                finished = true;
            }
        }
    }
    #endregion

    #region Static Methods
    public static Number operator +(Number c1, Number c2)
    {
        return Add(c1, c2);
    }
    public static Number operator *(Number c1, Number c2)
    {
        return Mul(c1, c2);
    }
    private static Number Add(Number value1, Number value2)
    {
        Number result = new Number();
        int count = Math.Max(value1._value.Count, value2._value.Count);
        long reminder = 0;
        long firstValue, secondValue;
        for (int i = 0; i < count; i++)
        {
            firstValue = 0;
            secondValue = 0;
            if (value1._value.Count > i)
            {
                firstValue = value1._value[i];
            }
            if (value2._value.Count > i)
            {
                secondValue = value2._value[i];
            }
            reminder += firstValue + secondValue;
            result._value.Add(reminder % _valueDivision);
            reminder /= _valueDivision;
        }
        while (reminder > 0)
        {
            result._value.Add(reminder % _valueDivision);
            reminder /= _valueDivision;
        }
        return result;
    }
    private static Number Mul(Number value1, Number value2)
    {
        List<List<long>> values = new List<List<long>>();
        for (int i = 0; i < value2._value.Count; i++)
        {
            values.Add(new List<long>());
            long lastremain = 0;
            for (int j = 0; j < value1._value.Count; j++)
            {
                values[i].Add(((value1._value[j] * value2._value[i] + lastremain) % _valueDivision));
                lastremain = ((value1._value[j] * value2._value[i] + lastremain) / _valueDivision);
                //result.Add(();
            }
            while (lastremain > 0)
            {
                values[i].Add((lastremain % _valueDivision));
                lastremain /= _valueDivision;
            }
        }
        List<long> result = new List<long>();
        for (int i = 0; i < values.Count; i++)
        {
            for (int j = 0; j < i; j++)
            {
                values[i].Insert(0, 0);
            }
        }
        int count = values.Select(list => list.Count).Max();
        int index = 0;
        long lastRemain = 0;
        while (count > 0)
        {
            for (int i = 0; i < values.Count; i++)
            {
                if (values[i].Count > index)
                    lastRemain += values[i][index];
            }
            result.Add((lastRemain % _valueDivision));
            lastRemain /= _valueDivision;
            count -= 1;
            index += 1;
        }
        while (lastRemain > 0)
        {
            result.Add((lastRemain % _valueDivision));
            lastRemain /= _valueDivision;
        }
        return new Number(result);
    }
    #endregion

    #region Overriden Methods Of Object
    public override string ToString()
    {
        string result = string.Empty;
        for (int i = 0; i < _value.Count; i++)
        {
            result = _value[i].ToString(_formatZeros) + result;
        }
        return result.TrimStart('0');
    }
    #endregion
}

class Program
{
    static void Main(string[] args)
    {
        Number number1 = new Number(5000);
        DateTime dateTime = DateTime.Now;
        string s = Factorial(number1).ToString();
        TimeSpan timeSpan = DateTime.Now - dateTime;
        long sum = s.Select(c => (long) (c - '0')).Sum();
    }
    static Number Factorial(Number value)
    {
        if( value.ValueCount==1 && value.ValueAsLong==2)
        {
            return value;
        }
        return Factorial(new Number(value.ValueAsLong - 1)) * value;
    }
}
psubsee2003
  • 8,563
  • 8
  • 61
  • 79
-1
using System;
//calculating factorial with recursion
namespace ConsoleApplication2
{
    class Program
    {
        long fun(long a)
        {
            if (a <= 1)
            {
                return 1;}
            else
            {
                long c = a * fun(a - 1);
                return c;
            }}

        static void Main(string[] args)
        {

            Console.WriteLine("enter the number");
            long num = Convert.ToInt64(Console.ReadLine());
            Console.WriteLine(new Program().fun(num));
            Console.ReadLine();
        }
    }
}
mmr
  • 14,781
  • 29
  • 95
  • 145
  • 2
    Please read the other answers. This is not a sufficient solution, because long does not have the capacity to handle the size of a number near 100!. As such, some other method (such as Stirling's approximation) need to be used instead. – mmr Feb 10 '11 at 22:30
  • Sorry, but recursive approach is not efficient in any way. – equintas Aug 02 '21 at 17:41
-4

hello everybody according to this solution i have my own solution where i calculate factorial of array 1D elements. the code is `int[] array = new int[5] { 4,3,4,3,8 };

        int fac = 1;

        int[] facs = new int[array.Length+1];

        for (int i = 0; i < array.Length; i++)
        {
            for (int j = array[i]; j > 0; j--)
            {
                fac *= j;
            }
            facs[i] = fac;
            textBox1.Text += facs[i].ToString() + " ";
            fac = 1;
        }`

copy and paste the code above ^ in the button , it solves factorial of elements of array 1D. best regards.

mehvan
  • 1
  • 1