P i = = 3,2

37

Zainspirowany tym filmem z serii Infinite .

Wprowadzenie

Pi jest zdefiniowane jako stosunek obwodu do średnicy koła. Ale jak zdefiniowano okrąg? Zwykle okrąg definiuje się jako punkty o stałej odległości od punktu środkowego (załóżmy, że środek znajduje się w punkcie (0,0)). Następne pytanie brzmi: jak zdefiniować odległość ? Poniżej rozważamy różne pojęcia odległości (wywołane przez Lp-normy):

Biorąc pod uwagę normę (= coś, co mierzy długość ), możemy łatwo skonstruować odległość (= odległość między dwoma punktami) w następujący sposób:

dist(A,B) := norm (A-B)

Normę euklidesową podaje:

norm((x,y)) = (x^2 + y^2)^(1/2)

Jest to również nazywane normą L2 . Pozostałe normy Lp są konstruowane przez zastąpienie 2powyższej formuły innymi wartościami od 1 do nieskończoności:

norm_p((x,y)) = (|x|^p + |y|^p)^(1/p)

Koła jednostek dla tych różnych norm mają dość wyraźne kształty:

Wyzwanie

Biorąc pod uwagę a p >= 1, obliczyć stosunek obwodu do średnicy koła Lp w odniesieniu do Lp-normy z dokładnością do czterech znaczących cyfr.

Przypadki testowe

Możemy użyć że p,qze 1 = 1/p + 1/qmamy ten sam stosunek dla Lpjak również Lqnormy. Co więcej, p = q = 2stosunek jest minimalny, a dla p = 1, q = infinitywspółczynnika 4, więc stosunki są zawsze pomiędzy pii 4.

p   or  q            ratio
1       infinity     4
2       2            3.141592
1.623   2.60513      3.200
1.5     3            3.25976
4       1.33333      3.39693
wada
źródło
2
Kształty są znane jako krzywe Lame'a lub supereliplipy i istnieją również dla 0 < p <1, nawet jeśli sama norma tego nie robi (ponieważ narusza to nierówność trójkąta). Artykuł w Wikipedii dotyczący superelipsu zawiera zamknięty formularz dla tego obszaru.
Neil
@ Neil Musimy jednak wziąć pod uwagę obwód, a nie obszar, który - o ile wiem - można obliczyć tylko za pomocą całki długości łuku.
flawr
7
Przepraszam, ale zanim skończyłem czytać o nich, zapomniałem, o co pytało pytanie.
Neil
2
Cudowne wyzwanie!
Luis Mendo,
1
Warto zauważyć, że formuła obszaru ( A = πr²) nie obowiązuje dlap ≠ 2
Mego

Odpowiedzi:

12

Python + scipy, 92 bajty

from scipy.integrate import*
lambda p:2/p*quad(lambda x:(x/x**p+(1-x)**(1-p))**(1/p),0,1)[0]

Formula pochodzi z tego pytania matematycznego .

orlp
źródło
Podczas testowania implementacji przy użyciu tej oceny miałem problem ze zbieżnością tego podejścia, ze względu na osobliwość w x=1, jak radzi sobie twoje poddanie?
flawr
Scipy nie jest częścią standardowej biblioteki Pythona. Może przestawisz się na Sage?
busukxuan
2
@busukxuan Nie ma wymagań dotyczących PPCG, które pozwalają na używanie tylko standardowych bibliotek. Ale i tak wspomnę o tym w tytule.
orlp
1
@ChristianSievers Zrobiłem własną integrację, aby uniknąć złego korzystania z zamkniętej formuły innej osoby :-P
Luis Mendo
1
@ChristianSievers Właściwie włączyłem też inną formułę do piaskownicy, na wypadek gdybyś był zainteresowany =)
flawr
10

MATL , 31 bajtów

0:1e-3:1lyG^-lG/^v!d|G^!slG/^sE

Wypróbuj online! Lub sprawdź wszystkie przypadki testowe .

Wyjaśnienie

To generuje współrzędne x , y jednej czwartej koła jednostki próbkowanej w 1001 punktach z krokiem 0,001 w x . Długość ćwiartki koła jest przybliżona przez długość linii wielokąta, która przechodzi przez te punkty; to jest suma długości 1000 segmentów. Długość jest oczywiście obliczana zgodnie z p-norm. Pomnożenie wyniku przez 2 daje przybliżoną długość półkola, czyli pi.

0:1e-3:1   % Push [0 0.001 0.002 ... 0.999 1]. These are the x coordinates of
           % the vertices of the polygonal line that will approximate a quarter
           % of the unit circle
l          % Push 1
y          % Duplicate [0 0.001 0.002 ... 0.999 1] onto the top of the stack.
G          % Push input, p
^          % Element-wise power: gives [0^p 0.001^p ... 1^p]
-          % Element-wise subtract from 1: gives [1-0^p 1-0.001^p ... 1-1^p]
lG/        % Push 1, push p, divide: gives 1/p
^          % Element-wise power: gives [(1-0^p)^(1/p) (1-0.001^p)^(1/p) ...
           % ... (1-1^p)^(1/p)]. These are the y coordinates of the vertices
           % of the polygonal line
