Este es otro desafío sobre los números de Fibonacci.

El objetivo es calcular el ^número 20'000'000 de Fibonacii lo más rápido posible. La salida decimal es aproximadamente 4 MiB grande; Empieza con:

28543982899108793710435526490684533031144309848579

La suma MD5 de la salida es

fa831ff5dd57a830792d8ded4c24c2cb

Tiene que enviar un programa que calcule el número mientras se ejecuta y coloca el resultado stdout. El programa más rápido, medido en mi propia máquina, gana.

Aquí hay algunas reglas adicionales:

• Debe enviar el código fuente y un binario ejecutable en un Linux x64
• El código fuente debe ser más corto que 1 MiB, en caso de ensamblaje también es aceptable si solo el binario es tan pequeño.
• No debe incluir el número que se calculará en su binario, incluso de forma encubierta. El número debe calcularse en tiempo de ejecución.
• Mi computadora tiene dos núcleos; puedes usar paralelismo

Tomé una pequeña implementación de Internet que se ejecuta en aproximadamente 4.5 segundos. No debería ser muy difícil superar esto, suponiendo que tenga un buen algoritmo.

C con GMP, 3.6s

Dioses, pero GMP hace que el código sea feo. Con un truco al estilo Karatsuba, logré reducir a 2 multiplicaciones por paso de duplicación. Ahora que estoy leyendo la solución de FUZxxl, no soy el primero en tener la idea. Tengo un par de trucos más bajo la manga ... tal vez los pruebe más tarde.

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

#define DBL mpz_mul_2exp(u,a,1);mpz_mul_2exp(v,b,1);mpz_add(u,u,b);mpz_sub(v,a,v);mpz_mul(b,u,b);mpz_mul(a,v,a);mpz_add(a,b,a);
#define ADD mpz_add(a,a,b);mpz_swap(a,b);

int main(){
    mpz_t a,b,u,v;
    mpz_init(a);mpz_set_ui(a,0);
    mpz_init(b);mpz_set_ui(b,1);
    mpz_init(u);
    mpz_init(v);

    DBL
    DBL
    DBL ADD
    DBL ADD
    DBL
    DBL
    DBL
    DBL ADD
    DBL
    DBL
    DBL ADD
    DBL
    DBL ADD
    DBL ADD
    DBL
    DBL ADD
    DBL
    DBL
    DBL
    DBL
    DBL
    DBL
    DBL
    DBL /*Comment this line out for F(10M)*/

    mpz_out_str(stdout,10,b);
    printf("\n");
}

Construido con gcc -O3 m.c -o m -lgmp.

Sabio

Hmm, parece suponer que lo más rápido será un programa compilado. ¡No hay binario para ti!

print fibonacci(2000000)

En mi máquina, toma 0,10 segundos de CPU, 0,15 segundos de pared.

editar: cronometrado en la consola, en lugar del cuaderno

Haskell

Este es mi propio intento, aunque no escribí el algoritmo por mí mismo. Lo copié bastante de haskell.org y lo adapté para usarlo Data.Vectorcon su famosa fusión de secuencias:

import Data.Vector as V
import Data.Bits

main :: IO ()
main = print $ fib 20000000

fib :: Int -> Integer
fib n = snd . V.foldl' fib' (1,0) . V.dropWhile not $ V.map (testBit n) $ V.enumFromStepN (s-1) (-1) s
    where
        s = bitSize n
        fib' (f,g) p
            | p         = (f*(f+2*g),ss)
            | otherwise = (ss,g*(2*f-g))
            where ss = f*f+g*g

Esto toma alrededor de 4.5 segundos cuando se compila con GHC 7.0.3 y los siguientes indicadores:

ghc -O3 -fllvm fib.hs

VACA

 MoO moO MoO mOo MOO OOM MMM moO moO
 MMM mOo mOo moO MMM mOo MMM moO moO
 MOO MOo mOo MoO moO moo mOo mOo moo

¡Mugir! (Toma un tiempo. Bebe un poco de leche ...)

Mathematica, interpretada:

First@Timing[Fibonacci[2 10^6]]

Cronometrado:

0.032 secs on my poor man's laptop.

Y, por supuesto, no binario.

Ocaml, 0.856s en mi laptop

Requiere la biblioteca zarith. Usé Big_int pero es lento para perros en comparación con zarith. ¡Tomó 10 minutos con el mismo código! ¡La mayor parte del tiempo la pasé imprimiendo el maldito número (aproximadamente 9½ minutos)!

module M = Map.Make
  (struct
    type t = int
    let compare = compare
   end)

let double b = Z.shift_left b 1
let ( +. ) b1 b2 = Z.add b1 b2
let ( *. ) b1 b2 = Z.mul b1 b2

let cache = ref M.empty 
let rec fib_log n =
  if n = 0
  then Z.zero
  else if n = 1
  then Z.one
  else if n mod 2 = 0
  then
    let f_n_half = fib_log_cached (n/2)
    and f_n_half_minus_one = fib_log_cached (n/2-1)
    in f_n_half *. (f_n_half +. double f_n_half_minus_one)
  else
    let f_n_half = fib_log_cached (n/2)
    and f_n_half_plus_one = fib_log_cached (n/2+1)
    in (f_n_half *. f_n_half) +.
    (f_n_half_plus_one *. f_n_half_plus_one)
and fib_log_cached n =
    try M.find n !cache
    with Not_found ->
      let res = fib_log n
      in cache := M.add n res !cache;
      res

let () =
  let res = fib_log 20_000_000 in
  Z.print res; print_newline ()

¡No puedo creer la diferencia que hizo la biblioteca!

Haskell

En mi sistema, esto funciona casi tan rápido como la respuesta de FUZxxl (~ 18 segundos en lugar de ~ 17 segundos).

main = print $ fst $ fib2 20000000

-- | fib2: Compute (fib n, fib (n+1)).
--
-- Having two adjacent Fibonacci numbers lets us
-- traverse up or down the series efficiently.
fib2 :: Int -> (Integer, Integer)

-- Guard against negative n.
fib2 n | n < 0 = error "fib2: negative index"

-- Start with a few base cases.
fib2 0 = (0, 1)
fib2 1 = (1, 1)
fib2 2 = (1, 2)
fib2 3 = (2, 3)

-- For larger numbers, derive fib2 n from fib2 (n `div` 2)
-- This takes advantage of the following identity:
--
--    fib(n) = fib(k)*fib(n-k-1) + fib(k+1)*fib(n-k)
--             where n > k
--               and k ≥ 0.
--
fib2 n =
    let (a, b) = fib2 (n `div` 2)
     in if even n
        then ((b-a)*a + a*b, a*a + b*b)
        else (a*a + b*b, a*b + b*(a+b))

C, algoritmo ingenuo

Tenía curiosidad, y no había usado gmp antes ... así que:

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

int main(int argc, char *argv[]){
    int n = (argc>1)?atoi(argv[1]):0;

    mpz_t temp,prev,result;
    mpz_init(temp);
    mpz_init_set_ui(prev, 0);
    mpz_init_set_ui(result, 1);

    for(int i = 2; i <= n; i++) {
        mpz_add(temp, result, prev);
        mpz_swap(temp, result);
        mpz_swap(temp, prev);
    }

    printf("fib(%d) = %s\n", n, mpz_get_str (NULL, 10, result));

    return 0;
}

fib (1 millón) toma alrededor de 7 segundos ... por lo que este algoritmo no ganará la carrera.

Implementé el método de multiplicación matricial (de sicp, http://sicp.org.ua/sicp/Exercise1-19 ) en SBCL, pero tarda unos 30 segundos en terminar. Lo porté a C usando GMP, y devuelve el resultado correcto en aproximadamente 1.36 segundos en mi máquina. Es casi tan rápido como la respuesta de Boothby.

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

int main()
{
  int n = 20000000;

  mpz_t a, b, p, q, psq, qsq, twopq, bq, aq, ap, bp;
  int count = n;

  mpz_init_set_si(a, 1);
  mpz_init_set_si(b, 0);
  mpz_init_set_si(p, 0);
  mpz_init_set_si(q, 1);
  mpz_init(psq);
  mpz_init(qsq);
  mpz_init(twopq);
  mpz_init(bq);
  mpz_init(aq);
  mpz_init(ap);
  mpz_init(bp);

  while(count > 0)
    {
      if ((count % 2) == 0)
        {
          mpz_mul(psq, p, p);
          mpz_mul(qsq, q, q);
          mpz_mul(twopq, p, q);
          mpz_mul_si(twopq, twopq, 2);

          mpz_add(p, psq, qsq);    // p -> (p * p) + (q * q)
          mpz_add(q, twopq, qsq);  // q -> (2 * p * q) + (q * q) 
          count/=2;
        }

      else
       {
          mpz_mul(bq, b, q);
          mpz_mul(aq, a, q);
          mpz_mul(ap, a, p);
          mpz_mul(bp, b, p);

          mpz_add(a, bq, aq);      // a -> (b * q) + (a * q)
          mpz_add(a, a, ap);       //              + (a * p)

          mpz_add(b, bp, aq);      // b -> (b * p) + (a * q)

          count--;
       }

    }

  gmp_printf("%Zd\n", b);
  return 0;
}

Java: 8 segundos para calcular, 18 segundos para escribir

public static BigInteger fibonacci1(int n) {
    if (n < 0) explode("non-negative please");
    short charPos = 32;
    boolean[] buf = new boolean[32];
    do {
        buf[--charPos] = (n & 1) == 1;
        n >>>= 1;
    } while (n != 0);
    BigInteger a = BigInteger.ZERO;
    BigInteger b = BigInteger.ONE;
    BigInteger temp;
    do {
        if (buf[charPos++]) {
            temp = b.multiply(b).add(a.multiply(a));
            b = b.multiply(a.shiftLeft(1).add(b));
            a = temp;
        } else {
            temp = b.multiply(b).add(a.multiply(a));
            a = a.multiply(b.shiftLeft(1).subtract(a));
            b = temp;
        }
    } while (charPos < 32);
    return a;
}

public static void main(String[] args) {
    BigInteger f;
    f = fibonacci1(20000000);
    // about 8 seconds
    System.out.println(f.toString());
    // about 18 seconds
}

Vamos

Es vergonzosamente lento. En mi computadora toma un poco menos de 3 minutos. Sin embargo, solo son 120 llamadas recursivas (después de agregar el caché). ¡Tenga en cuenta que esto puede usar mucha memoria (como 1.4 GiB)!

package main

import (
    "math/big"
    "fmt"
)

var cache = make(map[int64] *big.Int)

func fib_log_cache(n int64) *big.Int {
    if res, ok := cache[n]; ok {
        return res
    }
    res := fib_log(n)
    cache[n] = res
    return res
}

func fib_log(n int64) *big.Int {
    if n <= 1 {
        return big.NewInt(n)
    }

    if n % 2 == 0 {
        f_n_half := fib_log_cache(n/2)
        f_n_half_minus_one := fib_log_cache(n/2-1)
        res := new(big.Int).Lsh(f_n_half_minus_one, 1)
        res.Add(f_n_half, res)
        res.Mul(f_n_half, res)
        return res
    }
    f_n_half := fib_log_cache(n/2)
    f_n_half_plus_one := fib_log_cache(n/2+1)
    res := new(big.Int).Mul(f_n_half_plus_one, f_n_half_plus_one)
    tmp := new(big.Int).Mul(f_n_half, f_n_half)
    res.Add(res, tmp)
    return res
}

func main() {
    fmt.Println(fib_log(20000000))
}

pseudocódigo (no sé qué están usando)

product = 1
multiplier = 3 // 3 is fibonacci sequence, but this can be any number, 
      // generating an infinite amount of sequences
y = 28 // the 2^x-1 term, so 2^28-1=1,284,455,535th term
for (int i = 1; int < y; i++) {
  product= sum*multiplier-1
  multiplier= multiplier^2-2
}
multiplier=multiplier-product // 2^28+1 1,284,455,537th 

Le tomó a mi computadora 56 horas hacer esos dos términos. Mi computadora es un poco horrible. Tendré el número en un archivo de texto el 22 de octubre. 1.2 conciertos es un poco grande para compartir en mi conexión.