Какой самый быстрый способ получить значение π? - PullRequest
307 голосов
/ 01 августа 2008

Я ищу самый быстрый способ получить значение π, как личный вызов. Более конкретно, я использую способы, которые не включают использование #define констант, таких как M_PI, или жесткое кодирование числа в.

Программа ниже тестирует различные известные мне способы. Версия inline сборки, теоретически, является самым быстрым вариантом, хотя и явно не переносимым. Я включил его в качестве основы для сравнения с другими версиями. В моих тестах со встроенными модулями версия 4 * atan(1) была самой быстрой в GCC 4.2, потому что она автоматически складывает atan(1) в константу. С указанным -fno-builtin версия atan2(0, -1) самая быстрая.

Вот основная программа тестирования (pitimes.c):

#include <math.h>
#include <stdio.h>
#include <time.h>

#define ITERS 10000000
#define TESTWITH(x) {                                                       \
    diff = 0.0;                                                             \
    time1 = clock();                                                        \
    for (i = 0; i < ITERS; ++i)                                             \
        diff += (x) - M_PI;                                                 \
    time2 = clock();                                                        \
    printf("%s\t=> %e, time => %f\n", #x, diff, diffclock(time2, time1));   \
}

static inline double
diffclock(clock_t time1, clock_t time0)
{
    return (double) (time1 - time0) / CLOCKS_PER_SEC;
}

int
main()
{
    int i;
    clock_t time1, time2;
    double diff;

    /* Warmup. The atan2 case catches GCC's atan folding (which would
     * optimise the ``4 * atan(1) - M_PI'' to a no-op), if -fno-builtin
     * is not used. */
    TESTWITH(4 * atan(1))
    TESTWITH(4 * atan2(1, 1))

#if defined(__GNUC__) && (defined(__i386__) || defined(__amd64__))
    extern double fldpi();
    TESTWITH(fldpi())
#endif

    /* Actual tests start here. */
    TESTWITH(atan2(0, -1))
    TESTWITH(acos(-1))
    TESTWITH(2 * asin(1))
    TESTWITH(4 * atan2(1, 1))
    TESTWITH(4 * atan(1))

    return 0;
}

И встроенные компоненты сборки (fldpi.c), которые будут работать только для систем x86 и x64:

double
fldpi()
{
    double pi;
    asm("fldpi" : "=t" (pi));
    return pi;
}

И скрипт сборки, который собирает все конфигурации, которые я тестирую (build.sh):

#!/bin/sh
gcc -O3 -Wall -c           -m32 -o fldpi-32.o fldpi.c
gcc -O3 -Wall -c           -m64 -o fldpi-64.o fldpi.c

gcc -O3 -Wall -ffast-math  -m32 -o pitimes1-32 pitimes.c fldpi-32.o
gcc -O3 -Wall              -m32 -o pitimes2-32 pitimes.c fldpi-32.o -lm
gcc -O3 -Wall -fno-builtin -m32 -o pitimes3-32 pitimes.c fldpi-32.o -lm
gcc -O3 -Wall -ffast-math  -m64 -o pitimes1-64 pitimes.c fldpi-64.o -lm
gcc -O3 -Wall              -m64 -o pitimes2-64 pitimes.c fldpi-64.o -lm
gcc -O3 -Wall -fno-builtin -m64 -o pitimes3-64 pitimes.c fldpi-64.o -lm

Помимо тестирования между различными флагами компилятора (я также сравнил 32-битные с 64-битными, потому что оптимизации разные), я также попытался изменить порядок тестов. Но, тем не менее, версия atan2(0, -1) по-прежнему выходит на первое место каждый раз.

Ответы [ 23 ]

197 голосов
/ 02 августа 2008

В методе Монте-Карло , как уже упоминалось, применяются некоторые замечательные концепции, но он, очевидно, не самый быстрый, ни в дальнем плане, ни в какой-либо разумной мере. Кроме того, все зависит от того, какую точность вы ищете. Самый быстрый π, который я знаю, это тот, у которого цифры жестко запрограммированы. Глядя на Pi и Pi [PDF] , есть много формул.

Вот метод, который быстро сходится - около 14 цифр за итерацию. PiFast , самое быстрое текущее приложение, использует эту формулу с БПФ. Я просто напишу формулу, так как код прост. Эта формула была почти найдена Рамануджаном и обнаружена Чудновским . Это на самом деле, как он вычислил несколько миллиардов цифр числа, так что это не метод игнорировать. Формула быстро переполнится, и, поскольку мы делим факториалы, было бы целесообразно отложить такие вычисления, чтобы удалить термины.