v          % Concatenate vertically into a 2×1001 matrix. The first row contains
           % the x coordinates and the second row contains the y coordinates
!          % Transpose
d|         % Compute consecutive differences down each column. This gives a
           % 1000×2 matrix with the x and y increments of each segment. These
           % increments will be referred to as Δx, Δy
G          % Push p
^          % Element-wise power
!          % Transpose
s          % Sum of each column. This gives a 1×1000 vector containing
           % (Δx)^p+(Δy)^p for each segment
lG/        % Push 1/p
^          % Element-wise power. This gives a 1×1000 vector containing 
           % ((Δx)^p+(Δy)^p)^(1/p) for each segment, that is, the length of 
           % each segment according to p-norm
s          % Sum the lenghts of all segments. This approximates the length of
           % a quarter of the unit circle
E          % Multiply by 2. This gives the length of half unit circle, that is,
           % pi. Implicitly display
Luis Mendo
źródło
8

Mathematica, 49 46 bajtów

3 bajty zapisane z powodu alephalpha .

2NIntegrate[(1+(a^-#-1)^(1-#))^(1/#),{a,0,1}]&

Funkcja anonimowa. Pobiera liczbę jako dane wejściowe i zwraca liczbę jako dane wyjściowe.

LegionMammal978
źródło
1
2NIntegrate[(1+(a^-#-1)^(1-#))^(1/#),{a,0,1}]&
alephalpha
5

PARI / GP, 48 43 bajtów

Łatwo jest po znalezieniu formuły przez @orlp, a wersja @ alephalpha zapisuje 5 bajtów:

p->2*intnum(u=0,1,(1+(u^-p-1)^(1-p))^(1/p))

Aby dodać coś nieco przydatnego, obliczmy, pza co otrzymujemy 3.2:

? f=p->2*intnum(u=0,1,(1+(u^-p-1)^(1-p))^(1/p));
? solve(p=1,2,f(p)-3.2)
%2 = 1.623002382384469009676324702

Prawidłowa obsługa

Natomiast kod daje wyniki, które są znacznie bardziej dokładny niż żądań wyzwanie, może łatwo ulec poprawie wiele: czy możemy zastąpić górną granicę integracji 1z [1,1/p-1](co daje ręczne połączenia wykładnik osobliwość) to wszystko pokazane cyfry f(2)zgadzają się Pi. Jest to nadal prawdą, jeśli zwiększymy precyzję do 100 (typ \p100).

Jednak po tej zmianie solveobliczenia przestały działać. Zmieniłem wewnętrzny termin, aby jawnie obsłużyć sprawę, u=0a także zmieniłem na inny komputer z nowszą wersją PARI i 64-bitową (co oznacza wyższą domyślną precyzję).

Oto ulepszone obliczanie pwartości Pi=3.2i spójrzmy również na prawdziwe Pi:

? f=p->2*intnum(u=0,[1,1/p-1],if(u,(1+(u^-p-1)^(1-p))^(1/p),0));
? f(2)
%2 = 3.1415926535897932384626433832795028842
? Pi
%3 = 3.1415926535897932384626433832795028842
? solve(p=1,2,f(p)-3.2)
%4 = 1.6230023823844690096763253745604419761
Christian Sievers
źródło
p->2*intnum(u=0,1,(1+(u^-p-1)^(1-p))^(1/p))
alephalpha
0

JavaScript (ES7), 80 bajtów

Na podstawie odpowiedzi orlp . Ta implementacja JS jest dość powolna. Możesz spróbować i=1e-7(lub nawet wyżej), aby uzyskać szybsze zbliżenie.

Uwaga : jest to zasadniczo przeznaczone tylko dla Chrome i Edge. Równoważna wersja ES6 używająca Math.pow()przeglądarki Firefox 50.1 wydaje się być znacznie wolniejsza.

Edycja : Według Neila powinno to również działać dobrze w przeglądarce Firefox 52.

f=
p=>{for(i=5e-8,s=x=0;(x+=i)<1;)s+=i*(x**(1-p)+(1-x)**(1-p))**(1/p);return 2/p*s}

console.log(f(1).toFixed(3))
console.log(f(2).toFixed(3))
console.log(f(1.623).toFixed(3))

Arnauld
źródło
Wersja ES7 wydawała się dość dziwna, gdy wypróbowałem ją za pomocą Firefoksa 52 (nie zmierzyłem jej naukowo, ale czułem mniej więcej tę samą prędkość co Chrome; Edge zamarł na mnie).
Neil
@Neil Dziękujemy za opinię. Zaktualizowano odpowiednio.
Arnauld,