The divisor-based quantity in this problem can be reindexed as an ordered factor-pair sum:
$$S(N)=\sum_{n=1}^{N}\sum_{d\mid n}\gcd\!\left(d,\frac{n}{d}\right)=\sum_{ab\le N}\gcd(a,b).$$
So the task is to add \(\gcd(a,b)\) over all ordered pairs \((a,b)\) whose product is at most \(N\). At \(N=10^{15}\), direct enumeration is impossible, so the solution rewrites the sum using Euler's totient function and the divisor summatory function.
We work with the equivalent form
$$S(N)=\sum_{ab\le N}\gcd(a,b).$$
The efficient formula comes from separating the common-divisor contribution, then grouping repeated floor values so the large range can be handled in blocks.
For each \(n\) and each divisor \(d\mid n\), the pair \((d,n/d)\) is an ordered factor pair whose product is exactly \(n\). Summing over all \(n\le N\) and all divisors therefore gives exactly the same set of terms as summing over all ordered pairs \((a,b)\) with \(ab\le N\).
This is the first important simplification: the original divisor language becomes a two-variable summation problem with a clean product bound.
Use the classical identity
$$\sum_{t\mid m}\varphi(t)=m.$$
Applying it with \(m=\gcd(a,b)\) gives
$$\gcd(a,b)=\sum_{t\mid \gcd(a,b)}\varphi(t).$$
Substitute this into the sum and swap the order of summation:
$$S(N)=\sum_{ab\le N}\sum_{t\mid a,\ t\mid b}\varphi(t).$$
Now write \(a=tx\) and \(b=ty\). The condition becomes
$$t^2xy\le N.$$
Hence
$$S(N)=\sum_{t\le \sqrt N}\varphi(t)\,D\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right),$$
where \(D(m)\) counts ordered pairs \((x,y)\) with \(xy\le m\).
The pair-counting function is
$$D(m)=\#\{(x,y)\in\mathbb{Z}_{>0}^2:xy\le m\}.$$
For a fixed \(x\), there are exactly \(\left\lfloor m/x\right\rfloor\) valid choices for \(y\). Therefore
$$D(m)=\sum_{x=1}^{m}\left\lfloor\frac{m}{x}\right\rfloor.$$
Equivalently, every integer \(n\le m\) contributes once for each divisor pair \((x,y)\) with \(xy=n\), so
$$D(m)=\sum_{n=1}^{m}\tau(n),$$
where \(\tau(n)\) is the divisor-count function. The whole problem is now a single weighted sum of totients times divisor-summatory values.
A naive loop up to \(m\) is still too slow when \(m\) is large. The key observation is that the quotient
$$q=\left\lfloor\frac{m}{\ell}\right\rfloor$$
stays constant on an interval \([\ell,r]\), where
$$r=\left\lfloor\frac{m}{q}\right\rfloor.$$
So one whole block contributes
$$q\,(r-\ell+1).$$
Advancing from one block to the next yields
$$D(m)=\sum_{\text{blocks }[\ell,r]}\left\lfloor\frac{m}{\ell}\right\rfloor(r-\ell+1),$$
and the number of blocks is only \(O(\sqrt m)\). This is the hyperbola-style speedup used whenever \(D(m)\) has to be evaluated for a large argument.
Let
$$K=\left\lfloor N^{1/3}\right\rfloor,\qquad M=\left\lfloor\frac{N}{K^2}\right\rfloor.$$
Then the outer sum naturally splits into two regimes.
For \(t>K\), we have \(\left\lfloor N/t^2\right\rfloor\le M\), so all needed \(D(m)\) values are small. The implementation precomputes \(\tau(1),\dots,\tau(M)\), takes prefix sums, and obtains every small \(D(m)\) by table lookup.
For \(t\le K\), the argument \(\left\lfloor N/t^2\right\rfloor\) is large, but many consecutive \(t\) share the same value. If
$$v=\left\lfloor\frac{N}{t^2}\right\rfloor,$$
then the largest index with the same quotient is
$$r=\left\lfloor\sqrt{\frac{N}{v}}\right\rfloor.$$
That whole interval contributes
$$D(v)\sum_{u=t}^{r}\varphi(u).$$
Prefix sums of \(\varphi\) make the interval sum an \(O(1)\) lookup, so the expensive work is concentrated only in the distinct large \(D(v)\) calls.
Here \(\lfloor\sqrt{10}\rfloor=3\), so the transformed formula becomes
$$S(10)=\sum_{t=1}^{3}\varphi(t)\,D\!\left(\left\lfloor\frac{10}{t^2}\right\rfloor\right).$$
The required totients are
$$\varphi(1)=1,\qquad \varphi(2)=1,\qquad \varphi(3)=2.$$
The required divisor-summatory values are
$$D(10)=\sum_{x=1}^{10}\left\lfloor\frac{10}{x}\right\rfloor=27,\qquad D(2)=3,\qquad D(1)=1.$$
Therefore
$$S(10)=1\cdot 27+1\cdot 3+2\cdot 1=32,$$
which matches the checkpoint used by the implementation.
The C++, Python, and Java implementations all follow the same number-theoretic formula. They first compute the integer square-root limit \(\lfloor\sqrt N\rfloor\) and the cube-root split point \(\lfloor N^{1/3}\rfloor\). Next they build a totient sieve up to \(\lfloor\sqrt N\rfloor\) and store prefix sums so that \(\sum_{u=l}^{r}\varphi(u)\) can be recovered instantly.
For the small-argument region, the implementation precomputes divisor counts up to \(M=\lfloor N/K^2\rfloor\) by visiting multiples of each divisor, then turns those counts into prefix sums to obtain \(D(1),D(2),\dots,D(M)\). For the large-argument region, it groups equal values of \(\left\lfloor N/t^2\right\rfloor\), evaluates each large \(D(v)\) with quotient blocks, and multiplies by the corresponding totient-interval sum.
The C++ and Java implementations parallelize those independent large-group contributions before adding the small tail. The Python implementation is a thin wrapper: it builds the compiled solver when needed, runs it, and returns the printed result. The compiled version also checks the two checkpoints \(S(10)=32\) and \(S(1000)=12776\) before evaluating the full target.
Let \(L=\lfloor\sqrt N\rfloor\), \(K=\lfloor N^{1/3}\rfloor\), and \(M=\lfloor N/K^2\rfloor\). The totient sieve costs \(O(L\log\log L)\) time and \(O(L)\) memory. The small-table preparation for divisor counts costs
$$\sum_{d=1}^{M}\left\lfloor\frac{M}{d}\right\rfloor=O(M\log M)$$
time and \(O(M)\) memory.
In the large-value region there are only \(O(K)\) distinct quotient groups, and one call to the block method for \(D(v)\) costs \(O(\sqrt v)\). Summing those costs over \(t\le K\) gives
$$O\!\left(\sum_{t=1}^{K}\sqrt{\frac{N}{t^2}}\right)=O(\sqrt N\log K).$$
So the overall running time is \(O(\sqrt N\log N)\) in the worst case, while memory is dominated by the totient arrays and remains \(O(\sqrt N)\). Parallel execution improves wall-clock time for the grouped large-value part but does not change the asymptotic bound.
Die in der Aufgabe geforderte Divisorensumme kann als Summe über geordnete Faktorenpaare umgeschrieben werden:
$$S(N)=\sum_{n=1}^{N}\sum_{d\mid n}\gcd\!\left(d,\frac{n}{d}\right)=\sum_{ab\le N}\gcd(a,b).$$
Gesucht ist also die Summe von \(\gcd(a,b)\) über alle geordneten Paare \((a,b)\) mit \(ab\le N\). Für \(N=10^{15}\) ist direkte Enumeration unmöglich, daher formt die Lösung die Aufgabe mit der Eulerschen Phi-Funktion und der Divisor-Summenfunktion um.
Wir arbeiten mit der äquivalenten Darstellung
$$S(N)=\sum_{ab\le N}\gcd(a,b).$$
Die effiziente Formel entsteht, indem man zuerst den gemeinsamen Teilerbeitrag isoliert und danach gleiche Ganzzahlquotienten blockweise zusammenfasst.
Zu jedem \(n\) und jedem Teiler \(d\mid n\) gehört das geordnete Faktorenpaar \((d,n/d)\) mit Produkt \(n\). Wenn man also über alle \(n\le N\) und alle ihre Teiler summiert, erhält man exakt dieselben Terme wie bei der Summe über alle geordneten Paare \((a,b)\) mit \(ab\le N\).
Diese Umformung ist der erste wichtige Schritt, weil die saubere Produktbedingung \(ab\le N\) analytisch viel leichter zu behandeln ist als die ursprüngliche doppelte Divisorensumme.
Verwendet wird die klassische Identität
$$\sum_{t\mid m}\varphi(t)=m.$$
Für \(m=\gcd(a,b)\) folgt daraus
$$\gcd(a,b)=\sum_{t\mid \gcd(a,b)}\varphi(t).$$
Setzt man dies in \(S(N)\) ein und vertauscht die Summationen, so erhält man
$$S(N)=\sum_{ab\le N}\sum_{t\mid a,\ t\mid b}\varphi(t).$$
Schreibt man nun \(a=tx\) und \(b=ty\), dann lautet die Bedingung
$$t^2xy\le N.$$
Daraus wird
$$S(N)=\sum_{t\le \sqrt N}\varphi(t)\,D\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right),$$
wobei \(D(m)\) die Zahl geordneter Paare \((x,y)\) mit \(xy\le m\) bezeichnet.
Die Paarzählfunktion ist
$$D(m)=\#\{(x,y)\in\mathbb{Z}_{>0}^2:xy\le m\}.$$
Für festes \(x\) gibt es genau \(\left\lfloor m/x\right\rfloor\) zulässige Werte für \(y\). Also
$$D(m)=\sum_{x=1}^{m}\left\lfloor\frac{m}{x}\right\rfloor.$$
Äquivalent dazu trägt jede Zahl \(n\le m\) genau einmal pro Teilerpaar \((x,y)\) mit \(xy=n\) bei, also auch
$$D(m)=\sum_{n=1}^{m}\tau(n),$$
wobei \(\tau(n)\) die Anzahl der positiven Teiler von \(n\) ist. Das gesamte Problem ist damit auf eine gewichtete Summe aus Phi-Werten und Divisor-Summenwerten reduziert.
Ein naiver Lauf bis \(m\) bleibt für große \(m\) zu langsam. Die entscheidende Beobachtung ist, dass der Quotient
$$q=\left\lfloor\frac{m}{\ell}\right\rfloor$$
auf einem ganzen Intervall \([\ell,r]\) konstant bleibt, wobei
$$r=\left\lfloor\frac{m}{q}\right\rfloor.$$
Ein kompletter Block liefert daher den Beitrag
$$q\,(r-\ell+1).$$
Damit erhält man
$$D(m)=\sum_{\text{Blöcke }[\ell,r]}\left\lfloor\frac{m}{\ell}\right\rfloor(r-\ell+1),$$
und die Anzahl solcher Blöcke ist nur \(O(\sqrt m)\). Genau diese Hyperbeltechnik wird verwendet, wenn \(D(m)\) für große Argumente ausgewertet werden muss.
Setze
$$K=\left\lfloor N^{1/3}\right\rfloor,\qquad M=\left\lfloor\frac{N}{K^2}\right\rfloor.$$
Dann zerfällt die äußere Summe natürlich in zwei Bereiche.
Für \(t>K\) gilt \(\left\lfloor N/t^2\right\rfloor\le M\). Alle benötigten \(D(m)\)-Werte sind dann klein, also berechnet die Implementierung \(\tau(1),\dots,\tau(M)\) vorab, bildet Präfixsummen und liest jedes kleine \(D(m)\) direkt aus einer Tabelle ab.
Für \(t\le K\) ist das Argument \(\left\lfloor N/t^2\right\rfloor\) groß, aber viele aufeinanderfolgende \(t\) liefern denselben Wert. Ist
$$v=\left\lfloor\frac{N}{t^2}\right\rfloor,$$
dann ist der größte Index mit demselben Quotienten
$$r=\left\lfloor\sqrt{\frac{N}{v}}\right\rfloor.$$
Der gesamte Intervallbeitrag ist also
$$D(v)\sum_{u=t}^{r}\varphi(u).$$
Mit Präfixsummen der Phi-Funktion wird diese Intervallsumme zu einem \(O(1)\)-Lookup, sodass die teure Arbeit nur noch in den verschiedenen großen \(D(v)\)-Auswertungen steckt.
Hier ist \(\lfloor\sqrt{10}\rfloor=3\), also
$$S(10)=\sum_{t=1}^{3}\varphi(t)\,D\!\left(\left\lfloor\frac{10}{t^2}\right\rfloor\right).$$
Die benötigten Phi-Werte sind
$$\varphi(1)=1,\qquad \varphi(2)=1,\qquad \varphi(3)=2.$$
Und die benötigten Divisor-Summenwerte sind
$$D(10)=\sum_{x=1}^{10}\left\lfloor\frac{10}{x}\right\rfloor=27,\qquad D(2)=3,\qquad D(1)=1.$$
Damit folgt
$$S(10)=1\cdot 27+1\cdot 3+2\cdot 1=32,$$
genau der Kontrollwert der Implementierung.
Die C++-, Python- und Java-Implementierungen folgen alle derselben zahlentheoretischen Formel. Zuerst bestimmen sie die ganzzahligen Grenzen \(\lfloor\sqrt N\rfloor\) und \(\lfloor N^{1/3}\rfloor\). Danach wird bis \(\lfloor\sqrt N\rfloor\) eine Phi-Siebung aufgebaut und in Präfixsummen überführt, sodass \(\sum_{u=l}^{r}\varphi(u)\) sofort verfügbar ist.
Für den Bereich kleiner Argumente berechnet die Implementierung die Teileranzahlen bis \(M=\lfloor N/K^2\rfloor\), indem sie die Vielfachen jedes Divisors besucht, und bildet daraus Präfixsummen, um \(D(1),D(2),\dots,D(M)\) direkt nachzuschlagen. Für den Bereich großer Argumente gruppiert sie gleiche Werte von \(\left\lfloor N/t^2\right\rfloor\), berechnet jedes große \(D(v)\) blockweise und multipliziert es mit der passenden Phi-Intervallsumme.
Die C++- und Java-Implementierungen parallelisieren diese unabhängigen Großgruppen, bevor der kleine Restbereich addiert wird. Die Python-Implementierung ist ein dünner Wrapper: Sie baut den kompilierten Solver bei Bedarf, führt ihn aus und gibt das ausgegebene Ergebnis zurück. Die kompilierte Fassung prüft außerdem vor der Endauswertung die beiden Kontrollwerte \(S(10)=32\) und \(S(1000)=12776\).
Seien \(L=\lfloor\sqrt N\rfloor\), \(K=\lfloor N^{1/3}\rfloor\) und \(M=\lfloor N/K^2\rfloor\). Das Phi-Sieb benötigt \(O(L\log\log L)\) Zeit und \(O(L)\) Speicher. Die Vorbereitung der kleinen Divisortabelle kostet
$$\sum_{d=1}^{M}\left\lfloor\frac{M}{d}\right\rfloor=O(M\log M)$$
Zeit und \(O(M)\) Speicher.
Im Großwertbereich gibt es nur \(O(K)\) verschiedene Quotientengruppen, und ein Blockverfahren für \(D(v)\) kostet \(O(\sqrt v)\). Summiert über \(t\le K\) ergibt das
$$O\!\left(\sum_{t=1}^{K}\sqrt{\frac{N}{t^2}}\right)=O(\sqrt N\log K).$$
Damit liegt die Gesamtlaufzeit im Worst Case bei \(O(\sqrt N\log N)\), während der Speicherverbrauch von den Phi-Arrays dominiert wird und \(O(\sqrt N)\) bleibt. Parallelisierung verbessert die reale Laufzeit des großen gruppierten Teils, ändert aber nicht die asymptotische Schranke.
Bu problemdeki bölen-toplamı, sıralı çarpan çiftleri toplamı olarak yeniden yazılabilir:
$$S(N)=\sum_{n=1}^{N}\sum_{d\mid n}\gcd\!\left(d,\frac{n}{d}\right)=\sum_{ab\le N}\gcd(a,b).$$
Yani amaç, çarpımı en fazla \(N\) olan tüm sıralı \((a,b)\) çiftleri için \(\gcd(a,b)\) değerlerini toplamaktır. \(N=10^{15}\) düzeyinde doğrudan tarama imkansız olduğundan çözüm, toplamı Euler phi fonksiyonu ve bölen-toplam fonksiyonu cinsinden yeniden düzenler.
Eşdeğer biçim olarak
$$S(N)=\sum_{ab\le N}\gcd(a,b)$$
ifadesiyle çalışacağız. Verimli formül, önce ortak bölen katkısını ayırıp sonra tekrar eden taban bölüm değerlerini bloklar halinde toplamaktan gelir.
Her \(n\) ve her \(d\mid n\) için \((d,n/d)\) sıralı çarpan çifti vardır ve çarpımları tam olarak \(n\)'dir. Dolayısıyla tüm \(n\le N\) ve tüm bölenler üzerinden yapılan toplam, \(ab\le N\) koşulunu sağlayan tüm sıralı \((a,b)\) çiftleri üzerinden yapılan toplamla bire bir aynıdır.
Bu dönüşüm önemlidir; çünkü orijinal bölen-toplam yapısı yerine artık temiz bir \(ab\le N\) çarpım kısıtı elde ederiz.
Klasik özdeşlik şudur:
$$\sum_{t\mid m}\varphi(t)=m.$$
Bunu \(m=\gcd(a,b)\) için uygulayınca
$$\gcd(a,b)=\sum_{t\mid \gcd(a,b)}\varphi(t)$$
elde edilir. Bunu toplamın içine koyup toplam sırasını değiştirirsek
$$S(N)=\sum_{ab\le N}\sum_{t\mid a,\ t\mid b}\varphi(t)$$
olur. Şimdi \(a=tx\) ve \(b=ty\) yazalım. Yeni koşul
$$t^2xy\le N$$
şeklindedir. Böylece
$$S(N)=\sum_{t\le \sqrt N}\varphi(t)\,D\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right)$$
sonucuna ulaşırız; burada \(D(m)\), \(xy\le m\) koşulunu sağlayan sıralı \((x,y)\) çiftlerinin sayısıdır.
Bu çift sayma fonksiyonu
$$D(m)=\#\{(x,y)\in\mathbb{Z}_{>0}^2:xy\le m\}$$
şeklindedir. Sabit bir \(x\) için geçerli \(y\) sayısı tam olarak \(\left\lfloor m/x\right\rfloor\) olduğundan
$$D(m)=\sum_{x=1}^{m}\left\lfloor\frac{m}{x}\right\rfloor$$
olur. Aynı şey, her \(n\le m\) sayısının \(xy=n\) olacak her bölen çifti için bir kez sayılması biçiminde de görülebilir; dolayısıyla
$$D(m)=\sum_{n=1}^{m}\tau(n),$$
burada \(\tau(n)\) bölen sayısı fonksiyonudur. Artık tüm problem, phi değerleriyle bölen-toplam değerlerinin ağırlıklı toplamına indirgenmiştir.
\(m\)'ye kadar doğrudan döngü kurmak büyük değerlerde hâlâ çok yavaştır. Temel gözlem şudur: bölüm
$$q=\left\lfloor\frac{m}{\ell}\right\rfloor$$
bazı \([\ell,r]\) aralıklarında sabit kalır ve bu aralığın sonu
$$r=\left\lfloor\frac{m}{q}\right\rfloor$$
şeklinde bulunur. O halde tek bir blok katkısı
$$q\,(r-\ell+1)$$
olur. Böylece
$$D(m)=\sum_{\text{bloklar }[\ell,r]}\left\lfloor\frac{m}{\ell}\right\rfloor(r-\ell+1)$$
şeklinde hesaplanır ve blok sayısı sadece \(O(\sqrt m)\) olur. Büyük \(D(m)\) değerleri için kullanılan hızlandırma tam olarak budur.
Şunu tanımlayalım:
$$K=\left\lfloor N^{1/3}\right\rfloor,\qquad M=\left\lfloor\frac{N}{K^2}\right\rfloor.$$
Böylece dış toplam doğal olarak iki bölgeye ayrılır.
\(t>K\) için \(\left\lfloor N/t^2\right\rfloor\le M\) olur; yani gereken tüm \(D(m)\) değerleri küçüktür. Uygulama \(\tau(1),\dots,\tau(M)\) değerlerini önceden hesaplar, prefix toplam alır ve küçük \(D(m)\) değerlerini tablo bakışıyla elde eder.
\(t\le K\) için \(\left\lfloor N/t^2\right\rfloor\) büyüktür; fakat ardışık birçok \(t\) aynı değeri verir. Eğer
$$v=\left\lfloor\frac{N}{t^2}\right\rfloor$$
ise, aynı bölüm değerini veren en büyük indis
$$r=\left\lfloor\sqrt{\frac{N}{v}}\right\rfloor$$
olur. O zaman tüm aralığın katkısı
$$D(v)\sum_{u=t}^{r}\varphi(u)$$
şeklindedir. \(\varphi\) için prefix toplamlar bu aralık toplamını \(O(1)\) hale getirir.
Burada \(\lfloor\sqrt{10}\rfloor=3\) olduğundan dönüşmüş formül
$$S(10)=\sum_{t=1}^{3}\varphi(t)\,D\!\left(\left\lfloor\frac{10}{t^2}\right\rfloor\right)$$
olur. Gerekli phi değerleri
$$\varphi(1)=1,\qquad \varphi(2)=1,\qquad \varphi(3)=2$$
ve gerekli bölen-toplam değerleri
$$D(10)=\sum_{x=1}^{10}\left\lfloor\frac{10}{x}\right\rfloor=27,\qquad D(2)=3,\qquad D(1)=1$$
şeklindedir. Dolayısıyla
$$S(10)=1\cdot 27+1\cdot 3+2\cdot 1=32,$$
ki bu, uygulamanın kullandığı kontrol değeriyle aynıdır.
C++, Python ve Java uygulamaları aynı sayı kuramsal formülü izler. Önce tam sayı karekök sınırı \(\lfloor\sqrt N\rfloor\) ve küpkök ayırma noktası \(\lfloor N^{1/3}\rfloor\) hesaplanır. Sonra \(\lfloor\sqrt N\rfloor\)'ye kadar phi eleği kurulur ve prefix toplamlar saklanır; böylece \(\sum_{u=l}^{r}\varphi(u)\) anında bulunabilir.
Küçük argüman bölgesi için uygulama, \(M=\lfloor N/K^2\rfloor\)'ye kadar bölen sayılarını her bölenin katlarını gezerek önceden hesaplar ve bunları prefix toplamlarına dönüştürerek \(D(1),D(2),\dots,D(M)\) değerlerini tablo üzerinden alır. Büyük argüman bölgesi için \(\left\lfloor N/t^2\right\rfloor\) değerleri gruplanır, her büyük \(D(v)\) eşit bölüm bloklarıyla hesaplanır ve karşılık gelen phi aralık toplamıyla çarpılır.
C++ ve Java uygulamaları bu bağımsız büyük-grup katkılarını paralel toplar; ardından küçük kuyruk kısmı eklenir. Python uygulaması ise ince bir köprüdür: gerekirse derlenmiş çözücüyü oluşturur, onu çalıştırır ve yazdırılan sonucu döndürür. Derlenmiş sürüm ayrıca tam hedefi hesaplamadan önce \(S(10)=32\) ve \(S(1000)=12776\) denetimlerini yapar.
\(L=\lfloor\sqrt N\rfloor\), \(K=\lfloor N^{1/3}\rfloor\) ve \(M=\lfloor N/K^2\rfloor\) olsun. Phi eleği \(O(L\log\log L)\) zaman ve \(O(L)\) bellek gerektirir. Küçük bölen tablosunun hazırlanma maliyeti
$$\sum_{d=1}^{M}\left\lfloor\frac{M}{d}\right\rfloor=O(M\log M)$$
zamandır ve \(O(M)\) bellek kullanır.
Büyük değer bölgesinde sadece \(O(K)\) farklı bölüm grubu vardır ve \(D(v)\) için blok yöntemi \(O(\sqrt v)\) sürer. Bu maliyetler \(t\le K\) üzerinde toplandığında
$$O\!\left(\sum_{t=1}^{K}\sqrt{\frac{N}{t^2}}\right)=O(\sqrt N\log K)$$
elde edilir. Böylece toplam süre en kötü durumda \(O(\sqrt N\log N)\), bellek kullanımı ise phi dizilerinin baskın olması nedeniyle \(O(\sqrt N)\) olur. Paralellik, büyük gruplu bölümün gerçek çalışma süresini düşürür ama asimptotik sınırı değiştirmez.
La cantidad basada en divisores puede reindexarse como una suma sobre pares ordenados de factores:
$$S(N)=\sum_{n=1}^{N}\sum_{d\mid n}\gcd\!\left(d,\frac{n}{d}\right)=\sum_{ab\le N}\gcd(a,b).$$
Así, el objetivo es sumar \(\gcd(a,b)\) sobre todos los pares ordenados \((a,b)\) cuyo producto no supera \(N\). Para \(N=10^{15}\) no se puede enumerar directamente, de modo que la solución transforma la suma usando la función phi de Euler y la función sumatoria de divisores.
Trabajamos con la forma equivalente
$$S(N)=\sum_{ab\le N}\gcd(a,b).$$
La fórmula eficiente aparece al separar la contribución del divisor común y luego agrupar los mismos valores enteros de piso en bloques.
Para cada \(n\) y cada divisor \(d\mid n\), el par \((d,n/d)\) es un par ordenado de factores con producto exactamente igual a \(n\). Por tanto, sumar sobre todos los \(n\le N\) y sobre todos sus divisores es exactamente lo mismo que sumar sobre todos los pares ordenados \((a,b)\) con \(ab\le N\).
Esta reescritura es crucial porque la restricción \(ab\le N\) es mucho más manejable que la suma original sobre divisores.
Usamos la identidad clásica
$$\sum_{t\mid m}\varphi(t)=m.$$
Aplicada a \(m=\gcd(a,b)\), da
$$\gcd(a,b)=\sum_{t\mid \gcd(a,b)}\varphi(t).$$
Al sustituir esto en la suma e intercambiar el orden, obtenemos
$$S(N)=\sum_{ab\le N}\sum_{t\mid a,\ t\mid b}\varphi(t).$$
Ahora escribimos \(a=tx\) y \(b=ty\). La condición pasa a ser
$$t^2xy\le N.$$
Por consiguiente,
$$S(N)=\sum_{t\le \sqrt N}\varphi(t)\,D\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right),$$
donde \(D(m)\) cuenta los pares ordenados \((x,y)\) con \(xy\le m\).
La función de conteo de pares es
$$D(m)=\#\{(x,y)\in\mathbb{Z}_{>0}^2:xy\le m\}.$$
Para un \(x\) fijo, hay exactamente \(\left\lfloor m/x\right\rfloor\) valores posibles de \(y\). Por tanto,
$$D(m)=\sum_{x=1}^{m}\left\lfloor\frac{m}{x}\right\rfloor.$$
De manera equivalente, cada entero \(n\le m\) contribuye una vez por cada par de divisores \((x,y)\) con \(xy=n\), así que también
$$D(m)=\sum_{n=1}^{m}\tau(n),$$
donde \(\tau(n)\) es la función que cuenta divisores. Toda la tarea queda reducida a una sola suma ponderada de totientes por valores de la función sumatoria de divisores.
Un bucle ingenuo hasta \(m\) sigue siendo demasiado lento para \(m\) grande. La observación central es que el cociente
$$q=\left\lfloor\frac{m}{\ell}\right\rfloor$$
permanece constante en un intervalo \([\ell,r]\), donde
$$r=\left\lfloor\frac{m}{q}\right\rfloor.$$
Entonces un bloque completo aporta
$$q\,(r-\ell+1).$$
Así se obtiene
$$D(m)=\sum_{\text{bloques }[\ell,r]}\left\lfloor\frac{m}{\ell}\right\rfloor(r-\ell+1),$$
y el número de bloques baja a \(O(\sqrt m)\). Esta es la aceleración de tipo hipérbola que usa la implementación para los valores grandes de \(D(m)\).
Definimos
$$K=\left\lfloor N^{1/3}\right\rfloor,\qquad M=\left\lfloor\frac{N}{K^2}\right\rfloor.$$
Entonces aparecen dos regímenes naturales.
Si \(t>K\), se cumple \(\left\lfloor N/t^2\right\rfloor\le M\), así que todos los valores necesarios de \(D(m)\) son pequeños. La implementación precalcula \(\tau(1),\dots,\tau(M)\), forma sus prefijos y obtiene cada pequeño \(D(m)\) por consulta directa.
Si \(t\le K\), el argumento \(\left\lfloor N/t^2\right\rfloor\) es grande, pero muchos valores consecutivos de \(t\) producen el mismo cociente. Si
$$v=\left\lfloor\frac{N}{t^2}\right\rfloor,$$
entonces el mayor índice con ese mismo valor es
$$r=\left\lfloor\sqrt{\frac{N}{v}}\right\rfloor.$$
La contribución del intervalo completo es
$$D(v)\sum_{u=t}^{r}\varphi(u).$$
Los prefijos de \(\varphi\) convierten esa suma de intervalo en una consulta \(O(1)\), de modo que el trabajo costoso queda concentrado en los distintos valores grandes de \(D(v)\).
Aquí \(\lfloor\sqrt{10}\rfloor=3\), por lo que la fórmula transformada es
$$S(10)=\sum_{t=1}^{3}\varphi(t)\,D\!\left(\left\lfloor\frac{10}{t^2}\right\rfloor\right).$$
Los totientes necesarios son
$$\varphi(1)=1,\qquad \varphi(2)=1,\qquad \varphi(3)=2.$$
Y los valores necesarios de la función sumatoria son
$$D(10)=\sum_{x=1}^{10}\left\lfloor\frac{10}{x}\right\rfloor=27,\qquad D(2)=3,\qquad D(1)=1.$$
Entonces
$$S(10)=1\cdot 27+1\cdot 3+2\cdot 1=32,$$
que coincide con el valor de comprobación usado por la implementación.
Las implementaciones en C++, Python y Java siguen la misma fórmula de teoría de números. Primero calculan el límite entero \(\lfloor\sqrt N\rfloor\) y el punto de corte \(\lfloor N^{1/3}\rfloor\). Después construyen una criba de totientes hasta \(\lfloor\sqrt N\rfloor\) y guardan prefijos para recuperar \(\sum_{u=l}^{r}\varphi(u)\) de inmediato.
En la región de argumentos pequeños, la implementación precalcula el número de divisores hasta \(M=\lfloor N/K^2\rfloor\) recorriendo los múltiplos de cada divisor, y convierte esos conteos en prefijos para disponer de \(D(1),D(2),\dots,D(M)\). En la región de argumentos grandes, agrupa los valores iguales de \(\left\lfloor N/t^2\right\rfloor\), calcula cada \(D(v)\) grande por bloques de cociente constante y lo multiplica por la suma adecuada de totientes en el intervalo correspondiente.
Las implementaciones en C++ y Java paralelizan esas contribuciones independientes del tramo grande antes de añadir la cola pequeña. La implementación en Python es un envoltorio ligero: construye el solucionador compilado cuando hace falta, lo ejecuta y devuelve la salida. La versión compilada también comprueba previamente los dos puntos de control \(S(10)=32\) y \(S(1000)=12776\).
Sean \(L=\lfloor\sqrt N\rfloor\), \(K=\lfloor N^{1/3}\rfloor\) y \(M=\lfloor N/K^2\rfloor\). La criba de totientes cuesta \(O(L\log\log L)\) tiempo y \(O(L)\) memoria. La preparación de la tabla pequeña de divisores cuesta
$$\sum_{d=1}^{M}\left\lfloor\frac{M}{d}\right\rfloor=O(M\log M)$$
tiempo y \(O(M)\) memoria.
En la región de valores grandes solo hay \(O(K)\) grupos distintos, y una llamada al método por bloques para \(D(v)\) cuesta \(O(\sqrt v)\). Sumando sobre \(t\le K\) se obtiene
$$O\!\left(\sum_{t=1}^{K}\sqrt{\frac{N}{t^2}}\right)=O(\sqrt N\log K).$$
Así, el tiempo total en el peor caso es \(O(\sqrt N\log N)\), mientras que la memoria está dominada por los arreglos de totientes y permanece en \(O(\sqrt N)\). El paralelismo mejora el tiempo real del tramo agrupado grande, pero no modifica la cota asintótica.
这道题中的“约数 gcd 总和”可以改写成一个有序因子对求和：
$$S(N)=\sum_{n=1}^{N}\sum_{d\mid n}\gcd\!\left(d,\frac{n}{d}\right)=\sum_{ab\le N}\gcd(a,b).$$
也就是说，我们真正需要计算的是所有满足 \(ab\le N\) 的有序整数对 \((a,b)\) 的 \(\gcd(a,b)\) 之和。当 \(N=10^{15}\) 时，直接枚举这些二元组完全不可行，所以必须把问题改写成可以批量处理的数论求和。
以下都从等价形式
$$S(N)=\sum_{ab\le N}\gcd(a,b)$$
出发。高效算法的核心，是先把 gcd 拆成更容易求和的函数，再把重复出现的整除商整体打包处理。
对每个 \(n\) 以及每个约数 \(d\mid n\)，都对应一个有序因子对 \((d,n/d)\)，并且它们的乘积正好等于 \(n\)。因此，对所有 \(n\le N\) 及其全部约数求和，和对所有满足 \(ab\le N\) 的有序对 \((a,b)\) 求和，是完全同一件事。
这一步非常关键，因为原题中的“枚举每个数的约数”不容易直接优化，而条件 \(ab\le N\) 却非常适合做整除分块和前缀求和。
使用经典恒等式
$$\sum_{t\mid m}\varphi(t)=m.$$
把它应用到 \(m=\gcd(a,b)\) 上，就得到
$$\gcd(a,b)=\sum_{t\mid \gcd(a,b)}\varphi(t).$$
代回总和并交换求和顺序，有
$$S(N)=\sum_{ab\le N}\sum_{t\mid a,\ t\mid b}\varphi(t).$$
现在令 \(a=tx,\ b=ty\)，那么约束变成
$$t^2xy\le N.$$
于是得到核心公式
$$S(N)=\sum_{t\le \sqrt N}\varphi(t)\,D\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right),$$
其中 \(D(m)\) 表示满足 \(xy\le m\) 的有序正整数对 \((x,y)\) 的个数。
按定义，
$$D(m)=\#\{(x,y)\in\mathbb{Z}_{>0}^2:xy\le m\}.$$
如果固定 \(x\)，那么可以选择的 \(y\) 恰好有 \(\left\lfloor m/x\right\rfloor\) 个，因此
$$D(m)=\sum_{x=1}^{m}\left\lfloor\frac{m}{x}\right\rfloor.$$
另一方面，每个 \(n\le m\) 会因为每一组满足 \(xy=n\) 的约数对 \((x,y)\) 被计算一次，所以也可以写成
$$D(m)=\sum_{n=1}^{m}\tau(n),$$
其中 \(\tau(n)\) 是约数个数函数。这样一来，原问题就变成了“欧拉函数前缀信息”和“约数求和函数”之间的加权组合。
如果对每个 \(m\) 都直接做从 \(1\) 到 \(m\) 的循环，仍然会太慢。关键观察是，商
$$q=\left\lfloor\frac{m}{\ell}\right\rfloor$$
会在一个连续区间 \([\ell,r]\) 上保持不变，其中
$$r=\left\lfloor\frac{m}{q}\right\rfloor.$$
因此，这整个区间可以一次性贡献
$$q\,(r-\ell+1).$$
于是
$$D(m)=\sum_{\text{区间 }[\ell,r]}\left\lfloor\frac{m}{\ell}\right\rfloor(r-\ell+1).$$
这样的区间个数只有 \(O(\sqrt m)\)，这就是代码对大 \(D(m)\) 使用的双曲线式加速方法。
令
$$K=\left\lfloor N^{1/3}\right\rfloor,\qquad M=\left\lfloor\frac{N}{K^2}\right\rfloor.$$
外层求和会自然分成两段。
当 \(t>K\) 时，有 \(\left\lfloor N/t^2\right\rfloor\le M\)。也就是说，此时需要的 \(D(m)\) 全都是小值。实现会先预处理 \(\tau(1),\dots,\tau(M)\)，再做前缀和，从而一次性得到全部小范围的 \(D(m)\)。
当 \(t\le K\) 时，\(\left\lfloor N/t^2\right\rfloor\) 仍然很大，但连续很多个 \(t\) 会产生相同的值。若
$$v=\left\lfloor\frac{N}{t^2}\right\rfloor,$$
那么使得这个值保持不变的最大下标为
$$r=\left\lfloor\sqrt{\frac{N}{v}}\right\rfloor.$$
这一整段的贡献就可以写成
$$D(v)\sum_{u=t}^{r}\varphi(u).$$
而 \(\varphi\) 的前缀和可以把这个区间和变成 \(O(1)\) 查询，因此真正昂贵的部分只剩下那些彼此不同的大 \(D(v)\)。
此时 \(\lfloor\sqrt{10}\rfloor=3\)，于是
$$S(10)=\sum_{t=1}^{3}\varphi(t)\,D\!\left(\left\lfloor\frac{10}{t^2}\right\rfloor\right).$$
需要的欧拉函数值为
$$\varphi(1)=1,\qquad \varphi(2)=1,\qquad \varphi(3)=2.$$
需要的 \(D(m)\) 为
$$D(10)=\sum_{x=1}^{10}\left\lfloor\frac{10}{x}\right\rfloor=27,\qquad D(2)=3,\qquad D(1)=1.$$
所以
$$S(10)=1\cdot 27+1\cdot 3+2\cdot 1=32,$$
这正好与实现中使用的校验值一致。
C++、Python 和 Java 三个实现都遵循同一个数论公式。它们先求出整数范围内的 \(\lfloor\sqrt N\rfloor\) 和分界点 \(\lfloor N^{1/3}\rfloor\)，然后在线性数组上筛出直到 \(\lfloor\sqrt N\rfloor\) 为止的欧拉函数值，并构造前缀和，这样任意区间的 \(\sum_{u=l}^{r}\varphi(u)\) 都能立即得到。
对于小参数区间，实现先通过“枚举每个除数的倍数”的方式预处理到 \(M=\lfloor N/K^2\rfloor\) 为止的约数个数，再做前缀和，从而拿到全部 \(D(1),D(2),\dots,D(M)\)。对于大参数区间，实现把相同的 \(\left\lfloor N/t^2\right\rfloor\) 分成一组，对每个不同的大值 \(D(v)\) 用整除分块单独计算，再乘上对应区间的 totient 总和。
C++ 和 Java 实现会把这些彼此独立的大区间贡献并行累加，最后再加上小尾部。Python 实现则是一个很薄的包装层：它在需要时构建编译后的求解器，执行它，并解析输出结果。编译后的版本在计算最终目标之前，还会先验证 \(S(10)=32\) 和 \(S(1000)=12776\) 这两个检查点。
记 \(L=\lfloor\sqrt N\rfloor\)、\(K=\lfloor N^{1/3}\rfloor\)、\(M=\lfloor N/K^2\rfloor\)。欧拉函数筛法的时间复杂度为 \(O(L\log\log L)\)，空间复杂度为 \(O(L)\)。小范围约数表的预处理时间为
$$\sum_{d=1}^{M}\left\lfloor\frac{M}{d}\right\rfloor=O(M\log M),$$
空间复杂度为 \(O(M)\)。
在大值区域中，不同的商分组只有 \(O(K)\) 个，而一次整除分块求 \(D(v)\) 的代价是 \(O(\sqrt v)\)。把这些代价对 \(t\le K\) 累加，可得
$$O\!\left(\sum_{t=1}^{K}\sqrt{\frac{N}{t^2}}\right)=O(\sqrt N\log K).$$
因此，总时间复杂度在最坏情况下为 \(O(\sqrt N\log N)\)，空间则主要由 totient 数组主导，保持在 \(O(\sqrt N)\)。并行化只能改善大分组部分的实际运行时间，不改变渐近复杂度。
Искомая сумма по делителям может быть переиндексирована как сумма по упорядоченным парам множителей:
$$S(N)=\sum_{n=1}^{N}\sum_{d\mid n}\gcd\!\left(d,\frac{n}{d}\right)=\sum_{ab\le N}\gcd(a,b).$$
Иными словами, нужно просуммировать \(\gcd(a,b)\) по всем упорядоченным парам \((a,b)\), для которых \(ab\le N\). При \(N=10^{15}\) прямой перебор невозможен, поэтому решение переводит задачу в форму с функцией Эйлера и сумматорной функцией делителей.
Будем исходить из эквивалентной записи
$$S(N)=\sum_{ab\le N}\gcd(a,b).$$
Быстрая формула получается после выделения вклада общего делителя и последующей группировки одинаковых целых частных.
Для каждого \(n\) и каждого делителя \(d\mid n\) существует упорядоченная пара множителей \((d,n/d)\) с произведением \(n\). Поэтому суммирование по всем \(n\le N\) и по всем делителям эквивалентно суммированию по всем упорядоченным парам \((a,b)\), удовлетворяющим условию \(ab\le N\).
Это важная замена: исходная сумма по делителям превращается в задачу с простой границей на произведение.
Используем классическое тождество
$$\sum_{t\mid m}\varphi(t)=m.$$
Применяя его к \(m=\gcd(a,b)\), получаем
$$\gcd(a,b)=\sum_{t\mid \gcd(a,b)}\varphi(t).$$
Подставим это в сумму и поменяем порядок суммирования:
$$S(N)=\sum_{ab\le N}\sum_{t\mid a,\ t\mid b}\varphi(t).$$
Теперь положим \(a=tx\) и \(b=ty\). Тогда ограничение принимает вид
$$t^2xy\le N.$$
Отсюда следует
$$S(N)=\sum_{t\le \sqrt N}\varphi(t)\,D\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right),$$
где \(D(m)\) считает упорядоченные пары \((x,y)\) с условием \(xy\le m\).
Функция подсчета пар равна
$$D(m)=\#\{(x,y)\in\mathbb{Z}_{>0}^2:xy\le m\}.$$
При фиксированном \(x\) допустимых значений \(y\) ровно \(\left\lfloor m/x\right\rfloor\), поэтому
$$D(m)=\sum_{x=1}^{m}\left\lfloor\frac{m}{x}\right\rfloor.$$
Эквивалентно, каждое число \(n\le m\) учитывается один раз для каждой пары делителей \((x,y)\) с \(xy=n\), значит
$$D(m)=\sum_{n=1}^{m}\tau(n),$$
где \(\tau(n)\) обозначает число положительных делителей числа \(n\). Тем самым исходная задача сведена к одной взвешенной сумме значений \(\varphi\) и \(D\).
Прямой цикл до \(m\) остается слишком медленным для больших значений. Главная идея состоит в том, что частное
$$q=\left\lfloor\frac{m}{\ell}\right\rfloor$$
остается постоянным на некотором интервале \([\ell,r]\), где
$$r=\left\lfloor\frac{m}{q}\right\rfloor.$$
Целый такой блок дает вклад
$$q\,(r-\ell+1).$$
Следовательно,
$$D(m)=\sum_{\text{блоки }[\ell,r]}\left\lfloor\frac{m}{\ell}\right\rfloor(r-\ell+1),$$
а число блоков составляет лишь \(O(\sqrt m)\). Именно эта гиперболическая техника используется для больших значений \(D(m)\).
Положим
$$K=\left\lfloor N^{1/3}\right\rfloor,\qquad M=\left\lfloor\frac{N}{K^2}\right\rfloor.$$
Тогда внешняя сумма естественным образом разбивается на две области.
Если \(t>K\), то \(\left\lfloor N/t^2\right\rfloor\le M\), поэтому все нужные значения \(D(m)\) малы. Реализация заранее вычисляет \(\tau(1),\dots,\tau(M)\), строит префиксные суммы и получает каждое малое \(D(m)\) простым обращением к таблице.
Если \(t\le K\), то аргумент \(\left\lfloor N/t^2\right\rfloor\) велик, но многие последовательные \(t\) дают одно и то же значение. Если
$$v=\left\lfloor\frac{N}{t^2}\right\rfloor,$$
то наибольший индекс с тем же частным равен
$$r=\left\lfloor\sqrt{\frac{N}{v}}\right\rfloor.$$
Вклад всего интервала равен
$$D(v)\sum_{u=t}^{r}\varphi(u).$$
Префиксные суммы по \(\varphi\) делают эту интервальную сумму \(O(1)\)-запросом, так что дорогими остаются только различные большие вычисления \(D(v)\).
Здесь \(\lfloor\sqrt{10}\rfloor=3\), поэтому преобразованная формула выглядит так:
$$S(10)=\sum_{t=1}^{3}\varphi(t)\,D\!\left(\left\lfloor\frac{10}{t^2}\right\rfloor\right).$$
Нужные значения функции Эйлера:
$$\varphi(1)=1,\qquad \varphi(2)=1,\qquad \varphi(3)=2.$$
Нужные значения \(D(m)\):
$$D(10)=\sum_{x=1}^{10}\left\lfloor\frac{10}{x}\right\rfloor=27,\qquad D(2)=3,\qquad D(1)=1.$$
Следовательно,
$$S(10)=1\cdot 27+1\cdot 3+2\cdot 1=32,$$
что совпадает с контрольным значением в реализации.
Реализации на C++, Python и Java следуют одной и той же формуле из теории чисел. Сначала они вычисляют целочисленный предел \(\lfloor\sqrt N\rfloor\) и точку разбиения \(\lfloor N^{1/3}\rfloor\). Затем строится решето значений \(\varphi\) до \(\lfloor\sqrt N\rfloor\), а также их префиксные суммы, чтобы \(\sum_{u=l}^{r}\varphi(u)\) можно было получать мгновенно.
Для области малых аргументов реализация заранее вычисляет числа делителей до \(M=\lfloor N/K^2\rfloor\), проходя по кратным каждого делителя, а затем берет префиксные суммы, получая значения \(D(1),D(2),\dots,D(M)\). Для области больших аргументов одинаковые значения \(\left\lfloor N/t^2\right\rfloor\) группируются, каждый большой \(D(v)\) вычисляется блоками одинаковых частных и умножается на соответствующую сумму \(\varphi\) по интервалу.
Реализации на C++ и Java распараллеливают эти независимые крупные группы, а затем добавляют малый хвост. Реализация на Python представляет собой тонкую обертку: при необходимости она собирает скомпилированный решатель, запускает его и возвращает напечатанный результат. Скомпилированная версия также проверяет контрольные значения \(S(10)=32\) и \(S(1000)=12776\) перед вычислением полного ответа.
Пусть \(L=\lfloor\sqrt N\rfloor\), \(K=\lfloor N^{1/3}\rfloor\) и \(M=\lfloor N/K^2\rfloor\). Решето для функции Эйлера работает за \(O(L\log\log L)\) времени и требует \(O(L)\) памяти. Подготовка малой таблицы чисел делителей стоит
$$\sum_{d=1}^{M}\left\lfloor\frac{M}{d}\right\rfloor=O(M\log M)$$
времени и \(O(M)\) памяти.
В области больших значений существует лишь \(O(K)\) различных групп частных, а один блочный вызов для \(D(v)\) стоит \(O(\sqrt v)\). Суммирование этих затрат по \(t\le K\) дает
$$O\!\left(\sum_{t=1}^{K}\sqrt{\frac{N}{t^2}}\right)=O(\sqrt N\log K).$$
Следовательно, общая трудоемкость в худшем случае равна \(O(\sqrt N\log N)\), а память определяется в основном массивами \(\varphi\) и остается \(O(\sqrt N)\). Параллельное выполнение уменьшает реальное время работы на крупном сгруппированном участке, но не меняет асимптотику.
يمكن إعادة كتابة المجموع المعتمد على القواسم في هذه المسألة على صورة مجموع على أزواج مرتبة من العوامل:
$$S(N)=\sum_{n=1}^{N}\sum_{d\mid n}\gcd\!\left(d,\frac{n}{d}\right)=\sum_{ab\le N}\gcd(a,b).$$
أي إن المطلوب فعليًا هو جمع \(\gcd(a,b)\) على جميع الأزواج المرتبة \((a,b)\) التي تحقق \(ab\le N\). وعندما يكون \(N=10^{15}\) فإن التعداد المباشر مستحيل عمليًا، لذا يحول الحل هذا المجموع إلى صيغة تعتمد على دالة أويلر وعلى دالة مجموع القواسم.
سننطلق من الصيغة المكافئة
$$S(N)=\sum_{ab\le N}\gcd(a,b).$$
تأتي السرعة من فصل مساهمة القاسم المشترك، ثم تجميع القيم المتكررة لعمليات القسمة الصحيحة في كتل كبيرة.
لكل عدد \(n\) ولكل قاسم \(d\mid n\)، يوجد زوج مرتب من العوامل هو \((d,n/d)\) وحاصل ضربه يساوي \(n\). لذلك فإن الجمع على جميع \(n\le N\) وعلى جميع قواسمها يساوي تمامًا الجمع على جميع الأزواج المرتبة \((a,b)\) التي تحقق \(ab\le N\).
هذه هي أول نقلة مهمة، لأنها تستبدل صياغة القواسم الأصلية بقيد بسيط على حاصل الضرب.
نستخدم المتطابقة الكلاسيكية
$$\sum_{t\mid m}\varphi(t)=m.$$
وعند تطبيقها على \(m=\gcd(a,b)\) نحصل على
$$\gcd(a,b)=\sum_{t\mid \gcd(a,b)}\varphi(t).$$
بالتعويض ثم تبديل ترتيب المجاميع نحصل على
$$S(N)=\sum_{ab\le N}\sum_{t\mid a,\ t\mid b}\varphi(t).$$
الآن نكتب \(a=tx\) و\(b=ty\)، فتتحول القيود إلى
$$t^2xy\le N.$$
ومن ثم
$$S(N)=\sum_{t\le \sqrt N}\varphi(t)\,D\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right),$$
حيث إن \(D(m)\) يعد الأزواج المرتبة \((x,y)\) التي تحقق \(xy\le m\).
دالة عد الأزواج هي
$$D(m)=\#\{(x,y)\in\mathbb{Z}_{>0}^2:xy\le m\}.$$
إذا ثبتنا \(x\)، فإن عدد القيم الممكنة لـ \(y\) يساوي بالضبط \(\left\lfloor m/x\right\rfloor\). ولذلك
$$D(m)=\sum_{x=1}^{m}\left\lfloor\frac{m}{x}\right\rfloor.$$
وبصورة مكافئة، فإن كل عدد \(n\le m\) يُحسب مرة واحدة لكل زوج قواسم \((x,y)\) يحقق \(xy=n\)، لذا يمكن كتابة
$$D(m)=\sum_{n=1}^{m}\tau(n),$$
حيث \(\tau(n)\) هي دالة عدد القواسم. وهكذا تتحول المسألة كلها إلى مجموع موزون بين قيم \(\varphi\) وقيم \(D\).
المرور المباشر حتى \(m\) يظل بطيئًا جدًا عندما يكون \(m\) كبيرًا. الملاحظة الأساسية هي أن الخارج
$$q=\left\lfloor\frac{m}{\ell}\right\rfloor$$
يبقى ثابتًا على فترة كاملة \([\ell,r]\)، حيث
$$r=\left\lfloor\frac{m}{q}\right\rfloor.$$
وعليه فإن مساهمة الكتلة الواحدة هي
$$q\,(r-\ell+1).$$
ومن ثم نحصل على
$$D(m)=\sum_{\text{blocks }[\ell,r]}\left\lfloor\frac{m}{\ell}\right\rfloor(r-\ell+1),$$
ويصبح عدد هذه الكتل فقط \(O(\sqrt m)\). هذه هي السرعة الشبيهة بطريقة القطع الزائد التي يستخدمها التنفيذ مع القيم الكبيرة.
لنأخذ
$$K=\left\lfloor N^{1/3}\right\rfloor,\qquad M=\left\lfloor\frac{N}{K^2}\right\rfloor.$$
عندها ينقسم المجموع الخارجي طبيعيًا إلى جزأين.
إذا كان \(t>K\)، فإن \(\left\lfloor N/t^2\right\rfloor\le M\)، أي إن جميع قيم \(D(m)\) المطلوبة صغيرة. ولذلك يسبق التنفيذ حساب \(\tau(1),\dots,\tau(M)\)، ثم يبني مجاميع تراكمية للحصول على كل قيمة صغيرة من \(D(m)\) مباشرة من جدول.
أما إذا كان \(t\le K\)، فإن \(\left\lfloor N/t^2\right\rfloor\) كبير، لكن كثيرًا من القيم المتتالية لـ \(t\) تعطي النتيجة نفسها. فإذا كان
$$v=\left\lfloor\frac{N}{t^2}\right\rfloor,$$
فإن أكبر فهرس يحقق القيمة نفسها هو
$$r=\left\lfloor\sqrt{\frac{N}{v}}\right\rfloor.$$
وبذلك تصبح مساهمة الفترة كلها
$$D(v)\sum_{u=t}^{r}\varphi(u).$$
ومجاميع \(\varphi\) التراكمية تجعل مجموع الفترة هذا استعلامًا من رتبة \(O(1)\)، فلا يبقى العمل الثقيل إلا في حسابات \(D(v)\) الكبيرة المختلفة.
هنا لدينا \(\lfloor\sqrt{10}\rfloor=3\)، ومن ثم
$$S(10)=\sum_{t=1}^{3}\varphi(t)\,D\!\left(\left\lfloor\frac{10}{t^2}\right\rfloor\right).$$
وقيم \(\varphi\) المطلوبة هي
$$\varphi(1)=1,\qquad \varphi(2)=1,\qquad \varphi(3)=2.$$
أما قيم \(D(m)\) المطلوبة فهي
$$D(10)=\sum_{x=1}^{10}\left\lfloor\frac{10}{x}\right\rfloor=27,\qquad D(2)=3,\qquad D(1)=1.$$
إذًا
$$S(10)=1\cdot 27+1\cdot 3+2\cdot 1=32,$$
وهو نفس مقدار التحقق الذي يستخدمه التنفيذ.
تتبع تطبيقات C++ وPython وJava الصيغة العددية نفسها. فهي تبدأ بحساب الحد الصحيح \(\lfloor\sqrt N\rfloor\) ونقطة التقسيم \(\lfloor N^{1/3}\rfloor\). ثم تبني غربالًا لقيم \(\varphi\) حتى \(\lfloor\sqrt N\rfloor\)، وتخزن المجاميع التراكمية حتى يمكن الحصول على \(\sum_{u=l}^{r}\varphi(u)\) فورًا.
في منطقة القيم الصغيرة، يسبق التنفيذ حساب عدد القواسم حتى \(M=\lfloor N/K^2\rfloor\) عبر المرور على مضاعفات كل قاسم، ثم يحول هذه القيم إلى مجاميع تراكمية للحصول على \(D(1),D(2),\dots,D(M)\). وفي منطقة القيم الكبيرة، يجمع القيم المتساوية لـ \(\left\lfloor N/t^2\right\rfloor\)، ويحسب كل \(D(v)\) كبيرًا على شكل كتل ذات خارج ثابت، ثم يضربه في مجموع \(\varphi\) المناسب لذلك المجال.
يقوم تطبيقا C++ وJava بتوزيع هذه المساهمات الكبيرة المستقلة على مسارات متوازية قبل إضافة الذيل الصغير. أما تطبيق Python فهو طبقة خفيفة: يبني المحرك المترجم عند الحاجة، ثم يشغله ويعيد النتيجة المطبوعة. كما أن النسخة المترجمة تتحقق أولًا من نقطتي الفحص \(S(10)=32\) و\(S(1000)=12776\) قبل حساب الجواب الكامل.
لنضع \(L=\lfloor\sqrt N\rfloor\)، و\(K=\lfloor N^{1/3}\rfloor\)، و\(M=\lfloor N/K^2\rfloor\). يحتاج غربال دالة أويلر إلى زمن \(O(L\log\log L)\) وذاكرة \(O(L)\). أما تجهيز جدول القواسم الصغير فيكلف
$$\sum_{d=1}^{M}\left\lfloor\frac{M}{d}\right\rfloor=O(M\log M)$$
زمنًا و\(O(M)\) ذاكرة.
في منطقة القيم الكبيرة يوجد فقط \(O(K)\) من مجموعات القيم المختلفة، أما حساب \(D(v)\) بطريقة الكتل فيكلف \(O(\sqrt v)\). وبجمع هذه الكلفة على \(t\le K\) نحصل على
$$O\!\left(\sum_{t=1}^{K}\sqrt{\frac{N}{t^2}}\right)=O(\sqrt N\log K).$$
إذن فإن زمن التنفيذ الكلي في أسوأ الأحوال هو \(O(\sqrt N\log N)\)، بينما تظل الذاكرة من رتبة \(O(\sqrt N)\) لأن المصفوفات الخاصة بـ \(\varphi\) هي المهيمنة. ويُحسن التوازي زمن التنفيذ الفعلي في الجزء الكبير المجمع، لكنه لا يغير الرتبة التقاربية.