Recursive factorial function program returns 0 for large number input in C

Question

I wrote program using recursion for factorial computation as taught in the class. I observed this : factorial of a number greater than 65 consistently gives an output as 0.

//program to recursively calculate the value of factorial for a given integer!

#include <stdio.h>
#include <stdint.h>

unsigned long rfacta(uint32_t n){
    if(n>0)
        return rfacta(n-1)*n;
    else
        return 1;
}

unsigned long rfactd(uint32_t n){
    if(n>0)
        return n*rfactd(n-1);
    else
        return 1;
}

int main(void){
    uint32_t number, methd;
    unsigned long fctorial;
    printf(" Enter 0 for call-time recursion\n Enter 1 for return-time recursion\n Enter 2 to print the table of passible factorials on this system\n\n");
    scanf("%d",&methd);
    switch(methd){
    case 0:
        printf("Compute factorial for number = ? \n");
        scanf("%d",&number);
        fctorial = rfacta(number);
        printf("Factorial(%d) = %lu\n",number,fctorial);
    break;

    case 1:
        printf("Compute factorial for number = ? \n");
        scanf("%d",&number);
        fctorial = rfactd(number);
        printf("Factorial(%d) = %lu\n",number,fctorial);
    break;

    case 2:
      for(int i=0;rfacta(i)!=0;i++){
        printf("Factorial(%d) = %lu\n",i,rfacta(i));
      }
      break;

    default:
      printf("Incorrect entry! \n");
      break; 
    }
    return 0;
}

output Screenshot attachment: Recursion program for factorial returning 0

Question:

As we can see in the output[screenshot image attached], there are three regions:

Region 1: In this region, output for factorial is correct. The reason for this is that,

log2(20!) ~= 62, [here logm(n) means log of n for the base m] i.e. assuming a 64-bit machine, then the largest factorial number than can be represented correctly at single address is 20!, and after which memory overflow happens, because log2(21!) ~= 65. Therefore, such a large number cannot be represented by a 64 bit address.

Region 2: In this region, value of factorial is non-zero but incorrect. The reason for this, as described above is because the factorial value integer is so large that it cannot be represented by single 64-bit address. Thus the what we are instead seeing as the function output is actually random 64 bit long unsigned integers.

Region 3: for Factorial(n > 65), factorial(n)=0.

I am unable to understand why exactly for n>65, n! = 0 as returned by the program?

From what I understand that region 2 must go on and not stop. But it consistently stopping on all machines I have tried my code at n=65 has me puzzled. Thanks in advance for help.

Answer 1

If you print the values in hex, what is happening becomes more apparent:

Factorial(0) = 0000000000000001
Factorial(1) = 0000000000000001
Factorial(2) = 0000000000000002
Factorial(3) = 0000000000000006
Factorial(4) = 0000000000000018
Factorial(5) = 0000000000000078
Factorial(6) = 00000000000002d0
Factorial(7) = 00000000000013b0
Factorial(8) = 0000000000009d80
Factorial(9) = 0000000000058980
Factorial(10) = 0000000000375f00
Factorial(11) = 0000000002611500
Factorial(12) = 000000001c8cfc00
Factorial(13) = 000000017328cc00
Factorial(14) = 000000144c3b2800
Factorial(15) = 0000013077775800
Factorial(16) = 0000130777758000
Factorial(17) = 0001437eeecd8000
Factorial(18) = 0016beecca730000
Factorial(19) = 01b02b9306890000
Factorial(20) = 21c3677c82b40000
Factorial(21) = c5077d36b8c40000
Factorial(22) = eea4c2b3e0d80000
Factorial(23) = 70cd7e2933680000
Factorial(24) = 9343d3dcd1c00000
Factorial(25) = 619fb0907bc00000
Factorial(26) = ea37eeac91800000
Factorial(27) = b3e62c3358800000
Factorial(28) = ad2cd59dae000000
Factorial(29) = 9e1432dcb6000000
Factorial(30) = 865df5dd54000000
Factorial(31) = 4560c5cd2c000000
Factorial(32) = ac18b9a580000000
Factorial(33) = 2f2fee5580000000
Factorial(34) = 445da75b00000000
Factorial(35) = 58cde17100000000
Factorial(36) = 7cf3b3e400000000
Factorial(37) = 0f38fff400000000
Factorial(38) = 4275fe3800000000
Factorial(39) = 1ff9ba8800000000
Factorial(40) = ff05254000000000
Factorial(41) = d7d2f74000000000
Factorial(42) = 689c908000000000
Factorial(43) = 924c458000000000
Factorial(44) = 251bf20000000000
Factorial(45) = 85e98a0000000000
Factorial(46) = 0ff6cc0000000000
Factorial(47) = ee4f740000000000
Factorial(48) = aee5c00000000000
Factorial(49) = 79f9c00000000000
Factorial(50) = d2c7800000000000
Factorial(51) = fdbe800000000000
Factorial(52) = 8ab2000000000000
Factorial(53) = b6da000000000000
Factorial(54) = 91fc000000000000
Factorial(55) = 5d24000000000000
Factorial(56) = 5fe0000000000000
Factorial(57) = 58e0000000000000
Factorial(58) = 22c0000000000000
Factorial(59) = 0240000000000000
Factorial(60) = 8700000000000000
Factorial(61) = 2b00000000000000
Factorial(62) = 6a00000000000000
Factorial(63) = 1600000000000000
Factorial(64) = 8000000000000000
Factorial(65) = 8000000000000000

As you continue to multiply in values containing 2 as a factor, the number of trailing zeros continues to increase. By the time you multiply in 66, all nonzero bits have been pushed out to the left, so you're left with 0.

Also, the values from 21! to 65! are not actually random values but the low order 64 bits of the result. Unsigned integer arithmetic is carried out modulo 2^bitlen where "bitlen" is the bit length of the type in question which is 64 in this case.

Answer 2

As long as you multiply with two’s complement and wrapping (which is not guaranteed by the C standard but appears to be happening in your program), the computed result is the full mathematical result (without wrapping) modulo 2⁶⁴.

For n > 65, n! is a multiple of 2⁶⁴ (the factors from 1 to 66 have included 64 or more twos), so the value modulo 2⁶⁴ is zero:

• When you multiply by 2, one factor of two is provided.
• When you multiply by 4, two factors are provided.
• When you multiply by 6, one factor is provided.
• When you multiply by 8, three factors are provided.
• The total number of factors of two in n! is n/2 + n/4 + n/8 + …, using integer arithmetic (that is, truncating the division, so 3/2 = 1) and continuing until the term reaches zero. For 66!, this is 66/2 + 66/4 + 66/8 + 66/16 + 66/32 + 66/64 = 33 + 16 + 8 + 4 + 2 + 1 = 64.