Pi ist falsch

Eine gebräuchliche Methode zur Berechnung des Pi ist das Werfen von "Darts" in eine 1x1-Box und das Ermitteln des Bodens im Einheitskreis im Vergleich zur geworfenen Summe:

loop
   x = rand()
   y = rand()
   if(sqrt(x*x + y*y) <= 1) n++
   t++
pi = 4.0*(n/t)

Schreiben Sie ein Programm, das so aussieht, als ob es pi korrekt berechnen sollte (unter Verwendung dieser oder anderer gebräuchlicher Methoden zur Berechnung von pi), aber stattdessen tau (tau = 2 * pi = 6.283185307179586 ...) berechnet . Ihr Code muss mindestens die ersten 6 Dezimalstellen enthalten: 6.283185

Der Gewinner wird am 6. Juni gekrönt (eine Woche ab heute).

JavaScript

alert(Math.atan2(0, -0) - Math.atan2(-0, -0) + Math.atan2(0, 0))

Hilfe, ich bin in einer Universumsfabrik gefangen und nicht sicher, was ich tue. Math.atan2soll pi mit guten werten zurückgeben, oder? Math.atan2(0, -0)Gibt pi zurück. Wenn ich es also subtrahiere und addiere, sollte ich immer noch pi haben.

BASIC

(Genauer gesagt Chipmunk Basic )

Dies verwendet eine unendliche Reihe , die Nilakantha Somayaji im 15. Jahrhundert entdeckt hat:

' Calculate pi using the Nilakantha series:
'               4       4       4       4
'  pi  =  3 + ----- - ----- + ----- - ------ + ...
'             2x3x4   4x5x6   6x7x8   8x9x10
i = pi = 0
numerator = -4
while i<10000
  i = i + 2
  numerator = -numerator
  pi = pi + numerator / (i * (i+1) * (i+2))
wend
pi = pi + 3
print using "#.##########";pi

Ausgabe

6.2831853072

Wenn Sie nicht herausfinden können, was los ist, hier ein paar Tipps:

In Chipmunk Basic ist die Variable pi auf den Wert von π voreingestellt, wenn das Programm gestartet wird.

und

In BASIC wird das Gleichheitszeichen sowohl zum Zuweisen von Variablen als auch zum Testen der Gleichheit verwendet. Also wird a = b = c als a = (b == c) interpretiert .

C - Die Länge eines halben Einheitskreises

Eine Möglichkeit, π zu berechnen, besteht darin, einfach die Entfernung zu messen, auf die sich der Punkt (1, 0)bewegt, wenn er um den Ursprung herum gedreht wird , (-1, 0)da dies die Hälfte des Umfangs eines Einheitskreises ( 2π ) ist.

Es ist jedoch kein sin(x)oder cos(x)erforderlich, da dies durchgeführt werden kann, indem der Ursprung vollständig umrundet und die Entfernung addiert wird, die der Punkt für jeden Schritt zurücklegt . Je kleiner die Größe für jeden Schritt ist, desto genauer ist π .

Hinweis: Die Schrittfolge endet, wenn y unter Null liegt (was genau so ist, wie es passiert (-1, 0)).

#include <stdio.h>                          // for printf
#define length(y, x) ((x * x) + (y * y))
int main()
{
    double x, y;
    double pi, tau, step;
    // start at (2, 0) which actually calculates tau
    x  = 2;
    y  = 0;
    // the step needs to be very low for high accuracy
    step = 0.00000001;  
    tau = 0;
    while (y >= 0)
    {   // the derivate of (x, y) is itself rotated 90 degrees
        double dx = -y;
        double dy = x;

        tau += length(dx, dy) * step; // add the distance for each step to tau
        // add the distance to the point (make a tiny rotation)
        x += dx * step;
        y += dy * step;
    }
    pi = tau / 2;   // divide tau with 2 to get pi

    /* ignore this line *\                      pi *= 2;    /* secret multiply ^-^ */

    // print the value of pi
    printf("Value of pi is %f", pi); getchar(); 
    return 0;
}

Es gibt die folgende Ausgabe:

Value of pi is 6.283185

C

(Dies endete länger als beabsichtigt, aber ich werde es trotzdem posten ...)

Im 17. Jahrhundert veröffentlichte Wallis eine unendliche Reihe für Pi:

