# [Pharo-project] Integer>>#factorial

Nicolas Cellier nicolas.cellier.aka.nice at gmail.com
Thu Jun 16 16:21:33 EDT 2011

```Anyway, algorithm will be dominated by multiplication cost, so you
need at least a Karatsuba algorithm or Toom-Cook or Fourier transform
for monsters factorials.

For Karatsuba, the best is to multiply well equilibrated LargeIntegers
(some having aproximately same size).

It would be optimal to split n! at p!*(n!/p!) such that p! = (n!/p!)
(approximately)
That means n! = p! squared
Using a Stirling like  n! log = (n*n log - n), we find a rough
approximation for the ratio
p/n  = ( n log + 1 ) / (2 * n log)

Thus we can use this approximation:
factorial
self < 3 ifTrue: [^self].
x := (self highBit - 1).
p := (x + 1) * self // x >> 1 - 1.
^(p factorial)*(p+1 productUpTo: self)

For spliting (n!/p!) I suggest using this ratio:
productUpTo: n
n > self ifFalse: [^n].
x := (n highBit - 1)*self highBit.
p = x * (n + self) + (n - self)  // x >> 1.
^(self productUpTo: p)*(p+1 productUpTo: n)

But the highBit and // evaluations cost, so it's not that usefull...

John Brant strategy is more efficient (evaluate 1 out of 2 and recursively).

But the Prime Swing wins on Cog (not the case in VWNC 7.7).

Just play with http://www.squeaksource.com/FactorialContest.html if you want.

Nicolas

2011/6/16 Sven Van Caekenberghe <sven at beta9.be>:
> Hi,
>
>
> I had some old Lisp code hanging around, implementing an algorithm that I once found somewhere, and I ported it to Smalltalk. I am sure there are faster solutions, but this one is twice as fast as the original while still being easy to understand and it has a useful, safe helper method:
>
> Integer>>#svcFactorial
>
>        self = 0 ifTrue: [ ^ 1 ].
>        self > 0 ifTrue: [ ^ 1 productUpTo: self ].
>        self error: 'Not valid for negative integers'
>
> Integer>>#productUpTo: anInteger
>        "Answer the product of all integers from the receiver (non-inclusive) up to anInteger (inclusive)."
>
>        | difference split |
>        self assert: (self between: 0 and: anInteger).
>        "the idea is to multiply LargePositiveIntegers of approximately the same size"
>        difference := anInteger - self.
>        difference = 0 ifTrue: [ ^ 1 ].
>        difference = 1 ifTrue: [ ^ anInteger ].
>        difference = 2 ifTrue: [ ^ (anInteger - 1) * anInteger ].
>        difference = 3 ifTrue: [ ^ (anInteger - 2) * (anInteger - 1) * anInteger ].
>        difference = 4 ifTrue: [ ^ (anInteger - 3) * (anInteger - 2) * (anInteger - 1) * anInteger ].
>        split := (self + anInteger) bitShift: -1.
>        ^ (self productUpTo: split) * (split productUpTo: anInteger)
>
> Would it be a good idea to replace the current implementation with this one ?
>
> Sven
>
>
>

```