Merge pull request #700 from yledu/shootout

Shootout pidigits code, new shootout directory in examples and moved existing shootout code to it

Author: pao

examples/shootout/README

@@ -0,0 +1,9 @@
+This directory contains the Julia version of the "The
+Computer Language Benchmarks Game":
+http://shootout.alioth.debian.org/
+
+The source code for all the benchmarks are available there:
+http://alioth.debian.org/scm/viewvc.php/shootout/bench/?root=shootout
+
+See specific Julia discussion here:
+https://github.com/JuliaLang/julia/issues/660

examples/fasta.jl renamed to examples/shootout/fasta.jl

@@ -0,0 +1,77 @@
+load("timing.jl")
+load("bigint.jl")
+
+function pidigits(N::Int, printOut::Bool)
+    """
+    See http://shootout.alioth.debian.org/u64q/performance.php?test=pidigits#about
+
+    Transliterated from Mario Pernici Python's program
+
+
+    INPUT:
+
+    - N -- a positive integer giving the number of digits of pi to be computed
+
+    - printOut -- a boolean specifying if we want intermediate printouts of digits in packets of 10
+
+    OUTPUT:
+
+    - returns the last ten digits anyway
+    
+    - prints all the digits in packets of 10 iff printOut == true
+
+    """
+
+    i = k = ns = 0
+    k1 = 1
+    n,a,d,t,u = map(BigInt,(1,0,1,0,0))
+
+    while true
+        k += 1
+        t = lshift(n,uint(1))
+        n *= k
+        a += t
+        k1 += 2
+        a *= k1
+        d *= k1
+
+        if a >= n
+            t,u = divmod(n*3 +a, d)
+            u += n
+            if d > u
+                ns = ns*10 + t
+                i += 1
+                if mod(i,10) == 0
+                    if printOut
+                        show(ns)
+                        printf("\t:%d\n", i)
+                    end
+                    if i >= N
+                        return ns
+                    end
+                    ns = 0
+                end
+                a -= d*t
+                a *= 10
+                n *= 10
+
+            end
+        end
+    end
+end
+
+function pidigits(N::Int)
+    pidigits(N,false)
+end
+
+
+@assert pidigits(1000) == 9216420198
+
+if length(ARGS)==1
+    N = int(ARGS[1])
+else
+    N = 1000
+end
+@timeit pidigits(N) "pidigits"
+
+

examples/shootout/pidigits.jl

@@ -0,0 +1,12 @@
+macro timeit(ex,name)
+    t, i = gensym(2)
+    quote
+        $t = Inf
+        for $i=1:5
+            $t = min($t, @elapsed $ex)
+        end
+        printf("julia, %s: best of 5 took %.1f ms\n", $name, $t*1000)
+    end
+end
+
+

examples/shootout/timing.jl

@@ -45,6 +45,10 @@ extern void _jl_mpz_div(mpz_t* rop, mpz_t* op1, mpz_t* op2) {
   mpz_fdiv_q(*rop, *op1, *op2);
 }
 
+extern void _jl_mpz_divmod(mpz_t* rop1, mpz_t* rop2, mpz_t* op1, mpz_t* op2) {
+  mpz_divmod(*rop1, *rop2, *op1, *op2);
+}
+  
 extern void _jl_mpz_rem(mpz_t* rop, mpz_t* op1, mpz_t* op2) {
   mpz_fdiv_r(*rop, *op1, *op2);
 }
@@ -61,6 +65,10 @@ extern int _jl_mpz_cmp(mpz_t* op1, mpz_t* op2) {
   return mpz_cmp(*op1, *op2);
 }
 
+extern void _jl_mpz_lshift(mpz_t* rop, mpz_t* base, unsigned long int count) {
+  mpz_mul_2exp(*rop, *base, count);
+}
+
 extern void _jl_mpz_pow_ui(mpz_t* rop, mpz_t* base, unsigned long int exp) {
   mpz_pow_ui(*rop, *base, exp);
 }

external/gmp_wrapper.c

@@ -1,4 +1,4 @@
-_jl_libgmp_wrapper = dlopen("libgmp_wrapper")
+_jl_libgmp_wrapper = dlopen("libgmp_wrapper.so")
 
 type BigInt <: Integer
 	mpz::Ptr{Void}
@@ -57,6 +57,12 @@ function -(x::BigInt)
 	BigInt(z)
 end
 
+function lshift(x::BigInt, c::Uint)
+    z= _jl_bigint_init()
+    ccall(dlsym(_jl_libgmp_wrapper, :_jl_mpz_lshift), Void, (Ptr{Void}, Ptr{Void}, Uint), z, x.mpz, c)
+    BigInt(z)
+end
+
 function -(x::BigInt, y::BigInt)
 	z= _jl_bigint_init()
 	ccall(dlsym(_jl_libgmp_wrapper, :_jl_mpz_sub), Void, (Ptr{Void}, Ptr{Void}, Ptr{Void}),z,x.mpz,y.mpz)
@@ -75,6 +81,13 @@ function div (x::BigInt, y::BigInt)
 	BigInt(z)
 end
 
+function divmod(x::BigInt, y::BigInt)
+	z1= _jl_bigint_init()
+	z2= _jl_bigint_init()
+	ccall(dlsym(_jl_libgmp_wrapper, :_jl_mpz_divmod), Void, (Ptr{Void}, Ptr{Void}, Ptr{Void}, Ptr{Void}), z1, z2, x.mpz, y.mpz)
+	BigInt(z1),BigInt(z2)
+end
+
 function rem (x::BigInt, y::BigInt)
 	z= _jl_bigint_init()
 	ccall(dlsym(_jl_libgmp_wrapper, :_jl_mpz_rem), Void, (Ptr{Void}, Ptr{Void}, Ptr{Void}),z,x.mpz,y.mpz)