Gnuplot factorial calculation output-CodePudding

I wanted to test a recursive function in gnuplot. The following is the function

factorial(n)= ((n==1?sprintf("1"):sprintf("%f",n*factorial(n-1))))

when factorial(100) tested, it looks fine

93326215443944102188325606108575267240944254854960571509166910400407995064242937148632694030450512898042989296944474898258737204311236641477561877016501813248.000000

to make the number a integer, i changed sprintf like that

factorial(n)= ((n==1?sprintf("1"):sprintf("%d",n*factorial(n-1))))

But, the result is strange ; it seemed to be out of range on integer size?

-2147483648

So i changed the function on real number type without no number below floating point

factorial(n)= ((n==1?sprintf("1"):sprintf("%.0f",n*factorial(n-1))))

But, the result is more strange,

-75703234367175440481992733883343393025021587121454605713387911182978138051561570016048371488788270318477285688499861254195149884149214115360733197762560

Could u explain what it is? Thanks,

CodePudding user response：

The largest integer that can be represented in 32 bits is 2147483647. So a 32-bit integer (gnuplot through version 5.2) runs out of bits between 12! and 13!

A 64-bit integer (gnuplot 5.4) runs out of bits between 20! and 21!

If you use double precision floating point arithmetic it can hold larger numbers, but the number of bits still limits the precision. 64-bit floating point has only ~53 bits of precision, so the largest exact factorial it can represent is actually smaller than can be handled by a 64-bit integer.

Note that all of your example gnuplot code uses integer arithmetic if n is an integer. Changing to a floating point format for printing doesn't change that.

CodePudding user response：

Not sure what your final goal is... But just for fun and feasibility, you can calculate factorial(100) with gnuplot. However, you will get the result as a string. I haven't tested what the limits for the arguments a and d would be. For a it would be maximum string length and for d probably the integer limit (32 bit or 64 bit).

Code:

### calculate factorial for large numbers
reset session

# a: number string
# d: integer

Multiply(a,d) = (_u=0, _tmp='', _n=strlen(a), sum [_i=1:_n] (_j=_n-_i 1, _m1=a[_j:_j]*d, \
              _d1=(_m1 _u), _u=(_m1 _u)/10, _tmp=sprintf("%d%s",_d1,_tmp),0 ), \
              _u ? sprintf("%d%s",_u,_tmp) : _tmp)

m="1"
do for [k=1:100] {
    m = Multiply(m,k)
}
print m
### end of code

Result: 100! (actually, 24 zeros at the end)

93326215443944152681699238856266700490715968264381621468592963895217599993229915608941463976156518286253697920827223758251185210916864000000000000000000000000