W związku z ostatnim uderzeniem w Python , oto próba pokazania mocnych stron Pythona. Twoim wyzwaniem jest napisanie programu, który oblicza silnię tak dużej liczby, jak to możliwe, w ciągu 10 sekund.n

Twój wynik będzie (highest n for your program on your machine)/(highest n for my program on your machine)

Zasady

• Musisz obliczyć dokładne rozwiązanie liczb całkowitych. Ponieważ silnia byłaby znacznie wyższa niż to, co może zmieścić się w 64-bitowej liczbie całkowitej bez znaku, możesz użyć ciągów, jeśli twój język nie obsługuje dużych liczb całkowitych
• Standardowe luki są zabronione. W szczególności nie można używać żadnych zasobów zewnętrznych.
• Tylko część obliczeniowa (obejmuje czas na wszelkie obejścia wykorzystujące ciągi) dodaje do całkowitego czasu, który powinien być średnio poniżej 10 sekund.
• Tylko programy jednowątkowe.
• Musisz przechowywać dane wyjściowe w formie łatwej do wydrukowania (ponieważ drukowanie wymaga czasu) (patrz mój program poniżej), ciąg znaków, zmienna, tablica znaków itp.

EDYTOWAĆ:

• Twój program musi dawać prawidłowe dane wyjściowe dla wszystkich n:1 <= n <= (your highest n)

EDYCJA 2:

• Nienawidzę mówić tego wprost, ale korzystanie z wbudowanych funkcji czynnikowych twojego języka podlega standardowym lukom http://meta.codegolf.stackexchange.com/a/1078/8766 Przepraszamy Mathematica i Sage

Mój program

from __future__ import print_function
import time


def factorial( n ):
    return reduce( ( lambda x , y : x * y ) , xrange( 1 , n + 1 ) , 1 )

start = time.clock()
answer = factorial( 90000 )
end = time.clock()

print ( answer )
print ( "Time:" , end - start , "sec" )

Najwyższy wynik wygrywa. Dla przypomnienia, mój kod może zarządzać n = 90000w około 9.89sekundach na Pentium 4 3,0 GHz

EDYCJA: Czy każdy może dodać wynik, a nie tylko najwyższą liczbę n . Po prostu najwyższy nnie ma żadnego znaczenia, ponieważ zależy od twojego sprzętu. W przeciwnym razie nie można mieć obiektywnego kryterium wygranej. Odpowiedź ali0sha robi to poprawnie.

Mamy zwycięzcę. Nie zaakceptowałem odpowiedzi java /codegolf//a/26974/8766, ponieważ to rodzaj spódnic blisko http://meta.codegolf.stackexchange.com/a/1080/8766

C ++ z GMP, wynik = 55,55 (10 000 000/180 000)

#include <stdio.h>
#include <stdint.h>
#include <stdlib.h>
#include <vector>
#include <iostream>
#include <queue>
#include <gmpxx.h>

int main(int argc, char *argv[]) {
  uint64_t n = atoi(argv[1]);

  // Iterate through 1..n.  Strip off powers of 2.  Multiply
  // remainders together into <= 64 bit chunks.
  uint64_t twos = 0;
  std::vector<uint64_t> terms;
  uint64_t m = 1;
  for(uint64_t i = 1; i <= n; i++) {
    uint64_t j = __builtin_ctzll(i);
    twos += j;
    uint64_t k = i >> j;
    if(__builtin_clzll(m) + __builtin_clzll(k) >= 64) {
      m *= k;
    } else {
      terms.push_back(m);
      m = k;
    }
  }
  if(m != 1) terms.push_back(m);

  // convert to gmp
  // why isn't there a 64-bit constructor?
  std::queue<mpz_class> gmpterms;
  for(int i = 0; i < terms.size(); i++) {
    mpz_class x = (uint32_t)(terms[i] >> 32);
    x <<= 32;
    x += (uint32_t)terms[i];
    gmpterms.push(x);
  }

  // pop two from the bottom, multiply them, push on the end.
  while(gmpterms.size() > 1) {
    mpz_class a = gmpterms.front();
    gmpterms.pop();
    mpz_class b = gmpterms.front();
    gmpterms.pop();
    gmpterms.push(a * b);
  }

  mpz_class r = gmpterms.front();
  r <<= twos;
  //std::cout << r << std::endl;
}

Python 2.7

42,575 = (6,812 000/160 000) ok

Kod:

import gmpy2

def fac1(n):
    m=lambda(L):([]if len(L)%2==0 else[L.pop()])+map(lambda(l):l[0]*l[1],zip(L[1::2],L[-2::-2]))
    L=map(gmpy2.mpz,xrange(1,n+1))
    Number = (len(L)-1).bit_length()
    while Number:Number-=1;L=m(L)
    return L[0]

def fac2(n):
    global E; E=0
    def f(i):
        global E; E+=i//2
        return[]if i==1 else f(i//2)+range(3,i,2)+[[1,i][i%2]]
    m=lambda(L):([]if len(L)%2==0 else[L.pop()])+map(lambda(l):l[0]*l[1],zip(L[1::2],L[-2::-2]))
    L=map(gmpy2.mpz,f(n))
    N=(len(L)-1).bit_length()
    while N: N-=1;L=m(L)
    return L[0]<<E

Test:

import time

start = time.time()
baseline(160000)
print time.time()-start

start = time.time()
fac1(6811000)
print time.time()-start

start = time.time()
fac2(6812000)
print time.time()-start

start = time.time()
gmpy2.fac(26000000)
print time.time()-start

Wynik:

10.0069999695
10.0729999542
10.0360000134
9.98699998856

Jak to działa:

Większe zwielokrotnienia zajmują więcej czasu, dlatego chcemy wykonać jak najwięcej małych zwielokrotnień. Jest to szczególnie prawdziwe w Pythonie, gdzie dla liczb mniejszych niż 2^64używamy arytmetyki sprzętowej, a powyżej używamy oprogramowania. Tak więc, m(L)zaczynamy od listy L; jeśli ma nieparzystą długość, usuwamy jedną liczbę z rozważań, aby była jeszcze większa. Następnie mnożymy element 1przez element -2, element za 3pomocą -4itd., Aby

m([1,2,3,4,5,6,7,8]) = [2*7, 4*5, 6*3, 8*1] = [14, 20, 18, 8]
m([10,12,6]) = [360,112]
m([120,6]) = [40320]

Takie podejście zapewnia, że ​​używamy arytmetyki sprzętowej tak długo, jak to możliwe, po czym przełączamy się na wydajną bibliotekę arytmetyczną gmc.

W fac2, przyjmujemy również bardziej klasyczne podejście dziel i zwyciężaj, w którym dzielimy każdą wielokrotność 2 i przesuwamy je pod koniec bitów, aby uzyskać trywialny wzrost wydajności. Umieściłem go tutaj, ponieważ zwykle jest o około 0,5% szybszy niż fac1.

Wersja w golfa fac1(bo mogę), 220B

import gmpy2
def f(n):
    m=lambda(L):([]if len(L)%2==0 else[L.pop()])+map(lambda(l):l[0]*l[1],zip(L[1::2],L[-2::-2]))
    L=map(gmpy2.mpz,xrange(1,n+1));N=(len(L)-1).bit_length()
    while N:N-=1;L=m(L)
return L[0]

Java - 125,15 (21 400 000 / 171,000)

Również bezwstydnie skopiowane z repozytorium Github Petera Luschnego (dzięki @ semi-extrinsic) i licencjonowane na licencji MIT, wykorzystuje algorytm „pierwiastkowania zagnieżdżonego do kwadratu”, jak zaproponowali Albert Schönhage i in. (zgodnie ze stroną opisu algorytmów czynnikowych Luschnego ).

Lekko dostosowałem algorytm, aby używał BigInteger Javy i nie używał tabeli odnośników dla n <20.

Skompilowany z gcj, który używa GMP do implementacji BigInteger i działał na Linuksie 3.12.4 (Gentoo), na Core i7 4700MQ przy 2,40 GHz

import java.math.BigInteger;

public class PrimeSieveFactorialSchoenhage {

    private static int[] primeList, multiList;

    public static BigInteger factorial(int n) {
        int log2n = 31 - Integer.numberOfLeadingZeros(n);
        int piN = log2n < 2 ? 1 : 2 + (15 * n) / (8 * (log2n - 1));

        primeList = new int[piN];
        multiList = new int[piN];

        int len = primeFactors(n);
        return nestedSquare(len).shiftLeft(n - Integer.bitCount(n));
    }

    private static BigInteger nestedSquare(int len) {
        if (len == 0) {
            return BigInteger.ONE;
        }

        int i = 0, mult = multiList[0];

        while (mult > 1) {
            if ((mult & 1) == 1) { // is mult odd ?
                primeList[len++] = primeList[i];
            }

            multiList[i++] = mult / 2;
            mult = multiList[i];
        }
        BigInteger ns = nestedSquare(i);
        if (len <= i) {
            return ns.multiply(ns);
        }

        return product(primeList, i, len - i).multiply(ns.multiply(ns));
    }

    private static BigInteger product(int[] a, int start, int length) {
        if (length == 0) {
            return BigInteger.ONE;
        }

        int len = (length + 1) / 2;
        long[] b = new long[len];

        int i, j, k;

        for (k = 0, i = start, j = start + length - 1; i < j; i++, k++, j--) {
            b[k] = a[i] * (long) a[j];
        }

        if (i == j) {
            b[k++] = a[j];
        }

        return recProduct(b, 0, k - 1);
    }

    private static BigInteger recProduct(long[] s, int n, int m) {
        if (n > m) {
            return BigInteger.ONE;
        }
        if (n == m) {
            return BigInteger.valueOf(s[n]);
        }
        int k = (n + m) >> 1;
        return recProduct(s, n, k).multiply(recProduct(s, k + 1, m));
    }

    private static int primeFactors(int n) {
        int[] primes = new int[n < 17 ? 6 : (int) Math.floor(n / (Math.log(n) - 1.5))];
        int numPrimes = makePrimeList(n, primes);

        int maxBound = n / 2, count = 0;

        int start = indexOf(primes, 2, 0, numPrimes - 1);
        int end = indexOf(primes, n, start, numPrimes);

        for (int i = start; i < end; i++) {
            int prime = primes[i];
            int m = prime > maxBound ? 1 : 0;

            if (prime <= maxBound) {
                int q = n;
                while (q >= prime) {
                    m += q /= prime;
                }
            }

            primeList[count] = prime;
            multiList[count++] = m;
        }
        return count;
    }

    private static int indexOf(final int[] data, int value, int low, int high) {
        while (low < high) {
            int mid = (low + high) >>> 1;

            if (data[mid] < value) {
                low = mid + 1;
            } else {
                high = mid;
            }
        }

        if (low >= data.length) {
            return low;
        }

        if (data[low] == value) {
            low++;
        }

        return low;
    }

    private static int makePrimeList(int n, int[] prime) {
        boolean[] composite = new boolean[n / 3];

        sieveOfEratosthenes(composite);

        boolean toggle = false;
        int p = 5, i = 0, j = 2;

        prime[0] = 2;
        prime[1] = 3;

        while (p <= n) {
            if (!composite[i++]) {
                prime[j++] = p;
            }
            // -- never mind, it's ok.
            p += (toggle = !toggle) ? 2 : 4;
        }

        return j; // number of primes
    }

    private static void sieveOfEratosthenes(final boolean[] composite) {
        int d1 = 8;
        int d2 = 8;
        int p1 = 3;
        int p2 = 7;
        int s1 = 7;
        int s2 = 3;
        int n = 0;
        int len = composite.length;
        boolean toggle = false;

        while (s1 < len) { // -- scan sieve
            if (!composite[n++]) { // -- if a prime is found, cancel its multiples
                int inc = p1 + p2;

                for (int k = s1; k < len; k += inc) {
                    composite[k] = true;
                }

                for (int k = s1 + s2; k < len; k += inc) {
                    composite[k] = true;
                }
            }

            if (toggle = !toggle) { // Never mind, it's ok.
                s1 += d2;
                d1 += 16;
                p1 += 2;
                p2 += 2;
                s2 = p2;
            } else {
                s1 += d1;
                d2 += 8;
                p1 += 2;
                p2 += 6;
                s2 = p1;
            }
        }
    }

    public static void main(String[] args) {
        int n = Integer.parseInt(args[0]);
        long nanos = System.nanoTime();
        BigInteger fact = factorial(n);
        nanos = System.nanoTime() - nanos;
        // Commented out because it takes ages to print
        //System.out.println(fact);
        System.out.println(nanos / 1e9);
    }
}

Python 3, n = 100000

Wystarczyła prosta zmiana algorytmu, aby podnieść kod przykładowy o 10000.

import time

def factorial(n):
    result = 1
    while n > 0:
        result *= n
        n = n - 1
    return result

start = time.clock()
answer = factorial(100000)
end = time.clock()

print(answer)
print("Time:", end - start, "sec")

Oczywiście nie jest to najbardziej kreatywna odpowiedź, ale tak naprawdę jest tylko jeden sposób, aby zrobić czynnik.

Perl + C, n = około 3 milionów

Tutaj używam biblioteki Math :: BigInt :: GMP dostępnej na CPAN, która zapewnia ogromne przyspieszenie podstawowych obiektów Perla Math :: BigInt.

use v5.14;
use Time::HiRes 'time';
use Math::BigInt only => 'GMP';

sub factorial { Math::BigInt::->new(@_)->bfac }

my $start  = time;
my $answer = factorial( 3_000_000 );
my $end    = time;

say $answer;
say "Time: ", $end - $start, " sec";

Pamiętaj, że mój komputer jest prawdopodobnie nieco wolniejszy od twojego. Używając twojego oryginalnego skryptu Python, mogę obliczyć tylko factorial(40000)w 10 sekund; factorial(90000)zajmuje dużo dłużej. (Wcisnąłem Ctrl + C po minucie.) Na twoim sprzęcie, używając Math :: BigInt :: GMP, możesz być w stanie obliczyć silnię 5 milionów lub więcej w mniej niż 10 sekund.

Jedną z rzeczy, które możesz zauważyć, jest to, że chociaż silnia jest obliczana niewiarygodnie szybko, wydrukowanie wyniku jest bardzo wolne, zajmuje około trzy razy więcej niż pierwotne obliczenie. Wynika to z tego, że GMP wewnętrznie używa reprezentacji binarnej zamiast dziesiętnej, a wydrukowanie go wymaga konwersji binarnej na dziesiętną.

Python 2.7
5,94 = 1'200'000 / 202'000

def fast_fac(n):
    def prod(start, fin):
            if fin - start <= 50:
                    return reduce(lambda x,y: x*y, xrange(start, fin+1), 1)
            else:
                    mid = (start+fin) / 2
                    return prod(start, mid) * prod(mid+1, fin)
    return prod(1, n)

Wykorzystuje względną łatwość mnożenia wielu grup małych liczb, a następnie mnożenia ich w porównaniu do dużej liczby mnożeń z udziałem dużej liczby.

C #: 0,48 (77 000/160 000)

Nie jestem z tego zadowolony.

Czy C # jest tak wolny?

Ale i tak jest mój wpis.

static void Main(string[] args)
    {
        Console.WriteLine("Enter N for fatorial:");
        int n = Convert.ToInt32(Console.ReadLine());

        Stopwatch s = Stopwatch.StartNew();


        BigInteger result = 1;
        while (0 <-- n) result *= n;

        s.Stop();

        Console.WriteLine("Output: {0} ", result);

        Console.WriteLine("Completed in {0}", s.Elapsed);

    }

Gdy n = 77000 trzeba 00:00:09:8708952obliczyć.

Pracuję w trybie Release, poza Visual Studio, używając Core i3-2330M @ 2.2GHz.

Edycja: Ponieważ nie robię niczego inteligentnego, akceptuję ten wynik. Być może .NET Framework 4.5 dodaje trochę narzutu (lub BigInteger nie jest tak szybki).

bc, wynik = 0,19

Co do cholery, oto mój pretendent do „Ile możesz powoli pomnożyć?”

bc jest „językiem kalkulatora o dowolnej precyzji”, ale niestety raczej powolnym:

n=read()
for(f=i=1;i<=n;i++)f*=i
f
quit

W około 10 sekund na moim MacBooku Pro z połowy 2012 roku (2,3 GHz Intel Core i7) odpowiedź na python referencyjny może obliczyć 122000 !, ale ten skrypt bc może obliczyć tylko 23600 !.

Z drugiej strony 10000! zajmuje 1.5s ze skryptem referencyjnym Pythona, ale skrypt bc zajmuje 50s.

O jej.

Bash: wynik = 0,001206 (181/150000)

Ukradłem funkcje matematyczne Rosettacode - długie mnożenie, którego nie analizowałem ani nie próbowałem zoptymalizować.
Możesz zmienić algorytm lub wypróbować inną metodę podziału ciągów.

#!/bin/bash


add() { # arbitrary-precision addition
  if (( ${#1} < ${#2} )); then
    local a="$2" b="$1" sum= carry=0
  else
    local a="$1" b="$2" sum= carry=0
  fi

  while (( ${#a} )); do
    local -i d1="${a##${a%?}}" d2="10#0${b##${b%?}}" s=carry+d1+d2
    sum="${s##${s%?}}$sum"
    carry="10#0${s%?}"
    a="${a%?}" b="${b%?}"
  done
  echo "$sum"
}

multiply() { # arbitrary-precision multiplication
  if (( ${#1} < ${#2} )); then
    local a="$2" b="$1" product=0
  else
    local a="$1" b="$2" product=0
  fi

  local zeroes=
  while (( ${#b} )); do
    local m1="$a"
    local m2="${b##${b%?}}"
    local partial=$zeroes 
    local -i carry=0
    while (( ${#m1} )); do 
      local -i d="${m1##${m1%?}}"
      m1="${m1%?}"
      local -i p=d*m2+carry
      partial="${p##${p%?}}$partial"
      carry="10#0${p%?}"
    done
    partial="${carry#0}$partial"
    product="$(add "$product" "$partial")"
    zeroes=0$zeroes
    b="${b%?}"
  done
  echo "$product"
}

# 'timerun' function
trap 'echo $((i -1)) $f; exit'  USR1  
(sleep 9.9; kill -USR1 $$)&

declare -i i 
f=1
for ((i=1; i< 10000 ; i++ ))   # 10000 is verry optimistic
do
    f=$(multiply $f $i)
done 

Python 3, Advanced Algo autor: Peter Luschny: 8,25x (1 280 000/155 000)

Bezwstydnie skopiowane od Petera Luschnego,
http://www.luschny.de/math/factorial/FastFactorialFunctions.htm ,
który udostępnia ten kod na licencji „Creative Commons Uznanie autorstwa-Na tych samych warunkach 3.0”.

Jest to właściwie dość zaawansowany algorytm, wykorzystujący coś, co nazywa się „zmiennym silnikiem” i listą liczb pierwszych. Podejrzewam, że mogłoby być jeszcze szybciej, gdyby podobało się wielu innym odpowiedziom i wykonał większość multiplikacji z 32-bitowymi liczbami całkowitymi.

#! /usr/bin/python3
import time
import bisect 

def Primes(n) : 
  primes = [2, 3] 
  lim, tog = n // 3, False 
  composite = [False for i in range(lim)] 

  d1 = 8; d2 = 8; p1 = 3; p2 = 7; s = 7; s2 = 3; m = -1 

  while s < lim :             # --  scan the sieve 
      m += 1                  # --  if a prime is found 
      if not composite[m] :   # --  cancel its multiples 
          inc = p1 + p2 
          for k in range(s,      lim, inc) : composite[k] = True 
          for k in range(s + s2, lim, inc) : composite[k] = True 

          tog = not tog 
          if tog: s += d2; d1 += 16; p1 += 2; p2 += 2; s2 = p2 
          else:   s += d1; d2 +=  8; p1 += 2; p2 += 6; s2 = p1 

  k, p, tog = 0, 5, False 
  while p <= n : 
      if not composite[k] : primes.append(p) 
      k += 1; 
      tog = not tog 
      p += 2 if tog else 4 

  return primes 

def isqrt(x): 
  ''' 
  Writing your own square root function
  ''' 
  if x < 0: raise ValueError('square root not defined for negative numbers') 
  n = int(x) 
  if n == 0: return 0 
  a, b = divmod(n.bit_length(), 2) 
  x = 2**(a + b) 
  while True: 
      y = (x + n // x) // 2 
      if y >= x: return x 
      x = y 

def product(s, n, m): 
  if n > m: return 1 
  if n == m: return s[n] 
  k = (n + m) // 2 
  return product(s, n, k) * product(s, k + 1, m) 

def factorialPS(n): 

  small_swing = [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] 

  def swing(m, primes): 
      if m < 33: return small_swing[m] 

      s = bisect.bisect_left(primes, 1 + isqrt(m)) 
      d = bisect.bisect_left(primes, 1 + m // 3) 
      e = bisect.bisect_left(primes, 1 + m // 2) 
      g = bisect.bisect_left(primes, 1 + m) 

      factors = primes[e:g] 
      factors += filter(lambda x: (m // x) & 1 == 1, primes[s:d]) 
      for prime in primes[1:s]:   
          p, q = 1, m 
          while True: 
              q //= prime 
              if q == 0: break 
              if q & 1 == 1: 
                  p *= prime 
          if p > 1: factors.append(p) 

      return product(factors, 0, len(factors) - 1) 

  def odd_factorial(n, primes): 
      if n < 2: return 1 
      return (odd_factorial(n // 2, primes)**2) * swing(n, primes) 

  def eval(n): 
      if n < 0: 
          raise ValueError('factorial not defined for negative numbers') 

      if n == 0: return 1 
      if n < 20: return product(range(2, n + 1), 0, n-2) 

      N, bits = n, n 
      while N != 0: 
          bits -= N & 1 
          N >>= 1 

      primes = Primes(n) 
      return odd_factorial(n, primes) * 2**bits 

  return eval(n)

start = time.time()
answer = factorialPS(1280000) 
print(time.time()-start)

Java - 10.9

n = 885000

Mergesort-y.

import java.math.BigInteger;

public class Factorials {

    public static BigInteger fac;

    public static BigInteger two = BigInteger.valueOf(2);

    static BigInteger mul(BigInteger start, BigInteger end) {
        if(start.equals(end)) {
            return start;
        } else {
            BigInteger mid = start.add(end.subtract(start).divide(Factorials.two));
            return Factorials.mul(start, mid).multiply(Factorials.mul(mid.add(BigInteger.ONE), end));
        }
    }

    public static void main(String[] args) {
        Factorials.fac = BigInteger.valueOf(Integer.parseInt(args[0]));
        long t = System.nanoTime();
        BigInteger result = mul(BigInteger.ONE, fac);
        t = System.nanoTime() - t;
        System.out.print(String.valueOf(((float) t) / 1000000000)); //result.toString()+" @ "+
    }
}

BigIntegers są powolne.

Zalecenia dotyczące bibliotek liczb całkowitych Java o dowolnej precyzji? : P

C ++ (specyficzny dla x86_64) - 3.0 (390000/130000)

(łatwo przenośny na x86-32, przenoszenie na inne architektury oznacza znaczną utratę prędkości)

Oto moja własna mikro-implementacja długiej arytmetyki.
Samo obliczenie zajmuje 10 sekund, a mimo że wydruk jest łatwy do wydrukowania (zobacz operator<<przeciążenie), wydrukowanie go zajmuje więcej czasu.

#include <vector>
#include <iostream>
#include <stdint.h>
#include <ctime>

typedef uint64_t digit;
typedef std::vector<digit> number;

std::ostream &operator<<(std::ostream &s, const number &x)
{
    std::vector<char> o;
    size_t size = x.size() * 21;
    o.resize(size);
    size_t lud = 0;
    for(number::const_reverse_iterator i = x.rbegin(), end = x.rend(); i != end; i++)
    {
        digit carry = 0;
        int j;
        for(j = 0; j <= lud || carry; j++)
        {
            digit r = o[j] * (1LL << 32) + carry;
            o[j] = r % 10;
            carry = r / 10;
        }
        lud = j;
        carry = 0;
        for(j = 0; j <= lud || carry; j++)
        {
            digit r = o[j] * (1LL << 32) + carry;
            o[j] = r % 10;
            carry = r / 10;
        }
        lud = j;
        carry = *i;
        for(j = 0; carry; j++)
        {
            digit r = o[j] + (carry % 10);
            carry /= 10;
            carry += r / 10;
            o[j] = r % 10;
        }
        if(j > lud)
            lud = j;
    }
    for(int j = lud; j--;)
        s.put(o[j] + '0');
    return s;
}

inline uint64_t dmul(uint64_t x, uint64_t y, uint64_t &carry)
{
    asm("mulq %2" : "+a"(x), "=d"(carry) : "r"(y));
    return x;
}
inline digit dadd(digit x, digit y, digit &carry)
{
    asm("movq $0, %1; addq %2, %0; adcq %1, %1" : "+r"(x), "=r"(carry), "+r"(y));
    return x;
}

void multiply(number &x, digit y)
{
    x.resize(x.size() + 2);
    digit carry = 0;
    for(number::iterator i = x.begin(), end = x.end(); i != end; i++)
    {
        digit nc, res = dmul(*i, y, nc);
        *i = dadd(res, carry, carry);
        carry += nc;
    }
    size_t sz = x.size();
    for(number::const_reverse_iterator i = x.rbegin(), end = x.rend(); i != end; i++)
    {
        if(*i)
            break;
        sz--;
    }
    x.resize(sz);
}

int main()
{
    const int r = 390000;
    clock_t start = clock();
    number n;
    digit mult = 1;
    n.push_back(1);
    for(digit a = 2; a <= r; a++)
    {
        digit carry, m = dmul(mult, a, carry);
        if(carry)
        {
            multiply(n, mult);
            mult = a;
        }
        else
            mult = m;
    }
    multiply(n, mult);
    std::cout << "Took: " << (clock() - start)/((double)CLOCKS_PER_SEC) << std::endl;
    std::cout << n << std::endl;
}

Python 3: 280000/168000

Czas działania programu: pomiędzy 9.87585953253i 10.3046453994. Czas działania mojego programu: około 10.35296977897559.

import time

def factorial(n):
    f = 1
    while n > 1:
        hn = n >> 1
        f = f * 2**hn * double_factorial(n) #dfl[hn + (n & 1) - 1]
        n = hn
    return f
def double_factorial(n):
    #dfl = [1]
    p = 1
    l = 3
    mh = n
    while l <= n:
        p *= l
        l += 2
        #dfl.append(p)
    return p

start = time.clock()
factorial(280000)
end = time.clock()

print(end - start)

Przeczytałem tę odpowiedź na cs.SE i postanowiłem spróbować ją zaimplementować w Pythonie. Jednak przypadkowo odkryłem, że n! = (⌊n / 2⌋)! * 2**(⌊n / 2⌋) * n!!(uwaga: !!jest podwójnym silnikiem ). Przekształciłem to w formę nierekurencyjną.

Komentarze pokazują moją próbę uniknięcia ponownego obliczenia podwójnej silni, ale odkryłem, że przechowywanie każdej wartości było zbyt kosztowne dla pamięci, co spowodowało, że mój komputer działał jeszcze wolniej. Mogę to poprawić, przechowując tylko to, co jest potrzebne.

O dziwo zaimplementowałem naiwne proste mnożenie w Pythonie 3 i robi to lepiej niż twój program: n = 169000w 10 sekund .:

def factorial(n):
    p=1
    for i in range(n):
        p*=i+1
    return p

Ruby 2.1

wynik = 1,80 = 176_000 / 98_000

EDYCJA: poprawiona z 1.35 = 132_000 / 98_000

Wziąłem pomysły z algorytmu czynnikowego GMP . Ten program używa standardowej biblioteki do generowania liczb pierwszych. Ruby to zły wybór, ponieważ mnożenie wydaje się wolniejsze w Rubim niż w Pythonie.

1. Mój program w Ruby 2.1: wynik = 1,80 = 176_000 / 98_000
2. Trywialny algorytm w Pythonie 2.7: wynik = 1 = 98_000 / 98_000
3. Trywialny algorytm w Ruby 2.1: wynik = 0,878 = 86_000 / 98_000

Tak, mój plik binarny ruby 2.1.0p0 (2013-12-25 revision 44422) [x86_64-openbsd]linków przeciwko GMP. Ruby 2.1 dodał funkcję do korzystania z GMP do dużego mnożenia, ale nadal wydaje się wolniejszy niż Python 2.7.

require 'benchmark'
require 'optparse'
require 'prime'

def factorial(n)
  # calculate primes up to n, drop the 2
  @odd_primes = Prime.each(n).drop(1)

  # count prime factors of factorial(n)
  @factors = Hash.new(0)
  factorial_recurse(n)

  shift = @factors.delete(2) || 0
  @factors.inject(1) {|product, (base, exp)|
    product * base**exp
  } << shift
end

def factorial_recurse(n)
  return if n < 2

  # collect prime factors of 2 * 4 * 6 * .. * n
  #  = (2 * 2 * 2 * .. * 2) * (1 * 2 * 3 * .. * exp)
  #  = 2**exp * factorial(exp) where exp = floor(n/2)
  exp = n >> 1
  factorial_recurse(exp)
  @factors[2] += exp

  # collect prime factors 3 * 5 * 7 * ... * n
  for prime in @odd_primes
    break if prime > n
    exp = 0
    # count occurences of prime, prime**2, prime**3, .. n
    prime_power = prime
    until prime_power > n
      # floor(n / prime_power) occurences in 1 * 2 * .. * n,
      # but only ceil(count / 2) occurences in 3 * 5 * .. * n
      @factors[prime] += (n / prime_power + 1) >> 1
      prime_power *= prime
    end
  end
end

# usage: factorial.rb [-ct] [number]
cflag = tflag = false
OptionParser.new {|opts|
  opts.on('-c', 'Check for bugs') { cflag = true }
  opts.on('-t', 'Use trivial algorithm') { tflag = true }
  opts.parse!
}
$*[1] and fail 'too many arguments'
n = Integer($*[0] || 176_000)

if cflag
  factorial(n) == (1..n).reduce(1, :*) or
    fail "bad program: factorial(#{n}) is wrong"
  puts "ok"
  exit
end

# measure processor time to calculate factorial
f = nil
if tflag
  time = Benchmark.measure { f = (1..n).reduce(1, :*) }
else
  time = Benchmark.measure { f = factorial(n) }
end
puts f
puts "Time #{time.total} sec"

Julia - wynik = 15,194

Wykorzystując dokładnie to samo podejście, co w programie referencyjnym ...

f(n)=reduce(*,1:big(n))

Wykorzystuje więc redukcję, podstawową operację mnożenia binarnego i zakres (w tym przypadku użycie dużej (n) do wymuszenia wykonania obliczeń w BigInt zamiast Int64). Z tego rozumiem

julia> @time K=f(2340000);
elapsed time: 9.991324093 seconds (814552840 bytes allocated)

Na moim komputerze, z programem referencyjnym działającym z wejściem 154000, otrzymuję Time: 10.041181 secdane wyjściowe (uruchom przy użyciu python ./test.py, gdzie test.py to nazwa pliku zawierającego kod referencyjny)

tcl, 13757

Moja odpowiedź to sprawdzić limity tcl.

Pierwszy wiersz służy tylko do ustawienia parametru wejściowego:

set n 13757

Pozostałe to sam algorytm:

set r 2
for {set i 3} {$i <= $n} {incr i} {set r [expr {$r*$i}]}   
puts $r

Przetestowałem swój kod na http://rextester.com/live/WEL36956 ; Jeśli powiększę n, dostanę SIGKILL; may n może stać się większy na lokalnym tłumaczu tcl, którego nie mam.