enter image description here

enter image description here

, где

enter image description here

Ниже приведен алгоритм Брента – Саламина . В Википедии упоминается, что когда a и b"достаточно близки", то (a + b) ² / 4t будет приближением π. Я не уверен, что означает «достаточно близко», но из моих тестов одна итерация получила 2 цифры, две - 7, а три - 15, конечно это с двойными числами, поэтому может возникнуть ошибка, основанная на ее представлении и true расчет может быть более точным.

let pi_2 iters =
    let rec loop_ a b t p i =
        if i = 0 then a,b,t,p
        else
            let a_n = (a +. b) /. 2.0 
            and b_n = sqrt (a*.b)
            and p_n = 2.0 *. p in
            let t_n = t -. (p *. (a -. a_n) *. (a -. a_n)) in
            loop_ a_n b_n t_n p_n (i - 1)
    in 
    let a,b,t,p = loop_ (1.0) (1.0 /. (sqrt 2.0)) (1.0/.4.0) (1.0) iters in
    (a +. b) *. (a +. b) /. (4.0 *. t)

Наконец, как насчет пи-гольфа (800 цифр)? 160 символов!

int a=10000,b,c=2800,d,e,f[2801],g;main(){for(;b-c;)f[b++]=a/5;for(;d=0,g=c*2;c-=14,printf("%.4d",e+d/a),e=d%a)for(b=c;d+=f[b]*a,f[b]=d%--g,d/=g--,--b;d*=b);}
110 голосов
/ 02 сентября 2008

Мне очень нравится эта программа, потому что она аппроксимирует π, глядя на свою собственную область.

IOCCC 1988: westley.c

#define _ -F<00||--F-OO--;
int F=00,OO=00;main(){F_OO();printf("%1.3f\n",4.*-F/OO/OO);}F_OO()
{
            _-_-_-_
       _-_-_-_-_-_-_-_-_
    _-_-_-_-_-_-_-_-_-_-_-_
  _-_-_-_-_-_-_-_-_-_-_-_-_-_
 _-_-_-_-_-_-_-_-_-_-_-_-_-_-_
 _-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
 _-_-_-_-_-_-_-_-_-_-_-_-_-_-_
 _-_-_-_-_-_-_-_-_-_-_-_-_-_-_
  _-_-_-_-_-_-_-_-_-_-_-_-_-_
    _-_-_-_-_-_-_-_-_-_-_-_
        _-_-_-_-_-_-_-_
            _-_-_-_
}
75 голосов
/ 01 августа 2008

Вот общее описание методики вычисления пи, которую я выучил в старшей школе.

Я делюсь этим только потому, что думаю, что это достаточно просто, чтобы каждый мог помнить это, на неопределенный срок, плюс он учит вас концепции методов "Монте-Карло" - статистических методов получения ответов, которые не сразу кажется выводимым через случайные процессы.

Нарисуйте квадрат и вписайте квадрант (одну четверть полукруга) внутри этого квадрата (квадрант с радиусом, равным стороне квадрата, поэтому он заполняет максимально возможную площадь)

Теперь бросьте дротик в квадрат и запишите, где он приземлится, то есть выберите случайную точку в любом месте квадрата. Конечно, он приземлился внутри квадрата, но внутри полукруга? Запишите этот факт.

Повторите этот процесс много раз - и вы обнаружите, что есть отношение числа точек внутри полукруга к общему количеству брошенных, назовите это соотношение x.

Поскольку площадь квадрата равна r, умноженной на r, можно сделать вывод, что площадь полукруга равна x, умноженной на r, умноженной на r (то есть, x раз, умноженных на квадрат). Следовательно, х раз 4 даст вам пи.

Это не быстрый метод для использования. Но это хороший пример метода Монте-Карло. И если вы оглянетесь вокруг, вы можете обнаружить, что многие методы, не относящиеся к вашим вычислительным навыкам, могут быть решены такими методами.

54 голосов
/ 22 декабря 2009

В интересах полноты, версия шаблона C ++, которая для оптимизированной сборки будет вычислять приближение PI во время компиляции и будет встроена в одно значение.

#include <iostream>

template<int I>
struct sign
{
    enum {value = (I % 2) == 0 ? 1 : -1};
};

template<int I, int J>
struct pi_calc
{
    inline static double value ()
    {
        return (pi_calc<I-1, J>::value () + pi_calc<I-1, J+1>::value ()) / 2.0;
    }
};

