tło

Paradoks urodzin jest popularnym problemem w teorii, która przeczy prawdopodobieństwa intuicji matematycznej (większość ludzi). Opis problemu jest następujący:

Biorąc pod uwagę N osób, jakie jest prawdopodobieństwo, że co najmniej dwa z nich mają takie same urodziny (bez względu na rok).

Problem zwykle upraszcza się, całkowicie ignorując dni przestępne. W tym przypadku odpowiedź dla N = 23 to P (23) ≈ 0,5072972 (jako częsty przykład). Powiązany artykuł w Wikipedii wyjaśnia, jak dojść do tego prawdopodobieństwa. Alternatywnie, ten film Numberphile wykonuje naprawdę dobrą robotę.

Jednak w przypadku tego wyzwania chcemy to zrobić dobrze i nie ignorujemy lat przestępnych. Jest to nieco bardziej skomplikowane, ponieważ teraz należy dodać 29 lutego, ale te szczególne urodziny są mniej prawdopodobne niż wszystkie inne.

Wykorzystamy również zasady dotyczące pełnego roku przestępnego :

• Jeśli rok dzieli się przez 400, to jest to rok przestępny.
• W przeciwnym razie, jeśli rok dzieli się przez 100, to nie jest to rok przestępny.
• W przeciwnym razie, jeśli rok dzieli się przez 4, to jest to rok przestępny.
• W przeciwnym razie to nie jest rok przestępny.

Zmieszany? Oznacza to, że lata 1700, 1800, 1900, 2100, 2200, 2300 nie są latami przestępnymi, ale 1600, 2000, 2400 są (podobnie jak każdy inny rok podzielny przez 4). Ten kalendarz powtarza się co 400 lat, a my przyjmiemy jednolity rozkład urodzin w ciągu tych 400 lat.

Skorygowany wynik dla N = 23 wynosi teraz P (23) ≈ 0,5068761 .

Wyzwanie

Biorąc pod uwagę liczbę całkowitą 1 ≤ N < 100, określ prawdopodobieństwo, że wśród Nosób co najmniej dwie mają takie same urodziny, biorąc pod uwagę zasady dotyczące roku przestępnego. Wynik powinien być liczbą zmiennoprzecinkową lub stałą, z dokładnością do co najmniej 6 miejsc po przecinku. Dopuszczalne jest obcinanie końcowych zer.

Możesz napisać program lub funkcję, pobierając dane wejściowe przez STDIN (lub najbliższą alternatywę), argument wiersza poleceń lub argument funkcji i wyprowadzić wynik przez STDOUT (lub najbliższą alternatywę), wartość zwracaną funkcji lub parametr funkcji (wyjściowej).

Twoje rozwiązanie musi być w stanie wygenerować dane wyjściowe dla wszystkich 99 danych wejściowych w ciągu kilku sekund. Ma to głównie na celu wykluczenie metod Monte Carlo z mnóstwem próbek, więc jeśli używasz zasadniczo szybkiego i dokładnego algorytmu w zbyt wolnym języku ezoterycznym, jestem gotów dać swobodę tej zasadzie.

Przypadki testowe

Oto pełna tabela wyników:

 1 => 0.000000
 2 => 0.002737
 3 => 0.008195
 4 => 0.016337
 5 => 0.027104
 6 => 0.040416
 7 => 0.056171
 8 => 0.074251
 9 => 0.094518
10 => 0.116818
11 => 0.140987
12 => 0.166844
13 => 0.194203
14 => 0.222869
15 => 0.252642
16 => 0.283319
17 => 0.314698
18 => 0.346578
19 => 0.378764
20 => 0.411063
21 => 0.443296
22 => 0.475287
23 => 0.506876
24 => 0.537913
25 => 0.568260
26 => 0.597796
27 => 0.626412
28 => 0.654014
29 => 0.680524
30 => 0.705877
31 => 0.730022
32 => 0.752924
33 => 0.774560
34 => 0.794917
35 => 0.813998
36 => 0.831812
37 => 0.848381
38 => 0.863732
39 => 0.877901
40 => 0.890932
41 => 0.902870
42 => 0.913767
43 => 0.923678
44 => 0.932658
45 => 0.940766
46 => 0.948060
47 => 0.954598
48 => 0.960437
49 => 0.965634
50 => 0.970242
51 => 0.974313
52 => 0.977898
53 => 0.981043
54 => 0.983792
55 => 0.986187
56 => 0.988266
57 => 0.990064
58 => 0.991614
59 => 0.992945
60 => 0.994084
61 => 0.995055
62 => 0.995880
63 => 0.996579
64 => 0.997169
65 => 0.997665
66 => 0.998080
67 => 0.998427
68 => 0.998715
69 => 0.998954
70 => 0.999152
71 => 0.999314
72 => 0.999447
73 => 0.999556
74 => 0.999645
75 => 0.999717
76 => 0.999775
77 => 0.999822
78 => 0.999859
79 => 0.999889
80 => 0.999913
81 => 0.999932
82 => 0.999947
83 => 0.999959
84 => 0.999968
85 => 0.999976
86 => 0.999981
87 => 0.999986
88 => 0.999989
89 => 0.999992
90 => 0.999994
91 => 0.999995
92 => 0.999996
93 => 0.999997
94 => 0.999998
95 => 0.999999
96 => 0.999999
97 => 0.999999
98 => 0.999999
99 => 1.000000

(Oczywiście, P (99) wynosi tylko 1,0 z powodu zaokrąglenia. Prawdopodobieństwo nie osiągnie dokładnie 1,0, dopóki P (367) .)

Pyth, 31 34 bajtów

Jt.2425K366-1c.xX0rK-KQ*JQ^+KJQ

Demonstracja , uprząż testowa

Działa to podobnie do starej wersji, z tą różnicą, że zamiast osobnego generowania wartości (366 + n * (.2425 - 1)) i mnożenia jej, zaczyna się od wygenerowania listy rozciągającej się od 366 do 365 - n + 2, następnie modyfikuje 366 na miejsce (366 + n * (.2425 - 1)), a następnie bierze wynik z listy. Również zamiast 366 i .2425 używane są stałe 366 i -.7575.

Stara wersja:

Pyth, 34 bajty

J.2425K365-1c*+hK*QtJ.xrK-hKQ^+KJQ

Istnieją dwa sposoby, aby nie było pary osób z tymi samymi urodzinami: każdy ma inne urodziny, i nikt nie ma urodzin 29 lutego, i ktoś ma urodziny 29 lutego, a wszyscy inni mają inne urodziny urodziny w normalne dni.

Prawdopodobieństwo pierwszego wystąpienia to (365 * 364 * ... 365 - n + 1) / (365.2425 ^ n).

Prawdopodobieństwo drugiego wystąpienia to (365 * 364 * ... 365 - n + 2) * .2425 * n / (365.2425 ^ n)

Można je zapisać razem jako (365 * 364 * ... 365 - n + 2) * (365 - n + 1 + .2425 * n) / (365.2425 ^ n) = (365 * 364 * ... 365 - n + 2) * (365 + 1 + (.2425-1) * n) / (365.2425 ^ n).

Jest to prawdopodobieństwo braku par, więc prawdopodobieństwo co najmniej jednej pary wynosi jeden minus powyższa liczba.

J = .2425
K = 365
.xrK-hKQ = (365 * 364 * ... 365 - n + 2)
+hK*QtJ = (365 + 1 + n * (.2425 - 1))
^+KJQ = (365.2425 ^ n)

Python, 179 178 144 143 140 136 136 133 133

f=.2425
g=365+f
a=lambda n:(n and(365-n)*a(n-1)or 365)/g
b=lambda n:(n<2and f or(367-n)*b(n-1)+a(n-2)*f)/g
p=lambda n:1-a(n-1)-b(n)

p(n)daje wynik. Zmień, .2425aby fractions.Fraction(97,400)uzyskać dokładny wynik.

Python 155 153 151 142 140 bajtów

d=146097
b=d/400
c=97/d
e=lambda n:n<2and 1-97/d or e(n-1)*(366-n)/b
f=lambda n:n<2and c or f(n-1)*(367-n)/b+e(n-1)*c
a=lambda n:1-e(n)-f(n)

Zadzwoń a(n)po wynik. Aby uzyskać dokładne wyniki, zawiń dułamek.

Przetestuj tutaj

Używa tej samej techniki, co tutaj , ale ze zmodyfikowanymi stałymi.

Kod maszynowy 80386, 43 bajty

Hexdump kodu:

68 75 6e 33 3b 68 5a eb 07 3b 8b c4 49 d9 e8 d9
e8 d8 00 d9 e8 d9 40 04 de ea d8 c9 d9 00 de eb
e2 f3 d8 ca d9 e8 d8 e1 58 58 c3

Zacząłem od następującego wzoru (dla prawdopodobieństwa uzupełniającego):

(tutaj d = 366 * 400 - 303jest liczba dni w 400 lat)

Oto kod c ++, który go implementuje (jest już trochę zoptymalizowany):

double it_opt(int k)
{
    double d = 366 * 400 - 303; // number of days in 400 years
    double result = 1; // probability of no coincidences
    const float const1 = float(400 / d);
    const float const2 = float(303 / d);
    double v1 = 1 + const2;
    double v2 = 1;

    for (int i = 0; i < k - 1; ++i)
    {
        v1 -= const1;
        result *= v1;
        v2 -= const2;
    }
    result *= v2;
    return 1 - result;
}

Kod jest tak ułożony, że wymaga minimalnej liczby stałych - tylko dwóch ( 400 / di 303 / d). Używam tego floattypu do ich reprezentowania, ponieważ zajmuje mniej miejsca (4 bajty na stałą). Ponadto nie chciałem pomnożyć const2przez k - 1(ponieważ wymagałoby to konwersji k - 1nafloat ); const2zamiast tego kod odejmuje się wielokrotnie.

Oto lista języków asemblera:

    ; // fastcall convention - parameter k is in ecx
    ; // result must be returned in st
    push dword ptr 0x3b336e75; // const1 = [esp + 4]
    push dword ptr 0x3b07eb5a; // const2 = [esp]
    mov eax, esp;              // use [eax] instead of [esp] - 1 byte less
    dec ecx;                   // calculate k - 1
    fld1;                      // initiaze result = 1
    fld1;                      // initialize v1
    fadd [eax];
    fld1;                      // initialilze v2
myloop:
    fld dword ptr [eax + 4];
    fsubp st(2), st;            // update v1
    fmul st, st(1);             // update result
    fld dword ptr [eax];
    fsubp st(3), st;            // update v2
    loop myloop;                // loop
    fmul st, st(2);             // update result by v2
    fld1;
    fsub st, st(1);             // complement result
    pop eax;                    // restore stack
    pop eax;                    // restore stack
    ret;                        // return