Problem Summary

For \(N=10^{11}\) and the prime \(p=2017\), we must compute

$$S(N,p)=\sum_{\substack{1\le n\le N\\ p\mid \sigma(n)}} n,$$

where \(\sigma(n)\) is the sum of positive divisors of \(n\). A brute-force scan would require evaluating \(\sigma(n)\) for every \(n\le N\), which is far too expensive. The solution instead identifies the exact prime-power exponents that force divisibility by \(2017\), rewrites the task as a weighted union of exact-valuation events, and then sums that union with inclusion-exclusion.

Mathematical Approach

Write the prime factorization of \(n\) as

$$n=\prod_q q^{a_q}.$$

Because the divisor-sum function is multiplicative,

$$\sigma(n)=\prod_q \sigma\!\left(q^{a_q}\right).$$

Therefore \(2017\mid \sigma(n)\) if and only if at least one prime \(q\) satisfies \(2017\mid \sigma(q^{a_q})\).

Step 1: Reduce the problem to trigger prime powers

Define the set of trigger pairs

$$B=\left\{(q,e): e\ge 1,\ q^e\le N,\ 2017\mid \sigma(q^e)\right\},$$

where \(q\) is prime.

Then the target sum is the weighted union

$$S(N,2017)=\sum_{\substack{1\le n\le N\\ \exists (q,e)\in B:\ v_q(n)=e}} n,$$

where \(v_q(n)\) denotes the exponent of \(q\) in \(n\). Exact valuations matter: the condition depends on the precise exponent \(a_q\), not merely on divisibility by \(q^e\).

The prime \(q=2017\) never contributes, because

$$\sigma(2017^e)=1+2017+\cdots+2017^e\equiv 1 \pmod{2017}.$$

Step 2: Characterize the bad exponents with multiplicative order

For a prime \(q\neq 2017\),

$$\sigma(q^e)=1+q+\cdots+q^e=\frac{q^{e+1}-1}{q-1}.$$

If \(q\not\equiv 1 \pmod{2017}\), the denominator is invertible modulo \(2017\), so

$$2017\mid \sigma(q^e)\iff q^{e+1}\equiv 1 \pmod{2017}.$$

Let \(\operatorname{ord}_{2017}(q)\) be the multiplicative order of \(q\) modulo \(2017\). Then

$$2017\mid \sigma(q^e)\iff \operatorname{ord}_{2017}(q)\mid (e+1),$$

so the bad exponents are exactly

$$e=m\operatorname{ord}_{2017}(q)-1,\qquad m\ge 1.$$

Since \(2017-1=2016=2^5\cdot 3^2\cdot 7\), every such order divides \(2016\).

The residue class \(q\equiv 1 \pmod{2017}\) behaves differently, because then

$$\sigma(q^e)\equiv e+1 \pmod{2017}.$$

That would require \(e+1\) to be a multiple of \(2017\), impossible under \(q^e\le 10^{11}\). So primes congruent to \(1\) modulo \(2017\) contribute nothing.

Step 3: Separate the order-2 branch from the rest

If \(q\equiv -1 \pmod{2017}\), then \(\operatorname{ord}_{2017}(q)=2\), hence all odd exponents are formally bad.

For the target scale, however, only \(e=1\) can occur. The smallest prime with \(q\equiv -1 \pmod{2017}\) is \(12101\), and

$$12101^3>10^{11}.$$

So the order-2 branch contributes only prime powers of the form \(q^1\).

All remaining bad cases come from primes \(q\not\equiv \pm 1 \pmod{2017}\) with order at least \(3\), so their smallest bad exponent is at least \(2\). Consequently such primes matter only when \(q\le \sqrt{N}\), which is why this branch is searched only up to \(\lfloor \sqrt{N}\rfloor\).

Step 4: Turn each trigger pair into a closed-form weighted sum

For each \((q,e)\in B\), define the event

$$A_{q,e}=\{n\le N:\ v_q(n)=e\}.$$

Every integer in this event has the form \(n=q^e m\) with \(q\nmid m\) and

$$m\le M=\left\lfloor \frac{N}{q^e}\right\rfloor.$$

If \(T(x)=x(x+1)/2\) denotes the triangular-number sum, then

$$\sum_{n\in A_{q,e}} n =q^e\sum_{\substack{m\le M\\ q\nmid m}} m =q^e\left(T(M)-q\,T\!\left(\left\lfloor \frac{M}{q}\right\rfloor\right)\right).$$

This gives every single-event contribution in constant time once \(q^e\) is known.

Step 5: Use pairwise inclusion-exclusion, and stop there exactly

Events with the same prime \(q\) but different exponents are disjoint, so overlaps only arise from distinct primes. If \((q,e)\) and \((r,f)\) are trigger pairs with \(q\neq r\), then

$$n=q^e r^f m,\qquad q\nmid m,\qquad r\nmid m,$$

with

$$m\le X=\left\lfloor \frac{N}{q^e r^f}\right\rfloor.$$

Therefore

$$\sum_{n\in A_{q,e}\cap A_{r,f}} n =q^e r^f\left(T(X)-q\,T\!\left(\left\lfloor \frac{X}{q}\right\rfloor\right)-r\,T\!\left(\left\lfloor \frac{X}{r}\right\rfloor\right)+qr\,T\!\left(\left\lfloor \frac{X}{qr}\right\rfloor\right)\right).$$

Hence

$$S(N,2017)=\sum_{(q,e)\in B}\sum_{n\in A_{q,e}} n-\sum_{\substack{(q,e),(r,f)\in B\\ q<r}}\sum_{n\in A_{q,e}\cap A_{r,f}} n,$$

because triple intersections cannot occur for \(N=10^{11}\). The three smallest admissible trigger prime powers are

$$2311^2=5{,}340{,}721,\qquad 229^3=12{,}008{,}989,\qquad 3739^2=13{,}980{,}121,$$

and their product already exceeds \(10^{11}\) by many orders of magnitude. So second-order inclusion-exclusion is exact, not approximate.

Worked Example

Two concrete trigger prime powers show the mechanism.

First, \(12101\equiv -1 \pmod{2017}\), so order \(2\) gives the bad exponent \(e=1\), and indeed

$$\sigma(12101)=1+12101=12102=6\cdot 2017.$$

Second, \(2311\) has multiplicative order \(3\) modulo \(2017\), so \(e=2\) is bad, and

$$\sigma(2311^2)=1+2311+2311^2=5{,}343{,}033=2649\cdot 2017.$$

For \(N=10^{11}\), their overlap uses

$$P=12101\cdot 2311^2=64{,}628{,}064{,}821\le 10^{11},$$

so exactly one multiplier \(m=1\) is possible. The pair contribution is therefore precisely \(P\), and it must be subtracted once from the two single-event contributions.

How the Code Works

The implementation first generates all primes up to \(\lfloor \sqrt{N}\rfloor\). That prime list serves two purposes: it determines multiplicative orders modulo \(2017\) for the non-\(\pm 1\) branch, and it also drives a sieve on the arithmetic progression \(2017k-1\) to collect every prime congruent to \(-1\) modulo \(2017\) up to \(N\).

Next, the C++, Python, and Java implementations enumerate all trigger pairs \((q,e)\). For the order-2 branch they record only exponent \(1\). For the other branch they compute the multiplicative order as a divisor of \(2016\), list every exponent \(e=m\operatorname{ord}_{2017}(q)-1\) with \(q^e\le N\), and discard primes that cannot contribute.

Once the trigger list is known, the implementation accumulates all single-event formulas above. It then subtracts all admissible pair overlaps in three groups: order-2 with order-2, order-2 with other triggers, and other triggers with each other. The Python implementation delegates execution to the compiled solver, while the Java implementation mirrors the same arithmetic directly with arbitrary-precision integers.

The solver also verifies the published checkpoints

$$S(10^6,2017)=150850429,\qquad S(10^9,2017)=249652238344557,$$

before evaluating the final \(N=10^{11}\) case.

Complexity Analysis

Generating all primes up to \(\sqrt{N}\) costs \(O(\sqrt{N}\log\log N)\) time and \(O(\sqrt{N})\) memory. Sieving the progression \(2017k-1\) up to \(N\) uses an array of length \(\lfloor (N+1)/2017\rfloor\), so that stage needs \(O(N/2017)\) memory and about \(O((N/2017)\log\log N)\) marking work. The remaining order checks and inclusion-exclusion sums run only over the finite trigger list and its admissible pairs. For the target parameters, those later steps are smaller than the two sieve phases.

Footnotes and References

  1. Problem page: https://projecteuler.net/problem=565
  2. Divisor function: Wikipedia — Divisor function
  3. Multiplicative order: Wikipedia — Multiplicative order
  4. \(p\)-adic valuation: Wikipedia — \(p\)-adic valuation
  5. Inclusion-exclusion principle: Wikipedia — Inclusion-exclusion principle

Problemzusammenfassung

Für \(N=10^{11}\) und die Primzahl \(p=2017\) soll

$$S(N,p)=\sum_{\substack{1\le n\le N\\ p\mid \sigma(n)}} n$$

