For each non-square integer \(d \lt 100\), we approximate \(\pi\) by numbers of the form
$$a+b\sqrt d,\qquad |a|,|b|\le 10^{13}.$$
Among all such pairs \((a,b)\), we choose the one that minimizes
$$\left|\pi-\left(a+b\sqrt d\right)\right|.$$
The quantity called the integral part is the coefficient \(a\), not the floor of the real number \(a+b\sqrt d\). The required output is therefore
$$\sum_{d\in S} |a_d|,$$
where \(S\) is the set of non-square integers from \(2\) to \(99\), and \(a_d\) is the integer coefficient belonging to the best approximation for that specific \(d\).
Write \(\alpha=\sqrt d\) and \(\tau=\pi-3\). The central observation is that once \(b\) is fixed, the optimal \(a\) is forced to be the nearest integer to \(\pi-b\alpha\). So the real difficulty is to find the right \(b\) without scanning all \(2\cdot 10^{13}+1\) possibilities.
For fixed \(b\), define
$$x_b=\pi-b\alpha.$$
Then minimizing \(\left|\pi-(a+b\alpha)\right|\) over all admissible \(a\) is the same as minimizing
$$|x_b-a|$$
over integers \(a\) with \(|a|\le 10^{13}\). Therefore the best \(a\) is the nearest integer to \(x_b\), and the two-variable problem becomes
$$\min_{|b|\le 10^{13}} \left\| \pi-b\alpha \right\|_{\mathbb Z},$$
where \(\|y\|_{\mathbb Z}\) means the distance from \(y\) to the nearest integer.
There is also a strong bound on \(b\). The trivial candidate \(a=3\), \(b=0\) has error \(\tau=\pi-3\lt 1\), so any true optimum must have error below \(1\). If \(|b|\alpha \gt 10^{13}+\pi+1\), then for every admissible \(a\),
$$\left|a+b\alpha-\pi\right|\ge |b|\alpha-|a|-\pi \gt 1,$$
which cannot be optimal. The implementation therefore uses the safe bound
$$|b|\le B,\qquad B=\min\!\left(10^{13},\left\lfloor\frac{10^{13}+\pi+2}{\alpha}\right\rfloor\right).$$
Because \(\alpha=\sqrt d\) is a quadratic irrational, its continued fraction is periodic, and its convergents \(p/q\) give exceptionally good rational approximations. For each convergent we write
$$q\alpha=p+\delta,$$
so \(\delta=q\alpha-p\) is very small.
This matters because changing \(b\) by \(q\) changes \(b\alpha\) by almost the integer \(p\):
$$ (b+q)\alpha=b\alpha+p+\delta.$$
So within the arithmetic progression \(b=b_0+tq\), the fractional part of \(b\alpha\) drifts only by the tiny amount \(t\delta\). That turns a huge search over \(b\) into a very small local search around each convergent.
Since every convergent satisfies \(\gcd(p,q)=1\), the inverse of \(p\) modulo \(q\) exists. For each residue \(j\in\{0,1,\dots,q-1\}\), choose \(b_0\) so that
$$b_0p\equiv j \pmod q.$$
Then \(b_0p=j+kq\) for some integer \(k\). If we now write
$$b=b_0+tq,$$
we get
$$\begin{aligned} b\alpha &=b_0\alpha+tq\alpha \\ &=\frac{b_0(p+\delta)}{q}+t(p+\delta) \\ &=k+tp+\frac{j}{q}+\frac{b_0\delta}{q}+t\delta. \end{aligned}$$
Since \(\pi=3+\tau\), it follows that
$$\pi-b\alpha=\left(3-k-tp\right)+\left(\tau-\frac{j}{q}-\frac{b_0\delta}{q}-t\delta\right).$$
The first parenthesis is an integer, so the quality of the approximation is controlled entirely by the small residual
$$\operatorname{diff}(t)=\frac{j}{q}+\frac{b_0\delta}{q}-\tau+t\delta.$$
We want this quantity to be as close to \(0\) as possible.
For fixed \(q\), \(p\), and residue class \(j\), the expression \(\operatorname{diff}(t)\) is linear in \(t\). Therefore the best \(t\) is near the real minimizer
$$t\approx -\frac{\frac{j}{q}+\frac{b_0\delta}{q}-\tau}{\delta}.$$
Only residue classes with \(j/q\) near \(\tau\) can compete, because \(|b_0\delta/q|\le |\delta|\) is tiny for a useful convergent. The implementations therefore inspect only a small fixed window of \(j\)-values around \(\tau q\), and for each such class they inspect only a small fixed window of \(t\)-values around the estimate above.
The admissible \(t\)-range is also clipped by \(|b|\le B\). Its endpoints are checked explicitly so that a best solution on the boundary is never missed.
Once a candidate \(b\) is known, we compute
$$x_b=\pi-b\alpha$$
and take the nearest integer \(a_0\) to \(x_b\). Because the computation uses finite-precision arithmetic, the implementations test the three integers
$$a_0-1,\qquad a_0,\qquad a_0+1$$
instead of trusting one rounding step blindly. This guarantees that the true nearest admissible integer is captured even if \(x_b\) lies very close to a half-integer.
The best pair is the one with smallest absolute error. If two candidates are numerically tied, the implementations keep the one with smaller \(|a|\), and then smaller \(|b|\), to make the output deterministic.
Here \(\alpha=\sqrt 2\approx 1.414213562\) and \(\tau=\pi-3\approx 0.141592654\).
Take the convergent
$$\frac{p}{q}=\frac{3}{2},$$
for which
$$\delta=2\sqrt 2-3\approx -0.171572875.$$
Since \(\tau q\approx 0.283\), the most relevant residue is \(j=0\). Because \(p\equiv 1\pmod 2\), the corresponding base residue is \(b_0\equiv 0\pmod 2\), so we take \(b_0=0\).
Then
$$\operatorname{diff}(t)=-\tau+t\delta.$$
The linear estimate gives
$$t\approx -\frac{-\tau}{\delta}\approx -0.825,$$
so \(t=-1\) is the natural nearby candidate. This yields
$$b=b_0+tq=-2.$$
Now compute
$$\pi-b\sqrt 2=\pi+2\sqrt 2\approx 5.969019778,$$
whose nearest integer is \(a=6\). Hence the approximation is
$$6-2\sqrt 2\approx 3.171572875,$$
with error about
$$|\,\pi-(6-2\sqrt 2)\,|\approx 0.029980221.$$
This matches the hard-coded sanity check used by the implementations.
The C++, Python, and Java implementations all follow the same mathematical search. For each non-square \(d\lt 100\), they compute \(\sqrt d\) with high precision, derive the safe bound \(B\), and generate continued-fraction convergents until the denominator exceeds that bound.
For every convergent, the implementation forms the small error term \(\delta=q\sqrt d-p\), looks only at residue classes \(j\) near \((\pi-3)q\), converts each class into a base value \(b_0\) by modular inversion, and derives the permitted interval of \(t\)-values from the constraint \(|b|\le B\).
Instead of scanning every possible \(t\), the implementation evaluates only a short interval around the linear estimate for the best \(t\), plus the interval endpoints. Each resulting \(b\) is converted into a nearby integer coefficient \(a\), and the error is checked explicitly.
The C++ version distributes the independent \(d\)-values across worker threads. The Python version performs the same search directly with decimal arithmetic. The Java version delegates execution to the same compiled numeric core, so its output agrees with the C++ behavior.
For fixed \(d\), the denominators of the convergents of a quadratic irrational grow exponentially, so only \(O(\log B)\) convergents are needed before the denominator passes the search bound. The inspected windows in \(j\) and \(t\) are fixed-size constants, so each convergent contributes only \(O(1)\) candidate evaluations.
Therefore the practical running time per \(d\) is \(O(\log B)\), which here is effectively \(O(\log N)\), and the memory usage is \(O(\log B)\) if the convergents are stored explicitly. Since there are only finitely many non-square values \(d\lt 100\), the total cost is this per-\(d\) work multiplied by a small constant, with straightforward parallelization in C++.
Für jede nichtquadratische Zahl \(d \lt 100\) soll \(\pi\) durch Ausdrücke der Form
$$a+b\sqrt d,\qquad |a|,|b|\le 10^{13}$$
möglichst gut approximiert werden. Gesucht ist also das Paar \((a,b)\), das den Fehler
$$\left|\pi-\left(a+b\sqrt d\right)\right|$$
minimiert.
Der sogenannte Integralteil ist hier nicht der ganzzahlige Anteil des reellen Werts, sondern einfach der Koeffizient \(a\). Am Ende wird daher
$$\sum_{d\in S} |a_d|$$
über alle zulässigen \(d\) aufsummiert, wobei \(S\) die Menge aller nichtquadratischen Zahlen von \(2\) bis \(99\) bezeichnet.
Schreibe \(\alpha=\sqrt d\) und \(\tau=\pi-3\). Die Schlüsselerkenntnis ist: Sobald \(b\) feststeht, ist \(a\) praktisch erzwungen, denn der beste Wert ist einfach die ganze Zahl, die \(\pi-b\alpha\) am nächsten liegt. Damit muss nicht mehr in zwei Dimensionen gesucht werden, sondern nur noch nach den richtigen \(b\)-Werten.
Für festes \(b\) setze
$$x_b=\pi-b\alpha.$$
Dann ist die Minimierung von \(\left|\pi-(a+b\alpha)\right|\) über alle erlaubten \(a\) genau dasselbe wie die Minimierung von
$$|x_b-a|.$$
Also ist \(a\) die nächste ganze Zahl zu \(x_b\), und das ursprüngliche Problem wird zu
$$\min_{|b|\le 10^{13}} \left\| \pi-b\alpha \right\|_{\mathbb Z}.$$
Dabei bezeichnet \(\|y\|_{\mathbb Z}\) den Abstand von \(y\) zur nächstgelegenen ganzen Zahl.
Außerdem lässt sich der Suchraum für \(b\) stark verkleinern. Der triviale Kandidat \(a=3\), \(b=0\) hat den Fehler \(\tau=\pi-3\lt 1\). Ein optimales Paar muss also Fehler kleiner als \(1\) haben. Falls \(|b|\alpha \gt 10^{13}+\pi+1\), gilt für jedes \(|a|\le 10^{13}\),
$$\left|a+b\alpha-\pi\right|\ge |b|\alpha-|a|-\pi \gt 1,$$
also kann ein solcher \(b\) nicht optimal sein. Deshalb verwendet die Implementierung die sichere Schranke
$$|b|\le B,\qquad B=\min\!\left(10^{13},\left\lfloor\frac{10^{13}+\pi+2}{\alpha}\right\rfloor\right).$$
Da \(\alpha=\sqrt d\) eine quadratische Irrationalzahl ist, besitzt sie einen periodischen Kettenbruch. Ihre Konvergenten \(p/q\) liefern besonders gute rationale Approximationen. Für jeden Konvergenten schreiben wir
$$q\alpha=p+\delta,$$
wobei \(\delta=q\alpha-p\) sehr klein ist.
Das ist genau der nützliche Effekt: Wenn man \(b\) um \(q\) erhöht, ändert sich \(b\alpha\) fast um die ganze Zahl \(p\), nämlich
$$ (b+q)\alpha=b\alpha+p+\delta.$$
Innerhalb einer Progression \(b=b_0+tq\) wandert der gebrochene Anteil von \(b\alpha\) also nur langsam, mit Schrittweite \(\delta\). Dadurch wird die Suche lokal und beherrschbar.
Für Konvergenten gilt stets \(\gcd(p,q)=1\), daher existiert das Inverse von \(p\) modulo \(q\). Für jede Restklasse \(j\in\{0,1,\dots,q-1\}\) wählen wir \(b_0\) so, dass
$$b_0p\equiv j \pmod q.$$
Dann gibt es ein ganzzahliges \(k\) mit \(b_0p=j+kq\). Für
$$b=b_0+tq$$
folgt
$$\begin{aligned} b\alpha &=b_0\alpha+tq\alpha \\ &=\frac{b_0(p+\delta)}{q}+t(p+\delta) \\ &=k+tp+\frac{j}{q}+\frac{b_0\delta}{q}+t\delta. \end{aligned}$$
Mit \(\pi=3+\tau\) erhält man
$$\pi-b\alpha=\left(3-k-tp\right)+\left(\tau-\frac{j}{q}-\frac{b_0\delta}{q}-t\delta\right).$$
Der erste Teil ist ganzzahlig. Die Approximationsgüte wird daher vollständig durch
$$\operatorname{diff}(t)=\frac{j}{q}+\frac{b_0\delta}{q}-\tau+t\delta$$
bestimmt. Diese Größe soll möglichst nahe bei \(0\) liegen.
Für feste Daten \(q\), \(p\) und \(j\) ist \(\operatorname{diff}(t)\) linear in \(t\). Das optimale \(t\) liegt also nahe bei
$$t\approx -\frac{\frac{j}{q}+\frac{b_0\delta}{q}-\tau}{\delta}.$$
Nur Restklassen mit \(j/q\) in der Nähe von \(\tau\) können überhaupt konkurrenzfähig sein, denn \(|b_0\delta/q|\le |\delta|\) ist für gute Konvergenten sehr klein. Die Implementierungen prüfen deshalb nur ein kleines festes Fenster von \(j\)-Werten um \(\tau q\) und für jede solche Klasse nur ein kleines Fenster von \(t\)-Werten um die lineare Schätzung.
Zusätzlich wird das durch \(|b|\le B\) vorgegebene Intervall für \(t\) an beiden Enden explizit getestet, damit auch Randfälle korrekt erfasst werden.
Sobald ein Kandidat \(b\) vorliegt, berechnen wir wieder
$$x_b=\pi-b\alpha$$
und nehmen die nächstgelegene ganze Zahl \(a_0\). Weil mit endlicher Präzision gerechnet wird, verlassen sich die Implementierungen nicht blind auf eine einzige Rundung, sondern testen
$$a_0-1,\qquad a_0,\qquad a_0+1.$$
So wird die wahre beste ganze Zahl auch dann gefunden, wenn \(x_b\) sehr nahe an einer Halbzahl liegt.
Behalten wird der Kandidat mit dem kleinsten absoluten Fehler. Bei numerischem Gleichstand gewinnt zuerst das kleinere \(|a|\), danach das kleinere \(|b|\), damit die Ausgabe deterministisch bleibt.
Hier ist \(\alpha=\sqrt 2\approx 1.414213562\) und \(\tau=\pi-3\approx 0.141592654\).
Wir nehmen den Konvergenten
$$\frac{p}{q}=\frac{3}{2},$$
also
$$\delta=2\sqrt 2-3\approx -0.171572875.$$
Weil \(\tau q\approx 0.283\), ist \(j=0\) die naheliegendste Restklasse. Da \(p\equiv 1\pmod 2\), gehört dazu \(b_0\equiv 0\pmod 2\), also wählen wir \(b_0=0\).
Dann lautet der Restterm
$$\operatorname{diff}(t)=-\tau+t\delta.$$
Die lineare Schätzung ergibt
$$t\approx -\frac{-\tau}{\delta}\approx -0.825,$$
also ist \(t=-1\) der natürliche Kandidat in der Nähe. Damit erhält man
$$b=b_0+tq=-2.$$
Nun ist
$$\pi-b\sqrt 2=\pi+2\sqrt 2\approx 5.969019778,$$
also ist die nächstgelegene ganze Zahl \(a=6\). Die Approximation lautet damit
$$6-2\sqrt 2\approx 3.171572875,$$
und der Fehler ist ungefähr
$$|\,\pi-(6-2\sqrt 2)\,|\approx 0.029980221.$$
Genau dieses Paar erscheint auch in der eingebauten Plausibilitätsprüfung der Implementierungen.
Die C++-, Python- und Java-Implementierungen folgen alle derselben mathematischen Suche. Für jedes nichtquadratische \(d\lt 100\) berechnen sie \(\sqrt d\) mit hoher Präzision, bestimmen die sichere Schranke \(B\) und erzeugen Kettenbruch-Konvergenten, bis deren Nenner diese Schranke überschreitet.
Für jeden Konvergenten wird der kleine Fehlerterm \(\delta=q\sqrt d-p\) gebildet. Danach betrachtet die Implementierung nur Restklassen \(j\) in der Nähe von \((\pi-3)q\), bestimmt daraus per modularem Inversen einen Basiswert \(b_0\) und leitet aus \(|b|\le B\) das zulässige Intervall für \(t\) ab.
Anstatt alle \(t\)-Werte zu durchsuchen, werden nur wenige Werte um die lineare Schätzung sowie die beiden Intervallenden geprüft. Jeder so gefundene \(b\)-Kandidat liefert einen reellen Zielwert \(\pi-b\sqrt d\), aus dem die nächstgelegenen ganzzahligen \(a\)-Kandidaten explizit getestet werden.
Die C++-Version verteilt die unabhängigen \(d\)-Werte auf mehrere Threads. Die Python-Version führt dieselbe Suche direkt mit Dezimalarithmetik aus. Die Java-Version delegiert an denselben kompilierten numerischen Kern und reproduziert daher das Verhalten der C++-Version.
Für festes \(d\) wachsen die Nenner der Konvergenten einer quadratischen Irrationalzahl exponentiell. Bis die Nennergrenze \(B\) erreicht wird, benötigt man daher nur \(O(\log B)\) Konvergenten. Die Fenster für \(j\) und \(t\) haben feste Größe, also kostet jeder Konvergente nur \(O(1)\) Kandidatenprüfungen.
Damit liegt die praktische Laufzeit pro \(d\) bei \(O(\log B)\), also hier effektiv bei \(O(\log N)\). Der Speicherbedarf ist \(O(\log B)\), wenn die Konvergenten explizit gespeichert werden. Da es nur endlich viele nichtquadratische \(d\lt 100\) gibt, ist der Gesamtaufwand dieser per-\(d\)-Arbeit mal einem kleinen konstanten Faktor, wobei C++ die äußere Schleife parallel ausführen kann.
\(d \lt 100\) olan ve tam kare olmayan her tamsayı için, \(\pi\) sayısını
$$a+b\sqrt d,\qquad |a|,|b|\le 10^{13}$$
biçimindeki sayılarla en iyi şekilde yaklaştırmak istiyoruz. Yani amaç,
$$\left|\pi-\left(a+b\sqrt d\right)\right|$$
değerini en küçük yapan \((a,b)\) çiftini bulmaktır.
Buradaki integral kısım, gerçek sayının tabanı değil, doğrudan \(a\) katsayısıdır. Bu yüzden sonuç olarak
$$\sum_{d\in S} |a_d|$$
toplamı istenir; burada \(S\), \(2\) ile \(99\) arasındaki tam kare olmayan sayıların kümesidir.
\(\alpha=\sqrt d\) ve \(\tau=\pi-3\) yazalım. Temel fikir şudur: \(b\) sabitlenince en iyi \(a\), \(\pi-b\alpha\)'ya en yakın tam sayıdır. Dolayısıyla iki değişkenli kaba arama yerine, doğru \(b\) değerlerini akıllıca seçmek yeterlidir.
Sabit bir \(b\) için
$$x_b=\pi-b\alpha$$
tanımını yapalım. O zaman \(\left|\pi-(a+b\alpha)\right|\) ifadesini en küçük yapmak,
$$|x_b-a|$$
ifadesini en küçük yapmakla aynıdır. Bu yüzden en iyi \(a\), \(x_b\)'ye en yakın tam sayıdır ve problem
$$\min_{|b|\le 10^{13}} \left\| \pi-b\alpha \right\|_{\mathbb Z}$$
şekline dönüşür. Burada \(\|y\|_{\mathbb Z}\), \(y\)'nin en yakın tam sayıya uzaklığıdır.
\(b\) için güçlü bir üst sınır da vardır. Basit aday \(a=3\), \(b=0\), hata olarak \(\tau=\pi-3\lt 1\) verir. Dolayısıyla gerçek optimumun hatası mutlaka \(1\)'den küçüktür. Eğer \(|b|\alpha \gt 10^{13}+\pi+1\) ise her \(|a|\le 10^{13}\) için
$$\left|a+b\alpha-\pi\right|\ge |b|\alpha-|a|-\pi \gt 1$$
olur; yani böyle bir \(b\) optimum olamaz. Bu nedenle uygulama güvenli sınır olarak
$$|b|\le B,\qquad B=\min\!\left(10^{13},\left\lfloor\frac{10^{13}+\pi+2}{\alpha}\right\rfloor\right)$$
kullanır.
\(\alpha=\sqrt d\) kuadratik irrasyonel olduğu için sürekli kesri periyodiktir. Bu sürekli kesrin yakınsakları \(p/q\), \(\alpha\)'yı olağanüstü iyi rasyonel sayılarla yaklaştırır. Her yakınsak için
$$q\alpha=p+\delta$$
yazarız; burada \(\delta=q\alpha-p\) çok küçüktür.
Bu şu anlama gelir: \(b\) değerini \(q\) kadar artırmak, \(b\alpha\)'yı neredeyse tam sayı olan \(p\) kadar artırır:
$$ (b+q)\alpha=b\alpha+p+\delta.$$
Yani \(b=b_0+tq\) aritmetik dizisinde, \(b\alpha\)'nın kesir kısmı yalnızca küçük \(t\delta\) adımlarıyla kayar. Böylece devasa arama alanı, her yakınsak çevresinde küçük yerel aramalara dönüşür.
Yakınsaklarda her zaman \(\gcd(p,q)=1\) olduğu için, \(p\)'nin \(q\)'ya göre tersi vardır. \(j\in\{0,1,\dots,q-1\}\) için \(b_0\) öyle seçilir ki
$$b_0p\equiv j \pmod q.$$
Böylece bir tam sayı \(k\) için \(b_0p=j+kq\) yazılabilir. Şimdi
$$b=b_0+tq$$
yazarsak,
$$\begin{aligned} b\alpha &=b_0\alpha+tq\alpha \\ &=\frac{b_0(p+\delta)}{q}+t(p+\delta) \\ &=k+tp+\frac{j}{q}+\frac{b_0\delta}{q}+t\delta \end{aligned}$$
elde edilir. \(\pi=3+\tau\) olduğu için
$$\pi-b\alpha=\left(3-k-tp\right)+\left(\tau-\frac{j}{q}-\frac{b_0\delta}{q}-t\delta\right).$$
İlk parantez tam sayıdır. Dolayısıyla yaklaşımın kalitesini belirleyen asıl terim
$$\operatorname{diff}(t)=\frac{j}{q}+\frac{b_0\delta}{q}-\tau+t\delta$$
olur ve bu değeri \(0\)'a mümkün olduğunca yaklaştırmak isteriz.
\(\operatorname{diff}(t)\), sabit \(q\), \(p\) ve \(j\) için \(t\)'de doğrusal olduğundan en iyi \(t\) yaklaşık olarak
$$t\approx -\frac{\frac{j}{q}+\frac{b_0\delta}{q}-\tau}{\delta}$$
noktasındadır. Sadece \(j/q\) değeri \(\tau\)'ya yakın olan sınıflar yarışabilir; çünkü \(|b_0\delta/q|\le |\delta|\) iyi yakınsaklarda zaten çok küçüktür. Bu yüzden uygulamalar yalnızca \(\tau q\) çevresindeki küçük sabit bir \(j\) penceresini ve her sınıf için bu tahminin çevresindeki küçük sabit bir \(t\) penceresini inceler.
\(|b|\le B\) kısıtından gelen \(t\) aralığının iki ucu da ayrıca denenir; böylece optimum sınırda kalıyorsa kaçırılmaz.
Bir \(b\) adayı bulunduğunda tekrar
$$x_b=\pi-b\alpha$$
hesaplanır ve buna en yakın tam sayı \(a_0\) alınır. Sonlu hassasiyet nedeniyle uygulamalar tek bir yuvarlamaya güvenmez; onun yerine
$$a_0-1,\qquad a_0,\qquad a_0+1$$
değerlerinin hepsini test eder. Böylece \(x_b\) yarım tamsayıya çok yakın olsa bile gerçek en iyi tam sayı yakalanır.
En küçük mutlak hatayı veren çift tutulur. Hata eşitse önce \(|a|\) küçük olan, sonra \(|b|\) küçük olan aday seçilerek sonuç deterministik hale getirilir.
Bu durumda \(\alpha=\sqrt 2\approx 1.414213562\) ve \(\tau=\pi-3\approx 0.141592654\).
Yakınsak olarak
$$\frac{p}{q}=\frac{3}{2}$$
alalım. O zaman
$$\delta=2\sqrt 2-3\approx -0.171572875.$$
\(\tau q\approx 0.283\) olduğundan en anlamlı sınıf \(j=0\)'dır. Ayrıca \(p\equiv 1\pmod 2\) olduğu için ilgili taban kalıntı \(b_0\equiv 0\pmod 2\), yani \(b_0=0\) seçilir.
Böylece
$$\operatorname{diff}(t)=-\tau+t\delta$$
olur. Doğrusal tahmin
$$t\approx -\frac{-\tau}{\delta}\approx -0.825$$
verir; dolayısıyla yakın aday \(t=-1\)'dir. Buradan
$$b=b_0+tq=-2$$
çıkar.
Şimdi
$$\pi-b\sqrt 2=\pi+2\sqrt 2\approx 5.969019778$$
olur ve en yakın tam sayı \(a=6\)'dır. Sonuçta yaklaşım
$$6-2\sqrt 2\approx 3.171572875$$
ve hata
$$|\,\pi-(6-2\sqrt 2)\,|\approx 0.029980221$$
elde edilir. Bu çift, uygulamaların içine yerleştirilmiş sağlamlık kontrolüyle de aynıdır.
C++, Python ve Java uygulamalarının üçü de aynı matematiksel aramayı izler. Her tam kare olmayan \(d\lt 100\) için yüksek hassasiyetli \(\sqrt d\) hesaplanır, güvenli \(B\) sınırı kurulur ve paydası bu sınırı aşana kadar sürekli kesir yakınsakları üretilir.
Her yakınsak için küçük hata terimi \(\delta=q\sqrt d-p\) hesaplanır. Ardından uygulama yalnızca \((\pi-3)q\) civarındaki birkaç \(j\) sınıfına bakar, modüler ters yardımıyla bunları birer \(b_0\) tabanına dönüştürür ve \(|b|\le B\) kısıtından izinli \(t\) aralığını çıkarır.
Tüm \(t\) değerlerini gezmek yerine, yalnızca doğrusal tahminin çevresindeki kısa aralık ve uç noktalar test edilir. Elde edilen her \(b\) için uygun \(a\) değerleri açıkça denenir ve hata doğrudan ölçülür.
C++ sürümü bağımsız \(d\) değerlerini iş parçacıklarına dağıtır. Python sürümü aynı aramayı ondalık aritmetikle doğrudan yürütür. Java sürümü ise aynı derlenmiş sayısal çekirdeğe devreder; bu yüzden çıktısı C++ davranışıyla aynıdır.
Sabit bir \(d\) için, kuadratik irrasyonellerin yakınsak paydaları üstel büyür. Bu yüzden arama sınırı \(B\)'yi aşmak için yalnızca \(O(\log B)\) yakınsak gerekir. \(j\) ve \(t\) pencereleri sabit boyutlu olduğundan her yakınsak yalnızca \(O(1)\) kadar aday üretir.
Dolayısıyla pratik çalışma süresi her \(d\) için \(O(\log B)\), yani burada etkin olarak \(O(\log N)\) olur. Yakınsaklar açıkça saklanırsa bellek kullanımı \(O(\log B)\) düzeyindedir. \(d\lt 100\) için yalnızca sonlu sayıda tam kare olmayan değer bulunduğundan toplam maliyet, küçük bir sabitle çarpılmış bu per-\(d\) maliyetidir; C++ sürümü bunu kolayca paralelleştirir.
Para cada entero no cuadrado \(d \lt 100\), queremos aproximar \(\pi\) mediante números de la forma
$$a+b\sqrt d,\qquad |a|,|b|\le 10^{13}.$$
Entre todos esos pares \((a,b)\), debemos elegir el que minimiza
$$\left|\pi-\left(a+b\sqrt d\right)\right|.$$
Aquí la parte integral no significa el piso del número real \(a+b\sqrt d\), sino simplemente el coeficiente \(a\). Por lo tanto, la cantidad pedida es
$$\sum_{d\in S} |a_d|,$$
donde \(S\) es el conjunto de enteros no cuadrados entre \(2\) y \(99\), y \(a_d\) es el coeficiente entero de la mejor aproximación correspondiente a ese valor de \(d\).
Escribamos \(\alpha=\sqrt d\) y \(\tau=\pi-3\). La idea clave es que, una vez fijado \(b\), el mejor \(a\) ya no requiere búsqueda complicada: es simplemente el entero más cercano a \(\pi-b\alpha\). Así, el reto real consiste en identificar los valores prometedores de \(b\) sin recorrer todo el intervalo hasta \(10^{13}\).
Para un \(b\) fijo definimos
$$x_b=\pi-b\alpha.$$
Minimizar \(\left|\pi-(a+b\alpha)\right|\) sobre los \(a\) permitidos equivale a minimizar
$$|x_b-a|.$$
Por eso el mejor \(a\) es el entero más cercano a \(x_b\), y el problema bidimensional se convierte en
$$\min_{|b|\le 10^{13}} \left\| \pi-b\alpha \right\|_{\mathbb Z},$$
donde \(\|y\|_{\mathbb Z}\) representa la distancia de \(y\) al entero más cercano.
Además, el rango útil de \(b\) se puede recortar mucho. El candidato trivial \(a=3\), \(b=0\) tiene error \(\tau=\pi-3\lt 1\), así que el óptimo verdadero debe tener error menor que \(1\). Si \(|b|\alpha \gt 10^{13}+\pi+1\), entonces para cualquier \(|a|\le 10^{13}\),
$$\left|a+b\alpha-\pi\right|\ge |b|\alpha-|a|-\pi \gt 1,$$
y por tanto ese \(b\) no puede ser óptimo. La implementación usa entonces la cota segura
$$|b|\le B,\qquad B=\min\!\left(10^{13},\left\lfloor\frac{10^{13}+\pi+2}{\alpha}\right\rfloor\right).$$
Como \(\alpha=\sqrt d\) es una irracional cuadrática, su fracción continua es periódica. Sus convergentes \(p/q\) proporcionan aproximaciones racionales excepcionalmente buenas. Para cada convergente escribimos
$$q\alpha=p+\delta,$$
de modo que \(\delta=q\alpha-p\) es muy pequeño.
Esto es exactamente lo útil: al aumentar \(b\) en \(q\), el valor \(b\alpha\) cambia casi en un entero \(p\):
$$ (b+q)\alpha=b\alpha+p+\delta.$$
Por ello, dentro de una progresión aritmética \(b=b_0+tq\), la parte fraccionaria de \(b\alpha\) se mueve lentamente, con paso \(t\delta\). Esa lentitud convierte una búsqueda enorme en varias búsquedas locales muy pequeñas.
Todo convergente cumple \(\gcd(p,q)=1\), por lo que \(p\) tiene inverso módulo \(q\). Para cada residuo \(j\in\{0,1,\dots,q-1\}\), elegimos \(b_0\) tal que
$$b_0p\equiv j \pmod q.$$
Entonces existe un entero \(k\) con \(b_0p=j+kq\). Si escribimos
$$b=b_0+tq,$$
obtenemos
$$\begin{aligned} b\alpha &=b_0\alpha+tq\alpha \\ &=\frac{b_0(p+\delta)}{q}+t(p+\delta) \\ &=k+tp+\frac{j}{q}+\frac{b_0\delta}{q}+t\delta. \end{aligned}$$
Como \(\pi=3+\tau\), se sigue que
$$\pi-b\alpha=\left(3-k-tp\right)+\left(\tau-\frac{j}{q}-\frac{b_0\delta}{q}-t\delta\right).$$
El primer paréntesis es entero. Por tanto, la calidad de la aproximación depende por completo del residuo pequeño
$$\operatorname{diff}(t)=\frac{j}{q}+\frac{b_0\delta}{q}-\tau+t\delta,$$
que queremos acercar lo más posible a \(0\).
Para valores fijos de \(q\), \(p\) y \(j\), la función \(\operatorname{diff}(t)\) es lineal en \(t\). Por eso el mejor \(t\) está cerca de
$$t\approx -\frac{\frac{j}{q}+\frac{b_0\delta}{q}-\tau}{\delta}.$$
Solo las clases con \(j/q\) cercano a \(\tau\) pueden competir, porque \(|b_0\delta/q|\le |\delta|\) es diminuto para un convergente útil. Por eso las implementaciones inspeccionan únicamente una pequeña ventana fija de valores \(j\) alrededor de \(\tau q\), y para cada una solo una pequeña ventana fija de valores \(t\) alrededor de la estimación lineal.
Los extremos del intervalo permitido por \(|b|\le B\) también se comprueban explícitamente, para no perder una solución óptima que caiga justo en el borde.
Una vez elegido un candidato \(b\), se calcula
$$x_b=\pi-b\alpha$$
y se toma \(a_0\), el entero más cercano a ese valor. Como el cálculo es numérico y de precisión finita, las implementaciones no se fían de una sola operación de redondeo: prueban
$$a_0-1,\qquad a_0,\qquad a_0+1.$$
De ese modo se garantiza que el mejor entero admisible no se pierda aunque \(x_b\) quede muy cerca de una media unidad.
Se conserva el par con menor error absoluto. Si dos candidatos empatan numéricamente, se elige primero el de menor \(|a|\) y luego el de menor \(|b|\), para que el resultado sea determinista.
Aquí \(\alpha=\sqrt 2\approx 1.414213562\) y \(\tau=\pi-3\approx 0.141592654\).
Tomemos el convergente
$$\frac{p}{q}=\frac{3}{2},$$
para el cual
$$\delta=2\sqrt 2-3\approx -0.171572875.$$
Como \(\tau q\approx 0.283\), la clase residual más natural es \(j=0\). Además, como \(p\equiv 1\pmod 2\), el residuo base correspondiente es \(b_0\equiv 0\pmod 2\), así que tomamos \(b_0=0\).
Entonces
$$\operatorname{diff}(t)=-\tau+t\delta.$$
La estimación lineal da
$$t\approx -\frac{-\tau}{\delta}\approx -0.825,$$
por lo que \(t=-1\) es el candidato cercano más razonable. Eso produce
$$b=b_0+tq=-2.$$
Ahora
$$\pi-b\sqrt 2=\pi+2\sqrt 2\approx 5.969019778,$$
y el entero más cercano es \(a=6\). Así, la aproximación es
$$6-2\sqrt 2\approx 3.171572875,$$
con error aproximado
$$|\,\pi-(6-2\sqrt 2)\,|\approx 0.029980221.$$
Ese mismo par coincide con la comprobación interna de consistencia usada por las implementaciones.
Las implementaciones en C++, Python y Java siguen la misma búsqueda matemática. Para cada \(d\) no cuadrado menor que \(100\), calculan \(\sqrt d\) con alta precisión, fijan la cota segura \(B\) y generan convergentes de la fracción continua hasta que el denominador supera esa cota.
Para cada convergente, la implementación forma el término pequeño \(\delta=q\sqrt d-p\), mira solo las clases residuales \(j\) cercanas a \((\pi-3)q\), convierte cada clase en un valor base \(b_0\) mediante inversión modular y deduce el intervalo permitido de \(t\) a partir de \(|b|\le B\).
En lugar de recorrer todos los posibles \(t\), solo evalúa una franja corta alrededor de la estimación lineal, además de los extremos del intervalo. Cada \(b\) resultante se convierte en candidatos enteros cercanos para \(a\), y el error se comprueba de forma explícita.
La versión de C++ reparte los distintos valores de \(d\) entre varios hilos. La versión de Python ejecuta la misma búsqueda directamente con aritmética decimal. La versión de Java delega la ejecución en el mismo núcleo numérico compilado, de modo que reproduce el comportamiento de C++.
Para un \(d\) fijo, los denominadores de los convergentes de una irracional cuadrática crecen exponencialmente, así que solo se necesitan \(O(\log B)\) convergentes antes de superar la cota de búsqueda. Como las ventanas de \(j\) y \(t\) tienen tamaño constante, cada convergente aporta solo \(O(1)\) evaluaciones reales.
Por ello, el coste práctico por cada \(d\) es \(O(\log B)\), es decir, aquí esencialmente \(O(\log N)\). El uso de memoria es \(O(\log B)\) si se almacenan explícitamente los convergentes. Como solo hay un número pequeño y fijo de valores no cuadrados con \(d\lt 100\), el coste total es ese trabajo por \(d\) multiplicado por una constante pequeña, con paralelización natural en C++.
对每个满足 \(d \lt 100\) 且不是完全平方数的整数 \(d\),都要用形如
$$a+b\sqrt d,\qquad |a|,|b|\le 10^{13}$$
的数去逼近 \(\pi\)。我们要在所有这样的整数对 \((a,b)\) 中,找出使
$$\left|\pi-\left(a+b\sqrt d\right)\right|$$
最小的那个。
题目里的 integral part 不是指实数 \(a+b\sqrt d\) 的下取整,而是直接指系数 \(a\)。因此最终要求的是
$$\sum_{d\in S} |a_d|,$$
其中 \(S\) 表示从 \(2\) 到 \(99\) 的所有非平方整数集合,\(a_d\) 表示该 \(d\) 的最佳逼近所对应的整数系数。
记 \(\alpha=\sqrt d\),再记 \(\tau=\pi-3\)。整道题的关键在于:一旦 \(b\) 固定,最优的 \(a\) 几乎立刻就能确定,因为它必须是 \(\pi-b\alpha\) 的最近整数。真正困难的不是同时枚举 \(a\) 与 \(b\),而是在极大的范围内快速找到值得检查的那些 \(b\)。
对固定的 \(b\),定义
$$x_b=\pi-b\alpha.$$
此时最小化 \(\left|\pi-(a+b\alpha)\right|\) 与最小化
$$|x_b-a|$$
完全等价,其中 \(a\) 必须满足 \(|a|\le 10^{13}\)。所以最优的 \(a\) 就是 \(x_b\) 的最近整数,原问题变成
$$\min_{|b|\le 10^{13}} \left\| \pi-b\alpha \right\|_{\mathbb Z},$$
这里 \(\|y\|_{\mathbb Z}\) 表示 \(y\) 到最近整数的距离。
同时,\(b\) 的搜索范围可以明显缩小。最简单的候选是 \(a=3\)、\(b=0\),它的误差正好是 \(\tau=\pi-3\lt 1\)。因此真正的最优解误差一定不会超过 \(1\)。如果 \(|b|\alpha \gt 10^{13}+\pi+1\),那么任意满足 \(|a|\le 10^{13}\) 的整数 \(a\) 都会有
$$\left|a+b\alpha-\pi\right|\ge |b|\alpha-|a|-\pi \gt 1,$$
这样的 \(b\) 不可能是最优候选。于是实现中使用安全上界
$$|b|\le B,\qquad B=\min\!\left(10^{13},\left\lfloor\frac{10^{13}+\pi+2}{\alpha}\right\rfloor\right).$$
这一步已经把原本巨大的搜索区域压缩到了与 \(\sqrt d\) 成反比的范围。
\(\alpha=\sqrt d\) 是二次无理数,它的连分数展开是周期性的。连分数的收敛分数 \(p/q\) 会给出非常好的有理逼近。对每个收敛分数,都写成
$$q\alpha=p+\delta,$$
其中 \(\delta=q\alpha-p\) 很小。
这一点特别有用,因为把 \(b\) 增加 \(q\) 时,\(b\alpha\) 几乎只增加了一个整数 \(p\):
$$ (b+q)\alpha=b\alpha+p+\delta.$$
也就是说,在等差数列
$$b=b_0+tq$$
里面,\(b\alpha\) 的小数部分只会以很慢的速度漂移,漂移量由 \(t\delta\) 控制。原本需要在 \(10^{13}\) 级别范围内做的搜索,就被分解成围绕这些收敛分数的若干个很小的局部搜索。
由于收敛分数总满足 \(\gcd(p,q)=1\),所以 \(p\) 在模 \(q\) 意义下有逆元。对每个
$$j\in\{0,1,\dots,q-1\},$$
选择一个 \(b_0\),使得
$$b_0p\equiv j \pmod q.$$
于是存在整数 \(k\) 满足 \(b_0p=j+kq\)。若进一步写
$$b=b_0+tq,$$
那么
$$\begin{aligned} b\alpha &=b_0\alpha+tq\alpha \\ &=\frac{b_0(p+\delta)}{q}+t(p+\delta) \\ &=k+tp+\frac{j}{q}+\frac{b_0\delta}{q}+t\delta. \end{aligned}$$
因为 \(\pi=3+\tau\),所以
$$\pi-b\alpha=\left(3-k-tp\right)+\left(\tau-\frac{j}{q}-\frac{b_0\delta}{q}-t\delta\right).$$
前一部分是整数,真正决定逼近效果的是后一部分的残差,也就是
$$\operatorname{diff}(t)=\frac{j}{q}+\frac{b_0\delta}{q}-\tau+t\delta.$$
只要这个量接近 \(0\),就说明存在合适的整数 \(a\),使得 \(a+b\sqrt d\) 非常靠近 \(\pi\)。
在固定 \(q\)、\(p\)、\(j\) 的情况下,\(\operatorname{diff}(t)\) 是关于 \(t\) 的线性函数,因此最佳 \(t\) 会出现在
$$t\approx -\frac{\frac{j}{q}+\frac{b_0\delta}{q}-\tau}{\delta}$$
附近。
还可以进一步剪枝。因为 \(|b_0\delta/q|\le |\delta|\) 很小,所以只有那些 \(j/q\) 已经接近 \(\tau\) 的剩余类才有可能给出真正优秀的候选。实现因此只检查 \(\tau q\) 附近一个很小的固定窗口中的 \(j\) 值;对于每个这样的 \(j\),再只检查线性估计附近一个很小的固定窗口中的 \(t\) 值。
此外,\(|b|\le B\) 会限制允许的 \(t\) 区间。区间的两个端点也会被显式检查,这样即使最优解恰好落在边界上,也不会漏掉。
一旦选定某个候选 \(b\),就重新计算
$$x_b=\pi-b\alpha$$
并取其最近整数 \(a_0\)。由于数值计算使用的是有限精度,实现不会只依赖一次四舍五入,而是同时测试
$$a_0-1,\qquad a_0,\qquad a_0+1.$$
这样一来,即便 \(x_b\) 非常接近半整数,真正的最优整数系数也一定会被覆盖到。
最终保留绝对误差最小的候选。如果数值误差相同,则先比较 \(|a|\),再比较 \(|b|\),从而保证结果具有确定性。
此时 \(\alpha=\sqrt 2\approx 1.414213562\),而 \(\tau=\pi-3\approx 0.141592654\)。
取收敛分数
$$\frac{p}{q}=\frac{3}{2},$$
于是
$$\delta=2\sqrt 2-3\approx -0.171572875.$$
因为 \(\tau q\approx 0.283\),最值得看的剩余类是 \(j=0\)。又因为 \(p\equiv 1\pmod 2\),对应的基准剩余就是 \(b_0\equiv 0\pmod 2\),所以取 \(b_0=0\)。
这时残差函数变成
$$\operatorname{diff}(t)=-\tau+t\delta.$$
线性估计给出
$$t\approx -\frac{-\tau}{\delta}\approx -0.825,$$
因此最近的自然整数候选是 \(t=-1\)。于是
$$b=b_0+tq=-2.$$
再计算
$$\pi-b\sqrt 2=\pi+2\sqrt 2\approx 5.969019778,$$
其最近整数为 \(a=6\)。于是得到近似值
$$6-2\sqrt 2\approx 3.171572875,$$
对应误差约为
$$|\,\pi-(6-2\sqrt 2)\,|\approx 0.029980221.$$
这正好与实现中内置的那个小规模校验结果一致。
C++、Python 和 Java 三个实现遵循同一套数学搜索流程。对于每个不是完全平方数的 \(d\lt 100\),它们都会先用高精度计算 \(\sqrt d\),再求出安全边界 \(B\),然后不断生成连分数收敛分数,直到分母超过这个边界为止。
对于每个收敛分数,实现先构造小误差项 \(\delta=q\sqrt d-p\),只查看靠近 \((\pi-3)q\) 的少量剩余类 \(j\),再利用模逆把这些剩余类转换成基准值 \(b_0\),并由 \(|b|\le B\) 推出允许的 \(t\) 区间。
实现不会遍历所有 \(t\),而只检查线性估计附近的一小段区间,以及区间的两个端点。每个由此得到的 \(b\) 都会被转换成少量附近的整数 \(a\) 候选,再直接计算误差并更新当前最优解。
C++ 版本把不同的 \(d\) 分配给多个线程并行处理。Python 版本用高精度十进制数直接执行同样的搜索。Java 版本则把执行委托给同一个已编译的数值核心,因此输出与 C++ 版本保持一致。
对固定的 \(d\) 而言,二次无理数收敛分数的分母增长是指数级的,因此在分母超过搜索上界 \(B\) 之前,只需要生成 \(O(\log B)\) 个收敛分数。由于 \(j\) 与 \(t\) 的检查窗口大小都是常数,每个收敛分数只会带来 \(O(1)\) 次真实候选评估。
因此每个 \(d\) 的实际运行代价是 \(O(\log B)\),在本题中可以视为 \(O(\log N)\)。如果把收敛分数显式存储下来,空间复杂度是 \(O(\log B)\)。而 \(d\lt 100\) 的非平方数总数本身就是一个很小的常数,所以总时间就是这个单个 \(d\) 成本乘上一个小常数;C++ 版本还能进一步并行化外层循环。
Для каждого целого \(d \lt 100\), которое не является квадратом, нужно приблизить число \(\pi\) выражениями вида
$$a+b\sqrt d,\qquad |a|,|b|\le 10^{13}.$$
Из всех таких пар \((a,b)\) выбирается та, для которой величина
$$\left|\pi-\left(a+b\sqrt d\right)\right|$$
минимальна.
Под integral part здесь понимается не пол реального числа \(a+b\sqrt d\), а именно коэффициент \(a\). Поэтому требуется вычислить
$$\sum_{d\in S} |a_d|,$$
где \(S\) есть множество всех неквадратных целых чисел от \(2\) до \(99\), а \(a_d\) является целым коэффициентом у лучшего приближения для данного \(d\).
Обозначим \(\alpha=\sqrt d\) и \(\tau=\pi-3\). Главная идея состоит в том, что после фиксации \(b\) значение \(a\) фактически определяется автоматически: нужно взять целое число, ближайшее к \(\pi-b\alpha\). Значит, вместо полного перебора по двум переменным достаточно быстро находить хорошие кандидаты для \(b\).
Для фиксированного \(b\) введем
$$x_b=\pi-b\alpha.$$
Тогда минимизация \(\left|\pi-(a+b\alpha)\right|\) по допустимым \(a\) эквивалентна минимизации
$$|x_b-a|.$$
Следовательно, лучший \(a\) есть ближайшее целое к \(x_b\), а исходная задача превращается в
$$\min_{|b|\le 10^{13}} \left\| \pi-b\alpha \right\|_{\mathbb Z},$$
где \(\|y\|_{\mathbb Z}\) обозначает расстояние от \(y\) до ближайшего целого.
Диапазон для \(b\) можно резко сузить. Тривиальный кандидат \(a=3\), \(b=0\) имеет ошибку \(\tau=\pi-3\lt 1\), значит настоящий оптимум обязан иметь ошибку меньше \(1\). Если \(|b|\alpha \gt 10^{13}+\pi+1\), то для любого \(|a|\le 10^{13}\),
$$\left|a+b\alpha-\pi\right|\ge |b|\alpha-|a|-\pi \gt 1,$$
и такой \(b\) уже не может быть лучшим. Поэтому в реализации используется безопасная граница
$$|b|\le B,\qquad B=\min\!\left(10^{13},\left\lfloor\frac{10^{13}+\pi+2}{\alpha}\right\rfloor\right).$$
Так как \(\alpha=\sqrt d\) является квадратичной иррациональностью, ее цепная дробь периодична. Конвергенты \(p/q\) дают особенно хорошие рациональные приближения. Для каждого такого конвергента пишем
$$q\alpha=p+\delta,$$
где \(\delta=q\alpha-p\) мало.
Это ключевой эффект: если увеличить \(b\) на \(q\), то величина \(b\alpha\) почти увеличится на целое \(p\):
$$ (b+q)\alpha=b\alpha+p+\delta.$$
Значит, внутри арифметической прогрессии \(b=b_0+tq\) дробная часть \(b\alpha\) меняется медленно, с шагом \(t\delta\). Именно это превращает гигантский перебор в серию очень маленьких локальных поисков.
Для конвергентов всегда выполняется \(\gcd(p,q)=1\), поэтому у \(p\) существует обратный элемент по модулю \(q\). Для каждого остатка \(j\in\{0,1,\dots,q-1\}\) выбираем \(b_0\) так, чтобы
$$b_0p\equiv j \pmod q.$$
Тогда для некоторого целого \(k\) имеем \(b_0p=j+kq\). Если далее положить
$$b=b_0+tq,$$
то
$$\begin{aligned} b\alpha &=b_0\alpha+tq\alpha \\ &=\frac{b_0(p+\delta)}{q}+t(p+\delta) \\ &=k+tp+\frac{j}{q}+\frac{b_0\delta}{q}+t\delta. \end{aligned}$$
Поскольку \(\pi=3+\tau\), получаем
$$\pi-b\alpha=\left(3-k-tp\right)+\left(\tau-\frac{j}{q}-\frac{b_0\delta}{q}-t\delta\right).$$
Первая скобка является целым числом. Поэтому качество приближения полностью определяется малым остатком
$$\operatorname{diff}(t)=\frac{j}{q}+\frac{b_0\delta}{q}-\tau+t\delta,$$
который нужно сделать как можно ближе к нулю.
При фиксированных \(q\), \(p\) и \(j\) выражение \(\operatorname{diff}(t)\) линейно по \(t\). Значит, оптимальный \(t\) должен лежать около точки
$$t\approx -\frac{\frac{j}{q}+\frac{b_0\delta}{q}-\tau}{\delta}.$$
Конкурентоспособными могут быть только те классы, у которых \(j/q\) близко к \(\tau\), потому что \(|b_0\delta/q|\le |\delta|\) для хорошего конвергента очень мало. Поэтому реализации просматривают лишь небольшое фиксированное окно значений \(j\) вокруг \(\tau q\), а для каждого такого \(j\) только небольшое фиксированное окно значений \(t\) вокруг линейной оценки.
Кроме того, явно проверяются оба конца допустимого интервала, вытекающего из ограничения \(|b|\le B\), чтобы не потерять оптимум на границе.
Когда найден кандидат \(b\), вычисляется
$$x_b=\pi-b\alpha$$
и берется ближайшее целое \(a_0\). Поскольку вычисления выполняются с конечной точностью, реализации не полагаются на одно округление, а проверяют три соседних значения
$$a_0-1,\qquad a_0,\qquad a_0+1.$$
Это гарантирует, что истинное оптимальное целое не будет пропущено даже тогда, когда \(x_b\) находится очень близко к полуцелому значению.
Далее сохраняется кандидат с наименьшей абсолютной ошибкой. Если два варианта численно совпадают, выбирается сначала меньший \(|a|\), а затем меньший \(|b|\), чтобы результат был детерминированным.
Здесь \(\alpha=\sqrt 2\approx 1.414213562\) и \(\tau=\pi-3\approx 0.141592654\).
Возьмем конвергент
$$\frac{p}{q}=\frac{3}{2},$$
для которого
$$\delta=2\sqrt 2-3\approx -0.171572875.$$
Так как \(\tau q\approx 0.283\), наиболее естественный остаток равен \(j=0\). Поскольку \(p\equiv 1\pmod 2\), соответствующий базовый остаток равен \(b_0\equiv 0\pmod 2\), то есть можно взять \(b_0=0\).
Тогда
$$\operatorname{diff}(t)=-\tau+t\delta.$$
Линейная оценка дает
$$t\approx -\frac{-\tau}{\delta}\approx -0.825,$$
поэтому ближайший разумный кандидат есть \(t=-1\). Отсюда
$$b=b_0+tq=-2.$$
Теперь
$$\pi-b\sqrt 2=\pi+2\sqrt 2\approx 5.969019778,$$
и ближайшее целое равно \(a=6\). Значит, получаем приближение
$$6-2\sqrt 2\approx 3.171572875,$$
с ошибкой примерно
$$|\,\pi-(6-2\sqrt 2)\,|\approx 0.029980221.$$
Именно эта пара совпадает со встроенной проверкой корректности в реализациях.
Реализации на C++, Python и Java используют один и тот же математический поиск. Для каждого неквадратного \(d\lt 100\) они вычисляют \(\sqrt d\) с высокой точностью, определяют безопасную границу \(B\) и порождают конвергенты цепной дроби до тех пор, пока знаменатель не превысит эту границу.
Для каждого конвергента строится малый член \(\delta=q\sqrt d-p\). Затем реализация рассматривает только остатки \(j\), лежащие рядом с \((\pi-3)q\), переводит каждый остаток в базовое значение \(b_0\) с помощью модульного обратного элемента и выводит допустимый интервал для \(t\) из ограничения \(|b|\le B\).
Вместо полного перебора всех \(t\) код проверяет лишь короткий отрезок вокруг линейной оценки, а также оба конца интервала. Для каждого полученного \(b\) явно проверяются ближайшие целые значения \(a\), после чего ошибка вычисляется напрямую.
Версия на C++ распределяет независимые значения \(d\) по нескольким потокам. Версия на Python выполняет тот же поиск напрямую в десятичной арифметике повышенной точности. Версия на Java передает вычисление тому же скомпилированному численному ядру, поэтому дает такой же результат, как и C++.
Для фиксированного \(d\) знаменатели конвергентов квадратичной иррациональности растут экспоненциально, поэтому до достижения границы \(B\) требуется только \(O(\log B)\) конвергентов. Окна по \(j\) и \(t\) имеют постоянный размер, так что каждый конвергент дает лишь \(O(1)\) реальных проверок кандидатов.
Следовательно, практическая трудоемкость для одного \(d\) равна \(O(\log B)\), то есть здесь по существу \(O(\log N)\). Если конвергенты сохраняются явно, память составляет \(O(\log B)\). Поскольку допустимых неквадратных значений \(d\lt 100\) всего конечное малое число, суммарная стоимость есть работа на одно \(d\), умноженная на небольшую константу; в C++ внешний цикл дополнительно распараллеливается.
لكل عدد صحيح \(d \lt 100\) ليس مربعًا كاملًا، نريد تقريب \(\pi\) بأعداد من الشكل
$$a+b\sqrt d,\qquad |a|,|b|\le 10^{13}.$$
ومن بين جميع الأزواج \((a,b)\) المسموح بها نختار الزوج الذي يصغّر
$$\left|\pi-\left(a+b\sqrt d\right)\right|.$$
المقصود بـ integral part هنا ليس الجزء الصحيح للعدد الحقيقي \(a+b\sqrt d\)، بل المعامل \(a\) نفسه. لذلك فالكمية المطلوبة هي
$$\sum_{d\in S} |a_d|,$$
حيث تمثل \(S\) مجموعة الأعداد الصحيحة غير المربعة من \(2\) إلى \(99\)، ويمثّل \(a_d\) المعامل الصحيح في أفضل تقريب خاص بذلك المقدار \(d\).
لنكتب \(\alpha=\sqrt d\) و \(\tau=\pi-3\). الفكرة الأساسية هي أن قيمة \(a\) تصبح شبه محسومة ما إن نثبت \(b\)، لأن أفضل اختيار لها هو أقرب عدد صحيح إلى \(\pi-b\alpha\). ولهذا فالتحدي الحقيقي ليس البحث ثنائي المتغير، بل العثور على القيم الواعدة لـ \(b\) من دون مسح هائل حتى \(10^{13}\).
عند تثبيت \(b\) نعرّف
$$x_b=\pi-b\alpha.$$
وعندها يصبح تصغير \(\left|\pi-(a+b\alpha)\right|\) على جميع قيم \(a\) المسموح بها مكافئًا تمامًا لتصغير
$$|x_b-a|.$$
إذن أفضل \(a\) هو أقرب عدد صحيح إلى \(x_b\)، وتتحول المسألة إلى
$$\min_{|b|\le 10^{13}} \left\| \pi-b\alpha \right\|_{\mathbb Z},$$
حيث \(\|y\|_{\mathbb Z}\) تمثل المسافة بين \(y\) وأقرب عدد صحيح.
كما يمكن تقليص مجال البحث في \(b\) كثيرًا. فالمرشح البسيط \(a=3\)، \(b=0\) يعطي خطأً يساوي \(\tau=\pi-3\lt 1\)، ولذلك يجب أن يكون الحل الأمثل الحقيقي ذا خطأ أصغر من \(1\). وإذا تحقق \(|b|\alpha \gt 10^{13}+\pi+1\)، فإننا نحصل لكل \(|a|\le 10^{13}\) على
$$\left|a+b\alpha-\pi\right|\ge |b|\alpha-|a|-\pi \gt 1,$$
ومن ثم لا يمكن لمثل هذه القيمة من \(b\) أن تكون مثلى. لهذا تستعمل التطبيقات الحد الآمن
$$|b|\le B,\qquad B=\min\!\left(10^{13},\left\lfloor\frac{10^{13}+\pi+2}{\alpha}\right\rfloor\right).$$
بما أن \(\alpha=\sqrt d\) عدد غير نسبي تربيعي، فإن كسره المستمر دوري. وتعطي المقاربات \(p/q\) الناتجة من هذا الكسر تقريبات كسرية ممتازة جدًا. ولكل مقاربة نكتب
$$q\alpha=p+\delta,$$
بحيث يكون \(\delta=q\alpha-p\) صغيرًا جدًا.
وهذا مفيد لأن زيادة \(b\) بمقدار \(q\) تغيّر \(b\alpha\) بمقدار يقارب العدد الصحيح \(p\):
$$ (b+q)\alpha=b\alpha+p+\delta.$$
وبالتالي، داخل المتتالية الحسابية
$$b=b_0+tq,$$
فإن الجزء الكسري من \(b\alpha\) يتحرك ببطء شديد، وسرعة هذا الانزياح تتحكم فيها الكمية \(t\delta\). وبهذا تتحول عملية البحث الضخمة إلى عدة عمليات بحث محلية صغيرة جدًا.
كل مقاربة تحقق \(\gcd(p,q)=1\)، ولذلك يوجد معكوس لـ \(p\) بترديد \(q\). لكل
$$j\in\{0,1,\dots,q-1\}$$
نختار \(b_0\) بحيث
$$b_0p\equiv j \pmod q.$$
وعندئذ يوجد عدد صحيح \(k\) يحقق \(b_0p=j+kq\). وإذا كتبنا
$$b=b_0+tq,$$
فإن
$$\begin{aligned} b\alpha &=b_0\alpha+tq\alpha \\ &=\frac{b_0(p+\delta)}{q}+t(p+\delta) \\ &=k+tp+\frac{j}{q}+\frac{b_0\delta}{q}+t\delta. \end{aligned}$$
ومع \(\pi=3+\tau\) نحصل على
$$\pi-b\alpha=\left(3-k-tp\right)+\left(\tau-\frac{j}{q}-\frac{b_0\delta}{q}-t\delta\right).$$
القوس الأول عدد صحيح، لذا فإن جودة التقريب تحددها تمامًا البقية الصغيرة
$$\operatorname{diff}(t)=\frac{j}{q}+\frac{b_0\delta}{q}-\tau+t\delta,$$
ونريد جعلها أقرب ما يمكن إلى الصفر.
عندما تكون \(q\) و\(p\) و\(j\) ثابتة، فإن \(\operatorname{diff}(t)\) دالة خطية في \(t\). لذلك يكون أفضل \(t\) قريبًا من
$$t\approx -\frac{\frac{j}{q}+\frac{b_0\delta}{q}-\tau}{\delta}.$$
ولا يمكن أن تنافس إلا الفئات التي يكون فيها \(j/q\) قريبًا من \(\tau\)، لأن \(|b_0\delta/q|\le |\delta|\) يبقى صغيرًا جدًا عند المقاربات الجيدة. ولهذا تفحص التطبيقات نافذة ثابتة صغيرة فقط من قيم \(j\) حول \(\tau q\)، ثم تفحص لكل فئة نافذة ثابتة صغيرة من قيم \(t\) حول التقدير الخطي.
كما تُفحص نهايتا المجال المسموح به الناتج من الشرط \(|b|\le B\) بشكل صريح، حتى لا تضيع حالة مثلى تقع على الحد.
بعد اختيار مرشح \(b\)، نحسب
$$x_b=\pi-b\alpha$$
ثم نأخذ \(a_0\)، وهو أقرب عدد صحيح إلى هذه القيمة. وبسبب الحسابات ذات الدقة المنتهية، لا تعتمد التطبيقات على عملية تقريب واحدة فقط، بل تختبر الأعداد الثلاثة
$$a_0-1,\qquad a_0,\qquad a_0+1.$$
وهذا يضمن عدم فقدان الحل الصحيح حتى لو كانت القيمة الحقيقية قريبة جدًا من نصف عدد صحيح.
ويُحتفَظ بالمرشح ذي أصغر خطأ مطلق. وإذا تساوى خطآن عدديًا، فيُفضَّل أولًا الأصغر \(|a|\)، ثم الأصغر \(|b|\)، حتى تكون النتيجة حتمية.
في هذه الحالة \(\alpha=\sqrt 2\approx 1.414213562\)، و\(\tau=\pi-3\approx 0.141592654\).
لنأخذ المقاربة
$$\frac{p}{q}=\frac{3}{2},$$
وعندها
$$\delta=2\sqrt 2-3\approx -0.171572875.$$
ولأن \(\tau q\approx 0.283\)، فإن الفئة الأنسب هي \(j=0\). وبما أن \(p\equiv 1\pmod 2\)، فإن الباقي الأساسي الموافق هو \(b_0\equiv 0\pmod 2\)، ولذلك نأخذ \(b_0=0\).
عندئذ تصبح
$$\operatorname{diff}(t)=-\tau+t\delta.$$
ويعطي التقدير الخطي
$$t\approx -\frac{-\tau}{\delta}\approx -0.825,$$
ومن ثم يكون \(t=-1\) مرشحًا طبيعيًا قريبًا. وهذا يؤدي إلى
$$b=b_0+tq=-2.$$
ثم نحسب
$$\pi-b\sqrt 2=\pi+2\sqrt 2\approx 5.969019778,$$
وأقرب عدد صحيح هو \(a=6\). وبذلك يكون التقريب
$$6-2\sqrt 2\approx 3.171572875,$$
وبخطأ يقارب
$$|\,\pi-(6-2\sqrt 2)\,|\approx 0.029980221.$$
وهذا يطابق فحص الاتساق الصغير المضمَّن في التطبيقات.
تتبع تطبيقات C++ وPython وJava البحث الرياضي نفسه. فلكل قيمة غير مربعة من \(d\lt 100\)، تحسب \(\sqrt d\) بدقة عالية، وتستخرج الحد الآمن \(B\)، ثم تولّد مقاربات الكسر المستمر حتى يتجاوز المقام ذلك الحد.
وبالنسبة لكل مقاربة، تبني الشيفرة الحد الصغير \(\delta=q\sqrt d-p\)، وتفحص فقط فئات البواقي \(j\) القريبة من \((\pi-3)q\)، ثم تحوّل كل فئة إلى قيمة أساس \(b_0\) باستعمال المعكوس بترديد \(q\)، وتستنتج المجال المسموح به لـ \(t\) من القيد \(|b|\le B\).
ولا تقوم الشيفرة بمسح جميع قيم \(t\)، بل تكتفي بنافذة قصيرة حول التقدير الخطي، مع فحص طرفي المجال أيضًا. وكل قيمة \(b\) ناتجة تُحوَّل إلى عدد قليل من المرشحين المجاورين لـ \(a\)، ثم يُحسب الخطأ مباشرة لاختيار الأفضل.
نسخة C++ توزع قيم \(d\) المستقلة على عدة خيوط تنفيذ. ونسخة Python تنفذ البحث نفسه مباشرة باستعمال حساب عشري عالي الدقة. أما نسخة Java فتفوض التنفيذ إلى النواة العددية المترجمة نفسها، ولذلك يطابق خرجها سلوك C++.
لـ \(d\) ثابت، تنمو مقامات المقاربات الخاصة بعدد غير نسبي تربيعي نموًا أسيًا، ولذلك لا نحتاج إلا إلى \(O(\log B)\) من المقاربات قبل أن يتجاوز المقام حد البحث. وبما أن نافذتي \(j\) و\(t\) ثابتتا الحجم، فإن كل مقاربة تولد فقط \(O(1)\) من التقييمات الحقيقية للمرشحين.
لذلك يكون الزمن العملي لكل قيمة \(d\) هو \(O(\log B)\)، أي فعليًا هنا \(O(\log N)\). وإذا خُزنت المقاربات صراحة فالاستهلاك المكاني يساوي \(O(\log B)\). وبما أن عدد القيم غير المربعة مع \(d\lt 100\) صغير وثابت، فإن الكلفة الكلية هي كلفة كل \(d\) مضروبة في ثابت صغير، مع إمكانية التوازي المباشر في C++.