For a prime \(p\), let \(R(p)\) be the number of ordered pairs \((x,y)\) with \(1\le x,y\le p(p-1)\) such that
$$x^y \equiv y^x \pmod p.$$
The full task is to evaluate
$$S(M,N)\equiv \sum_{\substack{M \le p \le N \\ p\text{ prime}}} R(p)\pmod{993353399}.$$
A direct search over all pairs is far too large, so the solution converts the congruence into a group-theoretic count that depends only on the factorization of \(p-1\).
Fix one prime \(p\) and write \(n=p-1\). The multiplicative group \(\mathbb{F}_p^\times\) has order \(n\), and that is the structure the solution exploits.
The interval \(1,2,\dots,pn\) contains exactly \(n\) multiples of \(p\).
If exactly one of \(x\) or \(y\) is divisible by \(p\), then one side of \(x^y \equiv y^x \pmod p\) is \(0\) and the other is nonzero, so such pairs never work.
If both \(x\) and \(y\) are divisible by \(p\), then both sides are \(0\), so every such pair works. Therefore these pairs contribute
$$n^2$$
solutions immediately.
Everything else comes from pairs with \(p\nmid x\) and \(p\nmid y\).
Choose a primitive root \(g\) modulo \(p\). Then every nonzero residue class can be written uniquely as \(g^a\) for some \(a\in \mathbb{Z}/n\mathbb{Z}\).
So for \(p\nmid x\) and \(p\nmid y\), write
$$x\equiv g^a \pmod p,\qquad y\equiv g^b \pmod p$$
with \(a,b\in \mathbb{Z}/n\mathbb{Z}\).
Now look at all integers in \(1\le x\le pn\) having the same residue modulo \(p\). They are
$$r,\ r+p,\ r+2p,\ \dots,\ r+(n-1)p.$$
Because \(p\equiv 1 \pmod n\), these numbers are congruent modulo \(n\) to
$$r,\ r+1,\ r+2,\ \dots,\ r+(n-1),$$
so each residue class modulo \(n\) occurs exactly once. Hence every pair
$$\bigl(a,u\bigr)\in (\mathbb{Z}/n\mathbb{Z})^2$$
corresponds to a unique integer \(x\) with
$$x\equiv g^a \pmod p,\qquad x\equiv u \pmod n,$$
and similarly every \((b,v)\) determines a unique \(y\).
Since \(g\) has order \(n\), we have
$$x^y \equiv y^x \pmod p \iff g^{ay}\equiv g^{bx}\pmod p \iff ay\equiv bx\pmod n.$$
Only the classes of \(x\) and \(y\) modulo \(n\) matter in the exponents, so with \(u\equiv x\pmod n\) and \(v\equiv y\pmod n\), the condition becomes
$$av\equiv bu\pmod n.$$
For fixed \(a\) and \(b\), define
$$\Phi_{a,b}(u,v)=bu-av \pmod n.$$
This is a homomorphism from \((\mathbb{Z}/n\mathbb{Z})^2\) to \(\mathbb{Z}/n\mathbb{Z}\). If
$$d=\gcd(a,b,n),$$
then the image of \(\Phi_{a,b}\) is exactly the subgroup of multiples of \(d\), which has size \(n/d\). Therefore the kernel has size
$$\frac{n^2}{n/d}=nd.$$
So for each fixed pair \((a,b)\), the number of admissible pairs \((u,v)\) is
$$n\,\gcd(a,b,n).$$
Summing over all \(a,b\in \mathbb{Z}/n\mathbb{Z}\) gives the nonzero contribution
$$n\,H(n),\qquad H(n)=\sum_{a=0}^{n-1}\sum_{b=0}^{n-1}\gcd(a,b,n).$$
Use the standard identity
$$\gcd(a,b,n)=\sum_{d\mid \gcd(a,b,n)} \varphi(d),$$
where \(\varphi\) is Euler's totient function.
Substituting this into \(H(n)\) and swapping the order of summation yields
$$H(n)=\sum_{d\mid n}\varphi(d)\left(\frac{n}{d}\right)^2.$$
Indeed, for a fixed divisor \(d\mid n\), there are exactly \(n/d\) residue classes modulo \(n\) divisible by \(d\), both for \(a\) and for \(b\).
This formula already shows that \(H(n)\) is multiplicative in \(n\).
Let
$$n=\prod_{i=1}^t q_i^{e_i}.$$
Because \(H\) is multiplicative, it is enough to evaluate one prime power:
$$H(q^e)=\sum_{k=0}^{e}\varphi(q^k)\,q^{2e-2k}.$$
Now \(\varphi(q^0)=1\) and \(\varphi(q^k)=q^k-q^{k-1}\) for \(k\ge 1\), so
$$H(q^e)=q^{2e}+\sum_{k=1}^{e}(q^k-q^{k-1})q^{2e-2k}.$$
This simplifies to
$$H(q^e)=q^{2e}+q^{2e-1}-q^{e-1}=q^{e-1}(q^{e+1}+q^e-1).$$
Therefore
$$H(n)=\prod_{i=1}^{t} q_i^{e_i-1}(q_i^{e_i+1}+q_i^{e_i}-1),$$
and the final prime-local count is
$$\boxed{R(p)=n^2+n\prod_{q^e\parallel n} q^{e-1}(q^{e+1}+q^e-1),\qquad n=p-1.}$$
Here \(n=p-1=4\), so the search range is \(1\le x,y\le 20\).
There are \(4\) multiples of \(5\) in that range, namely \(5,10,15,20\). Pairs with both entries divisible by \(5\) contribute
$$4^2=16$$
solutions.
Next compute
$$H(4)=\sum_{d\mid 4}\varphi(d)\left(\frac{4}{d}\right)^2=\varphi(1)\cdot 4^2+\varphi(2)\cdot 2^2+\varphi(4)\cdot 1^2=16+4+2=22.$$
So the nonzero part contributes
$$4\cdot 22=88,$$
and hence
$$R(5)=16+88=104.$$
A direct brute-force check over \(1\le x,y\le 20\) gives the same value.
The C++, Python, and Java implementations scan the interval \([M,N]\), handle the even prime separately, and test odd candidates for primality. Only primes contribute to the final sum.
For each such prime \(p\), the implementation factors \(n=p-1\). Fast modular exponentiation is used throughout, and Pollard-Rho is used to split composite factors until the full prime factorization is known. Equal prime factors are then grouped to obtain the exponents \(e\).
Once the factorization is available, the implementation evaluates the local prime-power factors
$$q^{e-1}(q^{e+1}+q^e-1)\pmod{993353399}$$
and multiplies them to obtain \(H(n)\) modulo \(993353399\). Finally it adds
$$n^2+nH(n)\pmod{993353399}$$
to the running interval sum. The program never searches over \((x,y)\); all counting is done through the closed formula above.
If \(W=N-M+1\), the interval scan examines \(O(W)\) candidates, skipping even numbers after \(2\). For each prime \(p\), the dominant cost is factoring \(p-1\); once the factorization is known, evaluating the product formula only needs a number of modular multiplications proportional to the number of distinct prime factors. In practice the method is fast because it replaces a quadratic search over pairs by arithmetic on the much smaller factorization of \(p-1\). Memory usage per prime stays small, essentially the temporary factor list and recursion stack.
Für eine Primzahl \(p\) bezeichne \(R(p)\) die Anzahl geordneter Paare \((x,y)\) mit \(1\le x,y\le p(p-1)\), für die
$$x^y \equiv y^x \pmod p$$
gilt. Gesucht ist insgesamt
$$S(M,N)\equiv \sum_{\substack{M \le p \le N \\ p\text{ prime}}} R(p)\pmod{993353399}.$$
Eine direkte Enumeration aller Paare ist viel zu teuer. Die Lösung ersetzt sie durch eine Zählung in der multiplikativen Gruppe modulo \(p\), die nur noch von der Primfaktorzerlegung von \(p-1\) abhängt.
Fixiere eine Primzahl \(p\) und schreibe \(n=p-1\). Die Gruppe \(\mathbb{F}_p^\times\) hat genau \(n\) Elemente; auf dieser Struktur basiert die gesamte Herleitung.
Im Bereich \(1,2,\dots,pn\) liegen genau \(n\) Vielfache von \(p\).
Ist genau eine der beiden Zahlen \(x\) oder \(y\) durch \(p\) teilbar, dann ist eine Seite der Kongruenz \(0\), die andere aber ungleich \(0\). Solche Paare scheiden also aus.
Sind dagegen beide Zahlen durch \(p\) teilbar, dann sind beide Seiten \(0\). Diese Paare liefern sofort
$$n^2$$
Lösungen.
Alle übrigen Lösungen stammen von Paaren mit \(p\nmid x\) und \(p\nmid y\).
Wähle eine primitive Wurzel \(g\) modulo \(p\). Dann lässt sich jede von \(0\) verschiedene Restklasse eindeutig als \(g^a\) mit \(a\in \mathbb{Z}/n\mathbb{Z}\) schreiben.
Für \(p\nmid x\) und \(p\nmid y\) gilt also
$$x\equiv g^a \pmod p,\qquad y\equiv g^b \pmod p$$
mit \(a,b\in \mathbb{Z}/n\mathbb{Z}\).
Betrachte nun alle Zahlen \(x\) mit festem Rest \(r\not\equiv 0\pmod p\) im Intervall \(1\le x\le pn\):
$$r,\ r+p,\ r+2p,\ \dots,\ r+(n-1)p.$$
Da \(p\equiv 1 \pmod n\), sind diese Zahlen modulo \(n\) kongruent zu
$$r,\ r+1,\ r+2,\ \dots,\ r+(n-1).$$
Jede Restklasse modulo \(n\) kommt also genau einmal vor. Deshalb bestimmt jedes Paar
$$\bigl(a,u\bigr)\in (\mathbb{Z}/n\mathbb{Z})^2$$
genau ein \(x\) mit
$$x\equiv g^a \pmod p,\qquad x\equiv u \pmod n,$$
und entsprechend bestimmt \((b,v)\) genau ein \(y\).
Weil \(g\) Ordnung \(n\) hat, folgt
$$x^y \equiv y^x \pmod p \iff g^{ay}\equiv g^{bx}\pmod p \iff ay\equiv bx\pmod n.$$
In den Exponenten zählt also nur der Rest modulo \(n\). Mit \(u\equiv x\pmod n\) und \(v\equiv y\pmod n\) wird die Bedingung
$$av\equiv bu\pmod n.$$
Für feste \(a\) und \(b\) definiere
$$\Phi_{a,b}(u,v)=bu-av \pmod n.$$
Das ist ein Homomorphismus von \((\mathbb{Z}/n\mathbb{Z})^2\) nach \(\mathbb{Z}/n\mathbb{Z}\). Setze
$$d=\gcd(a,b,n).$$
Dann besteht das Bild von \(\Phi_{a,b}\) genau aus den Vielfachen von \(d\) und hat Größe \(n/d\). Folglich besitzt der Kern Größe
$$\frac{n^2}{n/d}=nd.$$
Für jedes feste \((a,b)\) gibt es also genau
$$n\,\gcd(a,b,n)$$
zulässige Paare \((u,v)\). Summiert über alle \(a,b\in \mathbb{Z}/n\mathbb{Z}\) ergibt das den Nichtnull-Anteil
$$n\,H(n),\qquad H(n)=\sum_{a=0}^{n-1}\sum_{b=0}^{n-1}\gcd(a,b,n).$$
Verwende die Standardidentität
$$\gcd(a,b,n)=\sum_{d\mid \gcd(a,b,n)} \varphi(d),$$
wobei \(\varphi\) die Eulersche Phi-Funktion ist.
Einsetzen in \(H(n)\) und Vertauschen der Summen liefert
$$H(n)=\sum_{d\mid n}\varphi(d)\left(\frac{n}{d}\right)^2.$$
Für jeden festen Teiler \(d\mid n\) gibt es nämlich genau \(n/d\) Restklassen modulo \(n\), die durch \(d\) teilbar sind, sowohl für \(a\) als auch für \(b\).
Damit ist sofort sichtbar, dass \(H(n)\) multiplikativ ist.
Schreibe
$$n=\prod_{i=1}^t q_i^{e_i}.$$
Wegen der Multiplikativität reicht es, \(H(q^e)\) zu berechnen:
$$H(q^e)=\sum_{k=0}^{e}\varphi(q^k)\,q^{2e-2k}.$$
Mit \(\varphi(q^0)=1\) und \(\varphi(q^k)=q^k-q^{k-1}\) für \(k\ge 1\) folgt
$$H(q^e)=q^{2e}+\sum_{k=1}^{e}(q^k-q^{k-1})q^{2e-2k}.$$
Das vereinfacht sich zu
$$H(q^e)=q^{2e}+q^{2e-1}-q^{e-1}=q^{e-1}(q^{e+1}+q^e-1).$$
Also gilt insgesamt
$$H(n)=\prod_{i=1}^{t} q_i^{e_i-1}(q_i^{e_i+1}+q_i^{e_i}-1),$$
und damit für die lokale Anzahl
$$\boxed{R(p)=n^2+n\prod_{q^e\parallel n} q^{e-1}(q^{e+1}+q^e-1),\qquad n=p-1.}$$
Hier ist \(n=p-1=4\), also reicht der Suchbereich von \(1\) bis \(20\).
Die Vielfachen von \(5\) in diesem Bereich sind \(5,10,15,20\). Paare, in denen beide Zahlen durch \(5\) teilbar sind, liefern
$$4^2=16$$
Lösungen.
Außerdem ist
$$H(4)=\sum_{d\mid 4}\varphi(d)\left(\frac{4}{d}\right)^2=\varphi(1)\cdot 4^2+\varphi(2)\cdot 2^2+\varphi(4)\cdot 1^2=16+4+2=22.$$
Der Nichtnull-Anteil ist also
$$4\cdot 22=88,$$
und damit insgesamt
$$R(5)=16+88=104.$$
Eine direkte Vollsuche über \(1\le x,y\le 20\) bestätigt denselben Wert.
Die C++-, Python- und Java-Implementierungen durchlaufen das Intervall \([M,N]\), behandeln die Primzahl \(2\) separat und unterziehen die ungeraden Kandidaten einem schnellen Primzahltest. Nur Primzahlen tragen zur Summe bei.
Für jede gefundene Primzahl \(p\) wird \(n=p-1\) faktorisiert. Dazu werden modulare Potenzen verwendet, und Pollard-Rho zerlegt zusammengesetzte Faktoren rekursiv, bis die vollständige Primfaktorzerlegung vorliegt. Gleiche Primfaktoren werden anschließend zu ihren Exponenten zusammengefasst.
Danach berechnet die Implementierung die lokalen Faktoren
$$q^{e-1}(q^{e+1}+q^e-1)\pmod{993353399}$$
und multipliziert sie zu \(H(n)\) modulo \(993353399\). Abschließend wird
$$n^2+nH(n)\pmod{993353399}$$
zur laufenden Summe addiert. Ein explizites Durchprobieren aller \((x,y)\) findet nirgends statt.
Mit \(W=N-M+1\) betrachtet der Intervall-Scan \(O(W)\) Kandidaten, nach \(2\) nur noch ungerade Zahlen. Für jede Primzahl \(p\) dominiert die Faktorisierung von \(p-1\); ist diese bekannt, dann benötigt die Auswertung der Produktformel nur noch eine Anzahl modularer Multiplikationen proportional zur Zahl der verschiedenen Primfaktoren. Praktisch ist die Methode schnell, weil sie die quadratische Paarsuche durch Arithmetik auf der viel kleineren Zerlegung von \(p-1\) ersetzt. Der Speicherbedarf pro Primzahl bleibt klein und besteht im Wesentlichen aus der temporären Faktorenliste und dem Rekursionsstack.
Bir asal \(p\) için \(R(p)\),
$$x^y \equiv y^x \pmod p$$
koşulunu sağlayan ve \(1\le x,y\le p(p-1)\) aralığında bulunan sıralı \((x,y)\) çiftlerinin sayısı olsun. İstenen toplam
$$S(M,N)\equiv \sum_{\substack{M \le p \le N \\ p\text{ prime}}} R(p)\pmod{993353399}$$
şeklindedir.
Tüm çiftleri doğrudan denemek çok pahalıdır. Çözüm bunun yerine sayımı mod \(p\) çarpımsal grubuna taşıyarak problemi yalnızca \(p-1\)'in asal çarpan yapısına indirger.
Tek bir asal \(p\) sabitleyelim ve \(n=p-1\) yazalım. \(\mathbb{F}_p^\times\) grubunun mertebesi \(n\)'dir; kullanılan tüm kapalı formül bu gözlemden çıkar.
\(1,2,\dots,pn\) aralığında tam olarak \(n\) tane \(p\) katı vardır.
Eğer \(x\) ve \(y\)'den yalnızca biri \(p\)'ye bölünüyorsa, kongruensin bir tarafı \(0\), diğer tarafı ise sıfır olmayan bir değer olur. Dolayısıyla bu tür çiftler çözüm değildir.
Her iki sayı da \(p\)'ye bölünüyorsa iki taraf da \(0\) olur ve tüm bu çiftler geçerlidir. Buradan doğrudan
$$n^2$$
adet çözüm gelir.
Geriye kalan tüm çözümler \(p\nmid x\) ve \(p\nmid y\) olan çiftlerden gelir.
Mod \(p\) için bir ilkel kök \(g\) seçelim. O zaman sıfır olmayan her kalıntı sınıfı, \(a\in \mathbb{Z}/n\mathbb{Z}\) olmak üzere tekil biçimde \(g^a\) şeklinde yazılır.
Bu nedenle \(p\nmid x\) ve \(p\nmid y\) iken
$$x\equiv g^a \pmod p,\qquad y\equiv g^b \pmod p$$
yazabiliriz.
Şimdi mod \(p\)'de aynı sıfır olmayan kalıntıya sahip ve \(1\le x\le pn\) aralığında yer alan tüm sayıları düşünelim:
$$r,\ r+p,\ r+2p,\ \dots,\ r+(n-1)p.$$
\(p\equiv 1 \pmod n\) olduğundan bunlar mod \(n\)'de
$$r,\ r+1,\ r+2,\ \dots,\ r+(n-1)$$
sınıflarına düşer. Yani mod \(n\)'deki her sınıf tam bir kez görülür. Böylece her
$$\bigl(a,u\bigr)\in (\mathbb{Z}/n\mathbb{Z})^2$$
çifti,
$$x\equiv g^a \pmod p,\qquad x\equiv u \pmod n$$
koşullarını sağlayan tek bir \(x\) belirler; aynı şekilde \((b,v)\) de tek bir \(y\) belirler.
\(g\)'nin mertebesi \(n\) olduğundan
$$x^y \equiv y^x \pmod p \iff g^{ay}\equiv g^{bx}\pmod p \iff ay\equiv bx\pmod n.$$
Üslerde yalnızca mod \(n\) sınıfları önemli olduğu için, \(u\equiv x\pmod n\) ve \(v\equiv y\pmod n\) ile koşul
$$av\equiv bu\pmod n$$
olur.
Sabit \(a\) ve \(b\) için
$$\Phi_{a,b}(u,v)=bu-av \pmod n$$
uygulamasını tanımlayalım. Bu, \((\mathbb{Z}/n\mathbb{Z})^2\)'den \(\mathbb{Z}/n\mathbb{Z}\)'ye bir homomorfizmadır. Eğer
$$d=\gcd(a,b,n)$$
ise, görüntü tam olarak \(d\)'nin katlarından oluşur ve boyutu \(n/d\)'dir. Dolayısıyla çekirdek boyutu
$$\frac{n^2}{n/d}=nd$$
olur.
Yani sabit bir \((a,b)\) çifti için uygun \((u,v)\) sayısı
$$n\,\gcd(a,b,n)$$
kadardır. Tüm \(a,b\in \mathbb{Z}/n\mathbb{Z}\) üzerinden toplarsak sıfır olmayan kalıntılardan gelen katkı
$$n\,H(n),\qquad H(n)=\sum_{a=0}^{n-1}\sum_{b=0}^{n-1}\gcd(a,b,n)$$
şeklinde yazılır.
Standart özdeşlik
$$\gcd(a,b,n)=\sum_{d\mid \gcd(a,b,n)} \varphi(d)$$
kullanılır; burada \(\varphi\) Euler totient fonksiyonudur.
Bunu \(H(n)\) içine koyup toplamların yerini değiştirince
$$H(n)=\sum_{d\mid n}\varphi(d)\left(\frac{n}{d}\right)^2$$
elde edilir.
Gerçekten de sabit bir \(d\mid n\) için mod \(n\)'de \(d\)'ye bölünen tam olarak \(n/d\) sınıf vardır; hem \(a\) hem \(b\) için sayı aynıdır.
Bu ifade aynı zamanda \(H(n)\)'nin çarpımsal olduğunu da hemen gösterir.
\(n\)'yi
$$n=\prod_{i=1}^t q_i^{e_i}$$
şeklinde ayıralım. Çarpımsallık nedeniyle yalnızca \(H(q^e)\)'yi hesaplamak yeterlidir:
$$H(q^e)=\sum_{k=0}^{e}\varphi(q^k)\,q^{2e-2k}.$$
\(\varphi(q^0)=1\) ve \(k\ge 1\) için \(\varphi(q^k)=q^k-q^{k-1}\) olduğundan
$$H(q^e)=q^{2e}+\sum_{k=1}^{e}(q^k-q^{k-1})q^{2e-2k}.$$
Bu ifade
$$H(q^e)=q^{2e}+q^{2e-1}-q^{e-1}=q^{e-1}(q^{e+1}+q^e-1)$$
şeklinde sadeleşir.
Dolayısıyla
$$H(n)=\prod_{i=1}^{t} q_i^{e_i-1}(q_i^{e_i+1}+q_i^{e_i}-1),$$
ve son yerel sayım formülü
$$\boxed{R(p)=n^2+n\prod_{q^e\parallel n} q^{e-1}(q^{e+1}+q^e-1),\qquad n=p-1.}$$
Burada \(n=p-1=4\) ve aralık \(1\le x,y\le 20\) olur.
Bu aralıktaki \(5\) katları \(5,10,15,20\) olduğundan her iki sayının da \(5\)'e bölündüğü çiftler
$$4^2=16$$
katkı verir.
Ayrıca
$$H(4)=\sum_{d\mid 4}\varphi(d)\left(\frac{4}{d}\right)^2=\varphi(1)\cdot 4^2+\varphi(2)\cdot 2^2+\varphi(4)\cdot 1^2=16+4+2=22$$
olur.
Dolayısıyla sıfır olmayan kısımdan gelen katkı
$$4\cdot 22=88,$$
toplam ise
$$R(5)=16+88=104$$
çıkar. Küçük aralıkta doğrudan sayım da aynı sonucu verir.
C++, Python ve Java uygulamaları \([M,N]\) aralığını tarar, \(2\)'yi ayrı ele alır ve tek adayları hızlı asal testinden geçirir. Yalnızca asal \(p\) değerleri toplamda yer alır.
Her asal için \(n=p-1\) çarpanlara ayrılır. Bileşik parçaları ayırmak için Pollard-Rho kullanılır; modüler üs alma işlemleri hızlı üs alma ile yapılır. Aynı asal çarpanlar gruplanarak üsler \(e\) elde edilir.
Ardından uygulama
$$q^{e-1}(q^{e+1}+q^e-1)\pmod{993353399}$$
yerel katsayılarını hesaplayıp çarparak \(H(n)\)'yi bulur ve son olarak
$$n^2+nH(n)\pmod{993353399}$$
değerini toplama ekler. Program hiçbir aşamada \((x,y)\) çiftlerini tek tek dolaşmaz.
\(W=N-M+1\) olmak üzere aralık taraması \(O(W)\) aday inceler; \(2\)'den sonra yalnızca tek sayılar kontrol edilir. Her asal \(p\) için baskın maliyet \(p-1\)'in çarpanlara ayrılmasıdır. Çarpanlar bulunduğunda kapalı formülün değerlendirilmesi yalnızca farklı asal çarpan sayısı kadar modüler çarpım gerektirir. Pratikte yöntem hızlıdır çünkü çiftler üzerinde ikinci dereceden bir arama yapmak yerine çok daha küçük olan \(p-1\)'in çarpan yapısı üzerinde çalışır. Bellek kullanımı da her asal için geçici çarpan listesi ve özyineleme yığını ile sınırlı kalır.
Para un primo \(p\), sea \(R(p)\) el número de pares ordenados \((x,y)\) con \(1\le x,y\le p(p-1)\) que satisfacen
$$x^y \equiv y^x \pmod p.$$
La cantidad total pedida es
$$S(M,N)\equiv \sum_{\substack{M \le p \le N \\ p\text{ prime}}} R(p)\pmod{993353399}.$$
Probar todos los pares directamente es inviable. La solución transforma la congruencia en un conteo dentro del grupo multiplicativo módulo \(p\), y al final todo depende solo de la factorización de \(p-1\).
Fijemos un primo \(p\) y escribamos \(n=p-1\). El grupo \(\mathbb{F}_p^\times\) tiene orden \(n\), y esa es la estructura que hace posible la fórmula cerrada.
En el intervalo \(1,2,\dots,pn\) hay exactamente \(n\) múltiplos de \(p\).
Si exactamente uno de los dos números \(x\) o \(y\) es divisible por \(p\), entonces un lado de la congruencia vale \(0\) y el otro no, así que esos pares no sirven.
Si ambos son divisibles por \(p\), entonces ambos lados valen \(0\). Eso aporta inmediatamente
$$n^2$$
soluciones.
Todo lo demás proviene de pares con \(p\nmid x\) y \(p\nmid y\).
Elige una raíz primitiva \(g\) módulo \(p\). Entonces cada clase no nula se escribe de forma única como \(g^a\) para algún \(a\in \mathbb{Z}/n\mathbb{Z}\).
Por tanto, cuando \(p\nmid x\) y \(p\nmid y\), podemos escribir
$$x\equiv g^a \pmod p,\qquad y\equiv g^b \pmod p$$
con \(a,b\in \mathbb{Z}/n\mathbb{Z}\).
Ahora fijemos un residuo no nulo \(r\) módulo \(p\). Los enteros del intervalo \(1\le x\le pn\) con ese residuo son
$$r,\ r+p,\ r+2p,\ \dots,\ r+(n-1)p.$$
Como \(p\equiv 1 \pmod n\), estos números son congruentes módulo \(n\) a
$$r,\ r+1,\ r+2,\ \dots,\ r+(n-1).$$
Eso significa que cada clase módulo \(n\) aparece exactamente una vez. En consecuencia, cada par
$$\bigl(a,u\bigr)\in (\mathbb{Z}/n\mathbb{Z})^2$$
determina un único entero \(x\) con
$$x\equiv g^a \pmod p,\qquad x\equiv u \pmod n,$$
y del mismo modo cada \((b,v)\) determina un único \(y\).
Como \(g\) tiene orden \(n\), se obtiene
$$x^y \equiv y^x \pmod p \iff g^{ay}\equiv g^{bx}\pmod p \iff ay\equiv bx\pmod n.$$
En los exponentes solo importa la clase módulo \(n\). Si \(u\equiv x\pmod n\) y \(v\equiv y\pmod n\), la condición pasa a ser
$$av\equiv bu\pmod n.$$
Para \(a\) y \(b\) fijos, definimos
$$\Phi_{a,b}(u,v)=bu-av \pmod n.$$
Es un homomorfismo de \((\mathbb{Z}/n\mathbb{Z})^2\) en \(\mathbb{Z}/n\mathbb{Z}\). Si
$$d=\gcd(a,b,n),$$
entonces la imagen de \(\Phi_{a,b}\) es exactamente el subgrupo de los múltiplos de \(d\), que tiene tamaño \(n/d\). Por tanto, el núcleo tiene tamaño
$$\frac{n^2}{n/d}=nd.$$
Así, para cada par fijo \((a,b)\), el número de pares \((u,v)\) admisibles es
$$n\,\gcd(a,b,n).$$
Al sumar sobre todos los \(a,b\in \mathbb{Z}/n\mathbb{Z}\), la contribución no nula queda
$$n\,H(n),\qquad H(n)=\sum_{a=0}^{n-1}\sum_{b=0}^{n-1}\gcd(a,b,n).$$
Usamos la identidad estándar
$$\gcd(a,b,n)=\sum_{d\mid \gcd(a,b,n)} \varphi(d),$$
donde \(\varphi\) es la función phi de Euler.
Al sustituir en \(H(n)\) e intercambiar el orden de las sumas, obtenemos
$$H(n)=\sum_{d\mid n}\varphi(d)\left(\frac{n}{d}\right)^2.$$
Para un divisor fijo \(d\mid n\), hay exactamente \(n/d\) clases módulo \(n\) divisibles por \(d\), tanto para \(a\) como para \(b\).
Esta fórmula muestra además que \(H(n)\) es multiplicativa.
Escribe
$$n=\prod_{i=1}^t q_i^{e_i}.$$
Como \(H\) es multiplicativa, basta calcular
$$H(q^e)=\sum_{k=0}^{e}\varphi(q^k)\,q^{2e-2k}.$$
Con \(\varphi(q^0)=1\) y \(\varphi(q^k)=q^k-q^{k-1}\) para \(k\ge 1\), resulta
$$H(q^e)=q^{2e}+\sum_{k=1}^{e}(q^k-q^{k-1})q^{2e-2k}.$$
Eso se simplifica a
$$H(q^e)=q^{2e}+q^{2e-1}-q^{e-1}=q^{e-1}(q^{e+1}+q^e-1).$$
Por lo tanto,
$$H(n)=\prod_{i=1}^{t} q_i^{e_i-1}(q_i^{e_i+1}+q_i^{e_i}-1),$$
y la cuenta local final es
$$\boxed{R(p)=n^2+n\prod_{q^e\parallel n} q^{e-1}(q^{e+1}+q^e-1),\qquad n=p-1.}$$
Aquí \(n=p-1=4\), así que el rango es \(1\le x,y\le 20\).
Los múltiplos de \(5\) en ese rango son \(5,10,15,20\). Los pares en los que ambos son divisibles por \(5\) aportan
$$4^2=16$$
soluciones.
Además,
$$H(4)=\sum_{d\mid 4}\varphi(d)\left(\frac{4}{d}\right)^2=\varphi(1)\cdot 4^2+\varphi(2)\cdot 2^2+\varphi(4)\cdot 1^2=16+4+2=22.$$
Entonces la contribución no nula vale
$$4\cdot 22=88,$$
y por consiguiente
$$R(5)=16+88=104.$$
Una comprobación directa sobre \(1\le x,y\le 20\) produce exactamente el mismo valor.
Las implementaciones en C++, Python y Java recorren el intervalo \([M,N]\), tratan el primo \(2\) por separado y someten los candidatos impares a una prueba rápida de primalidad. Solo los primos contribuyen a la suma final.
Para cada primo \(p\), la implementación factoriza \(n=p-1\). Se utiliza exponenciación modular rápida y Pollard-Rho para partir factores compuestos hasta obtener la factorización prima completa. Después se agrupan los factores iguales para recuperar cada exponente \(e\).
Con la factorización ya disponible, la implementación evalúa los factores locales
$$q^{e-1}(q^{e+1}+q^e-1)\pmod{993353399}$$
y los multiplica para obtener \(H(n)\) módulo \(993353399\). Finalmente añade
$$n^2+nH(n)\pmod{993353399}$$
al acumulador del intervalo. El programa no recorre pares \((x,y)\); toda la cuenta se hace mediante la fórmula cerrada.
Si \(W=N-M+1\), el barrido del intervalo examina \(O(W)\) candidatos, omitiendo los pares después de \(2\). Para cada primo \(p\), el coste dominante es factorizar \(p-1\); una vez conocida la factorización, evaluar el producto solo requiere un número de multiplicaciones modulares proporcional al número de primos distintos que aparecen. En la práctica el método es eficiente porque reemplaza una búsqueda cuadrática en \((x,y)\) por aritmética sobre la descomposición de \(p-1\), que es mucho más pequeña. La memoria por primo es reducida: esencialmente la lista temporal de factores y la pila de recursión.
对每个素数 \(p\)，记 \(R(p)\) 为满足下式的有序整数对 \((x,y)\) 的个数，其中 \(1\le x,y\le p(p-1)\)：
$$x^y \equiv y^x \pmod p.$$
整个题目要求计算
$$S(M,N)\equiv \sum_{\substack{M \le p \le N \\ p\text{ prime}}} R(p)\pmod{993353399}.$$
如果直接枚举所有 \((x,y)\)，规模会大得无法接受。真正可行的做法，是把同余条件改写成模 \(p\) 乘法群中的计数问题，最后只需要知道 \(p-1\) 的质因数分解。
固定一个素数 \(p\)，并记 \(n=p-1\)。有限域中的乘法群 \(\mathbb{F}_p^\times\) 恰好有 \(n\) 个元素，这正是整套推导的核心。
区间 \(1,2,\dots,pn\) 中恰好有 \(n\) 个数能被 \(p\) 整除。
如果 \(x\) 和 \(y\) 中只有一个能被 \(p\) 整除，那么同余式的一边是 \(0\)，另一边是非零值，因此这种情况不可能成立。
如果 \(x\) 和 \(y\) 都能被 \(p\) 整除，那么两边都是 \(0\)，所以所有这样的有序对都有效。它们一共贡献
$$n^2$$
个解。
因此，剩余部分只需要研究 \(p\nmid x\) 且 \(p\nmid y\) 的情况。
取模 \(p\) 的一个原根 \(g\)。那么每个非零剩余类都能唯一写成 \(g^a\) 的形式，其中 \(a\in \mathbb{Z}/n\mathbb{Z}\)。
于是当 \(p\nmid x\) 且 \(p\nmid y\) 时，可以写成
$$x\equiv g^a \pmod p,\qquad y\equiv g^b \pmod p$$
其中 \(a,b\in \mathbb{Z}/n\mathbb{Z}\)。
接下来固定模 \(p\) 的某个非零剩余 \(r\)。在区间 \(1\le x\le pn\) 内，与它同余的所有整数是
$$r,\ r+p,\ r+2p,\ \dots,\ r+(n-1)p.$$
因为 \(p\equiv 1 \pmod n\)，这些数在模 \(n\) 下分别同余于
$$r,\ r+1,\ r+2,\ \dots,\ r+(n-1).$$
也就是说，模 \(n\) 的每一个剩余类都恰好出现一次。于是任意一个
$$\bigl(a,u\bigr)\in (\mathbb{Z}/n\mathbb{Z})^2$$
都会唯一对应一个整数 \(x\)，满足
$$x\equiv g^a \pmod p,\qquad x\equiv u \pmod n,$$
同样地，任意 \((b,v)\) 都唯一对应一个 \(y\)。
由于 \(g\) 的阶为 \(n\)，有
$$x^y \equiv y^x \pmod p \iff g^{ay}\equiv g^{bx}\pmod p \iff ay\equiv bx\pmod n.$$
指数部分只依赖于模 \(n\) 的剩余类。因此若令 \(u\equiv x\pmod n\)、\(v\equiv y\pmod n\)，条件就变成
$$av\equiv bu\pmod n.$$
对固定的 \(a\) 和 \(b\)，定义映射
$$\Phi_{a,b}(u,v)=bu-av \pmod n.$$
这是从 \((\mathbb{Z}/n\mathbb{Z})^2\) 到 \(\mathbb{Z}/n\mathbb{Z}\) 的群同态。设
$$d=\gcd(a,b,n).$$
那么 \(\Phi_{a,b}\) 的像正好是模 \(n\) 下所有 \(d\) 的倍数，像的大小为 \(n/d\)。于是核的大小是
$$\frac{n^2}{n/d}=nd.$$
这说明：对每个固定的 \((a,b)\)，满足条件的 \((u,v)\) 恰有
$$n\,\gcd(a,b,n)$$
个。
把所有 \(a,b\in \mathbb{Z}/n\mathbb{Z}\) 累加起来，非零剩余部分的总贡献就是
$$n\,H(n),\qquad H(n)=\sum_{a=0}^{n-1}\sum_{b=0}^{n-1}\gcd(a,b,n).$$
使用标准恒等式
$$\gcd(a,b,n)=\sum_{d\mid \gcd(a,b,n)} \varphi(d),$$
其中 \(\varphi\) 是 Euler 的 totient 函数。
把它代入 \(H(n)\)，并交换求和顺序，可以得到
$$H(n)=\sum_{d\mid n}\varphi(d)\left(\frac{n}{d}\right)^2.$$
原因很直接：当 \(d\mid n\) 固定时，模 \(n\) 下能被 \(d\) 整除的剩余类恰好有 \(n/d\) 个，\(a\) 和 \(b\) 都是如此。
这个表达式同时说明 \(H(n)\) 是乘法性的。
把 \(n\) 分解为
$$n=\prod_{i=1}^t q_i^{e_i}.$$
因为 \(H\) 是乘法函数，只需先求
$$H(q^e)=\sum_{k=0}^{e}\varphi(q^k)\,q^{2e-2k}.$$
注意 \(\varphi(q^0)=1\)，而对 \(k\ge 1\) 有 \(\varphi(q^k)=q^k-q^{k-1}\)。因此
$$H(q^e)=q^{2e}+\sum_{k=1}^{e}(q^k-q^{k-1})q^{2e-2k}.$$
化简后得到
$$H(q^e)=q^{2e}+q^{2e-1}-q^{e-1}=q^{e-1}(q^{e+1}+q^e-1).$$
于是
$$H(n)=\prod_{i=1}^{t} q_i^{e_i-1}(q_i^{e_i+1}+q_i^{e_i}-1),$$
最后的素数局部计数公式就是
$$\boxed{R(p)=n^2+n\prod_{q^e\parallel n} q^{e-1}(q^{e+1}+q^e-1),\qquad n=p-1.}$$
此时 \(n=p-1=4\)，所以搜索范围是 \(1\le x,y\le 20\)。
区间中能被 \(5\) 整除的数有 \(5,10,15,20\)，因此两者都被 \(5\) 整除的有序对贡献
$$4^2=16$$
个解。
再计算
$$H(4)=\sum_{d\mid 4}\varphi(d)\left(\frac{4}{d}\right)^2=\varphi(1)\cdot 4^2+\varphi(2)\cdot 2^2+\varphi(4)\cdot 1^2=16+4+2=22.$$
于是非零部分贡献
$$4\cdot 22=88,$$
最终得到
$$R(5)=16+88=104.$$
对 \(1\le x,y\le 20\) 直接暴力验证，也正好得到 \(104\)。
C++、Python 和 Java 实现都会遍历区间 \([M,N]\)，把素数 \(2\) 单独处理，然后对奇数候选做快速素性测试。只有素数 \(p\) 才会产生贡献。
对每个素数，程序先分解 \(n=p-1\)。实现中使用快速模幂，并用 Pollard-Rho 递归拆分合数因子，直到得到完整的素因数分解；随后再把相同素因子合并成对应的指数 \(e\)。
分解完成后，程序计算各个局部因子
$$q^{e-1}(q^{e+1}+q^e-1)\pmod{993353399}$$
并把它们相乘得到 \(H(n)\) 在模 \(993353399\) 下的值，最后把
$$n^2+nH(n)\pmod{993353399}$$
累加到答案中。程序从头到尾都不会显式枚举 \((x,y)\) 对。
设 \(W=N-M+1\)。区间扫描会检查 \(O(W)\) 个候选数，除了 \(2\) 之外只检查奇数。对每个素数 \(p\) 来说，主要成本是分解 \(p-1\)；一旦分解完成，计算闭式乘积只需要与不同素因子个数成正比的模乘次数。实际运行中，这种方法之所以高效，是因为它避免了对 \((x,y)\) 的二次枚举，转而处理规模小得多的 \(p-1\) 的因子结构。每个素数的额外内存也很小，基本只是临时因子表和递归栈。
Для простого числа \(p\) обозначим через \(R(p)\) число упорядоченных пар \((x,y)\), где \(1\le x,y\le p(p-1)\), удовлетворяющих сравнению
$$x^y \equiv y^x \pmod p.$$
Нужно вычислить
$$S(M,N)\equiv \sum_{\substack{M \le p \le N \\ p\text{ prime}}} R(p)\pmod{993353399}.$$
Полный перебор всех пар слишком дорог. Решение переводит задачу в язык мультипликативной группы по модулю \(p\), после чего ответ выражается только через разложение числа \(p-1\).
Зафиксируем простое \(p\) и положим \(n=p-1\). Группа \(\mathbb{F}_p^\times\) имеет порядок \(n\), и именно это позволяет получить замкнутую формулу.
В интервале \(1,2,\dots,pn\) ровно \(n\) чисел делятся на \(p\).
Если только одно из чисел \(x\) или \(y\) делится на \(p\), то одна сторона сравнения равна \(0\), а другая нет, поэтому такие пары не подходят.
Если же на \(p\) делятся оба числа, то обе стороны равны \(0\). Это сразу даёт
$$n^2$$
решений.
Все остальные решения приходят из случая \(p\nmid x\) и \(p\nmid y\).
Выберем первообразный корень \(g\) по модулю \(p\). Тогда каждый ненулевой класс вычетов единственным образом записывается как \(g^a\), где \(a\in \mathbb{Z}/n\mathbb{Z}\).
Значит, для \(p\nmid x\) и \(p\nmid y\) можно написать
$$x\equiv g^a \pmod p,\qquad y\equiv g^b \pmod p$$
с некоторыми \(a,b\in \mathbb{Z}/n\mathbb{Z}\).
Рассмотрим все числа из диапазона \(1\le x\le pn\), имеющие фиксированный ненулевой остаток \(r\) по модулю \(p\):
$$r,\ r+p,\ r+2p,\ \dots,\ r+(n-1)p.$$
Так как \(p\equiv 1 \pmod n\), по модулю \(n\) они дают
$$r,\ r+1,\ r+2,\ \dots,\ r+(n-1).$$
То есть каждый класс по модулю \(n\) встречается ровно один раз. Следовательно, любая пара
$$\bigl(a,u\bigr)\in (\mathbb{Z}/n\mathbb{Z})^2$$
однозначно задаёт число \(x\), удовлетворяющее условиям
$$x\equiv g^a \pmod p,\qquad x\equiv u \pmod n,$$
и аналогично пара \((b,v)\) однозначно задаёт число \(y\).
Поскольку \(g\) имеет порядок \(n\), имеем
$$x^y \equiv y^x \pmod p \iff g^{ay}\equiv g^{bx}\pmod p \iff ay\equiv bx\pmod n.$$
В показателях важны только классы по модулю \(n\). Если \(u\equiv x\pmod n\) и \(v\equiv y\pmod n\), условие превращается в
$$av\equiv bu\pmod n.$$
Для фиксированных \(a\) и \(b\) рассмотрим отображение
$$\Phi_{a,b}(u,v)=bu-av \pmod n.$$
Это гомоморфизм из \((\mathbb{Z}/n\mathbb{Z})^2\) в \(\mathbb{Z}/n\mathbb{Z}\). Пусть
$$d=\gcd(a,b,n).$$
Тогда образ состоит ровно из кратных \(d\) и имеет размер \(n/d\). Поэтому ядро имеет размер
$$\frac{n^2}{n/d}=nd.$$
Значит, для каждой фиксированной пары \((a,b)\) число допустимых пар \((u,v)\) равно
$$n\,\gcd(a,b,n).$$
Суммируя по всем \(a,b\in \mathbb{Z}/n\mathbb{Z}\), получаем вклад ненулевых вычетов
$$n\,H(n),\qquad H(n)=\sum_{a=0}^{n-1}\sum_{b=0}^{n-1}\gcd(a,b,n).$$
Используем стандартное тождество
$$\gcd(a,b,n)=\sum_{d\mid \gcd(a,b,n)} \varphi(d),$$
где \(\varphi\) обозначает функцию Эйлера.
Подставляя его в \(H(n)\) и меняя порядок суммирования, получаем
$$H(n)=\sum_{d\mid n}\varphi(d)\left(\frac{n}{d}\right)^2.$$
Действительно, для фиксированного делителя \(d\mid n\) существует ровно \(n/d\) классов по модулю \(n\), делящихся на \(d\), и это верно и для \(a\), и для \(b\).
Из этой формулы сразу видно, что функция \(H(n)\) мультипликативна.
Запишем
$$n=\prod_{i=1}^t q_i^{e_i}.$$
Так как \(H\) мультипликативна, достаточно посчитать
$$H(q^e)=\sum_{k=0}^{e}\varphi(q^k)\,q^{2e-2k}.$$
При этом \(\varphi(q^0)=1\), а для \(k\ge 1\) имеем \(\varphi(q^k)=q^k-q^{k-1}\), поэтому
$$H(q^e)=q^{2e}+\sum_{k=1}^{e}(q^k-q^{k-1})q^{2e-2k}.$$
После упрощения получается
$$H(q^e)=q^{2e}+q^{2e-1}-q^{e-1}=q^{e-1}(q^{e+1}+q^e-1).$$
Следовательно,
$$H(n)=\prod_{i=1}^{t} q_i^{e_i-1}(q_i^{e_i+1}+q_i^{e_i}-1),$$
а итоговая локальная формула имеет вид
$$\boxed{R(p)=n^2+n\prod_{q^e\parallel n} q^{e-1}(q^{e+1}+q^e-1),\qquad n=p-1.}$$
Здесь \(n=p-1=4\), значит диапазон равен \(1\le x,y\le 20\).
Кратные \(5\) в этом диапазоне: \(5,10,15,20\). Пары, в которых оба числа делятся на \(5\), дают
$$4^2=16$$
решений.
Далее
$$H(4)=\sum_{d\mid 4}\varphi(d)\left(\frac{4}{d}\right)^2=\varphi(1)\cdot 4^2+\varphi(2)\cdot 2^2+\varphi(4)\cdot 1^2=16+4+2=22.$$
Следовательно, вклад ненулевой части равен
$$4\cdot 22=88,$$
и потому
$$R(5)=16+88=104.$$
Прямая проверка всех пар \((x,y)\) при \(1\le x,y\le 20\) даёт тот же результат.
Реализации на C++, Python и Java проходят по интервалу \([M,N]\), отдельно обрабатывают простое число \(2\) и проверяют нечётные кандидаты быстрым тестом на простоту. В сумму попадают только простые \(p\).
Для каждого такого простого программа факторизует \(n=p-1\). Используется быстрое модульное возведение в степень, а Pollard-Rho рекурсивно разбивает составные множители до полной факторизации. Затем одинаковые простые множители объединяются, чтобы получить показатели \(e\).
После этого реализация вычисляет локальные множители
$$q^{e-1}(q^{e+1}+q^e-1)\pmod{993353399}$$
и перемножает их, получая \(H(n)\) по модулю \(993353399\). Наконец, к ответу добавляется
$$n^2+nH(n)\pmod{993353399}.$$
Перебора всех \((x,y)\) в программе нет: всё сводится к замкнутой формуле.
Если \(W=N-M+1\), то просмотр интервала проверяет \(O(W)\) кандидатов, причём после \(2\) рассматриваются только нечётные числа. Для каждого простого \(p\) основная стоимость связана с факторизацией \(p-1\); когда разложение уже известно, вычисление произведения требует лишь числа модульных умножений, пропорционального количеству различных простых делителей. На практике метод работает быстро, потому что заменяет квадратичный перебор по \((x,y)\) арифметикой на существенно меньшей структуре разложения числа \(p-1\). Память на одно простое мала: в основном это временный список множителей и стек рекурсии.
لكل عدد أولي \(p\)، لنرمز بـ \(R(p)\) إلى عدد الأزواج المرتبة \((x,y)\) التي تحقق
$$x^y \equiv y^x \pmod p$$
مع الشرط \(1\le x,y\le p(p-1)\). والمطلوب في النهاية هو حساب
$$S(M,N)\equiv \sum_{\substack{M \le p \le N \\ p\text{ prime}}} R(p)\pmod{993353399}.$$
التعداد المباشر لكل الأزواج غير عملي إطلاقًا. لذلك تُحوِّل الحلول المسألة إلى مسألة عد داخل الزمرة الضربية بترديد \(p\)، ثم يصبح كل شيء معتمدًا فقط على تحليل \(p-1\) إلى عوامل أولية.
لنثبت عددًا أوليًا \(p\)، ونكتب \(n=p-1\). الزمرة \(\mathbb{F}_p^\times\) لها رتبة تساوي \(n\)، وهذه هي البنية التي تقود إلى الصيغة المغلقة.
في المجال \(1,2,\dots,pn\) يوجد بالضبط \(n\) عددًا يقبل القسمة على \(p\).
إذا كان أحد العددين \(x\) أو \(y\) فقط من مضاعفات \(p\)، فإن أحد طرفي الت合同 يكون \(0\) بينما يكون الطرف الآخر غير صفري، ولذلك لا نحصل على حل.
أما إذا كان كلاهما من مضاعفات \(p\)، فإن الطرفين يساويان \(0\)، وبالتالي تكون جميع هذه الأزواج صالحة. وهذه الحالة تعطي مباشرة
$$n^2$$
حلًا.
إذن تبقى فقط الأزواج التي تحقق \(p\nmid x\) و\(p\nmid y\).
اختر جذرًا بدائيًا \(g\) بترديد \(p\). عندئذ يمكن كتابة كل باقي غير صفري بشكل وحيد على صورة \(g^a\)، حيث \(a\in \mathbb{Z}/n\mathbb{Z}\).
لذلك عندما يكون \(p\nmid x\) و\(p\nmid y\)، نكتب
$$x\equiv g^a \pmod p,\qquad y\equiv g^b \pmod p$$
لبعض \(a,b\in \mathbb{Z}/n\mathbb{Z}\).
الآن ثبّت باقيًا غير صفري \(r\) modulo \(p\). جميع الأعداد في المجال \(1\le x\le pn\) التي تحقق هذا الباقي هي
$$r,\ r+p,\ r+2p,\ \dots,\ r+(n-1)p.$$
وبما أن \(p\equiv 1 \pmod n\)، فإن هذه القيم تعطي modulo \(n\)
$$r,\ r+1,\ r+2,\ \dots,\ r+(n-1).$$
أي إن كل فئة بترديد \(n\) تظهر مرة واحدة تمامًا. ومن ثم فإن كل زوج
$$\bigl(a,u\bigr)\in (\mathbb{Z}/n\mathbb{Z})^2$$
يحدد عددًا وحيدًا \(x\) يحقق
$$x\equiv g^a \pmod p,\qquad x\equiv u \pmod n,$$
وبالمثل فإن كل \((b,v)\) يحدد عددًا وحيدًا \(y\).
بما أن رتبة \(g\) تساوي \(n\)، فإن
$$x^y \equiv y^x \pmod p \iff g^{ay}\equiv g^{bx}\pmod p \iff ay\equiv bx\pmod n.$$
ولا يهم في الأسس إلا الباقي modulo \(n\). لذلك إذا كتبنا \(u\equiv x\pmod n\) و\(v\equiv y\pmod n\)، تصبح الشرطية
$$av\equiv bu\pmod n.$$
للقيمتين الثابتتين \(a\) و\(b\)، نعرّف
$$\Phi_{a,b}(u,v)=bu-av \pmod n.$$
وهذا تجانس من \((\mathbb{Z}/n\mathbb{Z})^2\) إلى \(\mathbb{Z}/n\mathbb{Z}\). إذا وضعنا
$$d=\gcd(a,b,n),$$
فإن صورة هذا التجانس هي بالضبط زمرة مضاعفات \(d\)، وحجمها \(n/d\). لذلك يكون حجم النواة
$$\frac{n^2}{n/d}=nd.$$
إذن، لكل زوج ثابت \((a,b)\)، يكون عدد الأزواج المسموح بها \((u,v)\) مساويًا لـ
$$n\,\gcd(a,b,n).$$
وبجمع ذلك على جميع \(a,b\in \mathbb{Z}/n\mathbb{Z}\)، نحصل على مساهمة الجزء غير الصفري:
$$n\,H(n),\qquad H(n)=\sum_{a=0}^{n-1}\sum_{b=0}^{n-1}\gcd(a,b,n).$$
نستخدم الهوية القياسية
$$\gcd(a,b,n)=\sum_{d\mid \gcd(a,b,n)} \varphi(d),$$
حيث \(\varphi\) هي دالة أويلر.
بالتعويض في \(H(n)\) ثم تبديل ترتيب الجمع نحصل على
$$H(n)=\sum_{d\mid n}\varphi(d)\left(\frac{n}{d}\right)^2.$$
فالسبب هو أنه عندما يكون \(d\mid n\) ثابتًا، يوجد بالضبط \(n/d\) فئة بترديد \(n\) قابلة للقسمة على \(d\)، سواءً لقيم \(a\) أو لقيم \(b\).
وهذا يبيّن مباشرة أن \(H(n)\) دالة ضربية.
اكتب
$$n=\prod_{i=1}^t q_i^{e_i}.$$
وبما أن \(H\) ضربية، يكفي حساب
$$H(q^e)=\sum_{k=0}^{e}\varphi(q^k)\,q^{2e-2k}.$$
ولأن \(\varphi(q^0)=1\)، ولأن \(\varphi(q^k)=q^k-q^{k-1}\) عندما \(k\ge 1\)، نحصل على
$$H(q^e)=q^{2e}+\sum_{k=1}^{e}(q^k-q^{k-1})q^{2e-2k}.$$
وبعد التبسيط يصبح
$$H(q^e)=q^{2e}+q^{2e-1}-q^{e-1}=q^{e-1}(q^{e+1}+q^e-1).$$
ومن ثم
$$H(n)=\prod_{i=1}^{t} q_i^{e_i-1}(q_i^{e_i+1}+q_i^{e_i}-1),$$
وتكون الصيغة المحلية النهائية
$$\boxed{R(p)=n^2+n\prod_{q^e\parallel n} q^{e-1}(q^{e+1}+q^e-1),\qquad n=p-1.}$$
هنا \(n=p-1=4\)، ولذلك يكون المجال \(1\le x,y\le 20\).
مضاعفات \(5\) في هذا المجال هي \(5,10,15,20\)، وبالتالي فإن الأزواج التي يكون فيها العددان كلاهما من مضاعفات \(5\) تعطي
$$4^2=16$$
حلًا.
كما أن
$$H(4)=\sum_{d\mid 4}\varphi(d)\left(\frac{4}{d}\right)^2=\varphi(1)\cdot 4^2+\varphi(2)\cdot 2^2+\varphi(4)\cdot 1^2=16+4+2=22.$$
إذًا مساهمة الجزء غير الصفري هي
$$4\cdot 22=88,$$
ومن ثم
$$R(5)=16+88=104.$$
والفحص المباشر لجميع الأزواج ضمن \(1\le x,y\le 20\) يعطي القيمة نفسها تمامًا.
تقوم تطبيقات C++ وPython وJava بمسح المجال \([M,N]\)، وتعالج العدد الأولي \(2\) على نحو منفصل، ثم تختبر المرشحين الفرديين باختبار أولية سريع. الأعداد الأولية فقط هي التي تدخل في المجموع النهائي.
لكل عدد أولي \(p\)، تقوم الشيفرة بتحليل \(n=p-1\). ويُستخدم الرفع السريع الأسّي modulo العدد، كما تُستخدم طريقة Pollard-Rho لتفكيك العوامل المركبة حتى نحصل على التحليل الأولي الكامل. بعد ذلك تُجمع العوامل المتساوية لاستخراج الأسس \(e\).
ثم تحسب الشيفرة العوامل المحلية
$$q^{e-1}(q^{e+1}+q^e-1)\pmod{993353399}$$
وتضربها للحصول على \(H(n)\) modulo \(993353399\)، ثم تضيف
$$n^2+nH(n)\pmod{993353399}$$
إلى المجموع التراكمي. لا يوجد في التنفيذ أي تعداد صريح للأزواج \((x,y)\)؛ فالعد كله يتم عبر الصيغة المغلقة.
إذا رمزنا إلى طول المجال بـ \(W=N-M+1\)، فإن المسح يفحص \(O(W)\) مرشحًا، مع تجاهل الأعداد الزوجية بعد \(2\). ولكل أولي \(p\)، تكون الكلفة الغالبة هي تحليل \(p-1\) إلى عوامله الأولية. وبعد معرفة هذا التحليل، لا يحتاج تقييم حاصل الضرب إلا إلى عدد من الضربات المعيارية يتناسب مع عدد العوامل الأولية المتميزة. عمليًا، هذه الطريقة سريعة لأنها تستبدل بحثًا تربيعيًا على الأزواج \((x,y)\) بعمليات حسابية على بنية أصغر بكثير، هي بنية عوامل \(p-1\). أما الذاكرة المطلوبة لكل أولي فتبقى صغيرة، وتكاد تقتصر على قائمة العوامل المؤقتة ومكدس الاستدعاء.