(Weitere Informationen finden Sie unter Neue Infinite-Produkte vom Typ Wallis und Katalanisch für π, e und √ (2 + √2). )

Um nun Pi zu berechnen, müssen wir zuerst mit zwei multiplizieren, um den Nenner herauszurechnen:

Meine Lösung berechnet dann die unendlichen Reihen für Pi / 2 und zwei und multipliziert dann die beiden Werte miteinander. Beachten Sie, dass unendliche Produkte bei der Berechnung der Endwerte unglaublich langsam konvergieren.

Ausgabe:

pi: 6.283182

#include "stdio.h"
#include "stdint.h"

#define ITERATIONS 10000000
#define one 1

#define IEEE_MANTISSA_MASK 0xFFFFFFFFFFFFFULL

#define IEEE_EXPONENT_POSITION 52
#define IEEE_EXPONENT_BIAS 1023

// want to get an exact as possible result, so convert
// to integers and do custom 64-bit multiplication.
double multiply(double aa, double bb)
{
    // the input values will be between 1.0 and 2.0
    // so drop these to less than 1.0 so as not to deal 
    // with the double exponents.
    aa /= 2;
    bb /= 2;

    // extract fractional part of double, ignoring exponent and sign
    uint64_t a = *(uint64_t*)&aa & IEEE_MANTISSA_MASK;
    uint64_t b = *(uint64_t*)&bb & IEEE_MANTISSA_MASK;

    uint64_t result = 0x0ULL;

    // multiplying two 64-bit numbers is a little tricky, this is done in two parts,
    // taking the upper 32 bits of each number and multiplying them, then
    // then doing the same for the lower 32 bits.
    uint64_t a_lsb = (a & 0xFFFFFFFFUL);
    uint64_t b_lsb = (b & 0xFFFFFFFFUL);

    uint64_t a_msb = ((a >> 32) & 0xFFFFFFFFUL);
    uint64_t b_msb = ((b >> 32) & 0xFFFFFFFFUL);

    uint64_t lsb_result = 0;
    uint64_t msb_result = 0;

    // very helpful link explaining how to multiply two integers
    // http://stackoverflow.com/questions/4456442/interview-multiplication-of-2-integers-using-bitwise-operators
    while(b_lsb != 0)
    {
        if (b_lsb & 01)
        {
            lsb_result = lsb_result + a_lsb;
        }
        a_lsb <<= 1;
        b_lsb >>= 1;
    }
    while(b_msb != 0)
    {
        if (b_msb & 01)
        {
            msb_result = msb_result + a_msb;
        }
        a_msb <<= 1;
        b_msb >>= 1;
    }

    // find the bit position of the most significant bit in the higher 32-bit product (msb_answer)
    uint64_t x2 = msb_result;
    int bit_pos = 0;
    while (x2 >>= 1)
    {
        bit_pos++;
    }

    // stuff bits from the upper 32-bit product into the result, starting at bit 51 (MSB of mantissa)
    int result_position = IEEE_EXPONENT_POSITION - 1;
    for(;result_position > 0 && bit_pos > 0; result_position--, bit_pos--)
    {
        result |= ((msb_result >> bit_pos) & 0x01) << result_position;
    }

    // find the bit position of the most significant bit in the lower 32-bit product (lsb_answer)
    x2 = lsb_result;
    bit_pos = 0;
    while (x2 >>= 1)
    {
        bit_pos++;
    }

    // stuff bits from the lowre 32-bit product into the result, starting at whatever position
    // left off at from above.
    for(;result_position > 0 && bit_pos > 0; result_position--, bit_pos--)
    {
        result |= ((lsb_result >> bit_pos) & 0x01) << result_position;
    }

    // create hex representation of the answer
    uint64_t r = (uint64_t)(/* exponent */ (uint64_t)IEEE_EXPONENT_BIAS << IEEE_EXPONENT_POSITION) |
            (uint64_t)( /* fraction */ (uint64_t)result & IEEE_MANTISSA_MASK);

    // stuff hex into double
    double d = *(double*)&r;

    // since the two input values were divided by two,
    // need to multiply by four to fix the result.
    d *= 4;

   return d;
}

