For each integer \(n \gt 1\), consider linear maps
$$f(x)\equiv ax+b \pmod{n},\qquad 0 \lt a \lt n,\quad 0 \le b \lt n.$$
The map is a retraction when it is idempotent on every residue class:
$$f(f(x))\equiv f(x)\pmod{n}\qquad \forall x.$$
Let \(R(n)\) be the number of such maps. Problem 447 asks for
$$F(N)=\sum_{n=2}^{N}R(n)\pmod{10^9+7},$$
with the checkpoint \(F(10^7)\equiv 638042271 \pmod{10^9+7}\) and target \(N=10^{14}\). Direct enumeration is impossible at this scale, so the implementation rewrites \(R(n)\) as an arithmetic function with a fast summatory formula.
Composing \(f(x)=ax+b\) once more gives
$$f(f(x))\equiv a^2x+ab+b \pmod{n}.$$
For this to equal \(ax+b\) for every residue class, both corrections must vanish modulo \(n\):
$$n \mid a(a-1),\qquad n \mid ab.$$
If we write \(d=\gcd(a,n)\), then \(d\mid n\) and \(\gcd(d,n/d)=1\), so \(d\) is a unitary divisor of \(n\). Conversely, each unitary divisor \(d \lt n\) determines one admissible residue class for \(a\), and once \(a\) is fixed, the number of valid \(b\) values is exactly \(d\). Therefore
$$R(n)=\sum_{\substack{d \mid n\\ \gcd(d,n/d)=1\\ d \lt n}} d=\sigma^*(n)-n,$$
where \(\sigma^*(n)\) is the sum of unitary divisors of \(n\). Since \(\sigma^*(1)=1\), the missing term at \(n=1\) contributes \(1-1=0\), so
$$F(N)=\sum_{n=1}^{N}\bigl(\sigma^*(n)-n\bigr)=\left(\sum_{n=1}^{N}\sigma^*(n)\right)-\frac{N(N+1)}{2}.$$
The unitary condition can be encoded by Möbius inversion:
$$\mathbf{1}_{\gcd(d,n/d)=1}=\sum_{t \mid \gcd(d,n/d)} \mu(t).$$
Insert this into the divisor sum:
$$\sigma^*(n)=\sum_{d \mid n} d \sum_{t \mid \gcd(d,n/d)} \mu(t).$$
If \(t \mid \gcd(d,n/d)\), then \(t^2 \mid n\), and we may write \(d=t\,m\) with \(m \mid n/t^2\). Hence
$$\sigma^*(n)=\sum_{t^2 \mid n} \mu(t)\, t \sum_{m \mid n/t^2} m=\sum_{t^2 \mid n} \mu(t)\, t\, \sigma\!\left(\frac{n}{t^2}\right).$$
This identity is the bridge from unitary divisors to the ordinary sum-of-divisors function \(\sigma\), which is much easier to sum over long intervals.
Define
$$T(x)=\sum_{m=1}^{x}\sigma(m).$$
Summing the previous identity for \(n \le N\) and swapping the order of summation gives
$$\sum_{n=1}^{N}\sigma^*(n)=\sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right).$$
Now use the classical divisor-summatory identity
$$T(x)=\sum_{m=1}^{x}\sum_{d \mid m} d=\sum_{d=1}^{x} d\left\lfloor\frac{x}{d}\right\rfloor.$$
So the original problem has been reduced to evaluating many instances of \(T(x)\) quickly.
The value \(\left\lfloor x/i \right\rfloor\) stays constant on intervals. If \(l\) is the first index of a block, define
$$q=\left\lfloor\frac{x}{l}\right\rfloor,\qquad r=\left\lfloor\frac{x}{q}\right\rfloor.$$
Then every \(i\in[l,r]\) has the same quotient \(q\), so the whole block contributes
$$q\sum_{i=l}^{r} i=q\cdot \frac{(l+r)(r-l+1)}{2}.$$
Advancing directly from \(l\) to \(r+1\) skips an entire constant-quotient range, reducing the cost of \(T(x)\) from \(O(x)\) to about \(O(\sqrt{x})\). This quotient-block decomposition is the central speedup in the implementation.
Putting everything together,
$$\boxed{F(N)\equiv \sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right)-\frac{N(N+1)}{2}\pmod{10^9+7}.}$$
For a small check, take \(N=10\). Quotient blocks give
$$T(10)=1\cdot 10+2\cdot 5+3\cdot 3+(4+5)\cdot 2+(6+7+8+9+10)\cdot 1=87.$$
Also \(T(2)=4\), \(T(1)=1\), and the nonzero Möbius values up to \(\sqrt{10}\) are \(\mu(1)=1\), \(\mu(2)=-1\), \(\mu(3)=-1\). Therefore
$$\sum_{n=1}^{10}\sigma^*(n)=1\cdot 87-2\cdot 4-3\cdot 1=76,$$
so
$$F(10)=76-\frac{10\cdot 11}{2}=21.$$
This agrees with direct evaluation of \(R(2),\dots,R(10)\), so the derivation is consistent before moving to \(10^{14}\).
The C++, Python, and Java implementations first sieve the Möbius function up to \(\lfloor\sqrt{N}\rfloor\). The outer sum then skips every \(t\) with \(\mu(t)=0\), because such terms contribute nothing. For each remaining \(t\), the implementation evaluates one inner query \(T\!\left(\lfloor N/t^2\rfloor\right)\) using quotient blocks instead of a linear loop.
After the summatory unitary-divisor part has been accumulated, the implementation subtracts the triangular number \(N(N+1)/2\) modulo \(10^9+7\). The C++ and Java implementations split the outer Möbius sum cyclically across worker threads so that expensive small \(t\) values are balanced more evenly; the Python implementation uses the same mathematics sequentially. In every language, modular subtraction is normalized back into the range \(0,\dots,10^9+6\).
The Möbius sieve up to \(\sqrt{N}\) needs \(O(\sqrt{N})\) memory and near-linear preprocessing in that range. The dominant work is
$$\sum_{t \le \sqrt{N}} O\!\left(\sqrt{\frac{N}{t^2}}\right)=\sum_{t \le \sqrt{N}} O\!\left(\frac{\sqrt{N}}{t}\right)=O(\sqrt{N}\log N).$$
So the overall method runs in \(O(\sqrt{N}\log N)\) time and \(O(\sqrt{N})\) memory. Parallel execution improves wall-clock time in the compiled implementations, but the asymptotic bound is unchanged.
Für jedes \(n \gt 1\) betrachten wir lineare Abbildungen
$$f(x)\equiv ax+b \pmod{n},\qquad 0 \lt a \lt n,\quad 0 \le b \lt n.$$
Eine solche Abbildung ist eine Retraktion, wenn sie auf jeder Restklasse idempotent ist:
$$f(f(x))\equiv f(x)\pmod{n}\qquad \forall x.$$
Sei \(R(n)\) die Anzahl dieser Abbildungen. Gesucht ist
$$F(N)=\sum_{n=2}^{N}R(n)\pmod{10^9+7},$$
mit dem Kontrollwert \(F(10^7)\equiv 638042271 \pmod{10^9+7}\) und dem Ziel \(N=10^{14}\). Direktes Aufsummieren ist aussichtslos, daher formt die Lösung \(R(n)\) in eine arithmetische Funktion mit schneller Summenformel um.
Durch nochmaliges Anwenden von \(f(x)=ax+b\) erhält man
$$f(f(x))\equiv a^2x+ab+b \pmod{n}.$$
Damit dies für alle Restklassen mit \(ax+b\) übereinstimmt, müssen beide Korrekturterme modulo \(n\) verschwinden:
$$n \mid a(a-1),\qquad n \mid ab.$$
Setzt man \(d=\gcd(a,n)\), dann gilt \(d\mid n\) und \(\gcd(d,n/d)=1\); also ist \(d\) ein unitärer Teiler von \(n\). Umgekehrt bestimmt jeder unitäre Teiler \(d \lt n\) genau eine zulässige Restklasse für \(a\), und für festes \(a\) gibt es genau \(d\) gültige Werte von \(b\). Daher
$$R(n)=\sum_{\substack{d \mid n\\ \gcd(d,n/d)=1\\ d \lt n}} d=\sigma^*(n)-n,$$
wobei \(\sigma^*(n)\) die Summe der unitären Teiler von \(n\) bezeichnet. Da \(\sigma^*(1)=1\) gilt, liefert der fehlende Term bei \(n=1\) genau \(1-1=0\). Also
$$F(N)=\sum_{n=1}^{N}\bigl(\sigma^*(n)-n\bigr)=\left(\sum_{n=1}^{N}\sigma^*(n)\right)-\frac{N(N+1)}{2}.$$
Die unitäre Bedingung lässt sich mit Möbius-Inversion kodieren:
$$\mathbf{1}_{\gcd(d,n/d)=1}=\sum_{t \mid \gcd(d,n/d)} \mu(t).$$
Eingesetzt in die Teilersumme ergibt das
$$\sigma^*(n)=\sum_{d \mid n} d \sum_{t \mid \gcd(d,n/d)} \mu(t).$$
Aus \(t \mid \gcd(d,n/d)\) folgt \(t^2 \mid n\), und man kann \(d=t\,m\) mit \(m \mid n/t^2\) schreiben. Daher
$$\sigma^*(n)=\sum_{t^2 \mid n} \mu(t)\, t \sum_{m \mid n/t^2} m=\sum_{t^2 \mid n} \mu(t)\, t\, \sigma\!\left(\frac{n}{t^2}\right).$$
Genau diese Identität verwenden die Implementierungen. Sie ersetzt die Summe über unitäre Teiler durch die gewöhnliche Teilersummenfunktion \(\sigma\), deren Summenformel deutlich günstiger ist.
Definiere
$$T(x)=\sum_{m=1}^{x}\sigma(m).$$
Summiert man die vorige Identität über alle \(n \le N\) und vertauscht die Reihenfolge der Summation, so erhält man
$$\sum_{n=1}^{N}\sigma^*(n)=\sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right).$$
Für \(T(x)\) nutzt man die klassische Identität
$$T(x)=\sum_{m=1}^{x}\sum_{d \mid m} d=\sum_{d=1}^{x} d\left\lfloor\frac{x}{d}\right\rfloor.$$
Damit ist das ursprüngliche Problem auf die schnelle Auswertung vieler Werte von \(T(x)\) reduziert.
Der Wert \(\left\lfloor x/i \right\rfloor\) bleibt auf ganzen Intervallen konstant. Ist \(l\) der erste Index eines Blocks, setze
$$q=\left\lfloor\frac{x}{l}\right\rfloor,\qquad r=\left\lfloor\frac{x}{q}\right\rfloor.$$
Dann gilt für jedes \(i\in[l,r]\) derselbe Quotient \(q\), und der ganze Block trägt bei mit
$$q\sum_{i=l}^{r} i=q\cdot \frac{(l+r)(r-l+1)}{2}.$$
Indem man direkt von \(l\) nach \(r+1\) springt, überspringt man jeweils einen kompletten Bereich mit konstantem Quotienten. So sinken die Kosten für \(T(x)\) von \(O(x)\) auf ungefähr \(O(\sqrt{x})\). Diese Blockzerlegung ist die entscheidende Beschleunigung der Lösung.
Zusammengesetzt ergibt sich
$$\boxed{F(N)\equiv \sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right)-\frac{N(N+1)}{2}\pmod{10^9+7}.}$$
Zur Kontrolle setzen wir \(N=10\). Mit Quotientenblöcken erhält man
$$T(10)=1\cdot 10+2\cdot 5+3\cdot 3+(4+5)\cdot 2+(6+7+8+9+10)\cdot 1=87.$$
Außerdem sind \(T(2)=4\), \(T(1)=1\), und die nichtverschwindenden Möbius-Werte bis \(\sqrt{10}\) sind \(\mu(1)=1\), \(\mu(2)=-1\), \(\mu(3)=-1\). Also
$$\sum_{n=1}^{10}\sigma^*(n)=1\cdot 87-2\cdot 4-3\cdot 1=76,$$
und damit
$$F(10)=76-\frac{10\cdot 11}{2}=21.$$
Das stimmt mit der direkten Auswertung von \(R(2),\dots,R(10)\) überein und bestätigt die Herleitung im kleinen Fall.
Die C++-, Python- und Java-Implementierungen sieben zunächst die Möbius-Funktion bis \(\lfloor\sqrt{N}\rfloor\). In der äußeren Summe werden dann alle \(t\) mit \(\mu(t)=0\) übersprungen, weil diese Beiträge verschwinden. Für jedes verbleibende \(t\) wird genau ein innerer Wert \(T\!\left(\lfloor N/t^2\rfloor\right)\) per Quotientenblock-Technik berechnet, nicht durch lineares Durchlaufen aller Indizes.
Nach dem Aufsummieren des unitären Teilersummenanteils wird die Dreieckszahl \(N(N+1)/2\) modulo \(10^9+7\) abgezogen. Die C++- und Java-Implementierungen verteilen die äußere Möbius-Summe zyklisch auf mehrere Threads, damit die teuren kleinen \(t\)-Werte besser ausbalanciert sind; die Python-Implementierung verwendet dieselbe Mathematik sequentiell. In allen drei Sprachen wird die modulare Subtraktion wieder in den Bereich \(0,\dots,10^9+6\) zurückgeführt.
Das Möbius-Sieb bis \(\sqrt{N}\) benötigt \(O(\sqrt{N})\) Speicher und nahezu lineare Vorverarbeitung in diesem Bereich. Der dominierende Aufwand ist
$$\sum_{t \le \sqrt{N}} O\!\left(\sqrt{\frac{N}{t^2}}\right)=\sum_{t \le \sqrt{N}} O\!\left(\frac{\sqrt{N}}{t}\right)=O(\sqrt{N}\log N).$$
Damit läuft das Verfahren in \(O(\sqrt{N}\log N)\) Zeit und \(O(\sqrt{N})\) Speicher. Parallelisierung verbessert die Laufzeit in der Praxis, ändert aber nicht die asymptotische Ordnung.
Her \(n \gt 1\) için şu doğrusal dönüşümleri ele alıyoruz:
$$f(x)\equiv ax+b \pmod{n},\qquad 0 \lt a \lt n,\quad 0 \le b \lt n.$$
Bir dönüşüm, her kalıntı sınıfında idempotent ise retraksiyon olur:
$$f(f(x))\equiv f(x)\pmod{n}\qquad \forall x.$$
\(R(n)\), bu tür dönüşümlerin sayısı olsun. Problem 447 şu toplamı ister:
$$F(N)=\sum_{n=2}^{N}R(n)\pmod{10^9+7},$$
ve verilen kontrol değeri \(F(10^7)\equiv 638042271 \pmod{10^9+7}\), hedef ise \(N=10^{14}\)'tür. Bu ölçekte doğrudan toplama mümkün olmadığından, çözüm \(R(n)\)'yi hızlı önek toplamı alınabilen bir aritmetik fonksiyona dönüştürür.
\(f(x)=ax+b\) dönüşümünü bir kez daha bileştirince
$$f(f(x))\equiv a^2x+ab+b \pmod{n}$$
elde edilir. Bunun her \(x\) için yeniden \(ax+b\)'ye eşit olması için iki koşul gerekir:
$$n \mid a(a-1),\qquad n \mid ab.$$
\(d=\gcd(a,n)\) yazalım. O zaman \(d\mid n\) ve \(\gcd(d,n/d)=1\) olur; yani \(d\), \(n\)'nin bir unitary bölenidir. Tersine, her \(d \lt n\) unitary böleni için tam bir uygun \(a\) kalıntı sınıfı vardır ve \(a\) sabitlenince geçerli \(b\) değerlerinin sayısı tam olarak \(d\)'dir. Dolayısıyla
$$R(n)=\sum_{\substack{d \mid n\\ \gcd(d,n/d)=1\\ d \lt n}} d=\sigma^*(n)-n,$$
burada \(\sigma^*(n)\), \(n\)'nin unitary bölenlerinin toplamıdır. Ayrıca \(\sigma^*(1)=1\) olduğundan, eksik bırakılan \(n=1\) terimi \(1-1=0\) verir. Bu yüzden
$$F(N)=\sum_{n=1}^{N}\bigl(\sigma^*(n)-n\bigr)=\left(\sum_{n=1}^{N}\sigma^*(n)\right)-\frac{N(N+1)}{2}.$$
Unitary koşulu Möbius terslemesiyle gösterilebilir:
$$\mathbf{1}_{\gcd(d,n/d)=1}=\sum_{t \mid \gcd(d,n/d)} \mu(t).$$
Bunu bölen toplamına yerleştirirsek
$$\sigma^*(n)=\sum_{d \mid n} d \sum_{t \mid \gcd(d,n/d)} \mu(t)$$
olur. Eğer \(t \mid \gcd(d,n/d)\) ise \(t^2 \mid n\) olur ve \(d=t\,m\) biçiminde, \(m \mid n/t^2\) olacak şekilde yazılabilir. Böylece
$$\sigma^*(n)=\sum_{t^2 \mid n} \mu(t)\, t \sum_{m \mid n/t^2} m=\sum_{t^2 \mid n} \mu(t)\, t\, \sigma\!\left(\frac{n}{t^2}\right)$$
kimliğini elde ederiz. Uygulamanın kullandığı temel dönüşüm budur: unitary bölen toplamı, normal bölen-toplam fonksiyonuna \(\sigma\) indirgenir.
Şunu tanımlayalım:
$$T(x)=\sum_{m=1}^{x}\sigma(m).$$
Bir önceki kimliği \(n \le N\) için toplayıp toplam sırasını değiştirince
$$\sum_{n=1}^{N}\sigma^*(n)=\sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right)$$
sonucu çıkar. Ayrıca klasik bölen-toplam özdeşliği
$$T(x)=\sum_{m=1}^{x}\sum_{d \mid m} d=\sum_{d=1}^{x} d\left\lfloor\frac{x}{d}\right\rfloor$$
olduğundan, problem artık birçok \(T(x)\) değerini hızlı hesaplama işine indirgenmiştir.
\(\left\lfloor x/i \right\rfloor\) ifadesi aralıklar boyunca sabit kalır. Bir bloğun ilk indeksi \(l\) olsun:
$$q=\left\lfloor\frac{x}{l}\right\rfloor,\qquad r=\left\lfloor\frac{x}{q}\right\rfloor.$$
O zaman \(i\in[l,r]\) için bölüm değeri hep \(q\)'dur ve blok katkısı
$$q\sum_{i=l}^{r} i=q\cdot \frac{(l+r)(r-l+1)}{2}$$
şeklinde tek seferde hesaplanır. \(l\)'den doğrudan \(r+1\)'e atlamak, sabit bölüm veren tüm aralığı bir anda geçmek demektir. Böylece \(T(x)\)'in maliyeti \(O(x)\)'ten yaklaşık \(O(\sqrt{x})\)'e düşer. Uygulamadaki ana hız kazanımı budur.
Tüm parçaları birleştirince
$$\boxed{F(N)\equiv \sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right)-\frac{N(N+1)}{2}\pmod{10^9+7}.}$$
Küçük bir kontrol için \(N=10\) alalım. Bölüm bloklarıyla
$$T(10)=1\cdot 10+2\cdot 5+3\cdot 3+(4+5)\cdot 2+(6+7+8+9+10)\cdot 1=87$$
bulunur. Ayrıca \(T(2)=4\), \(T(1)=1\) ve \(\sqrt{10}\)'a kadar sıfır olmayan Möbius değerleri \(\mu(1)=1\), \(\mu(2)=-1\), \(\mu(3)=-1\)'dir. Dolayısıyla
$$\sum_{n=1}^{10}\sigma^*(n)=1\cdot 87-2\cdot 4-3\cdot 1=76,$$
ve sonuç
$$F(10)=76-\frac{10\cdot 11}{2}=21$$
olur. Bu değer, \(R(2),\dots,R(10)\)'un doğrudan hesabıyla da aynıdır.
C++, Python ve Java uygulamaları önce Möbius fonksiyonunu \(\lfloor\sqrt{N}\rfloor\)'ye kadar bir elek ile hesaplar. Dış toplamda \(\mu(t)=0\) olan tüm \(t\) değerleri atlanır; çünkü bu terimlerin katkısı sıfırdır. Kalan her \(t\) için bir adet \(T\!\left(\lfloor N/t^2\rfloor\right)\) sorgusu yapılır ve bu sorgu doğrusal döngü yerine bölüm bloklarıyla çözülür.
Unitary bölenlerin kümülatif katkısı toplandıktan sonra \(N(N+1)/2\) üçgensel toplamı mod \(10^9+7\) altında çıkarılır. C++ ve Java uygulamaları, pahalı küçük \(t\) değerlerini daha dengeli dağıtmak için dış Möbius toplamını çevrimsel biçimde iş parçacıklarına böler; Python uygulaması ise aynı matematiği ardışık olarak yürütür. Her üç dilde de modüler çıkarma sonrası sonuç yeniden \(0,\dots,10^9+6\) aralığına normalize edilir.
\(\sqrt{N}\)'ye kadar Möbius eleği \(O(\sqrt{N})\) bellek ve bu aralıkta yaklaşık doğrusal ön işlem ister. Baskın maliyet
$$\sum_{t \le \sqrt{N}} O\!\left(\sqrt{\frac{N}{t^2}}\right)=\sum_{t \le \sqrt{N}} O\!\left(\frac{\sqrt{N}}{t}\right)=O(\sqrt{N}\log N)$$
şeklindedir. Sonuç olarak yöntem \(O(\sqrt{N}\log N)\) zamanda ve \(O(\sqrt{N})\) bellekte çalışır. Paralelleştirme pratik süreyi düşürür, fakat asimptotik sınırı değiştirmez.
Para cada entero \(n \gt 1\), consideramos las aplicaciones lineales
$$f(x)\equiv ax+b \pmod{n},\qquad 0 \lt a \lt n,\quad 0 \le b \lt n.$$
La aplicación es una retracción cuando es idempotente en todas las clases residuales:
$$f(f(x))\equiv f(x)\pmod{n}\qquad \forall x.$$
Sea \(R(n)\) el número de tales aplicaciones. El problema 447 pide
$$F(N)=\sum_{n=2}^{N}R(n)\pmod{10^9+7},$$
con el punto de control \(F(10^7)\equiv 638042271 \pmod{10^9+7}\) y objetivo \(N=10^{14}\). Una suma directa es inviable, así que la solución transforma \(R(n)\) en una función aritmética con fórmula acumulada rápida.
Al componer de nuevo \(f(x)=ax+b\), se obtiene
$$f(f(x))\equiv a^2x+ab+b \pmod{n}.$$
Para que esto coincida con \(ax+b\) para toda clase residual, deben anularse ambos términos de corrección módulo \(n\):
$$n \mid a(a-1),\qquad n \mid ab.$$
Si escribimos \(d=\gcd(a,n)\), entonces \(d\mid n\) y \(\gcd(d,n/d)=1\), de modo que \(d\) es un divisor unitario de \(n\). Recíprocamente, cada divisor unitario \(d \lt n\) determina exactamente una clase residual admisible para \(a\), y una vez fijado \(a\), el número de valores válidos de \(b\) es exactamente \(d\). Por tanto,
$$R(n)=\sum_{\substack{d \mid n\\ \gcd(d,n/d)=1\\ d \lt n}} d=\sigma^*(n)-n,$$
donde \(\sigma^*(n)\) es la suma de los divisores unitarios de \(n\). Como \(\sigma^*(1)=1\), el término faltante con \(n=1\) aporta \(1-1=0\). Así,
$$F(N)=\sum_{n=1}^{N}\bigl(\sigma^*(n)-n\bigr)=\left(\sum_{n=1}^{N}\sigma^*(n)\right)-\frac{N(N+1)}{2}.$$
La condición unitaria puede codificarse mediante inversión de Möbius:
$$\mathbf{1}_{\gcd(d,n/d)=1}=\sum_{t \mid \gcd(d,n/d)} \mu(t).$$
Al insertarla en la suma sobre divisores se obtiene
$$\sigma^*(n)=\sum_{d \mid n} d \sum_{t \mid \gcd(d,n/d)} \mu(t).$$
Si \(t \mid \gcd(d,n/d)\), entonces \(t^2 \mid n\), y podemos escribir \(d=t\,m\) con \(m \mid n/t^2\). Por ello
$$\sigma^*(n)=\sum_{t^2 \mid n} \mu(t)\, t \sum_{m \mid n/t^2} m=\sum_{t^2 \mid n} \mu(t)\, t\, \sigma\!\left(\frac{n}{t^2}\right).$$
Esta es la identidad usada por la implementación. Convierte la suma de divisores unitarios en la función habitual \(\sigma\), mucho más cómoda de acumular.
Definimos
$$T(x)=\sum_{m=1}^{x}\sigma(m).$$
Sumando la identidad anterior para \(n \le N\) e intercambiando el orden de sumación, resulta
$$\sum_{n=1}^{N}\sigma^*(n)=\sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right).$$
Además, la identidad clásica de sumas sobre divisores dice que
$$T(x)=\sum_{m=1}^{x}\sum_{d \mid m} d=\sum_{d=1}^{x} d\left\lfloor\frac{x}{d}\right\rfloor.$$
Así, el problema original se reduce a calcular muchas veces \(T(x)\) de forma eficiente.
El valor \(\left\lfloor x/i \right\rfloor\) permanece constante en intervalos. Si \(l\) es el primer índice de un bloque, definimos
$$q=\left\lfloor\frac{x}{l}\right\rfloor,\qquad r=\left\lfloor\frac{x}{q}\right\rfloor.$$
Entonces todo \(i\in[l,r]\) tiene el mismo cociente \(q\), y la contribución completa del bloque es
$$q\sum_{i=l}^{r} i=q\cdot \frac{(l+r)(r-l+1)}{2}.$$
Al saltar directamente de \(l\) a \(r+1\), se evita recorrer uno por uno todos los índices del intervalo. Con eso, el coste de \(T(x)\) baja de \(O(x)\) a aproximadamente \(O(\sqrt{x})\). Este es el acelerador principal de la solución.
Uniendo todas las piezas,
$$\boxed{F(N)\equiv \sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right)-\frac{N(N+1)}{2}\pmod{10^9+7}.}$$
Para una verificación rápida, tomemos \(N=10\). Los bloques de cociente dan
$$T(10)=1\cdot 10+2\cdot 5+3\cdot 3+(4+5)\cdot 2+(6+7+8+9+10)\cdot 1=87.$$
Además, \(T(2)=4\), \(T(1)=1\), y los valores no nulos de Möbius hasta \(\sqrt{10}\) son \(\mu(1)=1\), \(\mu(2)=-1\), \(\mu(3)=-1\). Por tanto,
$$\sum_{n=1}^{10}\sigma^*(n)=1\cdot 87-2\cdot 4-3\cdot 1=76,$$
y entonces
$$F(10)=76-\frac{10\cdot 11}{2}=21.$$
Ese valor coincide con la evaluación directa de \(R(2),\dots,R(10)\), así que la derivación queda validada en un caso pequeño.
Las implementaciones en C++, Python y Java primero calculan la función de Möbius hasta \(\lfloor\sqrt{N}\rfloor\) mediante una criba. En la suma exterior se omiten todos los \(t\) con \(\mu(t)=0\), porque no aportan nada. Para cada \(t\) restante se evalúa una consulta interna \(T\!\left(\lfloor N/t^2\rfloor\right)\) usando bloques de cociente en lugar de un recorrido lineal.
Después de acumular la parte correspondiente a los divisores unitarios, la implementación resta el número triangular \(N(N+1)/2\) módulo \(10^9+7\). Las versiones de C++ y Java reparten cíclicamente la suma exterior entre hilos de ejecución para equilibrar mejor los pequeños \(t\), que son los más costosos; la versión de Python sigue exactamente las mismas fórmulas de forma secuencial. En los tres lenguajes, toda resta modular se normaliza de nuevo al rango \(0,\dots,10^9+6\).
La criba de Möbius hasta \(\sqrt{N}\) usa \(O(\sqrt{N})\) memoria y un preprocesamiento casi lineal en ese rango. El trabajo dominante es
$$\sum_{t \le \sqrt{N}} O\!\left(\sqrt{\frac{N}{t^2}}\right)=\sum_{t \le \sqrt{N}} O\!\left(\frac{\sqrt{N}}{t}\right)=O(\sqrt{N}\log N).$$
Por lo tanto, el método completo requiere \(O(\sqrt{N}\log N)\) tiempo y \(O(\sqrt{N})\) memoria. La paralelización mejora el tiempo real en las implementaciones compiladas, pero no cambia el orden asintótico.
对每个 \(n \gt 1\)，考虑线性映射
$$f(x)\equiv ax+b \pmod{n},\qquad 0 \lt a \lt n,\quad 0 \le b \lt n.$$
如果它在每个剩余类上都是幂等的，就称它为一个 retraction：
$$f(f(x))\equiv f(x)\pmod{n}\qquad \forall x.$$
记 \(R(n)\) 为这类映射的个数。第 447 题要求计算
$$F(N)=\sum_{n=2}^{N}R(n)\pmod{10^9+7},$$
并给出校验点 \(F(10^7)\equiv 638042271 \pmod{10^9+7}\)，目标是 \(N=10^{14}\)。在这个规模下不可能逐个枚举，因此实现采用数论变形，把 \(R(n)\) 改写成可以快速求前缀和的形式。
把 \(f(x)=ax+b\) 再复合一次，有
$$f(f(x))\equiv a^2x+ab+b \pmod{n}.$$
要让它对所有 \(x\) 都等于 \(ax+b\)，必须满足
$$n \mid a(a-1),\qquad n \mid ab.$$
设 \(d=\gcd(a,n)\)。那么 \(d\mid n\) 且 \(\gcd(d,n/d)=1\)，所以 \(d\) 是 \(n\) 的一个 unitary 因子。反过来，每个 \(d \lt n\) 的 unitary 因子都唯一决定一个可行的 \(a\) 剩余类，而一旦 \(a\) 固定，可行的 \(b\) 恰好有 \(d\) 个。因此
$$R(n)=\sum_{\substack{d \mid n\\ \gcd(d,n/d)=1\\ d \lt n}} d=\sigma^*(n)-n,$$
其中 \(\sigma^*(n)\) 表示 \(n\) 的 unitary 因子之和。由于 \(\sigma^*(1)=1\)，缺少的 \(n=1\) 项贡献正好是 \(1-1=0\)。于是
$$F(N)=\sum_{n=1}^{N}\bigl(\sigma^*(n)-n\bigr)=\left(\sum_{n=1}^{N}\sigma^*(n)\right)-\frac{N(N+1)}{2}.$$
unitary 条件可以用 Möbius 反演表示：
$$\mathbf{1}_{\gcd(d,n/d)=1}=\sum_{t \mid \gcd(d,n/d)} \mu(t).$$
代回因子求和后得到
$$\sigma^*(n)=\sum_{d \mid n} d \sum_{t \mid \gcd(d,n/d)} \mu(t).$$
若 \(t \mid \gcd(d,n/d)\)，则必有 \(t^2 \mid n\)，并且可写 \(d=t\,m\)，其中 \(m \mid n/t^2\)。因此
$$\sigma^*(n)=\sum_{t^2 \mid n} \mu(t)\, t \sum_{m \mid n/t^2} m=\sum_{t^2 \mid n} \mu(t)\, t\, \sigma\!\left(\frac{n}{t^2}\right).$$
这正是实现所使用的恒等式。它把 unitary 因子和转化成普通约数和函数 \(\sigma\)，后者更容易做前缀求和。
定义
$$T(x)=\sum_{m=1}^{x}\sigma(m).$$
对上式在 \(n \le N\) 上求和并交换求和顺序，可得
$$\sum_{n=1}^{N}\sigma^*(n)=\sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right).$$
再利用经典恒等式
$$T(x)=\sum_{m=1}^{x}\sum_{d \mid m} d=\sum_{d=1}^{x} d\left\lfloor\frac{x}{d}\right\rfloor,$$
原问题就变成了高效计算许多次 \(T(x)\)。
\(\left\lfloor x/i \right\rfloor\) 在一段区间内保持不变。设某个块的起点为 \(l\)，定义
$$q=\left\lfloor\frac{x}{l}\right\rfloor,\qquad r=\left\lfloor\frac{x}{q}\right\rfloor.$$
那么对所有 \(i\in[l,r]\)，商都等于 \(q\)。整个块的贡献可以一次写成
$$q\sum_{i=l}^{r} i=q\cdot \frac{(l+r)(r-l+1)}{2}.$$
把指针直接从 \(l\) 跳到 \(r+1\)，就能整段跳过同商区间。这样 \(T(x)\) 的代价从 \(O(x)\) 降到大约 \(O(\sqrt{x})\)。这是三种语言实现共同的核心优化。
综合起来，得到
$$\boxed{F(N)\equiv \sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right)-\frac{N(N+1)}{2}\pmod{10^9+7}.}$$
取 \(N=10\) 做一个快速检验。整除分块给出
$$T(10)=1\cdot 10+2\cdot 5+3\cdot 3+(4+5)\cdot 2+(6+7+8+9+10)\cdot 1=87.$$
同时 \(T(2)=4\)、\(T(1)=1\)，而到 \(\sqrt{10}\) 为止非零的 Möbius 值是 \(\mu(1)=1\)、\(\mu(2)=-1\)、\(\mu(3)=-1\)。所以
$$\sum_{n=1}^{10}\sigma^*(n)=1\cdot 87-2\cdot 4-3\cdot 1=76,$$
从而
$$F(10)=76-\frac{10\cdot 11}{2}=21.$$
这与直接计算 \(R(2),\dots,R(10)\) 的结果一致，因此推导在小范围内是自洽的。
C++、Python 和 Java 实现都会先筛出到 \(\lfloor\sqrt{N}\rfloor\) 为止的 Möbius 函数值。外层求和中，凡是 \(\mu(t)=0\) 的项都可以直接跳过。对每个剩余的 \(t\)，实现只需求一次 \(T\!\left(\lfloor N/t^2\rfloor\right)\)，而这个值通过整除分块获得，不会线性扫到 \(x\)。
累计完 unitary 因子部分之后，再减去三角和 \(N(N+1)/2\)，并始终在模 \(10^9+7\) 下处理。C++ 和 Java 实现把外层的 Möbius 求和按循环方式分配给多个线程，以便更均衡地分担代价最高的小 \(t\)；Python 实现则按同样的数学公式顺序执行。三种语言都会把负的模结果重新规范到 \(0,\dots,10^9+6\)。
筛到 \(\sqrt{N}\) 的 Möbius 值需要 \(O(\sqrt{N})\) 空间，并在这个范围内进行近线性的预处理。主导工作量为
$$\sum_{t \le \sqrt{N}} O\!\left(\sqrt{\frac{N}{t^2}}\right)=\sum_{t \le \sqrt{N}} O\!\left(\frac{\sqrt{N}}{t}\right)=O(\sqrt{N}\log N).$$
因此整体复杂度是 \(O(\sqrt{N}\log N)\) 时间和 \(O(\sqrt{N})\) 空间。并行化能改善编译型实现的实际运行时间，但不会改变渐近阶。
Для каждого \(n \gt 1\) рассматриваются линейные отображения
$$f(x)\equiv ax+b \pmod{n},\qquad 0 \lt a \lt n,\quad 0 \le b \lt n.$$
Такое отображение является ретракцией, если оно идемпотентно на всех классах вычетов:
$$f(f(x))\equiv f(x)\pmod{n}\qquad \forall x.$$
Пусть \(R(n)\) обозначает число таких отображений. В задаче 447 требуется вычислить
$$F(N)=\sum_{n=2}^{N}R(n)\pmod{10^9+7},$$
причем дан контрольный результат \(F(10^7)\equiv 638042271 \pmod{10^9+7}\), а целевое значение соответствует \(N=10^{14}\). Прямой перебор невозможен, поэтому решение сводит \(R(n)\) к арифметической функции с быстрой формулой для префиксной суммы.
Если снова подставить \(f(x)=ax+b\) в себя, получится
$$f(f(x))\equiv a^2x+ab+b \pmod{n}.$$
Чтобы это совпадало с \(ax+b\) для любого \(x\), должны выполняться условия
$$n \mid a(a-1),\qquad n \mid ab.$$
Положим \(d=\gcd(a,n)\). Тогда \(d\mid n\) и \(\gcd(d,n/d)=1\), то есть \(d\) является унитарным делителем числа \(n\). Обратно, каждый унитарный делитель \(d \lt n\) задает ровно один допустимый класс вычетов для \(a\), а при фиксированном \(a\) число допустимых значений \(b\) равно ровно \(d\). Следовательно,
$$R(n)=\sum_{\substack{d \mid n\\ \gcd(d,n/d)=1\\ d \lt n}} d=\sigma^*(n)-n,$$
где \(\sigma^*(n)\) обозначает сумму унитарных делителей числа \(n\). Поскольку \(\sigma^*(1)=1\), пропущенный член при \(n=1\) дает \(1-1=0\). Поэтому
$$F(N)=\sum_{n=1}^{N}\bigl(\sigma^*(n)-n\bigr)=\left(\sum_{n=1}^{N}\sigma^*(n)\right)-\frac{N(N+1)}{2}.$$
Условие унитарности можно записать через обращение Möbius:
$$\mathbf{1}_{\gcd(d,n/d)=1}=\sum_{t \mid \gcd(d,n/d)} \mu(t).$$
Подставляя это в сумму по делителям, получаем
$$\sigma^*(n)=\sum_{d \mid n} d \sum_{t \mid \gcd(d,n/d)} \mu(t).$$
Если \(t \mid \gcd(d,n/d)\), то \(t^2 \mid n\), и можно написать \(d=t\,m\), где \(m \mid n/t^2\). Тогда
$$\sigma^*(n)=\sum_{t^2 \mid n} \mu(t)\, t \sum_{m \mid n/t^2} m=\sum_{t^2 \mid n} \mu(t)\, t\, \sigma\!\left(\frac{n}{t^2}\right).$$
Именно эта формула используется в реализации. Она переводит сумму по унитарным делителям в обычную функцию суммы делителей \(\sigma\), с которой работать существенно проще.
Определим
$$T(x)=\sum_{m=1}^{x}\sigma(m).$$
Просуммировав предыдущую формулу по всем \(n \le N\) и поменяв порядок суммирования, получаем
$$\sum_{n=1}^{N}\sigma^*(n)=\sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right).$$
Далее используется классическое тождество
$$T(x)=\sum_{m=1}^{x}\sum_{d \mid m} d=\sum_{d=1}^{x} d\left\lfloor\frac{x}{d}\right\rfloor.$$
Тем самым исходная задача сводится к быстрому вычислению многих значений \(T(x)\).
Величина \(\left\lfloor x/i \right\rfloor\) постоянна на целых интервалах. Если \(l\) — первый индекс блока, положим
$$q=\left\lfloor\frac{x}{l}\right\rfloor,\qquad r=\left\lfloor\frac{x}{q}\right\rfloor.$$
Тогда для всех \(i\in[l,r]\) частное равно \(q\), и вклад всего блока равен
$$q\sum_{i=l}^{r} i=q\cdot \frac{(l+r)(r-l+1)}{2}.$$
Переход сразу от \(l\) к \(r+1\) позволяет пропускать целый диапазон одинаковых значений. Благодаря этому стоимость вычисления \(T(x)\) уменьшается с \(O(x)\) до примерно \(O(\sqrt{x})\). Это и есть главный ускоряющий прием во всех трех реализациях.
Итак, получаем
$$\boxed{F(N)\equiv \sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right)-\frac{N(N+1)}{2}\pmod{10^9+7}.}$$
Для быстрой проверки возьмем \(N=10\). Блочное разбиение дает
$$T(10)=1\cdot 10+2\cdot 5+3\cdot 3+(4+5)\cdot 2+(6+7+8+9+10)\cdot 1=87.$$
Кроме того, \(T(2)=4\), \(T(1)=1\), а ненулевые значения Möbius до \(\sqrt{10}\) равны \(\mu(1)=1\), \(\mu(2)=-1\), \(\mu(3)=-1\). Следовательно,
$$\sum_{n=1}^{10}\sigma^*(n)=1\cdot 87-2\cdot 4-3\cdot 1=76,$$
и потому
$$F(10)=76-\frac{10\cdot 11}{2}=21.$$
Это совпадает с прямым подсчетом \(R(2),\dots,R(10)\), так что вывод формулы проверяется уже на малом примере.
Реализации на C++, Python и Java сначала вычисляют функцию Möbius до \(\lfloor\sqrt{N}\rfloor\) с помощью решета. Во внешней сумме все \(t\) с \(\mu(t)=0\) сразу пропускаются, так как их вклад равен нулю. Для каждого оставшегося \(t\) вычисляется одно значение \(T\!\left(\lfloor N/t^2\rfloor\right)\), причем не линейным перебором, а проходом по блокам постоянного частного.
После накопления вклада от суммы унитарных делителей программа вычитает треугольное число \(N(N+1)/2\) по модулю \(10^9+7\). Реализации на C++ и Java циклически распределяют внешнюю сумму по потокам, чтобы более дорогие малые \(t\) делились между ними равномернее; реализация на Python использует ту же математику последовательно. Во всех языках отрицательный остаток после вычитания переводится обратно в диапазон \(0,\dots,10^9+6\).
Решето Möbius до \(\sqrt{N}\) требует \(O(\sqrt{N})\) памяти и почти линейной предобработки на этом диапазоне. Основной объем работы равен
$$\sum_{t \le \sqrt{N}} O\!\left(\sqrt{\frac{N}{t^2}}\right)=\sum_{t \le \sqrt{N}} O\!\left(\frac{\sqrt{N}}{t}\right)=O(\sqrt{N}\log N).$$
Итак, полный метод работает за \(O(\sqrt{N}\log N)\) времени и использует \(O(\sqrt{N})\) памяти. Параллелизация улучшает практическое время работы компилируемых версий, но асимптотическую оценку не меняет.
لكل عدد صحيح \(n \gt 1\)، ننظر إلى التطبيقات الخطية
$$f(x)\equiv ax+b \pmod{n},\qquad 0 \lt a \lt n,\quad 0 \le b \lt n.$$
تكون الدالة ارتدادًا إذا كانت idempotent على جميع فئات البواقي:
$$f(f(x))\equiv f(x)\pmod{n}\qquad \forall x.$$
ولتكن \(R(n)\) عدد هذه الدوال. المطلوب في المسألة 447 هو حساب
$$F(N)=\sum_{n=2}^{N}R(n)\pmod{10^9+7},$$
مع نقطة تحقق معطاة هي \(F(10^7)\equiv 638042271 \pmod{10^9+7}\)، والهدف النهائي \(N=10^{14}\). العد المباشر غير ممكن بهذا الحجم، لذلك تعيد الحلول كتابة \(R(n)\) على صورة دالة عددية لها صيغة تراكمية سريعة.
إذا ركّبنا \(f(x)=ax+b\) مرة أخرى نحصل على
$$f(f(x))\equiv a^2x+ab+b \pmod{n}.$$
ولكي يساوي هذا \(ax+b\) لكل \(x\)، يجب أن يتحقق الشرطان
$$n \mid a(a-1),\qquad n \mid ab.$$
لنكتب \(d=\gcd(a,n)\). عندئذٍ \(d\mid n\) و\(\gcd(d,n/d)=1\)، أي إن \(d\) قاسم أحادي للعدد \(n\). وبالعكس، كل قاسم أحادي \(d \lt n\) يحدد فئة بواقي واحدة صالحة لـ \(a\)، وبعد تثبيت \(a\) يكون عدد قيم \(b\) الصالحة مساويًا تمامًا لـ \(d\). لذلك
$$R(n)=\sum_{\substack{d \mid n\\ \gcd(d,n/d)=1\\ d \lt n}} d=\sigma^*(n)-n,$$
حيث ترمز \(\sigma^*(n)\) إلى مجموع القواسم الأحادية للعدد \(n\). وبما أن \(\sigma^*(1)=1\)، فإن الحد الناقص عند \(n=1\) يساوي \(1-1=0\). ومن ثم
$$F(N)=\sum_{n=1}^{N}\bigl(\sigma^*(n)-n\bigr)=\left(\sum_{n=1}^{N}\sigma^*(n)\right)-\frac{N(N+1)}{2}.$$
يمكن ترميز شرط القاسم الأحادي بواسطة قلب Möbius:
$$\mathbf{1}_{\gcd(d,n/d)=1}=\sum_{t \mid \gcd(d,n/d)} \mu(t).$$
وبالتعويض داخل مجموع القواسم نحصل على
$$\sigma^*(n)=\sum_{d \mid n} d \sum_{t \mid \gcd(d,n/d)} \mu(t).$$
إذا كان \(t \mid \gcd(d,n/d)\)، فلابد أن \(t^2 \mid n\)، ويمكن كتابة \(d=t\,m\) بحيث \(m \mid n/t^2\). إذًا
$$\sigma^*(n)=\sum_{t^2 \mid n} \mu(t)\, t \sum_{m \mid n/t^2} m=\sum_{t^2 \mid n} \mu(t)\, t\, \sigma\!\left(\frac{n}{t^2}\right).$$
هذه هي الهوية التي تعتمد عليها التنفيذات. فهي تنقلنا من مجموع القواسم الأحادية إلى دالة مجموع القواسم العادية \(\sigma\)، وهي أسهل بكثير في الجمع التراكمي.
نعرّف
$$T(x)=\sum_{m=1}^{x}\sigma(m).$$
وبجمع الهوية السابقة على جميع \(n \le N\) وتبديل ترتيب الجمع نحصل على
$$\sum_{n=1}^{N}\sigma^*(n)=\sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right).$$
ثم نستخدم الهوية الكلاسيكية
$$T(x)=\sum_{m=1}^{x}\sum_{d \mid m} d=\sum_{d=1}^{x} d\left\lfloor\frac{x}{d}\right\rfloor.$$
وبذلك تصبح المسألة الأصلية هي حساب قيم كثيرة من \(T(x)\) بسرعة.
القيمة \(\left\lfloor x/i \right\rfloor\) تبقى ثابتة على فترات كاملة. إذا كان \(l\) هو أول فهرس في كتلة ما، فنعرف
$$q=\left\lfloor\frac{x}{l}\right\rfloor,\qquad r=\left\lfloor\frac{x}{q}\right\rfloor.$$
عندئذٍ لكل \(i\in[l,r]\) يكون خارج القسمة مساويًا لـ \(q\)، ومساهمة الكتلة كلها هي
$$q\sum_{i=l}^{r} i=q\cdot \frac{(l+r)(r-l+1)}{2}.$$
والانتقال مباشرة من \(l\) إلى \(r+1\) يعني تخطي مجال كامل دفعة واحدة. بهذه الفكرة تنخفض كلفة حساب \(T(x)\) من \(O(x)\) إلى نحو \(O(\sqrt{x})\). وهذا هو التسريع الرئيسي في التطبيق.
بعد جمع كل الأجزاء نحصل على
$$\boxed{F(N)\equiv \sum_{t \le \sqrt{N}} \mu(t)\, t\, T\!\left(\left\lfloor\frac{N}{t^2}\right\rfloor\right)-\frac{N(N+1)}{2}\pmod{10^9+7}.}$$
للتحقق السريع خذ \(N=10\). تقسيم خارج القسمة إلى كتل يعطي
$$T(10)=1\cdot 10+2\cdot 5+3\cdot 3+(4+5)\cdot 2+(6+7+8+9+10)\cdot 1=87.$$
كذلك \(T(2)=4\) و\(T(1)=1\)، والقيم غير الصفرية لدالة Möbius حتى \(\sqrt{10}\) هي \(\mu(1)=1\) و\(\mu(2)=-1\) و\(\mu(3)=-1\). إذًا
$$\sum_{n=1}^{10}\sigma^*(n)=1\cdot 87-2\cdot 4-3\cdot 1=76,$$
ومن ثم
$$F(10)=76-\frac{10\cdot 11}{2}=21.$$
وهذا يطابق الحساب المباشر للقيم \(R(2),\dots,R(10)\)، وبالتالي فالمشتقة الرياضية صحيحة في مثال صغير قبل الانتقال إلى \(10^{14}\).
تبدأ تطبيقات C++ وPython وJava ببناء غربال لدالة Möbius حتى \(\lfloor\sqrt{N}\rfloor\). بعد ذلك تتخطى الحلقة الخارجية كل \(t\) يحقق \(\mu(t)=0\)، لأن مساهمته معدومة. ولكل \(t\) غير صفري تُحسب قيمة واحدة فقط من \(T\!\left(\lfloor N/t^2\rfloor\right)\) باستخدام كتل خارج القسمة بدلًا من مسح خطي كامل.
بعد تجميع مساهمة مجموع القواسم الأحادية، تُطرح القيمة \(N(N+1)/2\) تحت modulo \(10^9+7\). تطبقا C++ وJava توزعان المجموع الخارجي دوريًا على عدة خيوط حتى تتوازن القيم الصغيرة من \(t\) لأنها الأغلى حسابًا؛ أما Python فتطبق الفكرة نفسها بصورة تسلسلية. وفي جميع اللغات يُعاد تطبيع أي ناتج سالب بعد الطرح إلى المجال \(0,\dots,10^9+6\).
غربال Möbius حتى \(\sqrt{N}\) يحتاج إلى \(O(\sqrt{N})\) ذاكرة، ومعالجة تمهيدية شبه خطية في هذا المجال. والعمل الغالب يساوي
$$\sum_{t \le \sqrt{N}} O\!\left(\sqrt{\frac{N}{t^2}}\right)=\sum_{t \le \sqrt{N}} O\!\left(\frac{\sqrt{N}}{t}\right)=O(\sqrt{N}\log N).$$
إذًا الطريقة كلها تعمل بزمن \(O(\sqrt{N}\log N)\) وذاكرة \(O(\sqrt{N})\). التوازي يحسن الزمن العملي في التنفيذات المترجمة، لكنه لا يغير الرتبة التقاربية.