berechnet werden, wobei \(\sigma(n)\) die Summe aller positiven Teiler von \(n\) ist. Ein direkter Test aller \(n\le N\) wäre viel zu langsam. Die Lösung identifiziert stattdessen genau jene Primzahlpotenzen, deren Teilerfunktion durch \(2017\) teilbar ist, formuliert daraus eine gewichtete Vereinigungsmenge von Exakt-Bewertungs-Ereignissen und summiert diese dann mit Inklusions-Exklusion.

Mathematischer Ansatz

Schreibe die Primfaktorzerlegung von \(n\) als

$$n=\prod_q q^{a_q}.$$

Da die Teilerfunktion multiplikativ ist, gilt

$$\sigma(n)=\prod_q \sigma\!\left(q^{a_q}\right).$$

Also ist \(2017\mid \sigma(n)\) genau dann erfüllt, wenn es mindestens eine Primzahl \(q\) mit \(2017\mid \sigma(q^{a_q})\) gibt.

Schritt 1: Auf auslösende Primzahlpotenzen reduzieren

Definiere die Menge der Auslöserpaare

$$B=\left\{(q,e): e\ge 1,\ q^e\le N,\ 2017\mid \sigma(q^e)\right\},$$

wobei \(q\) prim ist.

Dann wird die Zielsumme zur gewichteten Vereinigung

$$S(N,2017)=\sum_{\substack{1\le n\le N\\ \exists (q,e)\in B:\ v_q(n)=e}} n,$$

wobei \(v_q(n)\) den Exponenten von \(q\) in \(n\) bezeichnet. Die exakte Bewertung ist entscheidend: Es geht um den genauen Exponenten \(a_q\), nicht nur darum, ob \(q^e\) teilt.

Die Primzahl \(q=2017\) trägt nie bei, denn

$$\sigma(2017^e)=1+2017+\cdots+2017^e\equiv 1 \pmod{2017}.$$

Schritt 2: Schlechte Exponenten über die multiplikative Ordnung bestimmen

Für eine Primzahl \(q\neq 2017\) gilt

$$\sigma(q^e)=1+q+\cdots+q^e=\frac{q^{e+1}-1}{q-1}.$$

Falls \(q\not\equiv 1 \pmod{2017}\), ist der Nenner modulo \(2017\) invertierbar, also

$$2017\mid \sigma(q^e)\iff q^{e+1}\equiv 1 \pmod{2017}.$$

Mit \(\operatorname{ord}_{2017}(q)\) als multiplikativer Ordnung von \(q\) modulo \(2017\) folgt

$$2017\mid \sigma(q^e)\iff \operatorname{ord}_{2017}(q)\mid (e+1),$$

und damit sind die schlechten Exponenten genau

$$e=m\operatorname{ord}_{2017}(q)-1,\qquad m\ge 1.$$

Da \(2017-1=2016=2^5\cdot 3^2\cdot 7\) ist, teilt jede solche Ordnung die Zahl \(2016\).

Die Restklasse \(q\equiv 1 \pmod{2017}\) ist anders, denn dann gilt

$$\sigma(q^e)\equiv e+1 \pmod{2017}.$$

Dafür müsste \(e+1\) ein Vielfaches von \(2017\) sein, was unter \(q^e\le 10^{11}\) nicht vorkommen kann. Solche Primzahlen liefern also keinen Beitrag.

Schritt 3: Den Ordnungs-2-Zweig separat behandeln

Ist \(q\equiv -1 \pmod{2017}\), so gilt \(\operatorname{ord}_{2017}(q)=2\), also sind formal alle ungeraden Exponenten schlecht.

Für die Zielgröße kann jedoch nur \(e=1\) auftreten. Die kleinste Primzahl mit \(q\equiv -1 \pmod{2017}\) ist \(12101\), und

$$12101^3>10^{11}.$$

Der Ordnungs-2-Zweig besteht daher nur aus Primzahlpotenzen der Form \(q^1\).

Alle übrigen schlechten Fälle stammen von Primzahlen \(q\not\equiv \pm 1 \pmod{2017}\) mit Ordnung mindestens \(3\). Dort ist der kleinste schlechte Exponent mindestens \(2\), weshalb solche Primzahlen nur für \(q\le \sqrt{N}\) relevant sind. Genau darum durchsucht die Implementierung diesen Zweig nur bis \(\lfloor \sqrt{N}\rfloor\).

Schritt 4: Jedes Auslöserpaar in eine geschlossene gewichtete Summe umwandeln

Für jedes \((q,e)\in B\) definieren wir das Ereignis

$$A_{q,e}=\{n\le N:\ v_q(n)=e\}.$$

Jede Zahl in diesem Ereignis hat die Form \(n=q^e m\) mit \(q\nmid m\) und

$$m\le M=\left\lfloor \frac{N}{q^e}\right\rfloor.$$

Mit der Dreieckszahlensumme \(T(x)=x(x+1)/2\) erhält man

$$\sum_{n\in A_{q,e}} n =q^e\sum_{\substack{m\le M\\ q\nmid m}} m =q^e\left(T(M)-q\,T\!\left(\left\lfloor \frac{M}{q}\right\rfloor\right)\right).$$

Damit ist jeder Einzelbeitrag in konstanter Zeit berechenbar, sobald \(q^e\) bekannt ist.

Schritt 5: Paarweise Inklusions-Exklusion und warum sie exakt endet

Ereignisse mit derselben Primzahl \(q\), aber unterschiedlichen Exponenten, sind disjunkt. Überlappungen entstehen also nur bei verschiedenen Primzahlen. Sind \((q,e)\) und \((r,f)\) zwei Auslöserpaare mit \(q\neq r\), dann hat jede Zahl im Schnitt die Form

$$n=q^e r^f m,\qquad q\nmid m,\qquad r\nmid m,$$

wobei

$$m\le X=\left\lfloor \frac{N}{q^e r^f}\right\rfloor.$$

Damit folgt

$$\sum_{n\in A_{q,e}\cap A_{r,f}} n =q^e r^f\left(T(X)-q\,T\!\left(\left\lfloor \frac{X}{q}\right\rfloor\right)-r\,T\!\left(\left\lfloor \frac{X}{r}\right\rfloor\right)+qr\,T\!\left(\left\lfloor \frac{X}{qr}\right\rfloor\right)\right).$$

Somit ist

$$S(N,2017)=\sum_{(q,e)\in B}\sum_{n\in A_{q,e}} n-\sum_{\substack{(q,e),(r,f)\in B\\ q<r}}\sum_{n\in A_{q,e}\cap A_{r,f}} n,$$

denn Dreifachschnitte können für \(N=10^{11}\) nicht auftreten. Die drei kleinsten zulässigen Auslöser-Primzahlpotenzen sind

$$2311^2=5{,}340{,}721,\qquad 229^3=12{,}008{,}989,\qquad 3739^2=13{,}980{,}121,$$

und ihr Produkt liegt bereits weit über \(10^{11}\). Die Inklusions-Exklusion zweiter Ordnung ist hier also exakt und keine Näherung.

Durchgerechnetes Beispiel

Zwei konkrete Auslöser-Primzahlpotenzen zeigen die Methode.

Erstens gilt \(12101\equiv -1 \pmod{2017}\), also liefert Ordnung \(2\) den schlechten Exponenten \(e=1\), und tatsächlich

$$\sigma(12101)=1+12101=12102=6\cdot 2017.$$

Zweitens hat \(2311\) die multiplikative Ordnung \(3\) modulo \(2017\), also ist \(e=2\) schlecht, und

$$\sigma(2311^2)=1+2311+2311^2=5{,}343{,}033=2649\cdot 2017.$$

Für \(N=10^{11}\) benutzt ihr Schnitt

$$P=12101\cdot 2311^2=64{,}628{,}064{,}821\le 10^{11},$$

sodass genau ein Multiplikator \(m=1\) möglich ist. Der Paarbeitrag ist also genau \(P\) und wird einmal von den beiden Einzelbeiträgen abgezogen.

Wie der Code arbeitet

Die Implementierung erzeugt zunächst alle Primzahlen bis \(\lfloor \sqrt{N}\rfloor\). Diese Liste hat zwei Aufgaben: Sie bestimmt die multiplikativen Ordnungen modulo \(2017\) im nicht-\(\pm 1\)-Zweig und sie treibt außerdem ein Sieb auf der arithmetischen Progression \(2017k-1\), um alle Primzahlen mit Rest \(-1\) modulo \(2017\) bis \(N\) zu sammeln.

Anschließend enumerieren die C++-, Python- und Java-Implementierungen alle Auslöserpaare \((q,e)\). Im Ordnungs-2-Zweig wird nur der Exponent \(1\) gespeichert. Im anderen Zweig wird die Ordnung als Teiler von \(2016\) bestimmt, danach werden alle Exponenten \(e=m\operatorname{ord}_{2017}(q)-1\) mit \(q^e\le N\) aufgelistet, und nicht beitragende Primzahlen werden verworfen.

Ist diese Liste bekannt, summiert die Implementierung alle Einzelereignisse mit den obigen Formeln und zieht danach alle zulässigen Paarüberschneidungen ab: Ordnungs-2 mit Ordnungs-2, Ordnungs-2 mit sonstigen Auslösern und sonstige Auslöser untereinander. Die Python-Implementierung delegiert die Ausführung an den kompilierten Solver, während die Java-Implementierung dieselbe Arithmetik direkt mit Ganzzahlen beliebiger Länge nachbildet.

