For each integer \(n \gt 1\), consider the 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 if it is idempotent on every residue class:
$$f(f(x))\equiv f(x)\pmod{n}\qquad\text{for all }x.$$
Let \(R(n)\) be the number of such maps. The target of the problem is
$$S(N)=\sum_{k=1}^{N-1} R\left(\binom{N}{k}\right),\qquad N=10^7,$$
and the final result is required modulo \(M=10^9+7\). The main challenge is that the binomial coefficients are enormous, so the implementations never factor them from scratch.
Starting from \(f(x)=ax+b\), we compose once more:
$$f(f(x))\equiv a(ax+b)+b \equiv a^2x+ab+b \pmod{n}.$$
For this to equal \(ax+b\) for every \(x\), the coefficient of \(x\) and the constant correction must both vanish modulo \(n\). Therefore
$$n\mid a(a-1),\qquad n\mid ab.$$
So every retraction is determined by two simultaneous conditions: \(a\) must satisfy an idempotence congruence, and once \(a\) is fixed, \(b\) must solve a linear congruence.
Write the modulus as
$$n=\prod_{i=1}^{m} p_i^{e_i}.$$
Because consecutive integers are coprime, \(\gcd(a,a-1)=1\). Hence if \(p_i^{e_i}\mid a(a-1)\), the full prime power \(p_i^{e_i}\) must divide exactly one of the two factors. For each \(p_i^{e_i}\parallel n\), we must have
$$a\equiv 0 \pmod{p_i^{e_i}}\qquad\text{or}\qquad a\equiv 1 \pmod{p_i^{e_i}}.$$
These choices are independent for distinct prime powers, and the Chinese remainder theorem combines them into exactly one residue class modulo \(n\).
If we define
$$d=\gcd(a,n),$$
then \(d\) is a unitary divisor of \(n\): it satisfies \(d\mid n\) and \(\gcd(d,n/d)=1\). Conversely, every unitary divisor \(d\lt n\) determines one admissible residue \(a\). The excluded case \(d=n\) corresponds to \(a\equiv 0\pmod{n}\), but the problem requires \(0 \lt a \lt n\).
Fix one valid \(a\), and let \(d=\gcd(a,n)\). Write
$$a=d\,a_1,\qquad n=d\,n_1,\qquad \gcd(a_1,n_1)=1.$$
The second condition becomes
$$n\mid ab \iff d\,n_1 \mid d\,a_1 b \iff n_1 \mid a_1 b.$$
Since \(a_1\) is invertible modulo \(n_1\), this is equivalent to
$$n_1\mid b.$$
Among the residues \(0\le b\lt n=d\,n_1\), exactly \(d\) values are multiples of \(n_1\). Therefore each admissible \(a\) contributes exactly \(d=\gcd(a,n)\) choices for \(b\).
Summing over all unitary divisors gives
$$R(n)=\sum_{\substack{d\parallel n\\ d\lt n}} d = \left(\sum_{d\parallel n} d\right)-n.$$
Thus
$$\boxed{R(n)=\sigma^*(n)-n},$$
where \(\sigma^*(n)\) is the sum of unitary divisors. Because each prime power is either chosen completely or not chosen at all,
$$\sigma^*(n)=\sum_{d\parallel n} d=\prod_{p^e\parallel n}(1+p^e).$$
Let
$$B_k=\binom{N}{k}.$$
If
$$B_k=\prod_p p^{e_p(k)},$$
then the closed form above becomes
$$R(B_k)=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)-B_k.$$
So the entire problem reduces to maintaining the prime exponents of \(B_k\) while \(k\) moves across the row of Pascal's triangle.
The standard recurrence
$$B_k=B_{k-1}\cdot \frac{N-k+1}{k}$$
implies, for every prime \(p\),
$$e_p(k)=e_p(k-1)+v_p(N-k+1)-v_p(k).$$
This means the implementations only factor the two ordinary integers \(N-k+1\) and \(k\) at each step. They do not factor \(\binom{N}{k}\) itself.
At the same time they maintain
$$Q_k=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)\pmod{M},$$
and then compute
$$R(B_k)\equiv Q_k-B_k \pmod{M}.$$
Because \(M\) is prime and \(N \lt M\), the recurrence for \(B_k \bmod M\) can also use modular inverses:
$$B_k \equiv B_{k-1}(N-k+1)k^{-1}\pmod{M}.$$
Binomial coefficients satisfy
$$\binom{N}{k}=\binom{N}{N-k},$$
so the sum is symmetric. The implementations only iterate up to \(\lfloor N/2\rfloor\), using weight \(2\) except for the middle term when \(N\) is even.
For \(N=10\), the distinct values are
$$\binom{10}{1}=10,\quad \binom{10}{2}=45,\quad \binom{10}{3}=120,\quad \binom{10}{4}=210,\quad \binom{10}{5}=252.$$
Applying \(R(n)=\sigma^*(n)-n\):
$$\begin{aligned} R(10)&=(1+2)(1+5)-10=8,\\ R(45)&=(1+9)(1+5)-45=15,\\ R(120)&=(1+8)(1+3)(1+5)-120=96,\\ R(210)&=(1+2)(1+3)(1+5)(1+7)-210=366,\\ R(252)&=(1+4)(1+9)(1+7)-252=148. \end{aligned}$$
Hence
$$S(10)=2(8+15+96+366)+148=1118,$$
which is the small checkpoint matched by the implementation.
The C++, Python, and Java implementations first build a smallest-prime-factor table up to \(N\). They also precompute modular inverses \(1^{-1},2^{-1},\dots,N^{-1}\pmod M\) so that the current binomial coefficient can be updated in constant time once the numerator and denominator factors are known.
During the main loop, the implementation stores for each active prime its current exponent in \(B_k\) and the current value of \(p^{e_p(k)} \bmod M\). When a prime exponent changes, the old factor \((1+p^e)\) is removed from the running product \(Q_k\), the exponent is adjusted, and the new factor is inserted.
Division inside the modulus is handled by multiplying with inverses. The implementations cache inverses of previously seen nonzero factors, and they also keep track of how many factors are currently \(0 \pmod M\). If at least one such factor exists, then \(Q_k\equiv 0\pmod M\), so the running product remains correct even in that edge case.
The smallest-prime-factor preprocessing uses \(O(N)\) memory. The total factor work over all updates is governed by the number of prime factors, counted with multiplicity, appearing in the integers \(1,2,\dots,N\), which is near \(O(N\log\log N)\) in aggregate. The remaining arithmetic per step is constant. Therefore the full method runs in near \(O(N\log\log N)\) time and \(O(N)\) memory, which is practical for \(N=10^7\).
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 allen Restklassen idempotent ist:
$$f(f(x))\equiv f(x)\pmod{n}\qquad\text{für alle }x.$$
Sei \(R(n)\) die Anzahl dieser Abbildungen. Gesucht ist
$$S(N)=\sum_{k=1}^{N-1} R\left(\binom{N}{k}\right),\qquad N=10^7,$$
wobei das Endergebnis modulo \(M=10^9+7\) berechnet werden soll. Die binomialen Koeffizienten sind viel zu groß für direkte Faktorisierung, daher wird ihre Primfaktorstruktur nur inkrementell gepflegt.
Aus \(f(x)=ax+b\) folgt durch erneute Anwendung
$$f(f(x))\equiv a(ax+b)+b \equiv a^2x+ab+b \pmod{n}.$$
Damit dies für alle \(x\) mit \(ax+b\) übereinstimmt, müssen sowohl der lineare als auch der konstante Korrekturterm modulo \(n\) verschwinden. Also gilt
$$n\mid a(a-1),\qquad n\mid ab.$$
Jede Retraktion wird somit durch zwei Bedingungen charakterisiert: \(a\) muss eine idempotente Kongruenz lösen, und für dieses \(a\) muss \(b\) eine lineare Kongruenz erfüllen.
Schreibe
$$n=\prod_{i=1}^{m} p_i^{e_i}.$$
Da aufeinanderfolgende ganze Zahlen teilerfremd sind, gilt \(\gcd(a,a-1)=1\). Wenn also \(p_i^{e_i}\mid a(a-1)\), muss die gesamte Primzahlpotenz \(p_i^{e_i}\) genau einen der beiden Faktoren teilen. Für jedes \(p_i^{e_i}\parallel n\) folgt daher
$$a\equiv 0 \pmod{p_i^{e_i}}\qquad\text{oder}\qquad a\equiv 1 \pmod{p_i^{e_i}}.$$
Die Entscheidungen für verschiedene Primzahlpotenzen sind unabhängig, und der chinesische Restsatz setzt sie zu genau einer Restklasse modulo \(n\) zusammen.
Definiert man
$$d=\gcd(a,n),$$
dann ist \(d\) ein unitärer Teiler von \(n\), also \(d\mid n\) und \(\gcd(d,n/d)=1\). Umgekehrt liefert jeder unitäre Teiler \(d\lt n\) genau eine zulässige Restklasse \(a\). Der ausgeschlossene Fall \(d=n\) entspräche \(a\equiv 0\pmod n\), was wegen \(0 \lt a \lt n\) nicht erlaubt ist.
Fixiere ein gültiges \(a\) und setze \(d=\gcd(a,n)\). Schreibe
$$a=d\,a_1,\qquad n=d\,n_1,\qquad \gcd(a_1,n_1)=1.$$
Dann wird aus der zweiten Bedingung
$$n\mid ab \iff d\,n_1 \mid d\,a_1 b \iff n_1 \mid a_1 b.$$
Weil \(a_1\) modulo \(n_1\) invertierbar ist, ist dies äquivalent zu
$$n_1\mid b.$$
Unter den Resten \(0\le b\lt n=d\,n_1\) gibt es genau \(d\) Vielfache von \(n_1\). Also trägt jedes zulässige \(a\) genau \(d=\gcd(a,n)\) Werte von \(b\) bei.
Summiert über alle unitären Teiler ergibt das
$$R(n)=\sum_{\substack{d\parallel n\\ d\lt n}} d = \left(\sum_{d\parallel n} d\right)-n.$$
Damit erhält man die geschlossene Formel
$$\boxed{R(n)=\sigma^*(n)-n},$$
wobei \(\sigma^*(n)\) die Summe der unitären Teiler ist. Da jede Primzahlpotenz entweder vollständig gewählt wird oder gar nicht, gilt
$$\sigma^*(n)=\sum_{d\parallel n} d=\prod_{p^e\parallel n}(1+p^e).$$
Setze
$$B_k=\binom{N}{k}.$$
Mit der Primfaktorzerlegung
$$B_k=\prod_p p^{e_p(k)}$$
wird die Formel zu
$$R(B_k)=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)-B_k.$$
Damit reduziert sich die Aufgabe darauf, die Primexponenten von \(B_k\) beim Durchlaufen der Pascal-Zeile effizient zu aktualisieren.
Aus der Rekursion
$$B_k=B_{k-1}\cdot \frac{N-k+1}{k}$$
folgt für jede Primzahl \(p\)
$$e_p(k)=e_p(k-1)+v_p(N-k+1)-v_p(k).$$
Die Implementierungen faktorisieren daher in jedem Schritt nur die beiden gewöhnlichen Zahlen \(N-k+1\) und \(k\), nicht aber den riesigen Binomialkoeffizienten selbst.
Zusätzlich wird
$$Q_k=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)\pmod{M}$$
fortlaufend gepflegt. Danach ist
$$R(B_k)\equiv Q_k-B_k \pmod{M}.$$
Da \(M\) prim und \(N \lt M\) ist, kann auch \(B_k \bmod M\) selbst über modulare Inversen aktualisiert werden:
$$B_k \equiv B_{k-1}(N-k+1)k^{-1}\pmod{M}.$$
Wegen
$$\binom{N}{k}=\binom{N}{N-k}$$
ist die Summe symmetrisch. Die Implementierungen laufen deshalb nur bis \(\lfloor N/2\rfloor\) und verwenden Gewicht \(2\), außer beim Mittelterm für gerades \(N\).
Für \(N=10\) lauten die verschiedenen Werte
$$\binom{10}{1}=10,\quad \binom{10}{2}=45,\quad \binom{10}{3}=120,\quad \binom{10}{4}=210,\quad \binom{10}{5}=252.$$
Mit \(R(n)=\sigma^*(n)-n\) erhält man
$$\begin{aligned} R(10)&=(1+2)(1+5)-10=8,\\ R(45)&=(1+9)(1+5)-45=15,\\ R(120)&=(1+8)(1+3)(1+5)-120=96,\\ R(210)&=(1+2)(1+3)(1+5)(1+7)-210=366,\\ R(252)&=(1+4)(1+9)(1+7)-252=148. \end{aligned}$$
Daher
$$S(10)=2(8+15+96+366)+148=1118,$$
genau der kleine Kontrollwert, den die Implementierung reproduziert.
Die Implementierungen in C++, Python und Java bauen zunächst eine Tabelle des kleinsten Primteilers bis \(N\) auf. Außerdem werden die modularen Inversen \(1^{-1},2^{-1},\dots,N^{-1}\pmod M\) vorab berechnet, damit sich der aktuelle Binomialkoeffizient nach Kenntnis von Zähler- und Nennerfaktoren in konstanter Zeit fortschreiben lässt.
In der Hauptschleife speichert die Implementierung für jede aktive Primzahl ihren aktuellen Exponenten in \(B_k\) sowie den Wert \(p^{e_p(k)} \bmod M\). Wenn sich ein Exponent ändert, wird der alte Faktor \((1+p^e)\) aus dem laufenden Produkt \(Q_k\) entfernt, der Exponent angepasst und der neue Faktor eingefügt.
Die Division im Modul erfolgt über Multiplikation mit modularen Inversen. Nichtverschwindende Faktoren werden bei Bedarf invertiert und zwischengespeichert. Zusätzlich zählt die Implementierung, wie viele Faktoren momentan \(0 \pmod M\) sind. Falls mindestens einer vorhanden ist, gilt automatisch \(Q_k\equiv 0\pmod M\), sodass auch dieser Spezialfall korrekt behandelt wird.
Die SPF-Vorverarbeitung benötigt \(O(N)\) Speicher. Die gesamte Faktorisierungsarbeit über alle Updates wird durch die Anzahl der Primfaktoren mit Vielfachheit in den Zahlen \(1,2,\dots,N\) bestimmt und liegt insgesamt nahe bei \(O(N\log\log N)\). Pro Schritt bleibt nur konstante Zusatzarithmetik. Insgesamt arbeitet das Verfahren daher in nahe \(O(N\log\log N)\) Zeit und \(O(N)\) Speicher, was für \(N=10^7\) praktikabel ist.
Her \(n \gt 1\) için şu doğrusal fonksiyonları ele alalım:
$$f(x)\equiv ax+b \pmod{n},\qquad 0 \lt a \lt n,\quad 0\le b \lt n.$$
Fonksiyon tüm kalan sınıflarında idempotent ise bir geri çekimtir:
$$f(f(x))\equiv f(x)\pmod{n}\qquad\text{tüm }x\text{ için}.$$
\(R(n)\), bu koşulu sağlayan fonksiyonların sayısı olsun. Problemde istenen büyüklük
$$S(N)=\sum_{k=1}^{N-1} R\left(\binom{N}{k}\right),\qquad N=10^7,$$
ve sonuç \(M=10^9+7\) modunda alınır. Asıl zorluk, binom katsayılarının çok büyük olmasıdır; bu yüzden uygulamalar onları baştan çarpanlara ayırmaz, yalnızca asal üslerini adım adım günceller.
\(f(x)=ax+b\) için bir kez daha bileşke alırsak
$$f(f(x))\equiv a(ax+b)+b \equiv a^2x+ab+b \pmod{n}$$
olur. Bunun her \(x\) için \(ax+b\)'ye eşit olması için hem \(x\)'in katsayısı hem de sabit düzeltme modulo \(n\) sıfır olmalıdır. Dolayısıyla
$$n\mid a(a-1),\qquad n\mid ab.$$
Yani bir geri çekimi saymak için önce \(a\)'nın idempotent kongruansı sağlaması, sonra da bu \(a\) için \(b\)'nin doğrusal kongruansı çözmesi gerekir.
\(n\)'yi
$$n=\prod_{i=1}^{m} p_i^{e_i}$$
şeklinde yazalım. Ardışık iki sayı aralarında asal olduğundan \(\gcd(a,a-1)=1\). O halde \(p_i^{e_i}\mid a(a-1)\) ise bu asal kuvvetin tamamı iki çarpandan tam olarak birini bölmek zorundadır. Her \(p_i^{e_i}\parallel n\) için
$$a\equiv 0 \pmod{p_i^{e_i}}\qquad\text{veya}\qquad a\equiv 1 \pmod{p_i^{e_i}}$$
elde edilir. Farklı asal kuvvetler için seçimler bağımsızdır ve Çin kalan teoremi bunları modulo \(n\) altında tek bir kalanda birleştirir.
Şimdi
$$d=\gcd(a,n)$$
olsun. Bu \(d\), \(n\)'nin bir üniter bölenidir; yani \(d\mid n\) ve \(\gcd(d,n/d)=1\). Tersine, \(d\lt n\) olan her üniter bölen tam bir uygun \(a\) kalanı üretir. Hariç bırakılan \(d=n\) durumu \(a\equiv 0\pmod n\) demektir; oysa problem \(0 \lt a \lt n\) istiyor.
Geçerli bir \(a\) sabitleyelim ve \(d=\gcd(a,n)\) diyelim. Şöyle yazalım:
$$a=d\,a_1,\qquad n=d\,n_1,\qquad \gcd(a_1,n_1)=1.$$
İkinci koşul
$$n\mid ab \iff d\,n_1 \mid d\,a_1 b \iff n_1 \mid a_1 b$$
olur. \(a_1\), \(n_1\) modunda tersinir olduğu için bu koşul
$$n_1\mid b$$
ile eşdeğerdir. \(0\le b\lt n=d\,n_1\) aralığında \(n_1\)'in tam olarak \(d\) katı vardır. Demek ki her uygun \(a\), tam \(d=\gcd(a,n)\) farklı \(b\) verir.
Böylece tüm üniter bölenler üzerinden toplayınca
$$R(n)=\sum_{\substack{d\parallel n\\ d\lt n}} d = \left(\sum_{d\parallel n} d\right)-n$$
olur. Sonuç olarak
$$\boxed{R(n)=\sigma^*(n)-n},$$
burada \(\sigma^*(n)\) üniter bölenlerin toplamıdır. Her asal kuvvet ya tamamen seçilir ya da hiç seçilmez; dolayısıyla
$$\sigma^*(n)=\sum_{d\parallel n} d=\prod_{p^e\parallel n}(1+p^e).$$
$$B_k=\binom{N}{k}$$
olsun. Eğer
$$B_k=\prod_p p^{e_p(k)}$$
ise kapalı form
$$R(B_k)=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)-B_k$$
şeklini alır. Böylece sorun, Pascal üçgeninin tek bir satırında ilerlerken \(B_k\)'lerin asal üslerini verimli biçimde güncel tutmaya dönüşür.
Standart bağıntı
$$B_k=B_{k-1}\cdot \frac{N-k+1}{k}$$
her asal \(p\) için
$$e_p(k)=e_p(k-1)+v_p(N-k+1)-v_p(k)$$
sonucunu verir. Bu yüzden uygulamalar her adımda yalnızca \(N-k+1\) ve \(k\) sayılarını çarpanlara ayırır; dev binom katsayısını değil.
Aynı anda
$$Q_k=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)\pmod{M}$$
ürünü de tutulur ve sonra
$$R(B_k)\equiv Q_k-B_k \pmod{M}$$
hesaplanır. \(M\) asal ve \(N \lt M\) olduğundan, \(B_k \bmod M\) değeri de modüler terslerle doğrudan güncellenebilir:
$$B_k \equiv B_{k-1}(N-k+1)k^{-1}\pmod{M}.$$
Binom katsayıları için
$$\binom{N}{k}=\binom{N}{N-k}$$
olduğu için toplam simetriktir. Uygulamalar yalnızca \(\lfloor N/2\rfloor\)'ye kadar gider; \(N\) çiftse orta terim hariç her katkı 2 ile ağırlıklandırılır.
\(N=10\) için farklı değerler
$$\binom{10}{1}=10,\quad \binom{10}{2}=45,\quad \binom{10}{3}=120,\quad \binom{10}{4}=210,\quad \binom{10}{5}=252$$
olur. \(R(n)=\sigma^*(n)-n\) ile
$$\begin{aligned} R(10)&=(1+2)(1+5)-10=8,\\ R(45)&=(1+9)(1+5)-45=15,\\ R(120)&=(1+8)(1+3)(1+5)-120=96,\\ R(210)&=(1+2)(1+3)(1+5)(1+7)-210=366,\\ R(252)&=(1+4)(1+9)(1+7)-252=148. \end{aligned}$$
Buradan
$$S(10)=2(8+15+96+366)+148=1118$$
elde edilir; bu da uygulamanın küçük kontrol değeridir.
C++, Python ve Java uygulamaları önce \(N\)'ye kadar en küçük asal bölen tablosu kurar. Ayrıca \(1^{-1},2^{-1},\dots,N^{-1}\pmod M\) modüler tersleri önceden hazırlanır; böylece mevcut binom katsayısı, pay ve payda çarpanları bilindiğinde sabit zamanda güncellenebilir.
Ana döngüde her etkin asal için \(B_k\) içindeki güncel üs ve \(p^{e_p(k)} \bmod M\) değeri saklanır. Bir asalın üssü değiştiğinde, \((1+p^e)\) çarpanının eski hali çalışan çarpımdan çıkarılır, üs düzeltilir ve yeni hali eklenir.
Mod altında bölme, modüler ters ile çarpma olarak yapılır. Sıfır olmayan faktörlerin tersleri ilk gerektiğinde önbelleğe alınır. Ayrıca uygulamalar şu anda \(0 \pmod M\) olan kaç çarpan bulunduğunu ayrı sayar. En az bir tane varsa doğrudan \(Q_k\equiv 0\pmod M\) olur; böylece bu köşe durum da doğru yönetilir.
SPF ön işleme \(O(N)\) bellek kullanır. Tüm güncellemeler boyunca yapılan asal çarpan işi, \(1,2,\dots,N\) sayılarındaki asal çarpanların çoklukla sayısına bağlıdır ve toplamda yaklaşık \(O(N\log\log N)\) mertebesindedir. Adım başına kalan aritmetik sabittir. Sonuçta yöntem yaklaşık \(O(N\log\log N)\) zamanda ve \(O(N)\) bellekte çalışır; bu da \(N=10^7\) için uygundur.
Para cada \(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 si es idempotente sobre todas las clases residuales:
$$f(f(x))\equiv f(x)\pmod{n}\qquad\text{para todo }x.$$
Sea \(R(n)\) el número de tales aplicaciones. El objetivo del problema es
$$S(N)=\sum_{k=1}^{N-1} R\left(\binom{N}{k}\right),\qquad N=10^7,$$
y el resultado final debe darse módulo \(M=10^9+7\). La dificultad central es que los coeficientes binomiales son enormes, así que la implementación nunca los factoriza desde cero.
Si \(f(x)=ax+b\), entonces
$$f(f(x))\equiv a(ax+b)+b \equiv a^2x+ab+b \pmod{n}.$$
Para que esto coincida con \(ax+b\) para todo \(x\), deben anularse módulo \(n\) tanto el coeficiente de \(x\) como la corrección constante. Por tanto,
$$n\mid a(a-1),\qquad n\mid ab.$$
Así, toda retracción queda determinada por dos condiciones simultáneas: \(a\) debe satisfacer una congruencia idempotente y, una vez fijado \(a\), \(b\) debe resolver una congruencia lineal.
Escribamos
$$n=\prod_{i=1}^{m} p_i^{e_i}.$$
Como enteros consecutivos son coprimos, \(\gcd(a,a-1)=1\). Entonces, si \(p_i^{e_i}\mid a(a-1)\), toda la potencia prima \(p_i^{e_i}\) debe dividir exactamente a uno de los factores. Para cada \(p_i^{e_i}\parallel n\), necesariamente
$$a\equiv 0 \pmod{p_i^{e_i}}\qquad\text{o}\qquad a\equiv 1 \pmod{p_i^{e_i}}.$$
Las elecciones para distintas potencias primas son independientes, y el teorema chino del resto las combina en una única clase residual módulo \(n\).
Si definimos
$$d=\gcd(a,n),$$
entonces \(d\) es un divisor unitario de \(n\): satisface \(d\mid n\) y \(\gcd(d,n/d)=1\). A la inversa, cada divisor unitario \(d\lt n\) produce exactamente una clase admisible de \(a\). El caso excluido \(d=n\) corresponde a \(a\equiv 0\pmod n\), pero el problema exige \(0 \lt a \lt n\).
Fijemos un \(a\) válido y sea \(d=\gcd(a,n)\). Escribimos
$$a=d\,a_1,\qquad n=d\,n_1,\qquad \gcd(a_1,n_1)=1.$$
La segunda condición se transforma en
$$n\mid ab \iff d\,n_1 \mid d\,a_1 b \iff n_1 \mid a_1 b.$$
Como \(a_1\) es invertible módulo \(n_1\), esto equivale a
$$n_1\mid b.$$
Entre los residuos \(0\le b\lt n=d\,n_1\), hay exactamente \(d\) múltiplos de \(n_1\). Por lo tanto, cada \(a\) admisible aporta exactamente \(d=\gcd(a,n)\) valores posibles de \(b\).
Sumando sobre todos los divisores unitarios obtenemos
$$R(n)=\sum_{\substack{d\parallel n\\ d\lt n}} d = \left(\sum_{d\parallel n} d\right)-n.$$
En consecuencia,
$$\boxed{R(n)=\sigma^*(n)-n},$$
donde \(\sigma^*(n)\) es la suma de divisores unitarios. Como cada potencia prima se toma completa o no se toma,
$$\sigma^*(n)=\sum_{d\parallel n} d=\prod_{p^e\parallel n}(1+p^e).$$
Sea
$$B_k=\binom{N}{k}.$$
Si
$$B_k=\prod_p p^{e_p(k)},$$
entonces la fórmula cerrada se convierte en
$$R(B_k)=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)-B_k.$$
Por tanto, el problema completo se reduce a mantener actualizados los exponentes primos de \(B_k\) mientras recorremos una fila del triángulo de Pascal.
La recurrencia clásica
$$B_k=B_{k-1}\cdot \frac{N-k+1}{k}$$
implica, para cada primo \(p\),
$$e_p(k)=e_p(k-1)+v_p(N-k+1)-v_p(k).$$
Así, las implementaciones solo factorizan en cada paso los enteros corrientes \(N-k+1\) y \(k\). Nunca factorizan directamente el enorme coeficiente binomial.
Al mismo tiempo mantienen
$$Q_k=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)\pmod{M},$$
y después calculan
$$R(B_k)\equiv Q_k-B_k \pmod{M}.$$
Como \(M\) es primo y \(N \lt M\), también se puede actualizar \(B_k \bmod M\) mediante inversos modulares:
$$B_k \equiv B_{k-1}(N-k+1)k^{-1}\pmod{M}.$$
Los coeficientes binomiales cumplen
$$\binom{N}{k}=\binom{N}{N-k},$$
de modo que la suma es simétrica. Las implementaciones solo recorren hasta \(\lfloor N/2\rfloor\), con peso \(2\) salvo para el término central cuando \(N\) es par.
Para \(N=10\), los valores distintos son
$$\binom{10}{1}=10,\quad \binom{10}{2}=45,\quad \binom{10}{3}=120,\quad \binom{10}{4}=210,\quad \binom{10}{5}=252.$$
Aplicando \(R(n)=\sigma^*(n)-n\):
$$\begin{aligned} R(10)&=(1+2)(1+5)-10=8,\\ R(45)&=(1+9)(1+5)-45=15,\\ R(120)&=(1+8)(1+3)(1+5)-120=96,\\ R(210)&=(1+2)(1+3)(1+5)(1+7)-210=366,\\ R(252)&=(1+4)(1+9)(1+7)-252=148. \end{aligned}$$
Así,
$$S(10)=2(8+15+96+366)+148=1118,$$
que coincide con la comprobación pequeña usada por la implementación.
Las implementaciones en C++, Python y Java construyen primero una tabla del menor factor primo hasta \(N\). También precalculan los inversos modulares \(1^{-1},2^{-1},\dots,N^{-1}\pmod M\), de modo que el coeficiente binomial actual puede actualizarse en tiempo constante una vez conocidos sus factores de numerador y denominador.
En el bucle principal, la implementación guarda para cada primo activo su exponente actual en \(B_k\) y el valor actual de \(p^{e_p(k)} \bmod M\). Cuando cambia un exponente, se elimina del producto acumulado \(Q_k\) el factor antiguo \((1+p^e)\), se ajusta el exponente y se inserta el factor nuevo.
La división en el módulo se realiza multiplicando por inversos. Los inversos de factores no nulos se almacenan en caché cuando aparecen por primera vez. Además, la implementación lleva la cuenta de cuántos factores valen actualmente \(0 \pmod M\). Si existe al menos uno, entonces \(Q_k\equiv 0\pmod M\), y el producto acumulado sigue siendo correcto incluso en ese caso extremo.
La precomputación SPF usa \(O(N)\) memoria. El trabajo total de factorización sobre todas las actualizaciones está controlado por el número de factores primos, con multiplicidad, que aparecen en los enteros \(1,2,\dots,N\), lo cual en conjunto es cercano a \(O(N\log\log N)\). El resto de la aritmética por paso es constante. En consecuencia, el método completo funciona en tiempo cercano a \(O(N\log\log N)\) y memoria \(O(N)\), suficiente para \(N=10^7\).
对每个 \(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\text{对所有 }x\text{ 都成立},$$
那么称它为一个回缩。记 \(R(n)\) 为这样的映射个数。题目要求计算
$$S(N)=\sum_{k=1}^{N-1} R\left(\binom{N}{k}\right),\qquad N=10^7,$$
并把最终结果对 \(M=10^9+7\) 取模。真正困难的地方在于二项式系数极大，因此实现不可能在每一步重新分解 \(\binom{N}{k}\)。
由 \(f(x)=ax+b\) 再复合一次可得
$$f(f(x))\equiv a(ax+b)+b \equiv a^2x+ab+b \pmod{n}.$$
若它要对每个 \(x\) 都等于 \(ax+b\)，那么 \(x\) 的系数和常数修正项都必须在模 \(n\) 下消失。因此有
$$n\mid a(a-1),\qquad n\mid ab.$$
也就是说，计数问题被拆成两部分：先找满足幂等同余的 \(a\)，再对每个这样的 \(a\) 统计满足线性同余的 \(b\)。
把 \(n\) 写成
$$n=\prod_{i=1}^{m} p_i^{e_i}.$$
因为相邻整数互素，\(\gcd(a,a-1)=1\)。所以如果 \(p_i^{e_i}\mid a(a-1)\)，那么整个素数幂 \(p_i^{e_i}\) 只能完整地落在两个因子中的一个上。对每个 \(p_i^{e_i}\parallel n\)，必须满足
$$a\equiv 0 \pmod{p_i^{e_i}}\qquad\text{或}\qquad a\equiv 1 \pmod{p_i^{e_i}}.$$
不同素数幂上的选择彼此独立，再由中国剩余定理拼成模 \(n\) 下唯一的一个剩余类。
若定义
$$d=\gcd(a,n),$$
那么 \(d\) 正好是 \(n\) 的一个单位因子，也就是 \(d\mid n\) 且 \(\gcd(d,n/d)=1\)。反过来，每个满足 \(d\lt n\) 的单位因子都对应一个合法的 \(a\)。被排除的 \(d=n\) 对应 \(a\equiv 0\pmod n\)，但题目要求 \(0 \lt a \lt n\)。
固定一个合法的 \(a\)，令 \(d=\gcd(a,n)\)。写成
$$a=d\,a_1,\qquad n=d\,n_1,\qquad \gcd(a_1,n_1)=1.$$
第二个条件变为
$$n\mid ab \iff d\,n_1 \mid d\,a_1 b \iff n_1 \mid a_1 b.$$
由于 \(a_1\) 在模 \(n_1\) 下可逆，这等价于
$$n_1\mid b.$$
在 \(0\le b\lt n=d\,n_1\) 这些剩余中，恰好有 \(d\) 个是 \(n_1\) 的倍数。因此，每个合法的 \(a\) 恰好贡献 \(d=\gcd(a,n)\) 个 \(b\)。
对所有单位因子求和，就得到
$$R(n)=\sum_{\substack{d\parallel n\\ d\lt n}} d = \left(\sum_{d\parallel n} d\right)-n.$$
于是有闭式
$$\boxed{R(n)=\sigma^*(n)-n},$$
其中 \(\sigma^*(n)\) 是单位因子和。因为每个素数幂要么整块选中，要么完全不选，所以
$$\sigma^*(n)=\sum_{d\parallel n} d=\prod_{p^e\parallel n}(1+p^e).$$
记
$$B_k=\binom{N}{k}.$$
若
$$B_k=\prod_p p^{e_p(k)},$$
那么上面的闭式就变成
$$R(B_k)=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)-B_k.$$
因此，整道题等价于沿着帕斯卡三角形的一行移动时，持续维护 \(B_k\) 的素因子指数。
标准递推式
$$B_k=B_{k-1}\cdot \frac{N-k+1}{k}$$
对每个素数 \(p\) 都给出
$$e_p(k)=e_p(k-1)+v_p(N-k+1)-v_p(k).$$
这说明实现中每一步只需分解两个普通整数 \(N-k+1\) 和 \(k\)，而不必直接分解巨大的 \(\binom{N}{k}\)。
同时维护
$$Q_k=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)\pmod{M},$$
然后计算
$$R(B_k)\equiv Q_k-B_k \pmod{M}.$$
由于 \(M\) 是素数且 \(N \lt M\)，\(B_k \bmod M\) 也可以用模逆递推：
$$B_k \equiv B_{k-1}(N-k+1)k^{-1}\pmod{M}.$$
二项式系数满足
$$\binom{N}{k}=\binom{N}{N-k},$$
所以总和是对称的。实现只遍历到 \(\lfloor N/2\rfloor\)，除了 \(N\) 为偶数时的中间项外，其余贡献都乘以 \(2\)。
当 \(N=10\) 时，不同的值为
$$\binom{10}{1}=10,\quad \binom{10}{2}=45,\quad \binom{10}{3}=120,\quad \binom{10}{4}=210,\quad \binom{10}{5}=252.$$
代入 \(R(n)=\sigma^*(n)-n\)：
$$\begin{aligned} R(10)&=(1+2)(1+5)-10=8,\\ R(45)&=(1+9)(1+5)-45=15,\\ R(120)&=(1+8)(1+3)(1+5)-120=96,\\ R(210)&=(1+2)(1+3)(1+5)(1+7)-210=366,\\ R(252)&=(1+4)(1+9)(1+7)-252=148. \end{aligned}$$
因此
$$S(10)=2(8+15+96+366)+148=1118,$$
这正是实现使用的小规模校验值。
C++、Python 和 Java 实现都会先预处理到 \(N\) 的最小素因子表。同时还会预计算 \(1^{-1},2^{-1},\dots,N^{-1}\pmod M\)，这样在知道分子和分母新增了哪些素因子之后，就能用常数次模运算更新当前的二项式系数。
在主循环中，实现会为每个当前活跃的素数保存它在 \(B_k\) 中的指数，以及 \(p^{e_p(k)} \bmod M\) 的当前值。当某个素数的指数发生变化时，先把旧的 \((1+p^e)\) 从乘积 \(Q_k\) 中移除，再更新指数，最后插入新的因子。
模意义下的除法通过乘以模逆完成。实现会缓存已经见过的非零因子的逆元。另外，还会单独记录当前有多少个因子等于 \(0 \pmod M\)。只要这个计数大于零，就有 \(Q_k\equiv 0\pmod M\)，从而保证这种边界情况也能被正确处理。
SPF 预处理需要 \(O(N)\) 空间。所有更新中的分解总工作量由 \(1,2,\dots,N\) 中素因子个数（按重数计）控制，整体接近 \(O(N\log\log N)\)。每一步剩余的模运算都是常数级。因此总复杂度为近似 \(O(N\log\log N)\) 时间和 \(O(N)\) 空间，对 \(N=10^7\) 是可行的。
Для каждого \(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\text{для всех }x.$$
Пусть \(R(n)\) обозначает число таких отображений. Требуется вычислить
$$S(N)=\sum_{k=1}^{N-1} R\left(\binom{N}{k}\right),\qquad N=10^7,$$
а затем взять ответ по модулю \(M=10^9+7\). Главная трудность в том, что биномиальные коэффициенты огромны, поэтому реализация не раскладывает каждый из них заново.
Если \(f(x)=ax+b\), то после ещё одной подстановки получаем
$$f(f(x))\equiv a(ax+b)+b \equiv a^2x+ab+b \pmod{n}.$$
Чтобы это совпадало с \(ax+b\) для любого \(x\), должны обращаться в нуль по модулю \(n\) и коэффициент при \(x\), и постоянная поправка. Значит,
$$n\mid a(a-1),\qquad n\mid ab.$$
Итак, задача распадается на две части: сначала нужно найти все подходящие \(a\), удовлетворяющие идемпотентной конгруэнции, а затем для каждого такого \(a\) посчитать допустимые \(b\).
Пусть
$$n=\prod_{i=1}^{m} p_i^{e_i}.$$
Так как соседние целые числа взаимно просты, \(\gcd(a,a-1)=1\). Поэтому если \(p_i^{e_i}\mid a(a-1)\), то вся степень \(p_i^{e_i}\) должна целиком делить ровно один из множителей. Для каждого \(p_i^{e_i}\parallel n\) необходимо
$$a\equiv 0 \pmod{p_i^{e_i}}\qquad\text{или}\qquad a\equiv 1 \pmod{p_i^{e_i}}.$$
Выборы для различных простых степеней независимы, а китайская теорема об остатках склеивает их в единственный класс вычетов по модулю \(n\).
Если положить
$$d=\gcd(a,n),$$
то \(d\) оказывается унитарным делителем числа \(n\): выполняются \(d\mid n\) и \(\gcd(d,n/d)=1\). Обратно, каждый унитарный делитель \(d\lt n\) задаёт ровно один допустимый остаток \(a\). Исключённый случай \(d=n\) соответствует \(a\equiv 0\pmod n\), а по условию нужно \(0 \lt a \lt n\).
Зафиксируем корректное \(a\) и обозначим \(d=\gcd(a,n)\). Тогда
$$a=d\,a_1,\qquad n=d\,n_1,\qquad \gcd(a_1,n_1)=1.$$
Второе условие переписывается так:
$$n\mid ab \iff d\,n_1 \mid d\,a_1 b \iff n_1 \mid a_1 b.$$
Поскольку \(a_1\) обратим по модулю \(n_1\), это эквивалентно
$$n_1\mid b.$$
Среди вычетов \(0\le b\lt n=d\,n_1\) ровно \(d\) кратны \(n_1\). Значит, каждое допустимое \(a\) даёт ровно \(d=\gcd(a,n)\) возможных значений \(b\).
Суммируя по всем унитарным делителям, получаем
$$R(n)=\sum_{\substack{d\parallel n\\ d\lt n}} d = \left(\sum_{d\parallel n} d\right)-n.$$
Следовательно,
$$\boxed{R(n)=\sigma^*(n)-n},$$
где \(\sigma^*(n)\) обозначает сумму унитарных делителей. Поскольку каждая простая степень либо выбирается целиком, либо не выбирается вовсе, имеем
$$\sigma^*(n)=\sum_{d\parallel n} d=\prod_{p^e\parallel n}(1+p^e).$$
Обозначим
$$B_k=\binom{N}{k}.$$
Если
$$B_k=\prod_p p^{e_p(k)},$$
то закрытая формула принимает вид
$$R(B_k)=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)-B_k.$$
Значит, вся задача сводится к поддержанию простых показателей \(B_k\) при движении вдоль одной строки треугольника Паскаля.
Стандартная рекурсия
$$B_k=B_{k-1}\cdot \frac{N-k+1}{k}$$
даёт для любого простого \(p\)
$$e_p(k)=e_p(k-1)+v_p(N-k+1)-v_p(k).$$
Поэтому реализация в каждом шаге раскладывает только два обычных числа \(N-k+1\) и \(k\), а не сам огромный биномиальный коэффициент.
Параллельно поддерживается
$$Q_k=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)\pmod{M},$$
после чего вычисляется
$$R(B_k)\equiv Q_k-B_k \pmod{M}.$$
Так как \(M\) простое и \(N \lt M\), сам биномиальный коэффициент по модулю тоже обновляется через обратные элементы:
$$B_k \equiv B_{k-1}(N-k+1)k^{-1}\pmod{M}.$$
Биномиальные коэффициенты удовлетворяют равенству
$$\binom{N}{k}=\binom{N}{N-k},$$
поэтому сумма симметрична. Реализация проходит только до \(\lfloor N/2\rfloor\), а вклад берётся с весом \(2\), кроме центрального члена при чётном \(N\).
Для \(N=10\) различные значения таковы:
$$\binom{10}{1}=10,\quad \binom{10}{2}=45,\quad \binom{10}{3}=120,\quad \binom{10}{4}=210,\quad \binom{10}{5}=252.$$
Подставляя \(R(n)=\sigma^*(n)-n\), получаем
$$\begin{aligned} R(10)&=(1+2)(1+5)-10=8,\\ R(45)&=(1+9)(1+5)-45=15,\\ R(120)&=(1+8)(1+3)(1+5)-120=96,\\ R(210)&=(1+2)(1+3)(1+5)(1+7)-210=366,\\ R(252)&=(1+4)(1+9)(1+7)-252=148. \end{aligned}$$
Следовательно,
$$S(10)=2(8+15+96+366)+148=1118,$$
и это совпадает с малой контрольной проверкой реализации.
Реализации на C++, Python и Java сначала строят таблицу наименьшего простого делителя до \(N\). Также заранее вычисляются обратные по модулю числа \(1^{-1},2^{-1},\dots,N^{-1}\pmod M\), чтобы текущий биномиальный коэффициент можно было обновлять за константное время после обработки числителя и знаменателя.
В основном цикле для каждого активного простого числа хранятся текущий показатель в \(B_k\) и значение \(p^{e_p(k)} \bmod M\). Когда показатель меняется, старый множитель \((1+p^e)\) убирается из накопленного произведения \(Q_k\), показатель корректируется, а затем вставляется новый множитель.
Деление по модулю реализуется умножением на обратный элемент. Обратные к ненулевым множителям кэшируются при первом использовании. Кроме того, реализация отдельно считает, сколько множителей сейчас равны \(0 \pmod M\). Если такой множитель хотя бы один, то сразу \(Q_k\equiv 0\pmod M\), и произведение остаётся корректным даже в этом пограничном случае.
Предобработка SPF требует \(O(N)\) памяти. Общая стоимость всех факторизационных обновлений определяется числом простых множителей с кратностями в числах \(1,2,\dots,N\) и в сумме близка к \(O(N\log\log N)\). Остальная арифметика на шаге константна. В итоге получаем близкое к \(O(N\log\log N)\) время и \(O(N)\) память, что достаточно для \(N=10^7\).
لكل عدد صحيح \(n \gt 1\)، ننظر إلى الدوال الخطية
$$f(x)\equiv ax+b \pmod{n},\qquad 0 \lt a \lt n,\quad 0\le b \lt n.$$
وتكون الدالة ارتدادًا إذا كانت idempotent على جميع فئات البواقي، أي إذا تحقق الشرط الآتي لكل \(x\):
$$f(f(x))\equiv f(x)\pmod{n}.$$
ولنرمز بعدد هذه الدوال بـ \(R(n)\). الكمية المطلوبة في المسألة هي
$$S(N)=\sum_{k=1}^{N-1} R\left(\binom{N}{k}\right),\qquad N=10^7,$$
ثم يؤخذ الناتج النهائي بترديد \(M=10^9+7\). الصعوبة الأساسية هي أن معاملات ثنائية الحد ضخمة جدًا، لذلك لا تعيد التنفيذات تحليل كل \(\binom{N}{k}\) من الصفر.
إذا كانت \(f(x)=ax+b\)، فإن تركيبها مع نفسها يعطي
$$f(f(x))\equiv a(ax+b)+b \equiv a^2x+ab+b \pmod{n}.$$
ولكي يساوي هذا التعبير \(ax+b\) لكل \(x\)، يجب أن يختفي كل من معامل \(x\) والجزء الثابت بترديد \(n\). إذن لا بد من تحقق
$$n\mid a(a-1),\qquad n\mid ab.$$
وبذلك ينقسم العد إلى جزأين: إيجاد جميع قيم \(a\) التي تحقق علاقة idempotence، ثم عدُّ قيم \(b\) الموافقة لكل \(a\).
لنكتب
$$n=\prod_{i=1}^{m} p_i^{e_i}.$$
وبما أن عددين متتاليين أوليان نسبيًا، فإن \(\gcd(a,a-1)=1\). لذلك إذا كان \(p_i^{e_i}\mid a(a-1)\)، فلا بد أن تقع القوة الأولية كلها في أحد العاملين فقط. ومن ثم لكل \(p_i^{e_i}\parallel n\) يجب أن يتحقق أحد الشرطين الآتيين:
$$a\equiv 0 \pmod{p_i^{e_i}}\qquad\text{or}\qquad a\equiv 1 \pmod{p_i^{e_i}}.$$
هذه الاختيارات مستقلة بين القوى الأولية المختلفة، ثم تجمعها مبرهنة البواقي الصينية في باقي واحد وحيد modulo \(n\).
إذا عرّفنا
$$d=\gcd(a,n),$$
فإن \(d\) يكون قاسمًا وحدويًا لـ \(n\)، أي \(d\mid n\) و\(\gcd(d,n/d)=1\). وبالعكس، كل قاسم وحدوي يحقق \(d\lt n\) يعطي باقيًا مناسبًا واحدًا لـ \(a\). أمّا الحالة المستبعدة \(d=n\) فتعني \(a\equiv 0\pmod n\)، وهذا ممنوع لأن الشرط يفرض \(0 \lt a \lt n\).
لنثبت \(a\) صالحًا واحدًا، ولتكن \(d=\gcd(a,n)\). نكتب
$$a=d\,a_1,\qquad n=d\,n_1,\qquad \gcd(a_1,n_1)=1.$$
حينئذ يصبح الشرط الثاني
$$n\mid ab \iff d\,n_1 \mid d\,a_1 b \iff n_1 \mid a_1 b.$$
وبما أن \(a_1\) قابل للعكس modulo \(n_1\)، فهذا يكافئ
$$n_1\mid b.$$
ومن بين البواقي \(0\le b\lt n=d\,n_1\)، يوجد بالضبط \(d\) عددًا من مضاعفات \(n_1\). إذن كل \(a\) صالح يساهم بعدد \(d=\gcd(a,n)\) من القيم الممكنة لـ \(b\).
وبالجمع على جميع القواسم الوحدوية نحصل على
$$R(n)=\sum_{\substack{d\parallel n\\ d\lt n}} d = \left(\sum_{d\parallel n} d\right)-n.$$
ومن ثم
$$\boxed{R(n)=\sigma^*(n)-n},$$
حيث \(\sigma^*(n)\) هو مجموع القواسم الوحدوية. ولأن كل قوة أولية تُختار كاملة أو لا تُختار أصلًا، فإن
$$\sigma^*(n)=\sum_{d\parallel n} d=\prod_{p^e\parallel n}(1+p^e).$$
لنضع
$$B_k=\binom{N}{k}.$$
فإذا كان
$$B_k=\prod_p p^{e_p(k)},$$
فإن الصيغة المغلقة تصبح
$$R(B_k)=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)-B_k.$$
وهكذا تتحول المسألة إلى المحافظة على أسس العوامل الأولية لـ \(B_k\) أثناء التحرك عبر صف واحد من مثلث باسكال.
العلاقة القياسية
$$B_k=B_{k-1}\cdot \frac{N-k+1}{k}$$
تعطي لكل عدد أولي \(p\)
$$e_p(k)=e_p(k-1)+v_p(N-k+1)-v_p(k).$$
لذلك لا تحتاج التنفيذات في كل خطوة إلا إلى تحليل العددين العاديين \(N-k+1\) و\(k\)، وليس تحليل \(\binom{N}{k}\) الهائل مباشرة.
وفي الوقت نفسه يُحفَظ المقدار
$$Q_k=\prod_{e_p(k)\gt 0}\left(1+p^{e_p(k)}\right)\pmod{M},$$
ثم يُحسب
$$R(B_k)\equiv Q_k-B_k \pmod{M}.$$
وبما أن \(M\) أولي و\(N \lt M\)، فيمكن كذلك تحديث \(B_k \bmod M\) باستعمال المعكوسات modular:
$$B_k \equiv B_{k-1}(N-k+1)k^{-1}\pmod{M}.$$
معاملات ثنائية الحد تحقق
$$\binom{N}{k}=\binom{N}{N-k},$$
ولذلك يكون المجموع متناظرًا. فتكتفي التنفيذات بالدوران حتى \(\lfloor N/2\rfloor\)، مع وزن \(2\) لكل مساهمة ما عدا الحد الأوسط عندما يكون \(N\) زوجيًا.
عند \(N=10\) تكون القيم المختلفة
$$\binom{10}{1}=10,\quad \binom{10}{2}=45,\quad \binom{10}{3}=120,\quad \binom{10}{4}=210,\quad \binom{10}{5}=252.$$
وباستخدام \(R(n)=\sigma^*(n)-n\) نجد
$$\begin{aligned} R(10)&=(1+2)(1+5)-10=8,\\ R(45)&=(1+9)(1+5)-45=15,\\ R(120)&=(1+8)(1+3)(1+5)-120=96,\\ R(210)&=(1+2)(1+3)(1+5)(1+7)-210=366,\\ R(252)&=(1+4)(1+9)(1+7)-252=148. \end{aligned}$$
ومن ثم
$$S(10)=2(8+15+96+366)+148=1118,$$
وهو بالضبط التحقق الصغير الذي تطابقه التنفيذات.
تبدأ تنفيذات C++ وPython وJava ببناء جدول أصغر عامل أولي حتى \(N\). كما تُحسب مسبقًا المعكوسات modular للأعداد \(1,2,\dots,N\) بحيث يمكن تحديث معامل ثنائية الحد الحالي بزمن ثابت بعد معرفة عوامل البسط والمقام الجديدة.
داخل الحلقة الرئيسية تحتفظ التنفيذات، لكل عدد أولي نشط، بأسه الحالي داخل \(B_k\) وبالقيمة الحالية \(p^{e_p(k)} \bmod M\). وعندما يتغير الأس، يُزال العامل القديم \((1+p^e)\) من حاصل الضرب الجاري \(Q_k\)، ثم يُعدَّل الأس، ثم يُدرج العامل الجديد.
القسمة داخل modulo تتم بضرب المعكوس. وتُخزن معكوسات العوامل غير الصفرية عند الحاجة الأولى. كما تحتفظ التنفيذات بعدّاد مستقل لعدد العوامل التي تساوي \(0 \pmod M\). فإذا وُجد عامل صفري واحد على الأقل صار \(Q_k\equiv 0\pmod M\)، وبذلك تبقى الصيغة صحيحة حتى في هذه الحالة الحدية.
تستخدم مرحلة SPF ذاكرة \(O(N)\). أما العمل الكلي لتحليل العوامل عبر جميع التحديثات فيتحكم فيه عدد العوامل الأولية مع التكرار في الأعداد \(1,2,\dots,N\)، وهو في المجمل قريب من \(O(N\log\log N)\). وما تبقى من حسابات في كل خطوة ثابت الكلفة. لذلك تعمل الطريقة كلها بزمن قريب من \(O(N\log\log N)\) وذاكرة \(O(N)\)، وهو مناسب للقيمة \(N=10^7\).