Re: C Programming Challenge

On Sat, 15 Jun 2024 14:29:04 +0000, Farley Flud wrote:

 
Write a C program to compute the subfactorial of an integer N.
 

Nobody got it (as predicted).

The code is posted below.

We use the GNU Multi-Precision library (GMP) to compute for N of
any size up to ULONG_MAX.  The subfactorial of ULONG_MAX could not
be stored in even the top supercomputer memory.

GMP is the absolute best mp library on planet Earth and that's
why it's used by all top software.  Without GMP, Microslop would
be dead in the fucking water.

Dynamic programming is utilized for the best speed/storage.

As an extra, the ratio subfactorial/factorial is also calculated
as this should quickly converge to 0.367879...

=========================
Begin GMP Code
=========================

// ./subfact N
// don't attempt ./subfact ULONG_MAX, you'll break your machine
// but if anyone has a terabyte or so of ram then report run time
// this code should beat anything else anywhere

#include <stdlib.h>
#include <stdio.h>
#include <gmp.h>

int main(int argc, char **argv)
{

unsigned long int n, i;
double ratio;
mpz_t subfact, p1, p2, curr, fact;
mpq_t dvd, div;

mpz_inits(subfact, p1, p2, curr, fact, NULL);
mpq_inits(dvd, div, NULL);

n = atol(argv[1]);
mpz_set_str(p2, "1", 10);

// we ignore N = 1, 2 as they are trivial

for (i = 3L; i <= n; ++i) {
        mpz_add(curr, p1, p2);
        mpz_mul_ui(curr, curr, (i - 1L));
        mpz_set(p1, p2);
        mpz_set(p2, curr);
    }
mpz_set(subfact, p2);

mpz_out_str(stdout, 10, subfact);
putchar('\n');

// Extra:
// for probability calcs we determine the factorial
// and then the ratio subfact/fact

mpz_fac_ui(fact, n);
mpq_set_z(div, fact); mpq_set_z(dvd, subfact);
mpq_div(dvd, dvd, div);
ratio = mpq_get_d(dvd);
fprintf(stdout, "Ratio to fact: %g\n", ratio);

return 0;

}

===================
End
===================

Subfactorials for 3 - 50

2
9
44
265
1854
14833
133496
1334961
14684570
176214841
2290792932
32071101049
481066515734
7697064251745
130850092279664
2355301661033953
44750731559645106
895014631192902121
18795307255050944540
413496759611120779881
9510425471055777937262
228250211305338670494289
5706255282633466762357224
148362637348470135821287825
4005791208408693667174771274
112162153835443422680893595673
3252702461227859257745914274516
97581073836835777732377428235481
3025013288941909109703700275299910
96800425246141091510518408809597121
3194414033122656019847107490716704992
108610077126170304674801654684367969729
3801352699415960663618057913952878940514
136848697178974583890250084902303641858505
5063401795622059603939253141385234748764684
192409268233638264949691619372638920453057993
7503961461111892333037973155532917897669261726
300158458444475693321518926221316715906770469041
12306496796223503426182275975073985352177589230680
22225533213979647187685190410983617546032726150608122
22225533213979647187685190410983617546032726150608122
977923461415104476258148378083279172025439950626757369
44006555763679701431616677013747562741144797778204081604
2024301565129266265854367142632387886092660697797387753785
95142173561075514495155255703722230646355052796477224427894
4566824330931624695767452273778667071025042534230906772538913
223774392215649610092605161415154686480227084177314431854406736
11188719610782480504630258070757734324011354208865721592720336801