Zusätzlich prüft der Solver die veröffentlichten Kontrollwerte

$$S(10^6,2017)=150850429,\qquad S(10^9,2017)=249652238344557,$$

bevor der endgültige Fall \(N=10^{11}\) berechnet wird.

Komplexitätsanalyse

Das Primzahlsieb bis \(\sqrt{N}\) benötigt \(O(\sqrt{N}\log\log N)\) Zeit und \(O(\sqrt{N})\) Speicher. Das Sieb über die Progression \(2017k-1\) bis \(N\) arbeitet mit einem Array der Länge \(\lfloor (N+1)/2017\rfloor\) und benötigt daher \(O(N/2017)\) Speicher sowie ungefähr \(O((N/2017)\log\log N)\) Markierungsarbeit. Die verbleibenden Ordnungsprüfungen und Inklusions-Exklusionssummen laufen nur über die endliche Auslöserliste und deren zulässige Paare. Für die Zielparameter sind diese Schritte kleiner als die beiden Siebphasen.

Fußnoten und Referenzen

  1. Problemseite: https://projecteuler.net/problem=565
  2. Teilerfunktion: Wikipedia — Divisor function
  3. Multiplikative Ordnung: Wikipedia — Multiplicative order
  4. \(p\)-adische Bewertung: Wikipedia — \(p\)-adic valuation
  5. Inklusions-Exklusions-Prinzip: Wikipedia — Inclusion-exclusion principle

Problem Özeti

\(N=10^{11}\) ve asal \(p=2017\) için amaç

$$S(N,p)=\sum_{\substack{1\le n\le N\\ p\mid \sigma(n)}} n$$

toplamını hesaplamaktır; burada \(\sigma(n)\), \(n\)'nin pozitif bölenleri toplamıdır. \(n\le N\) olan tüm sayılar için \(\sigma(n)\) değerini tek tek hesaplamak pratik değildir. Çözüm bunun yerine, hangi asal kuvvet üslerinin \(2017\) ile bölünebilirlik yarattığını tam olarak saptar, problemi tam değerlik olaylarının ağırlıklı birleşimi haline getirir ve sonra bunu içerme-dışlama ile toplar.

Matematiksel Yaklaşım

\(n\)'nin asal çarpanlara ayrılışını

$$n=\prod_q q^{a_q}$$

şeklinde yazalım. Bölen toplamı fonksiyonu çarpımsal olduğundan

$$\sigma(n)=\prod_q \sigma\!\left(q^{a_q}\right)$$

olur. Dolayısıyla \(2017\mid \sigma(n)\) koşulu, en az bir asal \(q\) için \(2017\mid \sigma(q^{a_q})\) olmasıyla eşdeğerdir.

Adım 1: Problemi tetikleyici asal kuvvetlere indirgeme

Tetikleyici çiftler kümesini

$$B=\left\{(q,e): e\ge 1,\ q^e\le N,\ 2017\mid \sigma(q^e)\right\}$$

ve burada \(q\) asaldır.

olarak tanımlayalım. O zaman hedef toplam

$$S(N,2017)=\sum_{\substack{1\le n\le N\\ \exists (q,e)\in B:\ v_q(n)=e}} n$$

biçimini alır. Burada \(v_q(n)\), \(q\)'nun \(n\) içindeki üssüdür. Tam değerlik kullanmak zorunludur; önemli olan yalnızca \(q^e\) bölmesi değil, tam olarak hangi üsse sahip olduğudur.

\(q=2017\) asla katkı vermez, çünkü

$$\sigma(2017^e)=1+2017+\cdots+2017^e\equiv 1 \pmod{2017}.$$

Adım 2: Kötü üsleri çarpımsal mertebe ile karakterize etme

\(q\neq 2017\) olan bir asal için

$$\sigma(q^e)=1+q+\cdots+q^e=\frac{q^{e+1}-1}{q-1}$$

yazılır. Eğer \(q\not\equiv 1 \pmod{2017}\) ise payda modulo \(2017\) terslenebilir, dolayısıyla

$$2017\mid \sigma(q^e)\iff q^{e+1}\equiv 1 \pmod{2017}.$$

\(q\)'nun modulo \(2017\) altındaki çarpımsal mertebesi \(\operatorname{ord}_{2017}(q)\) olsun. O zaman

$$2017\mid \sigma(q^e)\iff \operatorname{ord}_{2017}(q)\mid (e+1),$$

ve kötü üsler tam olarak

$$e=m\operatorname{ord}_{2017}(q)-1,\qquad m\ge 1$$

şeklindedir. \(2017-1=2016=2^5\cdot 3^2\cdot 7\) olduğundan her mertebe \(2016\)'yı böler.

\(q\equiv 1 \pmod{2017}\) sınıfı farklı davranır; bu durumda

$$\sigma(q^e)\equiv e+1 \pmod{2017}.$$

Bunun \(2017\)'ye bölünebilmesi için \(e+1\)'in \(2017\)'nin katı olması gerekir; bu ise \(q^e\le 10^{11}\) altında mümkün değildir. Bu yüzden \(1\) kalan sınıfındaki asallar katkı vermez.

Adım 3: Mertebesi 2 olan kolu ayrı ele alma

Eğer \(q\equiv -1 \pmod{2017}\) ise \(\operatorname{ord}_{2017}(q)=2\) olur; dolayısıyla teorik olarak tüm tek üsler kötüdür.

Fakat hedef ölçekte yalnızca \(e=1\) gerçekleşebilir. Çünkü \(q\equiv -1 \pmod{2017}\) koşulunu sağlayan en küçük asal \(12101\)'dir ve

$$12101^3>10^{11}.$$

Böylece bu kol yalnızca \(q^1\) biçimindeki asal kuvvetlerden oluşur.

Kalan tüm kötü durumlar \(q\not\equiv \pm 1 \pmod{2017}\) ve mertebesi en az \(3\) olan asallardan gelir. Bu durumda en küçük kötü üs en az \(2\) olduğu için böyle asallar ancak \(q\le \sqrt{N}\) iken anlamlıdır. Uygulamanın bu kolu yalnızca \(\lfloor \sqrt{N}\rfloor\)'ye kadar taramasının nedeni budur.

Adım 4: Her tetikleyici çifti kapalı biçimli ağırlıklı toplam haline getirme

Her \((q,e)\in B\) için

$$A_{q,e}=\{n\le N:\ v_q(n)=e\}$$

olayını tanımlayalım. Bu olaydaki her sayı \(n=q^e m\) biçimindedir; burada \(q\nmid m\) ve

$$m\le M=\left\lfloor \frac{N}{q^e}\right\rfloor.$$

\(T(x)=x(x+1)/2\) üçgensel toplamı ile

$$\sum_{n\in A_{q,e}} n =q^e\sum_{\substack{m\le M\\ q\nmid m}} m =q^e\left(T(M)-q\,T\!\left(\left\lfloor \frac{M}{q}\right\rfloor\right)\right)$$

elde edilir. Böylece \(q^e\) bilindiğinde her tekil katkı sabit zamanda hesaplanır.

Adım 5: Çiftli içerme-dışlama ve neden burada tam olarak durduğu

Aynı asal için farklı üsler veren olaylar ayrık olduğundan, kesişimler ancak farklı asallardan gelir. \((q,e)\) ve \((r,f)\) iki tetikleyici çift ve \(q\neq r\) olsun. O zaman ortak olaydaki her sayı

$$n=q^e r^f m,\qquad q\nmid m,\qquad r\nmid m$$

şeklindedir ve

$$m\le X=\left\lfloor \frac{N}{q^e r^f}\right\rfloor.$$

Bu durumda

$$\sum_{n\in A_{q,e}\cap A_{r,f}} n =q^e r^f\left(T(X)-q\,T\!\left(\left\lfloor \frac{X}{q}\right\rfloor\right)-r\,T\!\left(\left\lfloor \frac{X}{r}\right\rfloor\right)+qr\,T\!\left(\left\lfloor \frac{X}{qr}\right\rfloor\right)\right)$$

olur. Dolayısıyla

$$S(N,2017)=\sum_{(q,e)\in B}\sum_{n\in A_{q,e}} n-\sum_{\substack{(q,e),(r,f)\in B\\ q<r}}\sum_{n\in A_{q,e}\cap A_{r,f}} n,$$

çünkü \(N=10^{11}\) için üçlü kesişimler oluşamaz. En küçük üç uygun tetikleyici asal kuvvet

$$2311^2=5{,}340{,}721,\qquad 229^3=12{,}008{,}989,\qquad 3739^2=13{,}980{,}121$$

olup bunların çarpımı zaten \(10^{11}\)'den çok büyüktür. Yani ikinci mertebe içerme-dışlama burada yaklaşım değil, tam sonuç verir.

Çalışılmış Örnek

Yöntem iki somut tetikleyici asal kuvvet üzerinde görülebilir.

İlk olarak \(12101\equiv -1 \pmod{2017}\) olduğundan mertebe \(2\) bize kötü üs \(e=1\)'i verir ve gerçekten

$$\sigma(12101)=1+12101=12102=6\cdot 2017.$$

