1000 factorial

Well, someone wrote that they could write a program to do 1000! in one line of Haskell. I said that it would take less than 100 lines of C to do the same. I wrote the program, and it is less than 100 lines of C:

#include <stdio.h>
#include <string.h>

#define LENGTH 10000

int max (int a, int b) { return (a > b) ? a : b; }
int min (int a, int b) { return (b > a) ? a : b; }

int length (unsigned char * a)
 {
   int i;
   i = LENGTH - 1;
   while ((a[i] == 0) && (i >= 0)) i--;
   return i + 1;
 }

void add (unsigned char * a, unsigned char * b)
 {
   int i, c, C;
   c = min(max(length(a), length(b)) + 1, LENGTH);
   for (i = 0, C = 0; i < c; i++)
    {
      a[i] += b[i] + C;
      C = a[i] > 99;
      if (C) a[i] -= 100;
    }
 }

   unsigned char x [LENGTH], y [LENGTH];
void mul (unsigned char * a, unsigned char * b)
 {
   int i, j, c, d, C;
   memset(y, '\0', LENGTH);
   c = min(length(a) + length(b) + 1, LENGTH);
   d = length(b);
   for (i = 0; i < d; i++)
    {
      memset(x, '\0', LENGTH);
      for (j = 0, C = 0; (i + j) < c; j++)
       {
         x[i + j] = (a[j] * b[i] + C) % 100;
         C = (a[j] * b[i] + C) / 100;
       }
      add(y, x);
    }
   memcpy(a, y, LENGTH);
 }

void numset (unsigned char * a, unsigned long b)
 {
   int i;
   for (i = 0; i < LENGTH; i++)
    {
      a[i] = b % 100;
      b /= 100;
    }
 }

int tz (unsigned char * a)
 {
   int i;
   for (i = 0; i < LENGTH; i++) if (a[i] != 0) return 0;
   return 1;
 }

   unsigned char minus1 [LENGTH], prod [LENGTH];
void fact (unsigned char * a)
 {
   numset(prod, 1);
   while (!tz(a))
    {
      mul(prod, a);
      add(a, minus1);
    }
   memcpy(a, prod, LENGTH);
 }

   unsigned char arg_res [LENGTH];
int main (void)
 {
   unsigned long a;
   int i;
   memset (minus1, 99, LENGTH);
   printf("Input the number to take the factorial of: ");
   scanf("%lu", &a);
   numset(arg_res, a);
   fact(arg_res);
   i = LENGTH - 1;
   while (arg_res[i] == 0) i--;
   while (i >= 0)
    {
      printf("%02d", (int) arg_res[i]);
      i--;
    }
   putchar('\n');
   return 0;
 }

>>1
... Wow. It's even worsebetter than mine^[1].

[1] http://dis.4chan.org/read/prog/1296932382/5

From my calculator:
4.02387260077093773543702433923e+2567

The thought that I would ever have use for a number this great confounds my little brain.

'\0'

>>3
402387260077093773543702433923003985719374864210714632543799910429938512398629020592044208486969404800479988610197196058631666872994808558901323829669944590997424504087073759918823627727188732519779505950995276120874975462497043601418278094646496291056393887437886487337119181045825783647849977012476632889835955735432513185323958463075557409114262417474349347553428646576611667797396668820291207379143853719588249808126867838374559731746136085379534524221586593201928090878297308431392844403281231558611036976801357304216168747609675871348312025478589320767169132448426236131412508780208000261683151027341827977704784635868170164365024153691398281264810213092761244896359928705114964975419909342221566832572080821333186116811553615836546984046708975602900950537616475847728421889679646244945160765353408198901385442487984959953319101723355556602139450399736280750137837615307127761926849034352625200015888535147331611702103968175921510907788019393178114194545257223865541461062892187960223838971476088506276862967146674697562911234082439208160153780889893964518263243671616762179168909779911903754031274622289988005195444414282012187361745992642956581746628302955570299024324153181617210465832036786906117260158783520751516284225540265170483304226143974286933061690897968482590125458327168226458066526769958652682272807075781391858178889652208164348344825993266043367660176999612831860788386150279465955131156552036093988180612138558600301435694527224206344631797460594682573103790084024432438465657245014402821885252470935190620929023136493273497565513958720559654228749774011413346962715422845862377387538230483865688976461927383814900140767310446640259899490222221765904339901886018566526485061799702356193897017860040811889729918311021171229845901641921068884387121855646124960798722908519296819372388642614839657382291123125024186649353143970137428531926649875337218940694281434118520158014123344828015051399694290153483077644569099073152433278288269864602789864321139083506217095002597389863554277196742822248757586765752344220207573630569498825087968928162753848863396909959826280956121450994871701244516461260379029309120889086942028510640182154399457156805941872748998094254742173582401063677404595741785160829230135358081840096996372524230560855903700624271243416909004153690105933983835777939410970027753472000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000