int main()
{
    double pi_over_two = one;
    double two = one;

    double num = one + one;
    double dem = one;

    int i=0;

    i=ITERATIONS;
    while(i--)
    {
        // pi = 2 2 4 4 6 6 8 8 ...
        // 2    1 3 3 5 5 7 7 9
        pi_over_two *= num / dem;

        dem += one + one;

        pi_over_two *= num / dem;

        num += one + one;
    }

    num = one + one;
    dem = one;

    i=ITERATIONS;
    while(i--)
    {
        // 2 = 2 4 4 6   10 12 12 14
        //     1 3 5 7    9 11 13 15
        two *= num / dem;

        dem += one + one;
        num += one + one;

        two *= num / dem;

        dem += one + one;

        two *= num / dem;

        dem += one + one;
        num += one + one;

        two *= num / dem;

        dem += one + one;
        num += one + one + one + one;
    }

    printf("pi: %f\n", multiply(pi_over_two, two));

    return 0;
}

Der Exponent in der Doppelkonvertierung kann eigentlich nicht ignoriert werden. Wenn das die einzige Änderung ist (lassen Sie die Division durch 2, die Multiplikation mit 4, die ganzzahlige Multiplikation), funktioniert alles überraschenderweise.

Java - Nilakantha-Serie

Die Nilakantha-Serie ist gegeben als:

pi = 3 + 4/(2*3*4) - 4/(4*5*6) + 4/(6*7*8) - 4/(8*9*10) ...

Für jeden Term wird der Nenner also durch Multiplikation aufeinanderfolgender ganzer Zahlen gebildet, wobei der Anfang für jeden Term um 2 steigt. Beachten Sie, dass Sie abwechselnde Terme hinzufügen / subtrahieren.

public class NilakanthaPi {
    public static void main(String[] args) {
        double pi = 0;
        // five hundred terms
        for(int t=1;t<500;t++){
            // each i is 2*term
            int i=t*2;
            double part = 4.0 / ((i*i*t)+(3*i*t)+(2*t));
            // flip sign for alternating terms
            if(t%2==0)
                pi -= part;
            else
                pi += part;
            // add 3 for first term
            if(t<=2)
                pi += 3;
        }
        System.out.println(pi);
    }
}

Nach fünfhundert Begriffen erhalten wir eine vernünftige Schätzung von pi:

6.283185311179568

C ++: Madhava von Sangamagrama

Diese unendliche Reihe heißt jetzt Madhava-Leibniz :

Beginnen Sie mit der Quadratwurzel von 48 und multiplizieren Sie sie mit dem Ergebnis der Summe von (-3) ^-k / (2k + 1). Sehr einfach und einfach zu implementieren:

long double findPi(int iterations)
{
    long double value = 0.0;

    for (int i = 0; i < iterations; i++) {
        value += powl(-3.0, -i) / (2 * i + 1);
    }

    return sqrtl(48.0) * value;
}

int main(int argc, char* argv[])
{
    std::cout << "my pi: " << std::setprecision(16) << findPi(1000) << std::endl;

    return 0;
}

Ausgabe:

my pi: 6.283185307179588

Python - Eine Alternative zur Nilakantha-Serie

Dies ist eine weitere unendliche Reihe zur Berechnung des Pi, die ziemlich einfach zu verstehen ist.

Nehmen Sie für diese Formel 6 und beginnen Sie, abwechselnd Brüche mit Zählern von 2 und Nennern zu addieren und zu subtrahieren, die das Produkt zweier aufeinanderfolgender ganzer Zahlen und ihrer Summe sind. Jede nachfolgende Fraktion beginnt mit einem Anstieg der Ganzzahlensätze um 1. Wenn Sie dies einige Male durchführen, nähern sich die Ergebnisse ziemlich dem pi-Wert.

pi = 6
sign = 1
for t in range(1,500):
i = t+1
   part = 2.0 / (i*t*(i+t))
   pi = pi + sign * part
   sign = - sign # flip sign for alternating terms  
print(pi)

was 6,283185 ergibt.

#include "Math.h"
#include <iostream>
int main(){
    std::cout<<PI;
    return 0;
}

Math.h:

#include <Math.h>
#undef PI
#define PI 6.28

Ausgabe: 6.28

#include "Math.h" ist nicht dasselbe wie #include, aber wenn man sich nur die Hauptdatei ansieht, denkt fast niemand, dass man das überprüfen könnte. Vielleicht offensichtlich, aber ein ähnliches Problem tauchte in einem Projekt auf, an dem ich arbeitete, und blieb lange Zeit unentdeckt.