什么是获得＆＃960;的价值的最快方法？

•
•  |
•

题

#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;
}


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


#!/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


解决方案

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)


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);}


其他提示

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()
{
_-_-_-_
_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_-_-_-_-_
_-_-_-_-_-_-_-_
_-_-_-_
}


#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;
}


1）选择您喜欢的编程语言 2）加载其数学库 3）并发现Pi已经在那里定义 - 准备好使用了！

SECOND FASTEST 方式（更通用的解决方案）是：

http://www.eveandersson.com/pi/digits/1000000 （100万位...你的浮点精度是多少？）

BBP公式允许您计算第n位数 - 在基数2中（或16） - 首先不必担心以前的n-1数字：）

......他已发布了他的在这里工作

from math import pi
from time import time

precision = 10**6 # higher value -> higher precision
# lower  value -> higher speed

t = time()

calc = 0
for k in xrange(0, precision):
calc += ((-1)**k) / (2*k+1.)
calc *= 4. # this is just a little optimization

t = time()-t

print "Calculated: %.40f" % calc
print "Costant pi: %.40f" % pi
print "Difference: %.40f" % abs(calc-pi)
print "Time elapsed: %s" % repr(t)


from gmpy import pi
print pi(3000) # the rule is the same as
# the precision on the previous code


static char doc_pi[]="\
pi(n): returns pi with n bits of precision in an mpf object\n\
";

/* This function was originally from netlib, package bmp, by
* Richard P. Brent. Paulo Cesar Pereira de Andrade converted
* it to C and used it in his LISP interpreter.
*
* Original comments:
*
*   sets mp pi = 3.14159... to the available precision.
*   uses the gauss-legendre algorithm.
*   this method requires time o(ln(t)m(t)), so it is slower
*   than mppi if m(t) = o(t**2), but would be faster for
*   large t if a faster multiplication algorithm were used
*   (see comments in mpmul).
*   for a description of the method, see - multiple-precision
*   zero-finding and the complexity of elementary function
*   evaluation (by r. p. brent), in analytic computational
*   complexity (edited by j. f. traub), academic press, 1976, 151-176.
*   rounding options not implemented, no guard digits used.
*/
static PyObject *
Pygmpy_pi(PyObject *self, PyObject *args)
{
PympfObject *pi;
int precision;
mpf_t r_i2, r_i3, r_i4;
mpf_t ix;

ONE_ARG("pi", "i", &precision);
if(!(pi = Pympf_new(precision))) {
return NULL;
}

mpf_set_si(pi->f, 1);

mpf_init(ix);
mpf_set_ui(ix, 1);

mpf_init2(r_i2, precision);

mpf_init2(r_i3, precision);
mpf_set_d(r_i3, 0.25);

mpf_init2(r_i4, precision);
mpf_set_d(r_i4, 0.5);
mpf_sqrt(r_i4, r_i4);

for (;;) {
mpf_set(r_i2, pi->f);
mpf_add(pi->f, pi->f, r_i4);
mpf_div_ui(pi->f, pi->f, 2);
mpf_mul(r_i4, r_i2, r_i4);
mpf_sub(r_i2, pi->f, r_i2);
mpf_mul(r_i2, r_i2, r_i2);
mpf_mul(r_i2, r_i2, ix);
mpf_sub(r_i3, r_i3, r_i2);
mpf_sqrt(r_i4, r_i4);
mpf_mul_ui(ix, ix, 2);
/* Check for convergence */
if (!(mpf_cmp_si(r_i2, 0) &&
mpf_get_prec(r_i2) >= (unsigned)precision)) {
mpf_mul(pi->f, pi->f, r_i4);
mpf_div(pi->f, pi->f, r_i3);
break;
}
}

mpf_clear(ix);
mpf_clear(r_i2);
mpf_clear(r_i3);
mpf_clear(r_i4);

return (PyObject*)pi;
}


;''6666,-2%{2+.2/@*\/10.3??2*+}*1000<~\;


176 * arctan (1/57) + 28 * arctan (1/239) - 48 * arctan (1/682) + 96 * arctan(1/12943)

[; \left( 176 \arctan \frac{1}{57} + 28 \arctan \frac{1}{239} - 48 \arctan \frac{1}{682} + 96 \arctan \frac{1}{12943}\right) ;], for you TeX the World people.


（+（ - （+（* 176（atan（/ 1 57）））（* 28（atan（/ 1 239））））（* 48（atan（/ 1 682）））） （* 96（atan（/ 1 12943））））

4.0 * (4.0 * Math.Atan(0.2) - Math.Atan(1.0 / 239.0))


Seth，它是3.14159265358979323846 3 ，而不是64。

使用D。

（复制自 DSource.org

/** Calculate pi at compile time
*
* Compile with dmd -c pi.d
*/
module calcpi;

import meta.math;
import meta.conv;

/** real evaluateSeries!(real x, real metafunction!(real y, int n) term)
*
* Evaluate a power series at compile time.
*
* Given a metafunction of the form
*  real term!(real y, int n),
* which gives the nth term of a convergent series at the point y
* (where the first term is n==1), and a real number x,
* this metafunction calculates the infinite sum at the point x
* by adding terms until the sum doesn't change any more.
*/
template evaluateSeries(real x, alias term, int n=1, real sumsofar=0.0)
{
static if (n>1 && sumsofar == sumsofar + term!(x, n+1)) {
const real evaluateSeries = sumsofar;
} else {
const real evaluateSeries = evaluateSeries!(x, term, n+1, sumsofar + term!(x, n));
}
}

/*** Calculate atan(x) at compile time.
*
* Uses the Maclaurin formula
*  atan(z) = z - z^3/3 + Z^5/5 - Z^7/7 + ...
*/
template atan(real z)
{
const real atan = evaluateSeries!(z, atanTerm);
}

template atanTerm(real x, int n)
{
const real atanTerm =  (n & 1 ? 1 : -1) * pow!(x, 2*n-1)/(2*n-1);
}

/// Machin's formula for pi
/// pi/4 = 4 atan(1/5) - atan(1/239).
pragma(msg, "PI = " ~ fcvt!(4.0 * (4*atan!(1/5.0) - atan!(1/239.0))) );


Pi正好是3！ [教授弗林克（辛普森一家）]

using System;
using System.Text;

class Program {
static void Main(string[] args) {
int Digits = 100;

BigNumber x = new BigNumber(Digits);
BigNumber y = new BigNumber(Digits);
x.ArcTan(16, 5);
y.ArcTan(4, 239);
x.Subtract(y);
string pi = x.ToString();
Console.WriteLine(pi);
}
}

public class BigNumber {
private UInt32[] number;
private int size;
private int maxDigits;

public BigNumber(int maxDigits) {
this.maxDigits = maxDigits;
this.size = (int)Math.Ceiling((float)maxDigits * 0.104) + 2;
number = new UInt32[size];
}
public BigNumber(int maxDigits, UInt32 intPart)
: this(maxDigits) {
number[0] = intPart;
for (int i = 1; i < size; i++) {
number[i] = 0;
}
}
private void VerifySameSize(BigNumber value) {
if (Object.ReferenceEquals(this, value))
throw new Exception("BigNumbers cannot operate on themselves");
if (value.size != this.size)
throw new Exception("BigNumbers must have the same size");
}

public void Add(BigNumber value) {
VerifySameSize(value);

int index = size - 1;
while (index >= 0 && value.number[index] == 0)
index--;

UInt32 carry = 0;
while (index >= 0) {
UInt64 result = (UInt64)number[index] +
value.number[index] + carry;
number[index] = (UInt32)result;
if (result >= 0x100000000U)
carry = 1;
else
carry = 0;
index--;
}
}
public void Subtract(BigNumber value) {
VerifySameSize(value);

int index = size - 1;
while (index >= 0 && value.number[index] == 0)
index--;

UInt32 borrow = 0;
while (index >= 0) {
UInt64 result = 0x100000000U + (UInt64)number[index] -
value.number[index] - borrow;
number[index] = (UInt32)result;
if (result >= 0x100000000U)
borrow = 0;
else
borrow = 1;
index--;
}
}
public void Multiply(UInt32 value) {
int index = size - 1;
while (index >= 0 && number[index] == 0)
index--;

UInt32 carry = 0;
while (index >= 0) {
UInt64 result = (UInt64)number[index] * value + carry;
number[index] = (UInt32)result;
carry = (UInt32)(result >> 32);
index--;
}
}
public void Divide(UInt32 value) {
int index = 0;
while (index < size && number[index] == 0)
index++;

UInt32 carry = 0;
while (index < size) {
UInt64 result = number[index] + ((UInt64)carry << 32);
number[index] = (UInt32)(result / (UInt64)value);
carry = (UInt32)(result % (UInt64)value);
index++;
}
}
public void Assign(BigNumber value) {
VerifySameSize(value);
for (int i = 0; i < size; i++) {
number[i] = value.number[i];
}
}

public override string ToString() {
BigNumber temp = new BigNumber(maxDigits);
temp.Assign(this);

StringBuilder sb = new StringBuilder();
sb.Append(temp.number[0]);
sb.Append(System.Globalization.CultureInfo.CurrentCulture.NumberFormat.CurrencyDecimalSeparator);

int digitCount = 0;
while (digitCount < maxDigits) {
temp.number[0] = 0;
temp.Multiply(100000);
sb.AppendFormat("{0:D5}", temp.number[0]);
digitCount += 5;
}

return sb.ToString();
}
public bool IsZero() {
foreach (UInt32 item in number) {
if (item != 0)
return false;
}
return true;
}

public void ArcTan(UInt32 multiplicand, UInt32 reciprocal) {
BigNumber X = new BigNumber(maxDigits, multiplicand);
X.Divide(reciprocal);
reciprocal *= reciprocal;

this.Assign(X);

BigNumber term = new BigNumber(maxDigits);
UInt32 divisor = 1;
bool subtractTerm = true;
while (true) {
X.Divide(reciprocal);
term.Assign(X);
divisor += 2;
term.Divide(divisor);
if (term.IsZero())
break;

if (subtractTerm)
this.Subtract(term);
else
this.Add(term);
subtractTerm = !subtractTerm;
}
}
}


program calcpi;

{\$APPTYPE CONSOLE}

uses
SysUtils;

var
start, finish: TDateTime;

function CalculatePi(iterations: integer): double;
var
numerator, denominator, i: integer;
sum: double;
begin
{
PI may be approximated with this formula:
4 * (1 - 1/3 + 1/5 - 1/7 + 1/9 - 1/11 .......)
//}
numerator := 1;
denominator := 1;
sum := 0;
for i := 1 to iterations do begin
sum := sum + (numerator/denominator);
denominator := denominator + 2;
numerator := -numerator;
end;
Result := 4 * sum;
end;

begin
try
start := Now;
WriteLn(FloatToStr(CalculatePi(StrToInt(ParamStr(1)))));
finish := Now;
WriteLn('Seconds:' + FormatDateTime('hh:mm:ss.zz',finish-start));
except
on E:Exception do
Writeln(E.Classname, ': ', E.Message);
end;
end.


/* Return approximation of n * PI; n is integer */
#define pi_times(n) (((n) * 22) / 7)


PI = 3.141592654


＆＃960; = - 我 ln（-1）

<input id="range" type="range" min="10" max="960" value="10" step="50" oninput="calcPi()">
<br>
<div id="cont"></div>

<script>
function generateCircle(width) {
var c = width/2;
var delta = 1.0;
var str = "";
var xCount = 0;
for (var x=0; x <= width; x++) {
for (var y = 0; y <= width; y++) {
var d = Math.sqrt((x-c)*(x-c) + (y-c)*(y-c));
if (d > (width-1)/2) {
str += '.';
}
else {
xCount++;
str += 'o';
}
str += "&nbsp;"
}
str += "\n";
}
var pi = (xCount * 4) / (width * width);
return [str, pi];
}

function calcPi() {
var e = document.getElementById("cont");
var width = document.getElementById("range").value;
e.innerHTML = "<h4>Generating circle...</h4>";
setTimeout(function() {
var circ = generateCircle(width);
e.innerHTML  = "<pre>" + "π = " + circ[1].toFixed(2) + "\n" + circ[0] +"</pre>";
}, 200);
}
calcPi();
</script>

更好的方法

• 数学库的pi变量。数学库将变量pi存储为常量。

math_pi.py

import math
print math.pi


Command being timed: "python math_pi.py"
User time (seconds): 0.01
System time (seconds): 0.01
Percent of CPU this job got: 91%
Elapsed (wall clock) time (h:mm:ss or m:ss): 0:00.03

• 使用数学cos算法

acos_pi.py

import math
print math.acos(-1)


Command being timed: "python acos_pi.py"
User time (seconds): 0.02
System time (seconds): 0.01
Percent of CPU this job got: 94%
Elapsed (wall clock) time (h:mm:ss or m:ss): 0:00.03


bbp_pi.py

from decimal import Decimal, getcontext
getcontext().prec=100
print sum(1/Decimal(16)**k *
(Decimal(4)/(8*k+1) -
Decimal(2)/(8*k+4) -
Decimal(1)/(8*k+5) -
Decimal(1)/(8*k+6)) for k in range(100))


Command being timed: "python c.py"
User time (seconds): 0.05
System time (seconds): 0.01
Percent of CPU this job got: 98%
Elapsed (wall clock) time (h:mm:ss or m:ss): 0:00.06
`