here it OP, also a one-liner

printf("402387260077093773543702433923003985719374864210714632543799910429938512398629020592044208486969404800479988610197196058631666872994808558901323829669944590997424504087073759918823627727188732519779505950995276120874975462497043601418278094646496291056393887437886487337119181045825783647849977012476632889835955735432513185323958463075557409114262417474349347553428646576611667797396668820291207379143853719588249808126867838374559731746136085379534524221586593201928090878297308431392844403281231558611036976801357304216168747609675871348312025478589320767169132448426236131412508780208000261683151027341827977704784635868170164365024153691398281264810213092761244896359928705114964975419909342221566832572080821333186116811553615836546984046708975602900950537616475847728421889679646244945160765353408198901385442487984959953319101723355556602139450399736280750137837615307127761926849034352625200015888535147331611702103968175921510907788019393178114194545257223865541461062892187960223838971476088506276862967146674697562911234082439208160153780889893964518263243671616762179168909779911903754031274622289988005195444414282012187361745992642956581746628302955570299024324153181617210465832036786906117260158783520751516284225540265170483304226143974286933061690897968482590125458327168226458066526769958652682272807075781391858178889652208164348344825993266043367660176999612831860788386150279465955131156552036093988180612138558600301435694527224206344631797460594682573103790084024432438465657245014402821885252470935190620929023136493273497565513958720559654228749774011413346962715422845862377387538230483865688976461927383814900140767310446640259899490222221765904339901886018566526485061799702356193897017860040811889729918311021171229845901641921068884387121855646124960798722908519296819372388642614839657382291123125024186649353143970137428531926649875337218940694281434118520158014123344828015051399694290153483077644569099073152433278288269864602789864321139083506217095002597389863554277196742822248757586765752344220207573630569498825087968928162753848863396909959826280956121450994871701244516461260379029309120889086942028510640182154399457156805941872748998094254742173582401063677404595741785160829230135358081840096996372524230560855903700624271243416909004153690105933983835777939410970027753472000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000000");

New challenge: Write the shortest C program that outputs the factorial of 1000.

Deadline: 2011-04-17

Prize: 10 Susscoins, to be delivered via Sussmail

>>7
#include <stdio.h>
#include <gmp.h>
main(){mpz_t a;int b=1001;mpz_init_set_ui(a,1);while(--b)mpz_mul_ui(a,a,b);mpz_out_str(stdout,10,a);}

>>8
main(){system("dc -e '1000[d1-d1<F*]dsFxp'");}

#! /usr/bin/python

__import__("sys").setrecursionlimit(1002)

def fac(n): return n * fac(n-1) if n else 1

print fac(1000)

>>10
>>7
shortest C program
Your shitty, long^[1] phyok program is not valid.

[1]
(define (f n) (if (zero? n) 1 (* n (f (sub1 n)))))
(f 1000)

>>11
Not tail recursive.

(define (f n) (let loop ([n n] [i 1]) (if (= n 0) i (loop (sub1 n) (* i n)))))
(f 1000)