template<int J>
struct pi_calc<0, J>
{
    inline static double value ()
    {
        return (sign<J>::value * 4.0) / (2.0 * J + 1.0) + pi_calc<0, J-1>::value ();
    }
};


template<>
struct pi_calc<0, 0>
{
    inline static double value ()
    {
        return 4.0;
    }
};

template<int I>
struct pi
{
    inline static double value ()
    {
        return pi_calc<I, I>::value ();
    }
};

int main ()
{
    std::cout.precision (12);

    const double pi_value = pi<10>::value ();

    std::cout << "pi ~ " << pi_value << std::endl;

    return 0;
}

Примечание: для I> 10 оптимизированные сборки могут выполняться медленно, также как и для неоптимизированных запусков. Я полагаю, что для 12 итераций существует около 80 тыс. Вызовов value () (при отсутствии запоминания).

42 голосов
/ 24 августа 2008

На самом деле есть целая книга, посвященная (помимо прочего) быстрым методам вычисления \ pi: 'Pi и AGM' Джонатана и Питера Боровейна (, доступная на Amazon *). 1004 *).

Я немного изучил AGM и связанные с ним алгоритмы: это довольно интересно (хотя иногда нетривиально).

Обратите внимание, что для реализации большинства современных алгоритмов для вычисления \ pi вам понадобится многоточная арифметическая библиотека ( GMP - неплохой выбор, хотя с тех пор, как я последний раз использовал его, прошло довольно много времени).

Сложность по времени лучших алгоритмов находится в O (M (n) log (n)), где M (n) - сложность по времени для умножения двух n-битных целых чисел (M (n) = O (n log (n) log (log (n))) с использованием алгоритмов на основе FFT, которые обычно необходимы при вычислении цифр \ pi, и такой алгоритм реализован в GMP).

Обратите внимание, что хотя математика, лежащая в основе алгоритмов, может и не быть тривиальной, сами алгоритмы обычно представляют собой несколько строк псевдокода, и их реализация обычно очень проста (если вы решили не писать свою собственную арифметику мультиточности: - )).

41 голосов
/ 28 октября 2011

Следующие ответы точно, как это сделать максимально быстро - с наименьшими вычислительными затратами . Даже если вам не нравится ответ, вы должны признать, что это действительно самый быстрый способ получить значение PI.

FASTEST способ получить значение Pi:

1) выбрал свой любимый язык программирования 2) загрузить свою математическую библиотеку 3) и найдите, что Pi там уже определен - готов к использованию!

Если у вас нет библиотеки математики под рукой ..

Способ ВТОРОЙ БЫСТРЫЙ (более универсальное решение):

поиск Пи в Интернете, например здесь:

http://www.eveandersson.com/pi/digits/1000000 (1 миллион цифр .. какова ваша точность с плавающей запятой?)

или здесь:

http://3.141592653589793238462643383279502884197169399375105820974944592.com/

или здесь:

http://en.wikipedia.org/wiki/Pi

Очень быстро найти нужные цифры для любой точной арифметики, которую вы хотели бы использовать, и, определив константу, вы можете быть уверены, что не тратите драгоценное время процессора.

Мало того, что это отчасти юмористический ответ, но на самом деле, если бы кто-нибудь пошел дальше и вычислил значение Pi в реальном приложении ... это было бы довольно большой тратой процессорного времени, не так ли? По крайней мере, я не вижу реального приложения для попытки пересчитать это.

Уважаемый Модератор: обратите внимание, что ОП спросил: «Самый быстрый способ получить значение ПИ»

26 голосов
/ 29 августа 2008

Формула BBP позволяет вычислять n-ую цифру - в базе 2 (или 16) - без необходимости сначала беспокоиться о предыдущих n-1 цифрах:)

22 голосов
/ 08 марта 2009

Вместо того, чтобы определять пи как константу, я всегда использую acos(-1).

21 голосов
/ 06 января 2010

Если эта статья верна, то алгоритм , созданный Беллардом , может быть одним из самых быстрых из доступных. Он создал пи до 2,7 триллиона цифр с помощью настольного ПК!

... и он опубликовал здесь работу

Хорошая работа, Беллард, ты пионер!

http://www.theregister.co.uk/2010/01/06/very_long_pi/

21 голосов
/ 12 января 2009

Только что наткнулся на этот, который должен быть здесь для полноты:

Рассчитать PI в пите

Он имеет довольно приятное свойство: точность можно повысить, увеличив программу.

Здесь немного о самом языке

...