İkinci olarak \(2311\)'in modulo \(2017\) altındaki mertebesi \(3\)'tür; dolayısıyla \(e=2\) kötüdür ve

$$\sigma(2311^2)=1+2311+2311^2=5{,}343{,}033=2649\cdot 2017.$$

\(N=10^{11}\) için bu iki olayın ortak çarpanı

$$P=12101\cdot 2311^2=64{,}628{,}064{,}821\le 10^{11}$$

olduğundan yalnızca \(m=1\) mümkündür. Bu yüzden çift katkısı tam olarak \(P\)'dir ve tekil katkıların toplamından bir kez çıkarılır.

Kod Nasıl Çalışır

Uygulama önce \(\lfloor \sqrt{N}\rfloor\)'ye kadar tüm asalları üretir. Bu liste iki iş yapar: bir yandan \(\pm 1\) dışında kalan dal için modulo \(2017\) mertebelerini belirler, diğer yandan \(2017k-1\) aritmetik dizisi üzerinde yürütülen bir elek ile \(-1\) kalan sınıfındaki tüm asalları \(N\)'ye kadar toplar.

Sonra C++, Python ve Java uygulamaları tüm tetikleyici çiftleri \((q,e)\) çıkarır. Mertebe-2 kolunda sadece üs \(1\) kaydedilir. Diğer kolda mertebe, \(2016\)'nın bir böleni olarak hesaplanır; ardından \(q^e\le N\) koşulunu sağlayan tüm \(e=m\operatorname{ord}_{2017}(q)-1\) değerleri listelenir ve katkı veremeyen asallar elenir.

Bu liste hazır olduktan sonra uygulama önce tüm tekil olay toplamlarını ekler, sonra da izin verilen tüm çift kesişimleri çıkarır: mertebe-2 ile mertebe-2, mertebe-2 ile diğerleri ve diğerleri kendi aralarında. Python uygulaması çalışmayı derlenmiş çözücüye devreder; Java uygulaması ise aynı aritmetiği doğrudan keyfi büyüklükte tamsayılarla yürütür.

Çözücü ayrıca yayımlanmış kontrol değerlerini de doğrular:

$$S(10^6,2017)=150850429,\qquad S(10^9,2017)=249652238344557,$$

ve ardından son \(N=10^{11}\) hesabını yapar.

Karmaşıklık Analizi

\(\sqrt{N}\)'ye kadar asal üretimi \(O(\sqrt{N}\log\log N)\) zaman ve \(O(\sqrt{N})\) bellek gerektirir. \(2017k-1\) dizisi üzerindeki elek, uzunluğu \(\lfloor (N+1)/2017\rfloor\) olan bir dizi kullanır; dolayısıyla bu aşama \(O(N/2017)\) bellek ve yaklaşık \(O((N/2017)\log\log N)\) iş ister. Geriye kalan mertebe testleri ve içerme-dışlama toplamları yalnızca sonlu tetikleyici listesi ve izin verilen çiftleri üzerinde çalışır. Hedef parametrelerde bu kısım iki elek aşamasından daha küçüktür.

Dipnotlar ve Kaynaklar

  1. Problem sayfası: https://projecteuler.net/problem=565
  2. Bölen toplamı fonksiyonu: Wikipedia — Divisor function
  3. Çarpımsal mertebe: Wikipedia — Multiplicative order
  4. \(p\)-adik değerlik: Wikipedia — \(p\)-adic valuation
  5. İçerme-dışlama ilkesi: Wikipedia — Inclusion-exclusion principle

Resumen del Problema

Para \(N=10^{11}\) y el primo \(p=2017\), hay que calcular

$$S(N,p)=\sum_{\substack{1\le n\le N\\ p\mid \sigma(n)}} n,$$

donde \(\sigma(n)\) es la suma de los divisores positivos de \(n\). Probar todos los \(n\le N\) sería inviable. La solución identifica exactamente qué exponentes de potencias primas obligan a que \(\sigma(n)\) sea divisible por \(2017\), reescribe el problema como una unión ponderada de eventos de valuación exacta y suma esa unión mediante inclusión-exclusión.

Enfoque Matemático

Escribamos la factorización prima de \(n\) como

$$n=\prod_q q^{a_q}.$$

Como la función suma de divisores es multiplicativa,

$$\sigma(n)=\prod_q \sigma\!\left(q^{a_q}\right).$$

Por tanto, \(2017\mid \sigma(n)\) si y solo si existe al menos un primo \(q\) tal que \(2017\mid \sigma(q^{a_q})\).

Paso 1: Reducir el problema a potencias primas disparadoras

Definimos el conjunto de pares disparadores

$$B=\left\{(q,e): e\ge 1,\ q^e\le N,\ 2017\mid \sigma(q^e)\right\},$$

donde \(q\) es primo.

Entonces la suma buscada se convierte en la unión ponderada

$$S(N,2017)=\sum_{\substack{1\le n\le N\\ \exists (q,e)\in B:\ v_q(n)=e}} n,$$

donde \(v_q(n)\) es el exponente de \(q\) en la descomposición de \(n\). La valuación exacta es la cantidad correcta, porque lo que importa es el exponente preciso \(a_q\), no solo que \(q^e\) divida a \(n\).

El primo \(q=2017\) nunca contribuye, ya que

$$\sigma(2017^e)=1+2017+\cdots+2017^e\equiv 1 \pmod{2017}.$$

Paso 2: Describir los exponentes malos con el orden multiplicativo

Para un primo \(q\neq 2017\),

$$\sigma(q^e)=1+q+\cdots+q^e=\frac{q^{e+1}-1}{q-1}.$$

Si \(q\not\equiv 1 \pmod{2017}\), el denominador es invertible módulo \(2017\), así que

$$2017\mid \sigma(q^e)\iff q^{e+1}\equiv 1 \pmod{2017}.$$

Sea \(\operatorname{ord}_{2017}(q)\) el orden multiplicativo de \(q\) módulo \(2017\). Entonces

$$2017\mid \sigma(q^e)\iff \operatorname{ord}_{2017}(q)\mid (e+1),$$

de modo que los exponentes malos son exactamente

$$e=m\operatorname{ord}_{2017}(q)-1,\qquad m\ge 1.$$

Como \(2017-1=2016=2^5\cdot 3^2\cdot 7\), todo orden posible divide a \(2016\).

La clase \(q\equiv 1 \pmod{2017}\) se trata aparte, porque entonces

$$\sigma(q^e)\equiv e+1 \pmod{2017}.$$

Eso exigiría que \(e+1\) fuera múltiplo de \(2017\), algo imposible bajo \(q^e\le 10^{11}\). Por lo tanto, esos primos no aportan nada.

Paso 3: Separar la rama de orden 2

Si \(q\equiv -1 \pmod{2017}\), entonces \(\operatorname{ord}_{2017}(q)=2\), y en principio todos los exponentes impares son malos.

Sin embargo, para la escala del problema solo puede aparecer \(e=1\). El menor primo con \(q\equiv -1 \pmod{2017}\) es \(12101\), y

$$12101^3>10^{11}.$$

Así, la rama de orden 2 solo aporta potencias primas de la forma \(q^1\).

Todos los demás casos malos provienen de primos \(q\not\equiv \pm 1 \pmod{2017}\) con orden al menos \(3\), por lo que su menor exponente malo es al menos \(2\). En consecuencia, esos primos solo importan cuando \(q\le \sqrt{N}\), y por eso la implementación busca esa rama solo hasta \(\lfloor \sqrt{N}\rfloor\).

Paso 4: Convertir cada par disparador en una suma ponderada cerrada

Para cada \((q,e)\in B\), definimos el evento

$$A_{q,e}=\{n\le N:\ v_q(n)=e\}.$$

Todo entero en este evento tiene la forma \(n=q^e m\), con \(q\nmid m\) y

$$m\le M=\left\lfloor \frac{N}{q^e}\right\rfloor.$$

Si \(T(x)=x(x+1)/2\) es la suma triangular, entonces

$$\sum_{n\in A_{q,e}} n =q^e\sum_{\substack{m\le M\\ q\nmid m}} m =q^e\left(T(M)-q\,T\!\left(\left\lfloor \frac{M}{q}\right\rfloor\right)\right).$$

Así se obtiene cada contribución individual en tiempo constante una vez conocido \(q^e\).

Paso 5: Inclusión-exclusión por pares, y por qué ahí termina exactamente

Los eventos con el mismo primo \(q\) pero distintos exponentes son disjuntos, de modo que las intersecciones solo aparecen entre primos distintos. Si \((q,e)\) y \((r,f)\) son pares disparadores con \(q\neq r\), entonces

$$n=q^e r^f m,\qquad q\nmid m,\qquad r\nmid m,$$

con

$$m\le X=\left\lfloor \frac{N}{q^e r^f}\right\rfloor.$$

Por consiguiente,

$$\sum_{n\in A_{q,e}\cap A_{r,f}} n =q^e r^f\left(T(X)-q\,T\!\left(\left\lfloor \frac{X}{q}\right\rfloor\right)-r\,T\!\left(\left\lfloor \frac{X}{r}\right\rfloor\right)+qr\,T\!\left(\left\lfloor \frac{X}{qr}\right\rfloor\right)\right).$$

Luego

$$S(N,2017)=\sum_{(q,e)\in B}\sum_{n\in A_{q,e}} n-\sum_{\substack{(q,e),(r,f)\in B\\ q<r}}\sum_{n\in A_{q,e}\cap A_{r,f}} n,$$

porque para \(N=10^{11}\) no pueden existir intersecciones triples. Las tres menores potencias primas disparadoras son

$$2311^2=5{,}340{,}721,\qquad 229^3=12{,}008{,}989,\qquad 3739^2=13{,}980{,}121,$$

y su producto ya supera \(10^{11}\) por un margen enorme. Por eso la inclusión-exclusión de segundo orden es exacta.

Ejemplo Trabajado

Dos potencias primas concretas muestran el mecanismo.

Primero, \(12101\equiv -1 \pmod{2017}\), así que el orden \(2\) produce el exponente malo \(e=1\), y en efecto

$$\sigma(12101)=1+12101=12102=6\cdot 2017.$$

Segundo, \(2311\) tiene orden multiplicativo \(3\) módulo \(2017\), por lo que \(e=2\) es malo, y

$$\sigma(2311^2)=1+2311+2311^2=5{,}343{,}033=2649\cdot 2017.$$

Para \(N=10^{11}\), su solapamiento usa

$$P=12101\cdot 2311^2=64{,}628{,}064{,}821\le 10^{11},$$

de modo que solo es posible el multiplicador \(m=1\). La contribución del par es exactamente \(P\), y por ello se resta una vez de la suma de los dos eventos simples.

Cómo Funciona el Código

La implementación genera primero todos los primos hasta \(\lfloor \sqrt{N}\rfloor\). Esa lista se usa con dos fines: calcular órdenes multiplicativos módulo \(2017\) en la rama distinta de \(\pm 1\), y alimentar una criba sobre la progresión aritmética \(2017k-1\) para reunir todos los primos congruentes con \(-1\) módulo \(2017\) hasta \(N\).

Después, las implementaciones en C++, Python y Java enumeran todos los pares disparadores \((q,e)\). En la rama de orden 2 solo se guarda el exponente \(1\). En la otra rama se calcula el orden como divisor de \(2016\), se listan todos los exponentes \(e=m\operatorname{ord}_{2017}(q)-1\) con \(q^e\le N\), y se descartan los primos que no pueden aportar.

Una vez construida esa lista, la implementación acumula todas las fórmulas de eventos simples y luego resta todos los solapamientos admisibles por pares: orden 2 con orden 2, orden 2 con otros disparadores y otros disparadores entre sí. La implementación en Python delega la ejecución al solucionador compilado, mientras que la de Java reproduce la misma aritmética con enteros de precisión arbitraria.

El solucionador también verifica los puntos de control publicados

$$S(10^6,2017)=150850429,\qquad S(10^9,2017)=249652238344557,$$

antes de evaluar el caso final \(N=10^{11}\).

Análisis de Complejidad

Generar todos los primos hasta \(\sqrt{N}\) cuesta \(O(\sqrt{N}\log\log N)\) tiempo y \(O(\sqrt{N})\) memoria. La criba sobre la progresión \(2017k-1\) hasta \(N\) usa un arreglo de longitud \(\lfloor (N+1)/2017\rfloor\), por lo que requiere \(O(N/2017)\) memoria y aproximadamente \(O((N/2017)\log\log N)\) trabajo de marcado. Las comprobaciones de orden y las sumas de inclusión-exclusión restantes recorren solo la lista finita de disparadores y sus pares admisibles. Para los parámetros objetivo, esas fases son menores que las dos cribas.

Notas y Referencias

  1. Página del problema: https://projecteuler.net/problem=565
  2. Función suma de divisores: Wikipedia — Divisor function
  3. Orden multiplicativo: Wikipedia — Multiplicative order
  4. Valuación \(p\)-ádica: Wikipedia — \(p\)-adic valuation
  5. Principio de inclusión-exclusión: Wikipedia — Inclusion-exclusion principle

问题概述

当 \(N=10^{11}\)、素数 \(p=2017\) 时,需要计算

$$S(N,p)=\sum_{\substack{1\le n\le N\\ p\mid \sigma(n)}} n,$$

其中 \(\sigma(n)\) 表示 \(n\) 的所有正因子之和。若对每个 \(n\le N\) 都单独求一次 \(\sigma(n)\),规模上完全不可行。高效做法不是从整数 \(n\) 出发,而是从它的质因数分解出发,找出哪些“质数幂指数”会强制 \(\sigma(n)\) 被 \(2017\) 整除,再把问题改写成若干精确赋值事件的加权并集,最后用容斥原理求和。

数学方法

把 \(n\) 的质因数分解写成

$$n=\prod_q q^{a_q}.$$

因为除数和函数是乘法函数,所以

$$\sigma(n)=\prod_q \sigma\!\left(q^{a_q}\right).$$

因此 \(2017\mid \sigma(n)\) 当且仅当至少存在一个质数 \(q\),使得对应的质数幂因子满足 \(2017\mid \sigma(q^{a_q})\)。这一步把原题从“看整个整数”变成了“看每个质数的指数”。

步骤 1:先把问题化成“触发质数幂”的并集

定义触发对的集合

$$B=\left\{(q,e): e\ge 1,\ q^e\le N,\ 2017\mid \sigma(q^e)\right\},$$

其中 \(q\) 为质数。

这样一来,目标和就可以写成

$$S(N,2017)=\sum_{\substack{1\le n\le N\\ \exists (q,e)\in B:\ v_q(n)=e}} n,$$

其中 \(v_q(n)\) 表示质数 \(q\) 在 \(n\) 的分解中出现的次数。这里必须使用“精确指数”而不是“至少被 \(q^e\) 整除”,因为 \(\sigma(q^{a_q})\) 是否被 \(2017\) 整除,取决于准确的指数 \(a_q\)。

特别地,\(q=2017\) 本身永远不会触发,因为

$$\sigma(2017^e)=1+2017+\cdots+2017^e\equiv 1 \pmod{2017}.$$

步骤 2:用乘法阶刻画所有“坏指数”

对任意质数 \(q\neq 2017\),有

$$\sigma(q^e)=1+q+\cdots+q^e=\frac{q^{e+1}-1}{q-1}.$$

如果 \(q\not\equiv 1 \pmod{2017}\),那么分母 \(q-1\) 在模 \(2017\) 意义下可逆,于是

$$2017\mid \sigma(q^e)\iff q^{e+1}\equiv 1 \pmod{2017}.$$

设 \(\operatorname{ord}_{2017}(q)\) 表示 \(q\) 在模 \(2017\) 下的乘法阶,则上式等价于

$$2017\mid \sigma(q^e)\iff \operatorname{ord}_{2017}(q)\mid (e+1).$$

因此所有坏指数恰好是

$$e=m\operatorname{ord}_{2017}(q)-1,\qquad m\ge 1.$$

又因为 \(2017\) 是素数,且

$$2017-1=2016=2^5\cdot 3^2\cdot 7,$$

所以每个乘法阶都必须整除 \(2016\)。

而当 \(q\equiv 1 \pmod{2017}\) 时,情况不同,因为

$$\sigma(q^e)\equiv e+1 \pmod{2017}.$$

要让它被 \(2017\) 整除,就必须有 \(e+1\) 是 \(2017\) 的倍数。但在 \(q^e\le 10^{11}\) 的约束下,不可能出现这么大的指数,所以这一类质数根本不会贡献答案。

步骤 3:把阶为 2 的分支单独拿出来

若 \(q\equiv -1 \pmod{2017}\),那么 \(\operatorname{ord}_{2017}(q)=2\),从形式上看,所有奇数指数都是坏指数。

不过在本题的目标规模里,真正可能出现的只有 \(e=1\)。原因是满足 \(q\equiv -1 \pmod{2017}\) 的最小质数是 \(12101\),而

$$12101^3>10^{11}.$$

因此在这一分支中,任何更大的奇数指数都会让 \(q^e\) 超过上界,所以只需要考虑 \(q^1\)。

其余坏情况都来自 \(q\not\equiv \pm 1 \pmod{2017}\) 且乘法阶至少为 \(3\) 的质数。此时最小坏指数至少为 \(2\),于是只有 \(q\le \sqrt{N}\) 时才可能满足 \(q^e\le N\)。这正是实现中只把这一分支搜索到 \(\lfloor \sqrt{N}\rfloor\) 的原因。

步骤 4:把每个触发事件写成闭式求和

对每个 \((q,e)\in B\),定义事件

$$A_{q,e}=\{n\le N:\ v_q(n)=e\}.$$

事件中的每个整数都可以唯一写成

$$n=q^e m,\qquad q\nmid m,\qquad m\le M=\left\lfloor \frac{N}{q^e}\right\rfloor.$$

如果记三角和

$$T(x)=\frac{x(x+1)}{2},$$

那么单个事件的加权贡献就是

$$\sum_{n\in A_{q,e}} n =q^e\sum_{\substack{m\le M\\ q\nmid m}} m =q^e\left(T(M)-q\,T\!\left(\left\lfloor \frac{M}{q}\right\rfloor\right)\right).$$

这里的思路很直接:先把 \(m\le M\) 的总和全部算出来,再减去那些仍然被 \(q\) 整除的 \(m\)。这样每个单事件都能在常数时间内求出。

步骤 5:对不同质数做二阶容斥,而且到二阶就精确结束

同一个质数 \(q\) 对应的不同指数事件彼此不相交,所以真正的重叠只会来自不同质数。若 \((q,e)\) 与 \((r,f)\) 是两个触发对,且 \(q\neq r\),那么同时满足这两个事件的整数只能写成

$$n=q^e r^f m,\qquad q\nmid m,\qquad r\nmid m,$$

并且

$$m\le X=\left\lfloor \frac{N}{q^e r^f}\right\rfloor.$$

因此二重交的加权和为

$$\sum_{n\in A_{q,e}\cap A_{r,f}} n =q^e r^f\left(T(X)-q\,T\!\left(\left\lfloor \frac{X}{q}\right\rfloor\right)-r\,T\!\left(\left\lfloor \frac{X}{r}\right\rfloor\right)+qr\,T\!\left(\left\lfloor \frac{X}{qr}\right\rfloor\right)\right).$$

于是可以直接写出

$$S(N,2017)=\sum_{(q,e)\in B}\sum_{n\in A_{q,e}} n-\sum_{\substack{(q,e),(r,f)\in B\\ q<r}}\sum_{n\in A_{q,e}\cap A_{r,f}} n.$$

之所以没有三重交项,是因为在 \(N=10^{11}\) 下它们根本不可能出现。最小的三个触发质数幂分别是

$$2311^2=5{,}340{,}721,\qquad 229^3=12{,}008{,}989,\qquad 3739^2=13{,}980{,}121,$$

它们的乘积远远大于 \(10^{11}\)。所以这里只做到二阶容斥不是近似,而是严格正确。

例子:两个实际触发项如何产生一项交集

先看 \(12101\)。因为

$$12101\equiv -1 \pmod{2017},$$

它属于阶为 \(2\) 的分支,于是 \(e=1\) 是坏指数,并且

$$\sigma(12101)=1+12101=12102=6\cdot 2017.$$

再看 \(2311\)。它在模 \(2017\) 下的乘法阶是 \(3\),因此 \(e=2\) 是坏指数,而且

$$\sigma(2311^2)=1+2311+2311^2=5{,}343{,}033=2649\cdot 2017.$$

把这两个事件同时放在 \(N=10^{11}\) 下考察,公共因子为

$$P=12101\cdot 2311^2=64{,}628{,}064{,}821\le 10^{11}.$$

这时只剩下唯一的乘子 \(m=1\),所以这一对事件的交集贡献恰好就是 \(P\) 本身。容斥时必须把这一项从两个单事件贡献之和中减掉一次。

代码如何工作

实现的第一步是筛出所有不超过 \(\lfloor \sqrt{N}\rfloor\) 的质数。这一份质数表承担两个角色:一方面,它用于处理 \(q\not\equiv \pm 1 \pmod{2017}\) 时的乘法阶计算;另一方面,它也被用来在等差数列 \(2017k-1\) 上做筛法,从而枚举出所有不超过 \(N\) 且满足 \(q\equiv -1 \pmod{2017}\) 的质数。

接着,C++、Python 和 Java 实现都会枚举全部触发对 \((q,e)\)。阶为 \(2\) 的分支只保存指数 \(1\)。另一分支则把乘法阶当作 \(2016\) 的一个因子来求,再列出所有满足 \(e=m\operatorname{ord}_{2017}(q)-1\) 且 \(q^e\le N\) 的指数,不能触发的质数直接跳过。

有了这张触发表以后,实现先把所有单事件公式累加起来,再按三类情况减去合法的二重交:阶为 \(2\) 与阶为 \(2\)、阶为 \(2\) 与其他触发项、以及其他触发项彼此之间。Python 实现负责调用编译后的求解器,Java 实现则用任意精度整数直接复现同样的算术过程。

求解器还会先检查题目给出的两个校验值:

$$S(10^6,2017)=150850429,\qquad S(10^9,2017)=249652238344557,$$

然后才计算最终的 \(N=10^{11}\) 情形。

复杂度分析

筛出所有不超过 \(\sqrt{N}\) 的质数需要 \(O(\sqrt{N}\log\log N)\) 时间和 \(O(\sqrt{N})\) 空间。对等差数列 \(2017k-1\) 的筛法会使用长度为 \(\lfloor (N+1)/2017\rfloor\) 的数组,因此这一阶段需要 \(O(N/2017)\) 空间,以及大约 \(O((N/2017)\log\log N)\) 的标记工作。后续的乘法阶判断与容斥求和只在有限个触发项及其可行的二重交上运行。对题目的目标规模来说,这些后续步骤都小于前面的两次筛法。

脚注与参考资料

  1. 题目页面:https://projecteuler.net/problem=565
  2. 除数和函数:Wikipedia — Divisor function
  3. 乘法阶:Wikipedia — Multiplicative order
  4. \(p\)-进赋值:Wikipedia — \(p\)-adic valuation
  5. 容斥原理:Wikipedia — Inclusion-exclusion principle

Краткое описание задачи

Для \(N=10^{11}\) и простого числа \(p=2017\) требуется вычислить

$$S(N,p)=\sum_{\substack{1\le n\le N\\ p\mid \sigma(n)}} n,$$

где \(\sigma(n)\) означает сумму положительных делителей числа \(n\). Перебирать все \(n\le N\) и каждый раз считать \(\sigma(n)\) слишком дорого. Поэтому решение ищет не сами числа, а те простые степени, которые гарантируют делимость \(\sigma(n)\) на \(2017\), затем переписывает задачу как взвешенное объединение событий точной \(q\)-адической оценки и суммирует его по принципу включения-исключения.

Математический подход

Пусть разложение числа \(n\) на простые множители имеет вид

$$n=\prod_q q^{a_q}.$$

Поскольку функция суммы делителей мультипликативна, имеем

$$\sigma(n)=\prod_q \sigma\!\left(q^{a_q}\right).$$

Следовательно, условие \(2017\mid \sigma(n)\) эквивалентно существованию хотя бы одного простого \(q\), для которого выполнено \(2017\mid \sigma(q^{a_q})\).

Шаг 1: Свести задачу к «запускающим» простым степеням

Введем множество запускающих пар

$$B=\left\{(q,e): e\ge 1,\ q^e\le N,\ 2017\mid \sigma(q^e)\right\},$$

где \(q\) — простое число.

Тогда искомая сумма равна взвешенному объединению

$$S(N,2017)=\sum_{\substack{1\le n\le N\\ \exists (q,e)\in B:\ v_q(n)=e}} n,$$

где \(v_q(n)\) обозначает показатель простого \(q\) в разложении числа \(n\). Здесь нужна именно точная оценка: важно не просто деление на \(q^e\), а точное значение показателя \(a_q\).

Само простое \(q=2017\) никогда не дает вклад, потому что

$$\sigma(2017^e)=1+2017+\cdots+2017^e\equiv 1 \pmod{2017}.$$

Шаг 2: Описать плохие показатели через мультипликативный порядок

Для простого \(q\neq 2017\) имеем

$$\sigma(q^e)=1+q+\cdots+q^e=\frac{q^{e+1}-1}{q-1}.$$

Если \(q\not\equiv 1 \pmod{2017}\), знаменатель обратим по модулю \(2017\), поэтому

$$2017\mid \sigma(q^e)\iff q^{e+1}\equiv 1 \pmod{2017}.$$

Обозначим через \(\operatorname{ord}_{2017}(q)\) мультипликативный порядок \(q\) по модулю \(2017\). Тогда

$$2017\mid \sigma(q^e)\iff \operatorname{ord}_{2017}(q)\mid (e+1),$$

а значит, плохие показатели имеют вид

$$e=m\operatorname{ord}_{2017}(q)-1,\qquad m\ge 1.$$

Так как \(2017-1=2016=2^5\cdot 3^2\cdot 7\), любой такой порядок делит \(2016\).

Класс \(q\equiv 1 \pmod{2017}\) надо рассматривать отдельно, поскольку в этом случае

$$\sigma(q^e)\equiv e+1 \pmod{2017}.$$

Для делимости на \(2017\) пришлось бы требовать, чтобы \(e+1\) было кратно \(2017\), но при условии \(q^e\le 10^{11}\) такое невозможно. Значит, простые из этого класса не участвуют.

Шаг 3: Отделить ветвь порядка 2

Если \(q\equiv -1 \pmod{2017}\), то \(\operatorname{ord}_{2017}(q)=2\), и формально все нечетные показатели являются плохими.

Однако при масштабе задачи реально возможен только случай \(e=1\). Наименьшее простое с условием \(q\equiv -1 \pmod{2017}\) равно \(12101\), и

$$12101^3>10^{11}.$$

Поэтому ветвь порядка \(2\) состоит только из простых степеней вида \(q^1\).

Все остальные плохие случаи исходят от простых \(q\not\equiv \pm 1 \pmod{2017}\) с порядком не меньше \(3\), а значит их минимальный плохой показатель не меньше \(2\). Следовательно, такие простые важны только при \(q\le \sqrt{N}\), и именно поэтому реализация ищет эту ветвь лишь до \(\lfloor \sqrt{N}\rfloor\).