>>12
Not OMG OPTIMIZED

(define small-odd-swing
  (vector 1 1 1 3 3 15 5 35 35 315 63 693 231 3003 429 6435 6435 109395
          12155 230945 46189 969969 88179 2028117 676039 16900975 1300075
          35102025 5014575 145422675 9694845 300540195 300540195))

(define (swing n)
  (if (< n 33) (vector-ref small-odd-swing n)
      (let* ((primes (sieve (range 3 n #:by 2)))
             (rootn (integer-sqrt n))
             (primesa (take-while (curryr <= rootn) primes))
             (primesb (filter (λ (x) (odd? (quotient n x)))
                                 (take-while (curryr <= (quotient n 3))
                                             (drop-while (curryr <= rootn) primes)))))
        (apply
         *
         (append
          (map
           (λ (prime)
             (do ((q (quotient n prime) (quotient q prime))
                  (p 1 (if (odd? q) (* p prime) p)))
               ((zero? q) p)))
           primesa)
          (take-while (curryr <= n)
                      (drop-while (curryr <= (quotient n 2)) primes))
          primesb)))))

(define (rec-fact n)
  (if (< n 2) 1
      (* (expt (rec-fact (quotient n 2)) 2) (swing n))))

(define (fact n)
  (if (< n 20) (Π values 1 n)
      (arithmetic-shift (rec-fact n) (- n (bitwise-bit-count n)))))
(fact 1000)

>>13
That's... optimized?

>>14
Yes it is.

; just for reference, Π is:
#;
(define (Π f n m #:threads (c (processor-count)))
  (let* ((t (lambda (a)
              (lambda ()
                (let/ec k
                  (do ((x (+ n a) (+ x c))
                       (r 1 (* r (f x))))
                    ((> x m) (k r)))))))
         (xs (do ((x (sub1 c) (sub1 x))
                  (xs '() (cons (future (t x)) xs)))
               ((zero? x) xs)))
         (x ((t 0))))
    (foldr * x (map touch xs))))

(define (naive-fact n)
  (Π values 1 n))

(collect-garbage)
(time (void (fact 100000)))
(collect-garbage)
(time (void (naive-fact 100000)))

; cpu time: 1030 real time: 1029 gc time: 520
; cpu time: 21552 real time: 21551 gc time: 7063


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

int digits(int n)
{
    int digs = 0;
    for(;n!=0;n/=100) ++digs;
    return digs;
}

int mul(unsigned int len, unsigned char **num, unsigned int n)
{
    unsigned int rem = 0, i;
    for(i=0; i<len; i++) {
        unsigned int tmp = (*num)[i]*n + rem;
        (*num)[i] = tmp % 100;
        rem = tmp / 100;
    }
    if(rem!=0) {
        int addLen = digits(rem)+10;
        if((*num=realloc(*num, len+addLen))==NULL)
            return -1;      
        memset(&(*num)[len], 0, addLen);
        for(i=len; i<len+addLen; i++) {
            unsigned int tmp = (*num)[i]*n + rem;
            (*num)[i] = tmp % 100;
            rem = tmp / 100;
        }
        len += addLen;
    }
    return len;
}

int main() {
    unsigned int i, len=1000, factOf=100000;
    unsigned char *num = (unsigned char*)malloc(len*sizeof(unsigned char));
    memset(num, 0, len);
    num[0]=1;
    for(i=1; i<=factOf; i++)
        if ((len = mul(len, &num, i)) == -1) {
            printf("Hurrrr\n");
            return -1;
        }
    i = len;
    while(i--) printf("%02d", num[i]);
    putchar('\n');
    return 0;
}
real    0m41.340s
user    0m41.264s
sys     0m0.014s

>>16
real    0m41.340s
in >>15, those are milliseconds. It's 1.030s and 21.552s.

>>17
My levels of autism are not high enough.

>>18
C confirmed for RUBY AS FUCK

>>8 is a cheapo for using a bignum library. His entry will be downgraded.

>>9 is cheating.

>>10 doesn't seem to know the difference between C and Python.

>>11-12 seems to be confusing programming languages as well.

>>13,15 is an expert Scheme programer, but doesn't seem to know the difference between C and Scheme. I'd give it the most interesting implementation award, however.

>>16 qualifies.

>>17 >>7 doesn't say anything about performance.

>>11,20

I'm calling it right now, you guys are fucking morons, it wasn't a response to your retarded little SEE challenge, it was a response to the original post. It even runs faster than your pathetic entries, if you're writing C at least do it right.

>>21
Calm down, bro.

>>20
>>1 is good enough, without using a bignum library. We'd just reimplement it over and over. I'd say this should be a ``make a fast/complete C bignum library'' challenge.

Reminds me of http://www.frontiernet.net/~fys/hugi/compoold.htm#compo26

>>16
OP here. Nice. The math is problem specific, and faster for it.
>>21
Also, nice. You took a language which you wrote factorial in one line, and wrote it in alot of lines, and showed that version to be superior. How the hell does it work?

You have to be smarter than the language. How many C programmers can do ad hoc bignum math that does more than addition and subtraction? How many Scheme programmers can write >>13 ?

>>25
>>21
I'm not the deluded pythonista! I'm _usually respectful, and knew that >>13 is not a valid entry.
I tried to do it in C, and wanted to make the bignum part faster, smaller (using a whole uintmax) and general, but it seems I can't write C code anymore. I give up, it would be much easier in Assembly.

What I've done so far:

#include <stdio.h>
#include <stdlib.h>
#include <inttypes.h>
#include <string.h>

#define SIZE(sz) (sizeof(uintmax_t)*(sz))

typedef struct {
     uintmax_t *val;
     size_t len;
} bignum;

static inline max (size_t a, size_t b) {return a>b?a:b;}
static inline min (size_t a, size_t b) {return a<b?a:b;}

bignum *bignum_alloc(size_t sz) {
     bignum *n = malloc(sizeof(bignum));
     n->len = sz;
     n->val = malloc(SIZE(sz));
     return n;
}

bignum *bignum_copy(bignum *n) {
     bignum *m = bignum_alloc(n->len);
     memcpy(m->val, n->val, n->len);
     return m;
}

void bignum_realloc(bignum *n, size_t sz, uintmax_t val) {
     size_t oldsz = n->len;
     n->val = realloc(n->val, SIZE(sz));
     n->len = sz;
     n->val[oldsz] = val;
     if (oldsz < sz-1)
      memset(n+oldsz, 0, oldsz+1);
}

void bignum_dealloc(bignum *n) {
     free(n->val);
     free(n);
}

bignum *bignum_from_uintmax(uintmax_t n) {
     bignum *m = bignum_alloc(1);
     m->val[0] = n;
     return m;
}

void bignum_addbi(bignum *n, uintmax_t m) {
     size_t i, j = n->len;
     uintmax_t a;

     for (a = n->val[i], i = 0; m && i < j; i++) {
      a = n->val[i];
      n->val[i] += m;
      m = a > n->val[i];
     }
     if (m && i >= n->len)
      bignum_realloc(n, n->len+1, 1);
}


void bignum_addbb(bignum *n, bignum *m) {
     size_t i, j = max(n->len, m->len);
     uintmax_t a; int c;
     if (n->len < m->len)
      bignum_realloc(n, m->len, 0);
     for (i = 0, c = 0; i < j; i++) {
      a = n->val[i];
      n->val[i] += m->val[i] + c;
      c = a > n->val[i];
     }
     if (c && i >= j++)
      bignum_realloc(n, j, 1);
}

void bignum_subbi(bignum *n, uintmax_t m) {
     size_t i, j = n->len;
     uintmax_t a;

     for (a = n->val[i], i = 0; m && i < j; i++) {
      a = n->val[i];
      n->val[i] -= m;
      m = a < n->val[i];
     }
     if (--j && !n->val[j])
      bignum_realloc(n, j, n->val[j]);
}

void bignum_mulbb(bignum *n, bignum *m) {
     if (bignum_zerop(m)) {
      bignum_dealloc(n);
      n = bignum_from_uintmax(0);
     }
     if (bignum_equalpi(m, 1)) return;
     bignum *i = bignum_copy(m);
     bignum *l = bignum_copy(n);
     for (; !bignum_equalpi(i, 1); bignum_subbi(i, 1))
      bignum_addbb(n, l);
     bignum_dealloc(i);
     bignum_dealloc(l);
}

void bignum_mulbi(bignum *n, uintmax_t m) {
     if (!m) {
      bignum_dealloc(n);
      n = bignum_from_uintmax(0);
     }
     if (m == 1) return;
     bignum *l = bignum_copy(n);
     while (--m)
      bignum_addbb(n, l);
}

int bignum_zerop(bignum *n) {
     return (!(!(n->len-1) && n->val[0]));
}

int bignum_equalp(bignum *n, bignum *m) {
     if (n->len != m->len) return 0;
     size_t i = n->len;
     while (n->val[i] == m->val[i]) i--;
     return !!i;
}

inline int bignum_equalpi(bignum *n, uintmax_t m) {
     return (!(n->len > 1)) && (n->val[0] == m);
}

void debug_print_bn(bignum *n) {
     printf("Memory Location: %#08x; %#08x\n", n, n->val);
     printf("Length: %u\n", n->len);
     size_t i;
     for (i = 0; i < n->len; i++)
      printf("Cell#%08u: %016llx\n", i, n->val[i]);
}

>>26
OVER 9000 REALLOC OPERATION
REWRITE THAT SHIT

Or the weenie can learn to write proper code using a functional programming language like Haskell.

>>27
I will not, I fucking hate writing C and managing memory.

Listen to >>27-san, >>29-kun.

>>26
Is it me, or does that look like GNU bc code?
>>21
You can call me a fucking retard when you post a portable multithreaded FFT multiplication implementation that's one shot for a factorial program. I'm not a retard just because Scheme holds your little hand and hides the machine from you.

>>31
Shit, am I that bad? I've never seen GNU bc's code.

>>21
Scheme
I am the Schemer, that was the deluded pythonista troll. I'm aware I suck at C.

Also バンパ·チャレンジ.

>>32
``バンパ''?

>>33
Typo. *バンプ.

>>32
suck my cock faggot

>>32
Well, sorry. I still need to try to wrap my head around what that code does. What is a "Bumper Challenge" (according to Google Translate)?

>>36
It does prime factorization, and basically is a faster version of

(define (factor x p)
  (let loop ((x (quotient x p))
             (r 0))
    (if (< x 2) (+ x r)
        (loop (quotient x p)
              (+ x r)))))

(define (prime-factorization x)
  (let ((primes (sieve (cons 2 (range 3 x #:by 2)))))
    (map expt primes (map (λ (p) (factor x p)) primes))))

(apply * (prime-factorization 5)) ; => 120

I'm not exactly sure of how swing exactly works, and the rec-fact thing seems to be there to make it run in logarithmic time. The bit-count thing is a mystery too.
I'm still trying to understand it too.

pump

    ▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄
    █     ___                        █
    █    //  7                       █
    █   (_,_/\       WORSHIP THIS    █
    █    \    \   YOUR THROBBING GOD █
    █     \    \                     █
    █     _\    \__                  █
    █    (   \     )                 █
    █     \___\___/                  █
    █                                █
    ▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀

    ▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄▄
    █     ___                        █
    █    //  7                       █
    █   (_,_/\       WORSHIP THIS    █
    █    \    \   YOUR THROBBING GOD █
    █     \    \                     █
    █     _\    \__                  █
    █    (   \     )                 █
    █     \___\___/                  █
    █                                █
    ▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀▀