We want the number of primitive Pythagorean triples \((a,b,c)\) with hypotenuse bound
$$c \le N.$$
Euclid's parametrization says that every primitive triple appears uniquely in the form
$$a=m^2-n^2,\qquad b=2mn,\qquad c=m^2+n^2,$$
with
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
So the problem is equivalent to counting parameter pairs \((m,n)\) satisfying those three conditions together with \(m^2+n^2\le N\).
Let \(P(X)\) denote the number of primitive parameter pairs \((m,n)\) with \(m^2+n^2\le X\), \(m>n>0\), coprime, and of opposite parity. The desired answer is \(P(N)\). The implementations compute it in two layers: first a fast raw counter without the gcd condition, then a recurrence that removes non-primitive pairs.
The classical parametrization gives a bijection between primitive triples and pairs \((m,n)\) such that
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
The hypotenuse becomes
$$c=m^2+n^2,$$
so counting primitive triples with \(c\le N\) is exactly the same as counting those admissible pairs inside the quarter-circle
$$m^2+n^2\le N.$$
This change of viewpoint is the foundation of the whole solution: geometry gives the bound, and arithmetic gives the primitivity conditions.
Define \(R(X)\) as the number of pairs \((m,n)\) with
$$m>n>0,\qquad m^2+n^2\le X,\qquad m\not\equiv n \pmod{2},$$
but without the condition \(\gcd(m,n)=1\).
For a fixed \(m\), the largest allowed \(n\) is
$$n_{\max}(m)=\min\left(m-1,\left\lfloor\sqrt{X-m^2}\right\rfloor\right).$$
The row contribution is therefore
$$r_X(m)=\begin{cases} \left\lfloor \dfrac{n_{\max}(m)}{2}\right\rfloor, & m\text{ odd},\\[6pt] \left\lfloor \dfrac{n_{\max}(m)+1}{2}\right\rfloor, & m\text{ even}. \end{cases}$$
Hence
$$R(X)=\sum_{m\ge 1} r_X(m).$$
The code evaluates this quickly by separating a full prefix of rows. If
$$m^2+(m-1)^2\le X,$$
then the entire interval \(1\le n\le m-1\) is inside the circle. Solving this inequality gives
$$m\le \frac{1+\sqrt{2X-1}}{2}.$$
Let
$$m_0=\left\lfloor\frac{1+\sqrt{2X-1}}{2}\right\rfloor,\qquad K=\left\lfloor\frac{m_0}{2}\right\rfloor.$$
The first \(2K\) rows contribute exactly
$$1+3+5+\cdots+(2K-1)=K^2,$$
because row \(2j-1\) contributes \(j-1\) and row \(2j\) contributes \(j\). After that prefix, the remaining upper bound \(\lfloor\sqrt{X-m^2}\rfloor\) decreases monotonically, so one moving square-root pointer is enough for the tail.
Now restore the coprimality condition. Opposite parity already rules out a factor \(2\), so any common divisor of \(m\) and \(n\) must be odd.
If \((m,n)\) has odd gcd \(d\), then
$$m=d\,m',\qquad n=d\,n',$$
and the reduced pair \((m',n')\) is primitive. The circle bound becomes
$$m'^2+n'^2\le \left\lfloor \frac{X}{d^2}\right\rfloor.$$
Therefore every raw pair is obtained uniquely from a primitive pair and an odd scaling factor, which gives
$$R(X)=\sum_{\substack{d\ge 1\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
Rearranging,
$$P(X)=R(X)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
This is the recurrence used directly by the implementations. It is equivalent to odd-only Möbius inversion, but the code never needs an explicit table of Möbius values.
Let
$$L=\left\lfloor N^{1/3}\right\rfloor.$$
First compute \(P(x)\) for every \(1\le x\le L\). For a fixed \(x\), odd divisors \(d\le x^{1/3}\) are handled one by one in
$$P(x)=R(x)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{x}{d^2}\right\rfloor\right).$$
For larger \(d\), the quotient \(\left\lfloor x/d^2\right\rfloor\) is small, and many different \(d\) produce the same quotient. Those repeated values are grouped together. For a fixed small \(z\), the number of odd \(d\) with
$$\left\lfloor \frac{x}{d^2}\right\rfloor=z$$
is
$$\left\lfloor\frac{\sqrt{x/z}-1}{2}\right\rfloor-\left\lfloor\frac{\sqrt{x/(z+1)}-1}{2}\right\rfloor.$$
That converts many repeated recursive terms into one multiplication by a multiplicity.
After the small table is known, the implementations evaluate the transformed values
$$Q(t)=P\left(\left\lfloor \frac{N}{t^2}\right\rfloor\right),\qquad t\text{ odd},$$
processing odd \(t\) from large to small. When the recurrence for \(Q(t)\) reaches a reduced argument at most \(L\), it is read from the small table. Otherwise it is another transformed value \(Q(tu)\) with a larger odd factor, which has already been computed. Finally, the required answer is simply
$$Q(1)=P(N).$$
The opposite-parity pairs with \(m^2+n^2\le 50\) are
$$ (2,1),\ (3,2),\ (4,1),\ (4,3),\ (5,2),\ (5,4),\ (6,1),\ (6,3). $$
So \(R(50)=8\).
The only possible odd common divisor larger than \(1\) is \(3\), because \(5^2>50\). Thus
$$P(50)=R(50)-P\left(\left\lfloor\frac{50}{9}\right\rfloor\right)=8-P(5).$$
Now \(P(5)=1\), coming from the single pair \((2,1)\). Therefore
$$P(50)=8-1=7.$$
The seven primitive triples are
$$ (3,4,5),\ (5,12,13),\ (8,15,17),\ (7,24,25),\ (20,21,29),\ (12,35,37),\ (9,40,41). $$
This small example shows exactly how the raw geometric count and the odd-gcd recurrence fit together.
The C++, Python, and Java implementations begin with exact integer square-root and cube-root routines so that every cutoff such as \(\lfloor \sqrt{x}\rfloor\) and \(\lfloor N^{1/3}\rfloor\) is computed safely despite floating-point starting estimates.
They then evaluate the raw counter \(R(X)\) by combining a closed-form prefix \(K^2\) with a monotone tail scan. The tail keeps one decreasing square-root boundary and alternates the parity formula for odd and even rows, which avoids recomputing the whole range of \(n\) for each \(m\).
Next the implementation fills a first table for all primitive counts \(P(x)\) with \(x\le L\). Each entry uses the recurrence above, split into individually handled large quotients and grouped small quotients.
After that it fills a second table for transformed arguments \(\left\lfloor N/t^2\right\rfloor\) with odd \(t\), processed in descending order. Small reduced arguments are read from the first table; large reduced arguments correspond to already computed transformed values with larger odd factors.
The entry for \(t=1\) is \(P(N)\), so it is exactly the number of primitive Pythagorean triples with hypotenuse at most \(N\).
A single evaluation of the raw counter \(R(X)\) costs \(O(\sqrt{X})\) time because the tail scan visits only the remaining boundary rows. The recurrence work around that raw counter contributes about \(O(X^{1/3})\) grouped-divisor operations.
Building the small table up to \(L=\lfloor N^{1/3}\rfloor\) costs \(O(N^{1/2})\) time. In the transformed stage, the dominant raw-count cost is
$$\sum_{\substack{t\le L\\ t\equiv 1 \pmod{2}}} O\left(\sqrt{\frac{N}{t^2}}\right)=O(\sqrt{N}\log N).$$
The extra grouped-divisor work is smaller than that dominant term. The memory usage is \(O(N^{1/3})\), because both stored tables have size proportional to the cube-root threshold.
Gesucht ist die Anzahl primitiver pythagoreischer Tripel \((a,b,c)\) mit
$$c \le N.$$
Nach der Parametrisierung von Euklid hat jedes primitive Tripel eindeutig die Form
$$a=m^2-n^2,\qquad b=2mn,\qquad c=m^2+n^2,$$
wobei
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
Damit reduziert sich die Aufgabe auf das Zählen solcher Paare \((m,n)\) unter der Kreisbedingung \(m^2+n^2\le N\).
Sei \(P(X)\) die Anzahl primitiver Parameterpaare \((m,n)\) mit \(m^2+n^2\le X\), \(m>n>0\), teilerfremd und unterschiedlicher Parität. Gesucht ist also \(P(N)\). Die Implementierungen berechnen diesen Wert in zwei Stufen: zuerst einen schnellen Rohzähler ohne ggT-Bedingung, danach eine Rekursion, die die nicht-primitiven Paare entfernt.
Die klassische Parametrisierung liefert eine Bijektion zwischen primitiven Tripeln und Paaren \((m,n)\) mit
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
Die Hypotenuse ist dabei
$$c=m^2+n^2.$$
Primitive Tripel mit \(c\le N\) zu zählen ist also genau dasselbe wie diese zulässigen Paare im Viertelkreis
$$m^2+n^2\le N$$
zu zählen. Diese Umformulierung trennt die geometrische Schranke von den arithmetischen Primtivitätsbedingungen.
Definiere \(R(X)\) als Anzahl der Paare \((m,n)\) mit
$$m>n>0,\qquad m^2+n^2\le X,\qquad m\not\equiv n \pmod{2},$$
aber noch ohne die Bedingung \(\gcd(m,n)=1\).
Für festes \(m\) ist das größte zulässige \(n\)
$$n_{\max}(m)=\min\left(m-1,\left\lfloor\sqrt{X-m^2}\right\rfloor\right).$$
Der Beitrag dieser Zeile ist dann
$$r_X(m)=\begin{cases} \left\lfloor \dfrac{n_{\max}(m)}{2}\right\rfloor, & m\text{ ungerade},\\[6pt] \left\lfloor \dfrac{n_{\max}(m)+1}{2}\right\rfloor, & m\text{ gerade}. \end{cases}$$
Also gilt
$$R(X)=\sum_{m\ge 1} r_X(m).$$
Der Code beschleunigt diese Summe, indem er einen vollständigen Anfangsblock separat behandelt. Falls
$$m^2+(m-1)^2\le X,$$
dann liegt das ganze Intervall \(1\le n\le m-1\) noch innerhalb des Kreises. Aufgelöst nach \(m\) ergibt das
$$m\le \frac{1+\sqrt{2X-1}}{2}.$$
Mit
$$m_0=\left\lfloor\frac{1+\sqrt{2X-1}}{2}\right\rfloor,\qquad K=\left\lfloor\frac{m_0}{2}\right\rfloor$$
liefern die ersten \(2K\) Zeilen genau den Beitrag
$$1+3+5+\cdots+(2K-1)=K^2,$$
denn Zeile \(2j-1\) trägt \(j-1\) bei und Zeile \(2j\) trägt \(j\) bei. Dahinter fällt \(\lfloor\sqrt{X-m^2}\rfloor\) monoton, daher genügt ein einzelner beweglicher Wurzelzeiger.
Nun kommt die ggT-Bedingung zurück. Verschiedene Parität schließt bereits einen gemeinsamen Faktor \(2\) aus, also ist jeder gemeinsame Teiler von \(m\) und \(n\) ungerade.
Hat ein Paar den ungeraden ggT \(d\), dann gilt
$$m=d\,m',\qquad n=d\,n',$$
und \((m',n')\) ist primitiv. Die Kreisbedingung wird zu
$$m'^2+n'^2\le \left\lfloor \frac{X}{d^2}\right\rfloor.$$
Somit erhält man jede Rohpaarung eindeutig aus einem primitiven Paar und einem ungeraden Streckfaktor:
$$R(X)=\sum_{\substack{d\ge 1\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
Daraus folgt die Rekursion
$$P(X)=R(X)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
Diese Form setzt das Programm direkt um. Sie ist äquivalent zu einer Möbius-Inversion nur über ungerade Teiler, braucht aber keine explizite Möbius-Tabelle.
Setze
$$L=\left\lfloor N^{1/3}\right\rfloor.$$
Zuerst werden alle Werte \(P(x)\) für \(1\le x\le L\) berechnet. Für festes \(x\) behandelt man ungerade Teiler \(d\le x^{1/3}\) einzeln in
$$P(x)=R(x)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{x}{d^2}\right\rfloor\right).$$
Für größere \(d\) ist der Quotient \(\left\lfloor x/d^2\right\rfloor\) klein, und viele verschiedene \(d\) liefern denselben Quotienten. Diese Wiederholungen werden gruppiert. Für ein festes kleines \(z\) ist die Anzahl ungerader \(d\) mit
$$\left\lfloor \frac{x}{d^2}\right\rfloor=z$$
gleich
$$\left\lfloor\frac{\sqrt{x/z}-1}{2}\right\rfloor-\left\lfloor\frac{\sqrt{x/(z+1)}-1}{2}\right\rfloor.$$
Dadurch werden viele identische Rekursionsterme zu einer einzigen Multiplikation mit ihrer Häufigkeit zusammengefasst.
Ist die kleine Tabelle bekannt, berechnet die Implementierung zusätzlich die transformierten Werte
$$Q(t)=P\left(\left\lfloor \frac{N}{t^2}\right\rfloor\right),\qquad t\text{ ungerade},$$
wobei ungerade \(t\) absteigend verarbeitet werden. Fällt ein reduziertes Argument auf höchstens \(L\), wird es aus der kleinen Tabelle gelesen. Andernfalls ist es ein anderer transformierter Wert \(Q(tu)\) mit größerem ungeradem Faktor und daher bereits vorhanden. Am Ende ist
$$Q(1)=P(N)$$
genau die gesuchte Antwort.
Die Paare verschiedener Parität mit \(m^2+n^2\le 50\) sind
$$ (2,1),\ (3,2),\ (4,1),\ (4,3),\ (5,2),\ (5,4),\ (6,1),\ (6,3). $$
Also ist \(R(50)=8\).
Der einzige mögliche ungerade gemeinsame Teiler größer als \(1\) ist \(3\), denn \(5^2>50\). Daher
$$P(50)=R(50)-P\left(\left\lfloor\frac{50}{9}\right\rfloor\right)=8-P(5).$$
Es gilt \(P(5)=1\), nämlich für das einzige Paar \((2,1)\). Folglich
$$P(50)=8-1=7.$$
Die sieben primitiven Tripel lauten
$$ (3,4,5),\ (5,12,13),\ (8,15,17),\ (7,24,25),\ (20,21,29),\ (12,35,37),\ (9,40,41). $$
Das Beispiel zeigt sehr konkret, wie geometrische Rohzählung und ungerade-ggT-Rekursion zusammenwirken.
Die C++-, Python- und Java-Implementierungen beginnen mit exakten ganzzahligen Quadratwurzel- und Kubikwurzelroutinen. So werden alle Grenzwerte wie \(\lfloor \sqrt{x}\rfloor\) und \(\lfloor N^{1/3}\rfloor\) trotz Gleitkomma-Startwerten korrekt getroffen.
Danach berechnen sie den Rohzähler \(R(X)\) mit einer Kombination aus geschlossenem Präfix \(K^2\) und monotonem Schwanzscan. Im Schwanz wird nur ein fallender Wurzelgrenzwert gepflegt, und die Paritätsformel wechselt ab zwischen ungeraden und geraden Zeilen.
Anschließend füllt die Implementierung eine erste Tabelle aller primitiven Werte \(P(x)\) für \(x\le L\). Jeder Eintrag verwendet die obige Rekursion, aufgeteilt in einzeln behandelte große Quotienten und gebündelte kleine Quotienten.
Danach folgt eine zweite Tabelle für transformierte Argumente \(\left\lfloor N/t^2\right\rfloor\) mit ungeradem \(t\), diesmal in absteigender Reihenfolge von \(t\). Kleine reduzierte Argumente kommen aus der ersten Tabelle; große reduzierte Argumente entsprechen bereits berechneten transformierten Werten mit größerem ungeradem Faktor.
Der Eintrag für \(t=1\) ist \(P(N)\) und damit genau die Anzahl primitiver pythagoreischer Tripel mit Hypotenuse höchstens \(N\).
Eine einzelne Auswertung des Rohzählers \(R(X)\) benötigt \(O(\sqrt{X})\) Zeit, weil der Schwanzscan nur die verbleibenden Randzeilen durchläuft. Die Rekursionsarbeit um diesen Rohzähler herum fügt ungefähr \(O(X^{1/3})\) gruppierte Divisoroperationen hinzu.
Der Aufbau der kleinen Tabelle bis \(L=\lfloor N^{1/3}\rfloor\) kostet \(O(N^{1/2})\) Zeit. In der transformierten Phase ist der dominante Rohzählungsanteil
$$\sum_{\substack{t\le L\\ t\equiv 1 \pmod{2}}} O\left(\sqrt{\frac{N}{t^2}}\right)=O(\sqrt{N}\log N).$$
Die zusätzliche Quotientenbündelung ist asymptotisch kleiner als dieser führende Term. Der Speicherverbrauch ist \(O(N^{1/3})\), weil beide Tabellen nur proportional zur Kubikwurzelschranke groß sind.
Aranan şey, hipotenüsü
$$c \le N$$
olan primitive Pisagor üçlülerinin sayısıdır.
Euclid parametrelemesine göre her primitive üçlü tek bir şekilde
$$a=m^2-n^2,\qquad b=2mn,\qquad c=m^2+n^2,$$
biçiminde yazılır ve burada
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
Dolayısıyla problem, bu koşulları sağlayan \((m,n)\) çiftlerini \(m^2+n^2\le N\) çeyrek dairesi içinde saymaya dönüşür.
\(P(X)\), \(m^2+n^2\le X\), \(m>n>0\), aralarında asal ve zıt parity'li primitive parametre çiftlerinin sayısı olsun. İstenen sonuç \(P(N)\)'dir. Uygulamalar bunu iki aşamada hesaplar: önce aralarında asallık şartı olmadan hızlı bir ham sayaç kurulur, sonra primitive olmayan çiftler bir yineleme ile çıkarılır.
Klasik parametreleme, primitive üçlülerle şu koşulları sağlayan \((m,n)\) çiftleri arasında bire bir eşleme kurar:
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
Bu durumda hipotenüs
$$c=m^2+n^2$$
olur. O halde \(c\le N\) şartındaki primitive üçlüleri saymak, tam olarak
$$m^2+n^2\le N$$
bölgesindeki uygun parametre çiftlerini saymaktır. Çözümün temelinde bu geometrik-aritmetik dönüşüm bulunur.
\(R(X)\) fonksiyonunu,
$$m>n>0,\qquad m^2+n^2\le X,\qquad m\not\equiv n \pmod{2}$$
koşullarını sağlayan ama \(\gcd(m,n)=1\) şartını henüz istemeyen çiftlerin sayısı olarak tanımlayalım.
Sabit bir \(m\) için en büyük izinli \(n\)
$$n_{\max}(m)=\min\left(m-1,\left\lfloor\sqrt{X-m^2}\right\rfloor\right)$$
olur. Buna göre satır katkısı
$$r_X(m)=\begin{cases} \left\lfloor \dfrac{n_{\max}(m)}{2}\right\rfloor, & m\text{ tek},\\[6pt] \left\lfloor \dfrac{n_{\max}(m)+1}{2}\right\rfloor, & m\text{ çift}. \end{cases}$$
ve dolayısıyla
$$R(X)=\sum_{m\ge 1} r_X(m)$$
elde edilir.
Kod bu toplamı hızlandırmak için başta tam dolu bir satır bloğunu ayırır. Eğer
$$m^2+(m-1)^2\le X,$$
ise \(1\le n\le m-1\) aralığının tamamı dairenin içindedir. Bu eşitsizlikten
$$m\le \frac{1+\sqrt{2X-1}}{2}$$
çıkar. Şimdi
$$m_0=\left\lfloor\frac{1+\sqrt{2X-1}}{2}\right\rfloor,\qquad K=\left\lfloor\frac{m_0}{2}\right\rfloor$$
olsun. İlk \(2K\) satırın toplam katkısı tam olarak
$$1+3+5+\cdots+(2K-1)=K^2$$
olur; çünkü \(2j-1\). satır \(j-1\), \(2j\). satır ise \(j\) katkı verir. Bu tam bloktan sonra \(\lfloor\sqrt{X-m^2}\rfloor\) monoton biçimde azalır; bu yüzden kuyruk kısmı tek bir hareketli karekök sınırıyla taranabilir.
Şimdi aralarında asallık şartını geri getirelim. Zıt parity zaten ortak \(2\) çarpanını imkansız kılar; dolayısıyla ortak bölen varsa mutlaka tektir.
Eğer bir çiftin tek ortak böleni \(d\) ise
$$m=d\,m',\qquad n=d\,n'$$
yazılır ve \((m',n')\) primitive olur. Daire koşulu da
$$m'^2+n'^2\le \left\lfloor \frac{X}{d^2}\right\rfloor$$
şekline iner. Böylece her ham çift, primitive bir çift ile tek bir tek ölçek çarpanından gelir:
$$R(X)=\sum_{\substack{d\ge 1\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
Buradan yineleme bağıntısı
$$P(X)=R(X)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right)$$
elde edilir. Uygulamanın doğrudan hesapladığı şey budur. Bu formül, yalnızca tek bölenler üzerinde çalışan Möbius terslemesine denktir; fakat kodun açıkça Möbius tablosu oluşturmasına gerek yoktur.
Şunu tanımlayalım:
$$L=\left\lfloor N^{1/3}\right\rfloor.$$
Önce \(1\le x\le L\) için tüm \(P(x)\) değerleri hesaplanır. Sabit bir \(x\) için \(d\le x^{1/3}\) olan tek bölenler aşağıdaki bağıntıda tek tek ele alınır:
$$P(x)=R(x)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{x}{d^2}\right\rfloor\right).$$
Daha büyük \(d\) değerlerinde \(\left\lfloor x/d^2\right\rfloor\) küçülür ve birçok farklı \(d\) aynı bölümü üretir. Bu tekrarlar gruplanır. Sabit küçük bir \(z\) için
$$\left\lfloor \frac{x}{d^2}\right\rfloor=z$$
eşitliğini sağlayan tek \(d\) sayısı
$$\left\lfloor\frac{\sqrt{x/z}-1}{2}\right\rfloor-\left\lfloor\frac{\sqrt{x/(z+1)}-1}{2}\right\rfloor$$
kadardır. Böylece aynı yineleme terimi tekrar tekrar çıkarılmak yerine, bir kez kendi çokluğu ile çarpılır.
Küçük tablo tamamlandıktan sonra uygulama ayrıca
$$Q(t)=P\left(\left\lfloor \frac{N}{t^2}\right\rfloor\right),\qquad t\text{ tek}$$
dönüştürülmüş değerlerini de hesaplar ve tek \(t\) değerlerini büyükten küçüğe işler. Azaltılmış argüman \(L\)'nin altına düşerse küçük tablodan okunur; düşmezse daha büyük tek çarpanlı başka bir dönüştürülmüş değere karşılık gelir ve o değer zaten hesaplanmıştır. Sonuç olarak
$$Q(1)=P(N)$$
istenen cevaptır.
\(m^2+n^2\le 50\) altında zıt parity'li çiftler şunlardır:
$$ (2,1),\ (3,2),\ (4,1),\ (4,3),\ (5,2),\ (5,4),\ (6,1),\ (6,3). $$
Dolayısıyla \(R(50)=8\).
\(1\)'den büyük olabilecek tek ortak tek bölen yalnızca \(3\)'tür; çünkü \(5^2>50\). O halde
$$P(50)=R(50)-P\left(\left\lfloor\frac{50}{9}\right\rfloor\right)=8-P(5).$$
\(P(5)=1\) olur; bu da tek \((2,1)\) çiftinden gelir. Sonuç:
$$P(50)=8-1=7.$$
Buna karşılık gelen yedi primitive üçlü
$$ (3,4,5),\ (5,12,13),\ (8,15,17),\ (7,24,25),\ (20,21,29),\ (12,35,37),\ (9,40,41) $$
olur. Bu örnek ham geometrik sayım ile tek-ortak-bölen yinelemesinin nasıl birleştiğini açıkça gösterir.
C++, Python ve Java uygulamaları önce tam sayı karekökü ve küpkökünü güvenli biçimde hesaplar; böylece \(\lfloor\sqrt{x}\rfloor\) ve \(\lfloor N^{1/3}\rfloor\) gibi sınırlar kayan nokta başlangıç değerlerine rağmen doğru kalır.
Daha sonra \(R(X)\) ham sayacı, kapalı biçimde hesaplanan \(K^2\) ön eki ile monotondan yararlanan kuyruk taramasını birleştirerek elde edilir. Kuyrukta yalnızca azalan bir karekök sınırı tutulur ve tek/çift satırlar için uygun parity formülü dönüşümlü uygulanır.
Ardından uygulama \(x\le L\) için tüm primitive sayıları içeren ilk tabloyu doldurur. Her giriş, büyük bölümlerin tek tek işlendiği ve küçük bölümlerin çokluklarıyla gruplanıp çıkarıldığı yinelemeyi kullanır.
Sonra \(\left\lfloor N/t^2\right\rfloor\) biçimindeki dönüştürülmüş argümanlar için ikinci tablo oluşturulur; burada tek \(t\) değerleri büyükten küçüğe işlenir. Küçük azaltılmış argümanlar ilk tablodan, büyük azaltılmış argümanlar ise daha önce hesaplanmış başka dönüştürülmüş değerlerden alınır.
\(t=1\) girişindeki değer \(P(N)\)'dir; dolayısıyla bu değer doğrudan istenen primitive Pisagor üçlüleri sayısını verir.
\(R(X)\) ham sayacının tek bir değerlendirmesi, kuyruk taraması yalnızca sınır satırlarını geçtiği için \(O(\sqrt{X})\) zaman alır. Bu ham sayacın çevresindeki yineleme işi yaklaşık \(O(X^{1/3})\) gruplanmış bölen işlemi ekler.
\(L=\lfloor N^{1/3}\rfloor\) değerine kadar olan küçük tablo \(O(N^{1/2})\) zamanda kurulur. Dönüştürülmüş aşamada baskın ham-sayım maliyeti
$$\sum_{\substack{t\le L\\ t\equiv 1 \pmod{2}}} O\left(\sqrt{\frac{N}{t^2}}\right)=O(\sqrt{N}\log N)$$
şeklindedir. Ek bölüm gruplama maliyeti bu baskın terimden daha küçüktür. Bellek kullanımı \(O(N^{1/3})\)'tür; çünkü iki tablo da küpkök eşiği büyüklüğündedir.
Queremos contar las ternas pitagóricas primitivas \((a,b,c)\) cuya hipotenusa satisface
$$c \le N.$$
La parametrización de Euclides afirma que cada terna primitiva aparece de manera única como
$$a=m^2-n^2,\qquad b=2mn,\qquad c=m^2+n^2,$$
con
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
Por tanto, el problema equivale a contar los pares \((m,n)\) que cumplen esas condiciones dentro de la región \(m^2+n^2\le N\).
Sea \(P(X)\) el número de pares primitivos \((m,n)\) con \(m^2+n^2\le X\), \(m>n>0\), coprimos y de paridad opuesta. La respuesta buscada es \(P(N)\). Las implementaciones lo obtienen en dos capas: primero cuentan rápidamente todos los pares de paridad opuesta sin imponer coprimalidad; después eliminan los no primitivos con una recurrencia.
La parametrización clásica produce una biyección entre las ternas primitivas y los pares \((m,n)\) que satisfacen
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
La hipotenusa vale
$$c=m^2+n^2,$$
así que contar ternas primitivas con \(c\le N\) es exactamente contar esos pares admisibles dentro del cuarto de círculo
$$m^2+n^2\le N.$$
Esta reformulación es la clave del método: la geometría aporta la cota y la aritmética aporta la condición de primitividad.
Definamos \(R(X)\) como el número de pares \((m,n)\) con
$$m>n>0,\qquad m^2+n^2\le X,\qquad m\not\equiv n \pmod{2},$$
pero todavía sin imponer \(\gcd(m,n)=1\).
Para un \(m\) fijo, el mayor \(n\) permitido es
$$n_{\max}(m)=\min\left(m-1,\left\lfloor\sqrt{X-m^2}\right\rfloor\right).$$
La contribución de esa fila es entonces
$$r_X(m)=\begin{cases} \left\lfloor \dfrac{n_{\max}(m)}{2}\right\rfloor, & m\text{ impar},\\[6pt] \left\lfloor \dfrac{n_{\max}(m)+1}{2}\right\rfloor, & m\text{ par}. \end{cases}$$
Por tanto,
$$R(X)=\sum_{m\ge 1} r_X(m).$$
El código acelera esta suma separando un prefijo completo de filas. Si
$$m^2+(m-1)^2\le X,$$
entonces todo el intervalo \(1\le n\le m-1\) cabe dentro del círculo. Al resolver esta desigualdad obtenemos
$$m\le \frac{1+\sqrt{2X-1}}{2}.$$
Sea
$$m_0=\left\lfloor\frac{1+\sqrt{2X-1}}{2}\right\rfloor,\qquad K=\left\lfloor\frac{m_0}{2}\right\rfloor.$$
Las primeras \(2K\) filas aportan exactamente
$$1+3+5+\cdots+(2K-1)=K^2,$$
porque la fila \(2j-1\) aporta \(j-1\) y la fila \(2j\) aporta \(j\). Después de ese prefijo, \(\lfloor\sqrt{X-m^2}\rfloor\) decrece de forma monótona, de modo que basta un único puntero móvil.
Ahora imponemos la coprimalidad. La paridad opuesta ya descarta un factor común \(2\), así que cualquier divisor común debe ser impar.
Si un par tiene divisor común impar \(d\), entonces
$$m=d\,m',\qquad n=d\,n',$$
y \((m',n')\) es primitivo. La restricción geométrica pasa a ser
$$m'^2+n'^2\le \left\lfloor \frac{X}{d^2}\right\rfloor.$$
En consecuencia, cada par bruto proviene de un par primitivo y de un factor de escala impar:
$$R(X)=\sum_{\substack{d\ge 1\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
Reordenando,
$$P(X)=R(X)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
Ésa es la recurrencia que usan directamente las implementaciones. Es equivalente a una inversión de Möbius restringida a divisores impares, pero no hace falta construir explícitamente los valores de \(\mu\).
Definimos
$$L=\left\lfloor N^{1/3}\right\rfloor.$$
Primero se calculan todos los valores \(P(x)\) para \(1\le x\le L\). Para un \(x\) fijo, los divisores impares \(d\le x^{1/3}\) se tratan individualmente en
$$P(x)=R(x)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{x}{d^2}\right\rfloor\right).$$
Para \(d\) mayores, el cociente \(\left\lfloor x/d^2\right\rfloor\) es pequeño y muchos \(d\) distintos producen el mismo valor. Esos términos se agrupan. Para un pequeño \(z\), la cantidad de \(d\) impares que cumplen
$$\left\lfloor \frac{x}{d^2}\right\rfloor=z$$
es
$$\left\lfloor\frac{\sqrt{x/z}-1}{2}\right\rfloor-\left\lfloor\frac{\sqrt{x/(z+1)}-1}{2}\right\rfloor.$$
Así, muchos términos repetidos se reemplazan por una sola resta multiplicada por su multiplicidad.
Una vez conocida la tabla pequeña, la implementación calcula además
$$Q(t)=P\left(\left\lfloor \frac{N}{t^2}\right\rfloor\right),\qquad t\text{ impar},$$
procesando los \(t\) impares de mayor a menor. Si un argumento reducido cae por debajo de \(L\), se toma de la tabla pequeña. Si no, coincide con otro valor transformado \(Q(tu)\) con un factor impar mayor, ya calculado. Al final,
$$Q(1)=P(N)$$
es exactamente la respuesta.
Los pares de paridad opuesta con \(m^2+n^2\le 50\) son
$$ (2,1),\ (3,2),\ (4,1),\ (4,3),\ (5,2),\ (5,4),\ (6,1),\ (6,3). $$
Por lo tanto, \(R(50)=8\).
El único divisor común impar posible mayor que \(1\) es \(3\), ya que \(5^2>50\). Luego
$$P(50)=R(50)-P\left(\left\lfloor\frac{50}{9}\right\rfloor\right)=8-P(5).$$
Como \(P(5)=1\), procedente del único par \((2,1)\), obtenemos
$$P(50)=8-1=7.$$
Las siete ternas primitivas correspondientes son
$$ (3,4,5),\ (5,12,13),\ (8,15,17),\ (7,24,25),\ (20,21,29),\ (12,35,37),\ (9,40,41). $$
Este ejemplo muestra con claridad cómo se combinan el conteo geométrico bruto y la sustracción por divisores impares comunes.
Las implementaciones en C++, Python y Java comienzan con rutinas exactas de raíz cuadrada y raíz cúbica entera, para que todos los cortes como \(\lfloor\sqrt{x}\rfloor\) y \(\lfloor N^{1/3}\rfloor\) sean fiables a pesar de usar aproximaciones de punto flotante como punto de partida.
Después calculan el contador bruto \(R(X)\) combinando un prefijo cerrado \(K^2\) con un barrido final monótono. En ese barrido sólo se mantiene un límite decreciente basado en raíces cuadradas y se alterna la fórmula de paridad para filas impares y pares.
A continuación, la implementación rellena una primera tabla con todos los valores \(P(x)\) para \(x\le L\). Cada entrada usa la recurrencia, separando los cocientes grandes tratados individualmente de los cocientes pequeños agrupados por multiplicidad.
Luego se construye una segunda tabla para los argumentos transformados \(\left\lfloor N/t^2\right\rfloor\) con \(t\) impar, procesados de mayor a menor. Los argumentos reducidos pequeños salen de la primera tabla; los grandes coinciden con otros valores transformados ya calculados.
La entrada correspondiente a \(t=1\) es \(P(N)\), que coincide exactamente con el número de ternas pitagóricas primitivas con hipotenusa a lo sumo \(N\).
Una evaluación individual del contador bruto \(R(X)\) cuesta \(O(\sqrt{X})\) porque el barrido final sólo recorre las filas cercanas a la frontera. El trabajo adicional de la recurrencia aporta alrededor de \(O(X^{1/3})\) operaciones agrupadas sobre divisores.
Construir la tabla pequeña hasta \(L=\lfloor N^{1/3}\rfloor\) requiere \(O(N^{1/2})\) tiempo. En la fase transformada, el coste dominante del contador bruto suma
$$\sum_{\substack{t\le L\\ t\equiv 1 \pmod{2}}} O\left(\sqrt{\frac{N}{t^2}}\right)=O(\sqrt{N}\log N).$$
El trabajo extra de agrupación de cocientes es menor que ese término principal. La memoria usada es \(O(N^{1/3})\), porque las dos tablas almacenadas tienen tamaño proporcional al umbral cúbico.
我们要求的是所有满足
$$c \le N$$
的 primitive 勾股三元组 \((a,b,c)\) 的个数。
按照 Euclid 参数化,每一个 primitive 三元组都唯一对应于
$$a=m^2-n^2,\qquad b=2mn,\qquad c=m^2+n^2,$$
其中
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
因此,这道题等价于统计满足这些条件且位于区域 \(m^2+n^2\le N\) 内的参数对 \((m,n)\)。
记 \(P(X)\) 为满足 \(m^2+n^2\le X\)、\(m>n>0\)、互素且奇偶性相反的 primitive 参数对数量。题目要求的就是 \(P(N)\)。实现采用两层结构:先快速计算不带互素条件的原始计数,再用递推把非 primitive 的情况扣除掉。
经典参数化给出了 primitive 勾股三元组与参数对 \((m,n)\) 之间的一一对应关系,条件是
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
对应的斜边就是
$$c=m^2+n^2.$$
所以,统计 \(c\le N\) 的 primitive 三元组,完全等价于统计四分之一圆区域
$$m^2+n^2\le N$$
中的这些合法参数对。整个算法的本质就是把几何边界和数论约束拆开处理。
定义 \(R(X)\) 为满足
$$m>n>0,\qquad m^2+n^2\le X,\qquad m\not\equiv n \pmod{2}$$
的参数对数量,但暂时不要求 \(\gcd(m,n)=1\)。
对固定的 \(m\),允许的最大 \(n\) 为
$$n_{\max}(m)=\min\left(m-1,\left\lfloor\sqrt{X-m^2}\right\rfloor\right).$$
因此这一行的贡献是
$$r_X(m)=\begin{cases} \left\lfloor \dfrac{n_{\max}(m)}{2}\right\rfloor, & m\text{ 为奇数},\\[6pt] \left\lfloor \dfrac{n_{\max}(m)+1}{2}\right\rfloor, & m\text{ 为偶数}. \end{cases}$$
于是
$$R(X)=\sum_{m\ge 1} r_X(m).$$
实现并不是逐行暴力枚举,而是先把一个“完整前缀”单独拿出来。如果
$$m^2+(m-1)^2\le X,$$
那么整段 \(1\le n\le m-1\) 都在圆内。解这个不等式可得
$$m\le \frac{1+\sqrt{2X-1}}{2}.$$
设
$$m_0=\left\lfloor\frac{1+\sqrt{2X-1}}{2}\right\rfloor,\qquad K=\left\lfloor\frac{m_0}{2}\right\rfloor.$$
前 \(2K\) 行的总贡献恰好是
$$1+3+5+\cdots+(2K-1)=K^2,$$
因为第 \(2j-1\) 行贡献 \(j-1\),第 \(2j\) 行贡献 \(j\)。超过这个完整前缀后,\(\lfloor\sqrt{X-m^2}\rfloor\) 会单调下降,所以只需要维护一个不断减小的平方根边界即可完成尾部扫描。
现在恢复互素条件。由于 \(m\) 与 \(n\) 奇偶性相反,它们不可能同时被 \(2\) 整除,因此任何公因子都只能是奇数。
如果某个参数对的奇公因子为 \(d\),那么可以写成
$$m=d\,m',\qquad n=d\,n',$$
其中 \((m',n')\) 就是 primitive 对。圆形约束变成
$$m'^2+n'^2\le \left\lfloor \frac{X}{d^2}\right\rfloor.$$
因此,每一个原始参数对都可以唯一地表示为“一个 primitive 对”乘上“一个奇数缩放因子”,从而得到
$$R(X)=\sum_{\substack{d\ge 1\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
整理后便有递推式
$$P(X)=R(X)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
代码实际计算的正是这个递推。它与只在奇数因子上进行的 Möbius 反演等价,但实现并不需要显式构造 Möbius 函数表。
令
$$L=\left\lfloor N^{1/3}\right\rfloor.$$
首先顺序计算所有 \(1\le x\le L\) 的 \(P(x)\)。对于固定的 \(x\),较小的奇数因子 \(d\le x^{1/3}\) 在下面这个式子里逐个处理:
$$P(x)=R(x)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{x}{d^2}\right\rfloor\right).$$
而对于更大的 \(d\),\(\left\lfloor x/d^2\right\rfloor\) 已经变小,很多不同的 \(d\) 会给出同一个商值。实现把这些重复项合并。对固定的小整数 \(z\),满足
$$\left\lfloor \frac{x}{d^2}\right\rfloor=z$$
的奇数 \(d\) 的个数是
$$\left\lfloor\frac{\sqrt{x/z}-1}{2}\right\rfloor-\left\lfloor\frac{\sqrt{x/(z+1)}-1}{2}\right\rfloor.$$
这样,同一个递推项不再重复减很多次,而是减一次再乘上对应的出现次数。
小表算完之后,再计算变换后的值
$$Q(t)=P\left(\left\lfloor \frac{N}{t^2}\right\rfloor\right),\qquad t\text{ 为奇数},$$
并按奇数 \(t\) 从大到小处理。若递推中出现的缩小参数不超过 \(L\),就直接查小表;若仍然大于 \(L\),它一定对应另一个更大奇因子的变换值 \(Q(tu)\),而那个值已经在更早一步算好了。最终答案就是
$$Q(1)=P(N).$$
满足 \(m^2+n^2\le 50\) 且奇偶性相反的参数对共有
$$ (2,1),\ (3,2),\ (4,1),\ (4,3),\ (5,2),\ (5,4),\ (6,1),\ (6,3). $$
所以 \(R(50)=8\)。
大于 \(1\) 的奇公因子只可能是 \(3\),因为 \(5^2>50\)。因此
$$P(50)=R(50)-P\left(\left\lfloor\frac{50}{9}\right\rfloor\right)=8-P(5).$$
而 \(P(5)=1\),对应唯一的参数对 \((2,1)\)。于是
$$P(50)=8-1=7.$$
对应的七个 primitive 勾股三元组是
$$ (3,4,5),\ (5,12,13),\ (8,15,17),\ (7,24,25),\ (20,21,29),\ (12,35,37),\ (9,40,41). $$
这个小例子清楚地展示了“先做几何计数,再扣掉有奇公因子的项”这一核心思路。
C++、Python 和 Java 实现首先使用精确的整数平方根与整数立方根过程。虽然一开始会借助浮点近似值,但随后都会向上或向下微调,因此所有边界判断都保持正确。
接着,程序利用上面的分解方法计算原始计数 \(R(X)\):前面完整填满的行直接用 \(K^2\) 给出,后面的边界部分则通过一个单调下降的平方根指针完成扫描,从而避免重复工作。
然后实现按从小到大的顺序填充第一张表,得到所有 \(x\le L\) 的 primitive 计数 \(P(x)\)。在每个表项中,较大的商值单独处理,较小的商值则按相同结果分组后统一扣除。
之后再建立第二张表,存放形如 \(\left\lfloor N/t^2\right\rfloor\) 的变换参数,其中 \(t\) 为奇数,并且按 \(t\) 从大到小处理。小参数直接引用第一张表,大参数则对应第二张表中已经算出的后续项。
当 \(t=1\) 时,得到的就是 \(P(N)\),也就是题目要求的 primitive 勾股三元组总数。
单次计算原始计数 \(R(X)\) 的时间复杂度为 \(O(\sqrt{X})\),因为尾部扫描只经过靠近圆弧边界的那些行。围绕它的递推部分大约还需要 \(O(X^{1/3})\) 次按商分组的操作。
构造到 \(L=\lfloor N^{1/3}\rfloor\) 为止的小表,总耗时为 \(O(N^{1/2})\)。在变换阶段,主导项来自原始计数器的总成本:
$$\sum_{\substack{t\le L\\ t\equiv 1 \pmod{2}}} O\left(\sqrt{\frac{N}{t^2}}\right)=O(\sqrt{N}\log N).$$
额外的分组递推开销比这个主导项更小。内存复杂度为 \(O(N^{1/3})\),因为两张表的大小都与立方根阈值成正比。
Нужно посчитать число примитивных пифагоровых троек \((a,b,c)\), для которых
$$c \le N.$$
Параметризация Евклида говорит, что каждая примитивная тройка единственным образом задаётся формулами
$$a=m^2-n^2,\qquad b=2mn,\qquad c=m^2+n^2,$$
где
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
Значит, задача сводится к подсчёту таких пар \((m,n)\) внутри области \(m^2+n^2\le N\).
Обозначим через \(P(X)\) число примитивных пар \((m,n)\), удовлетворяющих условиям \(m^2+n^2\le X\), \(m>n>0\), взаимная простота и разная чётность. Тогда искомый ответ равен \(P(N)\). Реализация строится в два слоя: сначала быстро считается грубое число пар без условия взаимной простоты, а затем рекуррентно вычитаются непримитивные случаи.
Классическая параметризация устанавливает биекцию между примитивными тройками и парами \((m,n)\), для которых
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
Гипотенуза при этом равна
$$c=m^2+n^2.$$
Следовательно, подсчёт примитивных троек с \(c\le N\) полностью эквивалентен подсчёту допустимых пар в четверти круга
$$m^2+n^2\le N.$$
Именно поэтому программа работает не с самими тройками, а с параметрами \((m,n)\): геометрия задаёт границу, а арифметика отвечает за примитивность.
Пусть \(R(X)\) обозначает число пар \((m,n)\), для которых
$$m>n>0,\qquad m^2+n^2\le X,\qquad m\not\equiv n \pmod{2},$$
но взаимная простота пока не требуется.
Для фиксированного \(m\) максимальное допустимое \(n\) равно
$$n_{\max}(m)=\min\left(m-1,\left\lfloor\sqrt{X-m^2}\right\rfloor\right).$$
Тогда вклад строки имеет вид
$$r_X(m)=\begin{cases} \left\lfloor \dfrac{n_{\max}(m)}{2}\right\rfloor, & m\text{ нечётно},\\[6pt] \left\lfloor \dfrac{n_{\max}(m)+1}{2}\right\rfloor, & m\text{ чётно}. \end{cases}$$
Следовательно,
$$R(X)=\sum_{m\ge 1} r_X(m).$$
Код ускоряет эту сумму, выделяя полный начальный блок строк. Если
$$m^2+(m-1)^2\le X,$$
то весь интервал \(1\le n\le m-1\) целиком лежит внутри круга. Решая неравенство, получаем
$$m\le \frac{1+\sqrt{2X-1}}{2}.$$
Положим
$$m_0=\left\lfloor\frac{1+\sqrt{2X-1}}{2}\right\rfloor,\qquad K=\left\lfloor\frac{m_0}{2}\right\rfloor.$$
Первые \(2K\) строк дают ровно вклад
$$1+3+5+\cdots+(2K-1)=K^2,$$
потому что строка \(2j-1\) даёт \(j-1\), а строка \(2j\) даёт \(j\). После этого префикса величина \(\lfloor\sqrt{X-m^2}\rfloor\) убывает монотонно, и хвост можно пройти одним подвижным указателем.
Теперь возвращаем условие взаимной простоты. Разная чётность уже исключает общий множитель \(2\), значит любой общий делитель обязан быть нечётным.
Если у пары общий нечётный делитель равен \(d\), то
$$m=d\,m',\qquad n=d\,n',$$
и редуцированная пара \((m',n')\) уже примитивна. Ограничение по окружности превращается в
$$m'^2+n'^2\le \left\lfloor \frac{X}{d^2}\right\rfloor.$$
Следовательно, каждая грубая пара однозначно получается из примитивной пары и нечётного коэффициента растяжения:
$$R(X)=\sum_{\substack{d\ge 1\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
Отсюда следует рекурсия
$$P(X)=R(X)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
Именно эту форму используют реализации. Она эквивалентна нечётной версии обращения Мёбиуса, но явная таблица значений функции Мёбиуса не требуется.
Положим
$$L=\left\lfloor N^{1/3}\right\rfloor.$$
Сначала последовательно вычисляются все \(P(x)\) для \(1\le x\le L\). Для фиксированного \(x\) нечётные делители \(d\le x^{1/3}\) обрабатываются по одному в формуле
$$P(x)=R(x)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{x}{d^2}\right\rfloor\right).$$
Для больших \(d\) величина \(\left\lfloor x/d^2\right\rfloor\) становится маленькой, и многие разные \(d\) дают один и тот же результат. Такие повторяющиеся частные объединяются. Для фиксированного малого \(z\) число нечётных \(d\), удовлетворяющих
$$\left\lfloor \frac{x}{d^2}\right\rfloor=z,$$
равно
$$\left\lfloor\frac{\sqrt{x/z}-1}{2}\right\rfloor-\left\lfloor\frac{\sqrt{x/(z+1)}-1}{2}\right\rfloor.$$
Это заменяет множество одинаковых рекурсивных вычитаний одним вычитанием с коэффициентом кратности.
После построения малой таблицы реализации вычисляют ещё и преобразованные значения
$$Q(t)=P\left(\left\lfloor \frac{N}{t^2}\right\rfloor\right),\qquad t\text{ нечётно},$$
обрабатывая нечётные \(t\) в порядке убывания. Если уменьшенный аргумент не превосходит \(L\), он берётся из малой таблицы. Иначе это другой преобразованный элемент \(Q(tu)\) с большим нечётным множителем, уже посчитанный ранее. В итоге
$$Q(1)=P(N)$$
даёт окончательный ответ.
Пары разной чётности, удовлетворяющие \(m^2+n^2\le 50\), таковы:
$$ (2,1),\ (3,2),\ (4,1),\ (4,3),\ (5,2),\ (5,4),\ (6,1),\ (6,3). $$
Значит, \(R(50)=8\).
Единственный возможный нечётный общий делитель больше \(1\) здесь равен \(3\), потому что \(5^2>50\). Поэтому
$$P(50)=R(50)-P\left(\left\lfloor\frac{50}{9}\right\rfloor\right)=8-P(5).$$
А \(P(5)=1\), это единственная пара \((2,1)\). Следовательно,
$$P(50)=8-1=7.$$
Соответствующие семь примитивных троек равны
$$ (3,4,5),\ (5,12,13),\ (8,15,17),\ (7,24,25),\ (20,21,29),\ (12,35,37),\ (9,40,41). $$
Этот пример наглядно показывает, как грубый геометрический подсчёт соединяется с рекурсивным вычитанием по нечётным общим делителям.
Реализации на C++, Python и Java начинают с точных процедур вычисления целочисленных квадратных и кубических корней. Это гарантирует корректность всех пороговых сравнений, даже если первоначальная оценка берётся из вещественной арифметики.
Затем вычисляется грубый счётчик \(R(X)\): заполненный начальный префикс даётся формулой \(K^2\), а оставшаяся часть границы обходится монотонным сканированием с одним убывающим корневым указателем.
После этого заполняется первая таблица со всеми значениями \(P(x)\) при \(x\le L\). Для каждого элемента большие частные обрабатываются отдельно, а малые частные группируются по совпадающему значению.
Далее строится вторая таблица для аргументов вида \(\left\lfloor N/t^2\right\rfloor\), где \(t\) нечётно, причём обработка идёт от больших \(t\) к меньшим. Малые редуцированные аргументы берутся из первой таблицы, а большие совпадают с уже вычисленными преобразованными значениями.
Элемент при \(t=1\) равен \(P(N)\), то есть именно искомому числу примитивных пифагоровых троек с гипотенузой не больше \(N\).
Одна оценка грубого счётчика \(R(X)\) занимает \(O(\sqrt{X})\) времени, поскольку хвостовой проход затрагивает только строки у границы круга. Рекурсивная часть добавляет примерно \(O(X^{1/3})\) операций, связанных с группировкой одинаковых частных.
Построение малой таблицы до \(L=\lfloor N^{1/3}\rfloor\) требует \(O(N^{1/2})\) времени. На преобразованной стадии доминирующий вклад от грубого счётчика равен
$$\sum_{\substack{t\le L\\ t\equiv 1 \pmod{2}}} O\left(\sqrt{\frac{N}{t^2}}\right)=O(\sqrt{N}\log N).$$
Дополнительная работа по группировке частных меньше этого ведущего члена. Память имеет порядок \(O(N^{1/3})\), так как обе таблицы пропорциональны кубическому корню от \(N\).
نريد حساب عدد الثلاثيات الفيثاغورية البدائية \((a,b,c)\) التي تحقق
$$c \le N.$$
وفق تمثيل إقليدس، كل ثلاثية بدائية تُكتب بصورة وحيدة على الشكل
$$a=m^2-n^2,\qquad b=2mn,\qquad c=m^2+n^2,$$
مع الشروط
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
لذلك تصبح المسألة مساوية تمامًا لعدّ الأزواج \((m,n)\) التي تحقق هذه الشروط داخل المنطقة \(m^2+n^2\le N\).
لنرمز بـ \(P(X)\) إلى عدد أزواج المعاملات البدائية \((m,n)\) التي تحقق \(m^2+n^2\le X\)، و\(m>n>0\)، وهي متباينة الزوجية ومتباينة القاسم المشترك. المطلوب هو \(P(N)\). تعتمد التطبيقات على طبقتين: أولًا عدّ خام سريع من دون شرط البدائية، ثم علاقة تراجعية تزيل الأزواج غير البدائية.
التمثيل الكلاسيكي يعطي تقابلًا واحدًا لواحد بين الثلاثيات البدائية والأزواج \((m,n)\) التي تحقق
$$m>n>0,\qquad \gcd(m,n)=1,\qquad m\not\equiv n \pmod{2}.$$
وفي هذا التمثيل تكون الوتر مساويًا لـ
$$c=m^2+n^2.$$
إذًا عدّ الثلاثيات البدائية التي تحقق \(c\le N\) هو بالضبط عدّ هذه الأزواج المسموح بها داخل ربع الدائرة
$$m^2+n^2\le N.$$
هذه الصياغة هي أساس الحل كله: الحدّ الهندسي يأتي من الدائرة، وشرط البدائية يأتي من الحسابيات.
عرّف \(R(X)\) بأنه عدد الأزواج \((m,n)\) التي تحقق
$$m>n>0,\qquad m^2+n^2\le X,\qquad m\not\equiv n \pmod{2},$$
ولكن من دون أن نفرض بعدُ الشرط \(\gcd(m,n)=1\).
إذا ثبتنا \(m\)، فإن أكبر قيمة مسموحة لـ \(n\) هي
$$n_{\max}(m)=\min\left(m-1,\left\lfloor\sqrt{X-m^2}\right\rfloor\right).$$
وعليه يكون إسهام هذا الصف
$$r_X(m)=\begin{cases} \left\lfloor \dfrac{n_{\max}(m)}{2}\right\rfloor, & m\equiv 1 \pmod{2},\\[6pt] \left\lfloor \dfrac{n_{\max}(m)+1}{2}\right\rfloor, & m\equiv 0 \pmod{2}. \end{cases}$$
ومن ثم
$$R(X)=\sum_{m\ge 1} r_X(m).$$
تسرّع الشيفرة هذا المجموع بفصل بادئة كاملة من الصفوف. فإذا تحقق
$$m^2+(m-1)^2\le X,$$
فإن كامل المجال \(1\le n\le m-1\) يقع داخل الدائرة. بحل هذه المتباينة نحصل على
$$m\le \frac{1+\sqrt{2X-1}}{2}.$$
ولنضع
$$m_0=\left\lfloor\frac{1+\sqrt{2X-1}}{2}\right\rfloor,\qquad K=\left\lfloor\frac{m_0}{2}\right\rfloor.$$
عندئذ يكون مجموع الصفوف الأولى حتى \(2K\) مساويًا تمامًا لـ
$$1+3+5+\cdots+(2K-1)=K^2,$$
لأن الصف \(2j-1\) يضيف \(j-1\)، والصف \(2j\) يضيف \(j\). وبعد هذه البادئة يبدأ الحد الأعلى \(\lfloor\sqrt{X-m^2}\rfloor\) في الانخفاض الرتيب، ولذلك يكفي مؤشر جذر تربيعي متحرّك واحد لمسح الذيل.
الآن نعيد شرط التباين في القاسم المشترك. اختلاف الزوجية يمنع أصلًا وجود عامل مشترك \(2\)، ولذلك فإن أي قاسم مشترك يجب أن يكون فرديًا.
إذا كان للأزوجة قاسم مشترك فردي \(d\)، فإننا نكتب
$$m=d\,m',\qquad n=d\,n',$$
ويصبح الزوج المختزل \((m',n')\) بدائيًا. كما يتحول قيد الدائرة إلى
$$m'^2+n'^2\le \left\lfloor \frac{X}{d^2}\right\rfloor.$$
وبالتالي فإن كل زوج خام ينتج بصورة وحيدة من زوج بدائي ومعامل تمديد فردي، وهذا يعطي
$$R(X)=\sum_{\substack{d\ge 1\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
ومنها نستنتج العلاقة
$$P(X)=R(X)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{X}{d^2}\right\rfloor\right).$$
هذه هي الصيغة التي تنفذها التطبيقات مباشرة. وهي مكافئة لعكس Möbius على القواسم الفردية فقط، لكن التنفيذ لا يحتاج إلى جدول صريح لقيم هذه الدالة.
لنعرّف
$$L=\left\lfloor N^{1/3}\right\rfloor.$$
أولًا تُحسب جميع القيم \(P(x)\) عندما \(1\le x\le L\). ولعدد ثابت \(x\)، تُعالج القواسم الفردية الصغيرة \(d\le x^{1/3}\) واحدةً واحدة في العلاقة
$$P(x)=R(x)-\sum_{\substack{d\ge 3\\ d\equiv 1 \pmod{2}}} P\left(\left\lfloor \frac{x}{d^2}\right\rfloor\right).$$
أما القيم الأكبر من \(d\)، فإن \(\left\lfloor x/d^2\right\rfloor\) فيها يصبح صغيرًا، وتنتج كثير من القواسم المختلفة القيمة نفسها. لذلك تُجمَّع الحدود المتكررة. ولعدد صغير ثابت \(z\)، فإن عدد القواسم الفردية \(d\) التي تحقق
$$\left\lfloor \frac{x}{d^2}\right\rfloor=z$$
يساوي
$$\left\lfloor\frac{\sqrt{x/z}-1}{2}\right\rfloor-\left\lfloor\frac{\sqrt{x/(z+1)}-1}{2}\right\rfloor.$$
وبهذا يُستبدل طرح الحد نفسه مرات كثيرة بطرح واحد مضروبًا في عدد مرات ظهوره.
بعد اكتمال الجدول الصغير، تُحسب أيضًا القيم المحوّلة
$$Q(t)=P\left(\left\lfloor \frac{N}{t^2}\right\rfloor\right),\qquad t\equiv 1 \pmod{2},$$
مع معالجة القيم الفردية \(t\) من الأكبر إلى الأصغر. فإذا هبط الوسيط المختزل إلى ما دون \(L\)، نقرأه من الجدول الصغير. وإلا فإنه يساوي قيمة محوّلة أخرى \(Q(tu)\) ذات عامل فردي أكبر، وهي تكون محسوبة مسبقًا. وفي النهاية نحصل على
$$Q(1)=P(N),$$
وهو الجواب المطلوب.
الأزواج ذات الزوجية المختلفة التي تحقق \(m^2+n^2\le 50\) هي
$$ (2,1),\ (3,2),\ (4,1),\ (4,3),\ (5,2),\ (5,4),\ (6,1),\ (6,3). $$
إذن \(R(50)=8\).
القاسم المشترك الفردي الوحيد الممكن الأكبر من \(1\) هو \(3\)، لأن \(5^2>50\). لذلك
$$P(50)=R(50)-P\left(\left\lfloor\frac{50}{9}\right\rfloor\right)=8-P(5).$$
ولدينا \(P(5)=1\)، وهو الزوج الوحيد \((2,1)\). وعليه
$$P(50)=8-1=7.$$
والثلاثيات البدائية السبع المقابلة هي
$$ (3,4,5),\ (5,12,13),\ (8,15,17),\ (7,24,25),\ (20,21,29),\ (12,35,37),\ (9,40,41). $$
يبين هذا المثال بدقة كيف يتكامل العدّ الهندسي الخام مع طرح الأزواج ذات القواسم الفردية المشتركة.
تبدأ تطبيقات C++ وPython وJava بحسابات دقيقة للجذر التربيعي والجذر التكعيبي الصحيحين، بحيث تبقى كل حدود القطع مثل \(\lfloor\sqrt{x}\rfloor\) و\(\lfloor N^{1/3}\rfloor\) صحيحة حتى لو استُخدم تقريب عشري أولي.
بعد ذلك تُحسب الدالة الخام \(R(X)\) عبر جمع بادئة مغلقة مقدارها \(K^2\) مع مسح أحادي الاتجاه للذيل. وفي هذا الذيل يُحافَظ فقط على حد جذر تربيعي متناقص، وتُطبَّق صيغة التناوب بين الصفوف الفردية والزوجية.
ثم تملأ الشيفرة جدولًا أولًا يحتوي جميع القيم \(P(x)\) حتى \(x\le L\). كل عنصر من هذا الجدول يستخدم العلاقة التراجعية نفسها، مع فصل القيم الكبيرة للقسمة عن القيم الصغيرة المجمعة حسب التكرار.
بعد ذلك تُبنى جدولية ثانية للوسائط المحوّلة من الشكل \(\left\lfloor N/t^2\right\rfloor\) حيث \(t\) فردي، وتُعالَج هذه القيم من الأكبر إلى الأصغر. الوسائط الصغيرة تُقرأ من الجدول الأول، أما الوسائط الكبيرة فتطابق قيمًا محوّلة حُسبت سابقًا.
العنصر الموافق لـ \(t=1\) يساوي \(P(N)\)، وهو بالضبط عدد الثلاثيات الفيثاغورية البدائية ذات الوتر الذي لا يتجاوز \(N\).
الحساب الواحد للدالة الخام \(R(X)\) يحتاج إلى زمن \(O(\sqrt{X})\)، لأن مسح الذيل يزور فقط الصفوف الواقعة قرب حد الدائرة. وتضيف العلاقة التراجعية حوله نحو \(O(X^{1/3})\) من عمليات التجميع على القواسم.
بناء الجدول الصغير حتى \(L=\lfloor N^{1/3}\rfloor\) يكلّف \(O(N^{1/2})\) زمنًا. وفي المرحلة المحوّلة يكون الحد المسيطر القادم من العدّ الخام هو
$$\sum_{\substack{t\le L\\ t\equiv 1 \pmod{2}}} O\left(\sqrt{\frac{N}{t^2}}\right)=O(\sqrt{N}\log N).$$
أما عمل التجميع الإضافي فهو أصغر من هذا الحد المسيطر. واستهلاك الذاكرة هو \(O(N^{1/3})\)، لأن الجدولين المخزنين كلاهما بحجم متناسب مع عتبة الجذر التكعيبي.