Шаг 4: Превратить каждую запускающую пару в замкнутую формулу для суммы

Для каждой пары \((q,e)\in B\) введем событие

$$A_{q,e}=\{n\le N:\ v_q(n)=e\}.$$

Любое число из этого события имеет вид

$$n=q^e m,\qquad q\nmid m,\qquad m\le M=\left\lfloor \frac{N}{q^e}\right\rfloor.$$

Если обозначить через \(T(x)=x(x+1)/2\) сумму первых \(x\) натуральных чисел, то

$$\sum_{n\in A_{q,e}} n =q^e\sum_{\substack{m\le M\\ q\nmid m}} m =q^e\left(T(M)-q\,T\!\left(\left\lfloor \frac{M}{q}\right\rfloor\right)\right).$$

Итак, вклад каждого отдельного события вычисляется за константное время после того, как известна величина \(q^e\).

Шаг 5: Попарное включение-исключение и точное обрывание на этом месте

События с одним и тем же простым \(q\), но с разными показателями, попарно несовместимы, поэтому настоящие пересечения возникают только для разных простых. Если \((q,e)\) и \((r,f)\) — две запускающие пары с \(q\neq r\), то числа в их пересечении имеют вид

$$n=q^e r^f m,\qquad q\nmid m,\qquad r\nmid m,$$

где

$$m\le X=\left\lfloor \frac{N}{q^e r^f}\right\rfloor.$$

Следовательно,

$$\sum_{n\in A_{q,e}\cap A_{r,f}} n =q^e r^f\left(T(X)-q\,T\!\left(\left\lfloor \frac{X}{q}\right\rfloor\right)-r\,T\!\left(\left\lfloor \frac{X}{r}\right\rfloor\right)+qr\,T\!\left(\left\lfloor \frac{X}{qr}\right\rfloor\right)\right).$$

Поэтому

$$S(N,2017)=\sum_{(q,e)\in B}\sum_{n\in A_{q,e}} n-\sum_{\substack{(q,e),(r,f)\in B\\ q<r}}\sum_{n\in A_{q,e}\cap A_{r,f}} n,$$

поскольку тройные пересечения при \(N=10^{11}\) невозможны. Три наименьшие допустимые запускающие простые степени равны

$$2311^2=5{,}340{,}721,\qquad 229^3=12{,}008{,}989,\qquad 3739^2=13{,}980{,}121,$$

и их произведение уже намного больше \(10^{11}\). Значит, учет только одиночных и парных членов здесь дает точную формулу.

Разобранный пример

Рассмотрим два конкретных запускающих случая.

Во-первых, \(12101\equiv -1 \pmod{2017}\), значит порядок равен \(2\), плохой показатель равен \(e=1\), и действительно

$$\sigma(12101)=1+12101=12102=6\cdot 2017.$$

Во-вторых, число \(2311\) имеет мультипликативный порядок \(3\) по модулю \(2017\), поэтому \(e=2\) — плохой показатель, и

$$\sigma(2311^2)=1+2311+2311^2=5{,}343{,}033=2649\cdot 2017.$$

При \(N=10^{11}\) их совместный множитель равен

$$P=12101\cdot 2311^2=64{,}628{,}064{,}821\le 10^{11},$$

поэтому возможен ровно один множитель \(m=1\). Значит, вклад пересечения равен точно \(P\), и его нужно один раз вычесть из суммы двух одиночных вкладов.

Как работает код

Сначала реализация генерирует все простые числа до \(\lfloor \sqrt{N}\rfloor\). Этот список используется в двух местах: он помогает вычислять мультипликативные порядки по модулю \(2017\) в ветви \(q\not\equiv \pm 1 \pmod{2017}\), а также запускает решето по арифметической прогрессии \(2017k-1\), чтобы перечислить все простые числа, сравнимые с \(-1\) по модулю \(2017\), вплоть до \(N\).

Затем реализации на C++, Python и Java перечисляют все запускающие пары \((q,e)\). Для ветви порядка \(2\) сохраняется только показатель \(1\). Для остальных простых порядок вычисляется как делитель числа \(2016\), затем перечисляются все показатели \(e=m\operatorname{ord}_{2017}(q)-1\) с условием \(q^e\le N\), а простые, которые не могут дать вклад, отбрасываются.

После построения этого списка код суммирует все одиночные события по формулам выше и затем вычитает все допустимые парные пересечения: ветвь порядка \(2\) с самой собой, ветвь порядка \(2\) с остальными триггерами и остальные триггеры между собой. Python-реализация передает выполнение компилируемому решателю, тогда как Java-реализация повторяет те же вычисления напрямую с целыми числами произвольной длины.

Перед окончательным вычислением также проверяются опубликованные контрольные значения

$$S(10^6,2017)=150850429,\qquad S(10^9,2017)=249652238344557,$$

после чего считается конечный случай \(N=10^{11}\).

Анализ сложности

Построение решета простых до \(\sqrt{N}\) требует \(O(\sqrt{N}\log\log N)\) времени и \(O(\sqrt{N})\) памяти. Решето по прогрессии \(2017k-1\) до \(N\) использует массив длины \(\lfloor (N+1)/2017\rfloor\), поэтому эта стадия требует \(O(N/2017)\) памяти и примерно \(O((N/2017)\log\log N)\) работы по отметке составных. Оставшиеся проверки порядков и суммы по включению-исключению выполняются только по конечному списку триггеров и допустимых пар. Для целевых параметров эти этапы меньше двух ситовых проходов.

Сноски и ссылки

  1. Страница задачи: https://projecteuler.net/problem=565
  2. Функция суммы делителей: Wikipedia — Divisor function
  3. Мультипликативный порядок: Wikipedia — Multiplicative order
  4. \(p\)-адическая оценка: Wikipedia — \(p\)-adic valuation
  5. Принцип включения-исключения: Wikipedia — Inclusion-exclusion principle

ملخص المسألة

عندما يكون \(N=10^{11}\) والعدد الأولي \(p=2017\)، نريد حساب

$$S(N,p)=\sum_{\substack{1\le n\le N\\ p\mid \sigma(n)}} n,$$

حيث تمثل \(\sigma(n)\) مجموع القواسم الموجبة للعدد \(n\). فحص كل \(n\le N\) مباشرة وحساب \(\sigma(n)\) لكل واحد منها غير عملي إطلاقًا. الفكرة الفعالة هي تحديد أسس القوى الأولية التي تجعل \(\sigma(n)\) قابلة للقسمة على \(2017\)، ثم إعادة صياغة المسألة كاتحاد موزون لأحداث مرتبطة بالقيم الدقيقة للتقييمات الأولية، وبعد ذلك جمع هذا الاتحاد باستخدام مبدأ الشمول والاستبعاد.

المنهج الرياضي

لنكتب التحليل الأولي للعدد \(n\) على الصورة

$$n=\prod_q q^{a_q}.$$

وبما أن دالة مجموع القواسم دالة ضربية، فإن

$$\sigma(n)=\prod_q \sigma\!\left(q^{a_q}\right).$$

إذن الشرط \(2017\mid \sigma(n)\) يتحقق إذا وفقط إذا وُجد عدد أولي واحد على الأقل \(q\) بحيث يكون \(2017\mid \sigma(q^{a_q})\).

الخطوة 1: اختزال المسألة إلى قوى أولية مُحفِّزة

نعرّف مجموعة الأزواج المحفِّزة بالصيغة

$$B=\left\{(q,e): e\ge 1,\ q^e\le N,\ 2017\mid \sigma(q^e)\right\},$$

حيث إن \(q\) عدد أولي.

وبذلك تصبح القيمة المطلوبة هي الاتحاد الموزون

$$S(N,2017)=\sum_{\substack{1\le n\le N\\ \exists (q,e)\in B:\ v_q(n)=e}} n,$$

حيث يرمز \(v_q(n)\) إلى أس \(q\) في تحليل \(n\). من المهم هنا استعمال التقييم الدقيق، لأن قابلية \(\sigma(q^{a_q})\) للقسمة تعتمد على الأس الحقيقي \(a_q\)، لا على مجرد كون \(q^e\) يقسم \(n\).

أما \(q=2017\) نفسه فلا يعطي أي مساهمة، لأن

$$\sigma(2017^e)=1+2017+\cdots+2017^e\equiv 1 \pmod{2017}.$$

الخطوة 2: توصيف الأسس السيئة بواسطة الرتبة الضربية

لكل عدد أولي \(q\neq 2017\) لدينا

$$\sigma(q^e)=1+q+\cdots+q^e=\frac{q^{e+1}-1}{q-1}.$$

إذا كان \(q\not\equiv 1 \pmod{2017}\)، فإن المقام قابل للعكس بترديد \(2017\)، ومن ثم

$$2017\mid \sigma(q^e)\iff q^{e+1}\equiv 1 \pmod{2017}.$$

ولنرمز إلى الرتبة الضربية لـ \(q\) بترديد \(2017\) بالرمز \(\operatorname{ord}_{2017}(q)\). عندئذ نحصل على

$$2017\mid \sigma(q^e)\iff \operatorname{ord}_{2017}(q)\mid (e+1),$$

أي إن جميع الأسس السيئة تأخذ الشكل

$$e=m\operatorname{ord}_{2017}(q)-1,\qquad m\ge 1.$$

وبما أن

$$2017-1=2016=2^5\cdot 3^2\cdot 7,$$

فكل رتبة ممكنة لا بد أن تقسم \(2016\).

أما الفئة \(q\equiv 1 \pmod{2017}\) فلها سلوك مختلف، لأن

$$\sigma(q^e)\equiv e+1 \pmod{2017}.$$

وهذا يفرض أن يكون \(e+1\) مضاعفًا لـ \(2017\)، وهو أمر مستحيل تحت القيد \(q^e\le 10^{11}\). لذلك لا تدخل هذه الأعداد الأولية في الحل.

الخطوة 3: فصل فرع الرتبة 2 عن بقية الحالات

إذا كان \(q\equiv -1 \pmod{2017}\)، فإن \(\operatorname{ord}_{2017}(q)=2\)، وبالتالي تكون جميع الأسس الفردية سيئة من الناحية النظرية.

لكن على مقياس هذه المسألة لا يمكن أن يظهر إلا \(e=1\). فأصغر عدد أولي يحقق \(q\equiv -1 \pmod{2017}\) هو \(12101\)، ولدينا

$$12101^3>10^{11}.$$

إذن هذا الفرع يساهم فقط بقوى أولية من الشكل \(q^1\).

أما بقية الحالات السيئة فتأتي من أعداد أولية تحقق \(q\not\equiv \pm 1 \pmod{2017}\) ورتبتها على الأقل \(3\)، وبالتالي فإن أصغر أس سيئ فيها لا يقل عن \(2\). لذلك لا تصبح هذه الأعداد مهمة إلا عندما يكون \(q\le \sqrt{N}\)، ولهذا السبب يبحث التنفيذ في هذا الفرع فقط حتى \(\lfloor \sqrt{N}\rfloor\).

الخطوة 4: تحويل كل زوج محفِّز إلى مجموع مغلق

لكل \((q,e)\in B\) نعرّف الحدث

$$A_{q,e}=\{n\le N:\ v_q(n)=e\}.$$

كل عدد ضمن هذا الحدث يمكن كتابته على الصورة

$$n=q^e m,\qquad q\nmid m,\qquad m\le M=\left\lfloor \frac{N}{q^e}\right\rfloor.$$

وإذا رمزنا إلى مجموع الأعداد من \(1\) إلى \(x\) بالرمز

$$T(x)=\frac{x(x+1)}{2},$$

فإن مساهمة الحدث الواحد تساوي

$$\sum_{n\in A_{q,e}} n =q^e\sum_{\substack{m\le M\\ q\nmid m}} m =q^e\left(T(M)-q\,T\!\left(\left\lfloor \frac{M}{q}\right\rfloor\right)\right).$$

المعنى واضح: نحسب مجموع كل القيم الممكنة لـ \(m\)، ثم نطرح القيم التي لا تزال قابلة للقسمة على \(q\). وهكذا تصبح مساهمة كل حدث أحادي قابلة للحساب مباشرة.

الخطوة 5: شمول واستبعاد ثنائي فقط، وهو كافٍ بدقة

الأحداث المرتبطة بالعدد الأولي نفسه ولكن بأسس مختلفة تكون متباينة، لذا فإن التداخل الحقيقي لا يظهر إلا بين عددين أوليين مختلفين. إذا كان \((q,e)\) و\((r,f)\) زوجين محفزين مع \(q\neq r\)، فإن عناصر التقاطع تأخذ الصورة

$$n=q^e r^f m,\qquad q\nmid m,\qquad r\nmid m,$$

حيث

$$m\le X=\left\lfloor \frac{N}{q^e r^f}\right\rfloor.$$

ومن ثم يكون مجموع التقاطع الموزون

$$\sum_{n\in A_{q,e}\cap A_{r,f}} n =q^e r^f\left(T(X)-q\,T\!\left(\left\lfloor \frac{X}{q}\right\rfloor\right)-r\,T\!\left(\left\lfloor \frac{X}{r}\right\rfloor\right)+qr\,T\!\left(\left\lfloor \frac{X}{qr}\right\rfloor\right)\right).$$

إذًا نحصل على الصيغة

$$S(N,2017)=\sum_{(q,e)\in B}\sum_{n\in A_{q,e}} n-\sum_{\substack{(q,e),(r,f)\in B\\ q<r}}\sum_{n\in A_{q,e}\cap A_{r,f}} n.$$

ولا تظهر حدود ثلاثية هنا لأن ذلك مستحيل أصلًا عندما \(N=10^{11}\). فأصغر ثلاث قوى أولية محفِّزة هي

$$2311^2=5{,}340{,}721,\qquad 229^3=12{,}008{,}989,\qquad 3739^2=13{,}980{,}121,$$

وحاصل ضربها أكبر بكثير من \(10^{11}\). لذلك التوقف عند الرتبة الثانية ليس تقريبًا، بل صيغة دقيقة تمامًا.

مثال تطبيقي

لنأخذ مثالين فعليين من الأزواج المحفزة.

أولًا، بما أن

$$12101\equiv -1 \pmod{2017},$$

فإنه يقع في فرع الرتبة \(2\)، ومن ثم يكون \(e=1\) أسًا سيئًا، ونجد فعلًا أن

$$\sigma(12101)=1+12101=12102=6\cdot 2017.$$

ثانيًا، للعدد \(2311\) رتبة ضربية تساوي \(3\) بترديد \(2017\)، ولهذا يكون \(e=2\) سيئًا، وكذلك

$$\sigma(2311^2)=1+2311+2311^2=5{,}343{,}033=2649\cdot 2017.$$

وعند \(N=10^{11}\) يكون عامل التقاطع بين الحدثين

$$P=12101\cdot 2311^2=64{,}628{,}064{,}821\le 10^{11}.$$

وهذا يعني أن المضاعف الممكن الوحيد هو \(m=1\)، فتكون مساهمة هذا التقاطع مساوية بالضبط لـ \(P\)، ويجب طرحها مرة واحدة من مجموع المساهمتين الأحاديتين.

كيف تعمل الشيفرة

تبدأ الشيفرة بتوليد جميع الأعداد الأولية حتى \(\lfloor \sqrt{N}\rfloor\). وتُستعمل هذه القائمة في أمرين: حساب الرتب الضربية بترديد \(2017\) في الفرع الذي لا يساوي فيه \(q\) أيًا من \(\pm 1\)، وتشغيل غربال على المتتالية الحسابية \(2017k-1\) لجمع كل الأعداد الأولية الموافقة لـ \(-1\) بترديد \(2017\) حتى \(N\).

بعد ذلك تقوم تطبيقات C++ وPython وJava بحصر جميع الأزواج المحفِّزة \((q,e)\). في فرع الرتبة \(2\) يُحتفظ فقط بالأس \(1\). وفي الفرع الآخر تُحسب الرتبة على أنها قاسم للعدد \(2016\)، ثم تُسرد كل الأسس من الشكل \(e=m\operatorname{ord}_{2017}(q)-1\) التي تحقق \(q^e\le N\)، وتُهمَل الأعداد الأولية التي لا يمكن أن تعطي مساهمة.

بعد معرفة هذه القائمة، تجمع الشيفرة جميع مساهمات الأحداث الأحادية، ثم تطرح كل التقاطعات الثنائية المسموح بها: فرع الرتبة \(2\) مع نفسه، وفرع الرتبة \(2\) مع بقية الأزواج، وبقية الأزواج فيما بينها. تنفيذ Python يمرر العمل إلى الحل المترجم، أما تنفيذ Java فيكرر الحسابات نفسها مباشرة باستعمال أعداد صحيحة ذات دقة غير محدودة.

كما يجري التحقق أولًا من نقطتي الفحص المنشورتين

$$S(10^6,2017)=150850429,\qquad S(10^9,2017)=249652238344557,$$

ثم تُحسب الحالة النهائية عند \(N=10^{11}\).

تحليل التعقيد

توليد جميع الأعداد الأولية حتى \(\sqrt{N}\) يحتاج إلى زمن \(O(\sqrt{N}\log\log N)\) وذاكرة \(O(\sqrt{N})\). أما الغربال على المتتالية \(2017k-1\) حتى \(N\) فيستعمل مصفوفة طولها \(\lfloor (N+1)/2017\rfloor\)، ولذلك يحتاج إلى ذاكرة \(O(N/2017)\) ونحو \(O((N/2017)\log\log N)\) من أعمال الوسم. أما حساب الرتب ومجاميع الشمول والاستبعاد بعد ذلك فيجري فقط على القائمة المنتهية من الأزواج المحفزة وتقاطعها الثنائي الممكن. وعند معاملات المسألة تكون هذه المرحلة أخف من مرحلتي الغربلة.

حواشٍ ومراجع

  1. صفحة المسألة: https://projecteuler.net/problem=565
  2. دالة مجموع القواسم: Wikipedia — Divisor function
  3. الرتبة الضربية: Wikipedia — Multiplicative order
  4. التقييم \(p\)-أدي: Wikipedia — \(p\)-adic valuation
  5. مبدأ الشمول والاستبعاد: Wikipedia — Inclusion-exclusion principle