Problem Summary

Fix a bound \(B\) and look only at primes \(p \le B\). The relevant primes split into the two nonzero residue classes modulo 3:

$$\mathcal P_1=\{p \le B : p \text{ prime},\ p \equiv 1 \pmod{3}\}, \qquad \mathcal P_2=\{p \le B : p \text{ prime},\ p \equiv 2 \pmod{3}\}.$$

The prime \(3\) is intentionally excluded, because it belongs to neither nonzero class and does not change the congruence that drives the count. For any integer \(n\), define

$$\omega_1(n)=\#\{p\in \mathcal P_1 : p\mid n\},\qquad \omega_2(n)=\#\{p\in \mathcal P_2 : p\mid n\}.$$

The computation counts

$$F(N,B)=\#\{1\le n\le N : \omega_1(n)+2\omega_2(n)\equiv0\pmod3\}.$$

Equivalently, it counts the integers for which \(\omega_1(n)\equiv\omega_2(n)\pmod3\). Only the presence or absence of relevant prime divisors matters, so prime exponents never enter the formula. The challenge is to evaluate this count for \(N=10^{18}\) and \(B=120\) without scanning all integers individually.

Mathematical Approach

The solution is a two-dimensional inclusion-exclusion over subsets of the relevant primes. The goal is to turn the modular condition on \(\omega_1(n)\) and \(\omega_2(n)\) into a weighted sum over squarefree products.

The relevant prime signature

Let

$$\mathcal P_B=\mathcal P_1\cup\mathcal P_2.$$

For a subset \(S\subseteq \mathcal P_B\), write

$$d(S)=\prod_{p\in S} p,\qquad a(S)=\#(S\cap\mathcal P_1),\qquad b(S)=\#(S\cap\mathcal P_2).$$

If \(d(S)\mid n\), then every prime in \(S\) divides \(n\). Since the definition depends only on which relevant primes divide \(n\), every integer is summarized by the signature \((\omega_1(n),\omega_2(n))\). Factors of \(3\), powers of primes already present, and prime divisors larger than \(B\) are all invisible to this signature.

Detecting the congruence class

Define the indicator

$$G(x,y)=\begin{cases}1,& x+2y\equiv0\pmod3,\\0,& \text{otherwise}.\end{cases}$$

A clean way to isolate this condition is with a cubic root-of-unity filter. If \(\zeta\neq1\) and \(\zeta^3=1\), then

$$G(x,y)=\frac{1+\zeta^{x+2y}+\zeta^{2x+y}}{3}.$$

This identity shows that only the counts modulo 3 matter. The implementations never use complex arithmetic directly, but this filter explains why the final weight table depends only on the two residue-class counts.

Binomial inversion turns the indicator into integer weights

Suppose an integer \(n\) has \(x\) relevant prime divisors from \(\mathcal P_1\) and \(y\) from \(\mathcal P_2\). We want weights \(w(a,b)\) such that

$$G(x,y)=\sum_{a=0}^{x}\sum_{b=0}^{y}\binom{x}{a}\binom{y}{b}w(a,b).$$

The combinatorial meaning is that we choose any \(a\) of the \(x\) class-\(1\) primes and any \(b\) of the \(y\) class-\(2\) primes. Each choice corresponds to one subset of the relevant prime divisors of \(n\).

Applying binomial inversion in both variables gives

$$w(a,b)=\sum_{i=0}^{a}\sum_{j=0}^{b}(-1)^{a-i+b-j}\binom{a}{i}\binom{b}{j}G(i,j).$$

Because \(G(i,j)\) vanishes unless \(i+2j\equiv0\pmod3\), this becomes

$$w(a,b)=(-1)^{a+b}\sum_{\substack{0\le i\le a\\0\le j\le b\\i+2j\equiv0\!\!\!\pmod3}}(-1)^{i+j}\binom{a}{i}\binom{b}{j}.$$

This is exactly the formula precomputed by the C++, Python, and Java implementations.

Collapsing the count to squarefree products

For each integer \(n\), let \(R_B(n)\subseteq\mathcal P_B\) be the set of relevant primes dividing \(n\), and let \(\mathbf{1}_B(n)\) be the indicator of \(B\)-trivisibility. The inversion formula implies

$$\mathbf{1}_B(n)=\sum_{S\subseteq R_B(n)} w\big(a(S),b(S)\big).$$

Now sum over all \(n\le N\) and swap the order of summation. A fixed subset \(S\) contributes once for every multiple of \(d(S)\), so its total contribution is

$$\left\lfloor\frac{N}{d(S)}\right\rfloor.$$

Therefore

$$\boxed{F(N,B)=\sum_{S\subseteq\mathcal P_B} w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor.}$$

This is the core reduction. The original search over \(1,\dots,N\) is replaced by a sum over subsets of primes up to \(B\), and only squarefree products appear because each prime is either selected once or not selected at all.

Worked example: \(B=10\)

For \(B=10\), the relevant primes are

$$\mathcal P_1=\{7\},\qquad \mathcal P_2=\{2,5\}.$$

The condition is \(\omega_1(n)+2\omega_2(n)\equiv0\pmod3\). Among the integers \(1\le n\le10\), the numbers \(1,3,9\) have signature \((0,0)\) and are counted. The numbers \(2,4,5,6,8\) have signature \((0,1)\), so they fail the congruence. The number \(7\) has signature \((1,0)\), and \(10\) has signature \((0,2)\); both are rejected as well. Hence

$$F(10,10)=3,$$

which matches the small checkpoint used by the implementations. The same reasoning explains \(F(10,4)=5\): if the only relevant prime is \(2\), then the counted integers are exactly the odd ones.

How the Code Works

The C++, Python, and Java implementations follow the derivation almost literally. They preprocess the prime classes, build the integer weight table, and then enumerate prime subsets recursively.

Generating the relevant primes

First, the implementation runs a sieve up to \(B\), removes the prime \(3\), and counts how many remaining primes lie in \(\mathcal P_1\) and \(\mathcal P_2\). These counts determine the dimensions of the weight table.

Building Pascal's triangle and the weight table

Next, it constructs binomial coefficients \(\binom{n}{k}\) by Pascal's rule. Using those coefficients, it fills the table \(w(a,b)\) with the inversion formula above. The root-of-unity viewpoint remains only a proof device; the actual computation stays entirely in integer arithmetic.

Recursive enumeration of squarefree products

The main recursion decides, for each relevant prime, whether that prime is excluded from the current subset or included in it. Along the recursion the implementation tracks the current product \(d(S)\), the number of selected primes from \(\mathcal P_1\), and the number selected from \(\mathcal P_2\). At a leaf it adds

$$w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor$$

to the answer. If including the next prime would make the product exceed \(N\), the branch is cut immediately, because every deeper extension would contribute zero. The C++ and Java implementations also verify a few small checkpoint values before printing the final result.

Complexity Analysis

Let \(P=|\mathcal P_B|\), \(M_1=|\mathcal P_1|\), \(M_2=|\mathcal P_2|\), and \(D=\max(M_1,M_2)\). Building the binomial table costs \(O(D^2)\) time and memory. Filling the weight table by the direct inversion formula costs \(O(M_1^2M_2^2)\) time and \(O(M_1M_2)\) memory, which is tiny for \(B=120\).

The recursive subset traversal has worst-case size \(O(2^P)\), since every relevant prime can be absent or present. In practice it is much smaller because branches are pruned as soon as their product exceeds \(N\). For the actual problem, \(P=29\), so the search tree is manageable and the algorithm avoids any dependence on the full scale of \(N=10^{18}\).

Footnotes and References

  1. Problem page: Project Euler 967 - B-Trivisible Numbers
  2. Inclusion-exclusion principle: Wikipedia - Inclusion-exclusion principle
  3. Binomial transform and binomial inversion: Wikipedia - Binomial transform
  4. Modular arithmetic: Wikipedia - Modular arithmetic
  5. Roots of unity: Wikipedia - Root of unity

Problemzusammenfassung

Für ein festes \(B\) sind nur Primzahlen \(p \le B\) relevant. Diese zerfallen in die beiden von null verschiedenen Restklassen modulo 3:

$$\mathcal P_1=\{p \le B : p \text{ prim},\ p \equiv 1 \pmod{3}\}, \qquad \mathcal P_2=\{p \le B : p \text{ prim},\ p \equiv 2 \pmod{3}\}.$$

Die Primzahl \(3\) wird bewusst ausgeschlossen, weil sie zu keiner dieser beiden Klassen gehört und die entscheidende Kongruenz nicht verändert. Für eine ganze Zahl \(n\) definieren wir

$$\omega_1(n)=\#\{p\in \mathcal P_1 : p\mid n\},\qquad \omega_2(n)=\#\{p\in \mathcal P_2 : p\mid n\}.$$

Gesucht ist damit

$$F(N,B)=\#\{1\le n\le N : \omega_1(n)+2\omega_2(n)\equiv0\pmod3\}.$$

Äquivalent dazu werden genau die Zahlen gezählt, für die \(\omega_1(n)\equiv\omega_2(n)\pmod3\) gilt. Wichtig ist: Es zählt nur, ob ein relevanter Primfaktor vorkommt, nicht mit welcher Potenz. Für \(N=10^{18}\) und \(B=120\) ist daher eine direkte Prüfung aller Zahlen ausgeschlossen.

Mathematischer Ansatz

Die Lösung ist eine zweidimensionale Inklusions-Exklusions-Rechnung über Teilmengen der relevanten Primzahlen. Die Kongruenzbedingung auf \(\omega_1(n)\) und \(\omega_2(n)\) wird dabei in eine gewichtete Summe über quadratfreie Produkte übersetzt.

Die relevante Primsignatur

Setze

$$\mathcal P_B=\mathcal P_1\cup\mathcal P_2.$$

Für eine Teilmenge \(S\subseteq \mathcal P_B\) schreiben wir

$$d(S)=\prod_{p\in S} p,\qquad a(S)=\#(S\cap\mathcal P_1),\qquad b(S)=\#(S\cap\mathcal P_2).$$

Aus \(d(S)\mid n\) folgt, dass alle Primzahlen aus \(S\) den Wert \(n\) teilen. Da die Eigenschaft nur davon abhängt, welche relevanten Primfaktoren vorkommen, wird jede Zahl vollständig durch die Signatur \((\omega_1(n),\omega_2(n))\) beschrieben. Faktoren \(3\), höhere Potenzen derselben Primzahlen und Primteiler größer als \(B\) spielen für diese Signatur keine Rolle.

Die Kongruenzklasse erkennen

Definiere die Indikatorfunktion

$$G(x,y)=\begin{cases}1,& x+2y\equiv0\pmod3,\\0,& \text{sonst}.\end{cases}$$

Eine elegante Darstellung benutzt einen Filter mit dritten Einheitswurzeln. Ist \(\zeta\neq1\) und \(\zeta^3=1\), dann gilt

$$G(x,y)=\frac{1+\zeta^{x+2y}+\zeta^{2x+y}}{3}.$$

Damit sieht man sofort, dass nur die Restklassen der Zähler modulo 3 relevant sind. Die Implementierungen rechnen später nur noch mit ganzen Zahlen, aber diese Identität erklärt die Form der Gewichtstabelle.

Binomialinversion liefert ganzzahlige Gewichte

Hat eine Zahl \(n\) genau \(x\) relevante Primteiler aus \(\mathcal P_1\) und \(y\) aus \(\mathcal P_2\), so suchen wir Gewichte \(w(a,b)\) mit

$$G(x,y)=\sum_{a=0}^{x}\sum_{b=0}^{y}\binom{x}{a}\binom{y}{b}w(a,b).$$

Die kombinatorische Bedeutung ist direkt: Man wählt beliebige \(a\) der \(x\) Primzahlen aus Klasse \(1\) und beliebige \(b\) der \(y\) Primzahlen aus Klasse \(2\). Jede solche Wahl entspricht einer Teilmenge der relevanten Primteiler von \(n\).

Die Binomialinversion in beiden Variablen ergibt

$$w(a,b)=\sum_{i=0}^{a}\sum_{j=0}^{b}(-1)^{a-i+b-j}\binom{a}{i}\binom{b}{j}G(i,j).$$

Da \(G(i,j)\) nur dann ungleich null ist, wenn \(i+2j\equiv0\pmod3\), folgt daraus

$$w(a,b)=(-1)^{a+b}\sum_{\substack{0\le i\le a\\0\le j\le b\\i+2j\equiv0\!\!\!\pmod3}}(-1)^{i+j}\binom{a}{i}\binom{b}{j}.$$

Genau diese Formel wird in den Implementierungen vorab ausgewertet.

Die Zählung auf quadratfreie Produkte reduzieren

Für jede Zahl \(n\) sei \(R_B(n)\subseteq\mathcal P_B\) die Menge der relevanten Primteiler von \(n\), und \(\mathbf{1}_B(n)\) sei der Indikator dafür, dass \(n\) \(B\)-trivisibel ist. Aus der Inversion folgt

$$\mathbf{1}_B(n)=\sum_{S\subseteq R_B(n)} w\big(a(S),b(S)\big).$$

Summiert man nun über alle \(n\le N\) und vertauscht die Summationsreihenfolge, dann trägt eine feste Teilmenge \(S\) genau einmal für jedes Vielfache von \(d(S)\) bei. Ihr Gesamtbeitrag ist also

$$\left\lfloor\frac{N}{d(S)}\right\rfloor.$$

Damit erhält man

$$\boxed{F(N,B)=\sum_{S\subseteq\mathcal P_B} w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor.}$$

Das ist die eigentliche Reduktion. Die Suche über \(1,\dots,N\) wird durch eine Summe über Primzahl-Teilmengen bis \(B\) ersetzt, und es treten nur quadratfreie Produkte auf, weil jede Primzahl höchstens einmal gewählt wird.

Durchgerechnetes Beispiel: \(B=10\)

Für \(B=10\) sind die relevanten Primzahlen

$$\mathcal P_1=\{7\},\qquad \mathcal P_2=\{2,5\}.$$

Die Bedingung lautet \(\omega_1(n)+2\omega_2(n)\equiv0\pmod3\). Unter den Zahlen \(1\le n\le10\) haben \(1,3,9\) die Signatur \((0,0)\) und werden gezählt. Die Zahlen \(2,4,5,6,8\) haben die Signatur \((0,1)\) und scheitern an der Kongruenz. Die Zahl \(7\) hat \((1,0)\), die Zahl \(10\) hat \((0,2)\); beide werden ebenfalls nicht gezählt. Also gilt

$$F(10,10)=3,$$

genau wie in der kleinen Kontrolle der Implementierungen. Derselbe Gedanke erklärt auch \(F(10,4)=5\): Wenn nur die Primzahl \(2\) relevant ist, werden genau die ungeraden Zahlen gezählt.

Wie der Code arbeitet

Die C++-, Python- und Java-Implementierungen folgen dieser Herleitung fast wörtlich. Sie bereiten die Primklassen vor, bauen die Gewichtstabelle auf und enumerieren dann rekursiv die Primzahl-Teilmengen.

Erzeugung der relevanten Primzahlen

Zunächst läuft ein Sieb bis \(B\), danach wird die Primzahl \(3\) entfernt. Anschließend wird gezählt, wie viele verbleibende Primzahlen in \(\mathcal P_1\) bzw. \(\mathcal P_2\) liegen. Diese beiden Zahlen bestimmen die Größe der Gewichtstabelle.

Pascal-Dreieck und Gewichtstabelle

Danach werden die Binomialkoeffizienten \(\binom{n}{k}\) über Pascalsche Rekursion aufgebaut. Mit ihnen füllt die Implementierung die Tabelle \(w(a,b)\) gemäß der Inversionsformel. Die Einheitswurzeln dienen also nur der mathematischen Motivation; die eigentliche Rechnung bleibt vollständig ganzzahlig.

Rekursive Enumeration quadratfreier Produkte

Die Hauptrekursion entscheidet für jede relevante Primzahl, ob sie in der aktuellen Teilmenge vorkommt oder nicht. Dabei werden das aktuelle Produkt \(d(S)\), die Anzahl gewählter Primzahlen aus \(\mathcal P_1\) und die Anzahl aus \(\mathcal P_2\) mitgeführt. An einem Blatt addiert die Implementierung

$$w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor$$

zum Ergebnis. Würde das Hinzufügen der nächsten Primzahl das Produkt größer als \(N\) machen, wird der Ast sofort abgeschnitten, denn jede Fortsetzung hätte Beitrag null. Die C++- und Java-Versionen prüfen zusätzlich einige kleine Kontrollwerte, bevor sie den Endwert ausgeben.

Komplexitätsanalyse

Seien \(P=|\mathcal P_B|\), \(M_1=|\mathcal P_1|\), \(M_2=|\mathcal P_2|\) und \(D=\max(M_1,M_2)\). Das Aufbauen der Binomialtabelle benötigt \(O(D^2)\) Zeit und Speicher. Das Füllen der Gewichtstabelle mit der direkten Inversionsformel kostet \(O(M_1^2M_2^2)\) Zeit und \(O(M_1M_2)\) Speicher. Für \(B=120\) sind diese Tabellen sehr klein.

Die rekursive Teilmengensuche hat im Worst Case Größe \(O(2^P)\), weil jede relevante Primzahl fehlen oder vorkommen kann. Praktisch ist sie deutlich kleiner, da Äste abgeschnitten werden, sobald ihr Produkt \(N\) überschreitet. Im konkreten Problem gilt \(P=29\), sodass die Rechnung durch einen handhabbaren Tiefensuchbaum dominiert wird und nicht durch die Größe von \(N=10^{18}\).

Fußnoten und Referenzen

  1. Problemseite: Project Euler 967 - B-Trivisible Numbers
  2. Prinzip von Inklusion und Exklusion: Wikipedia - Prinzip von Inklusion und Exklusion
  3. Binomialtransformation und Inversion: Wikipedia - Binomial transform
  4. Modulare Arithmetik: Wikipedia - Modulare Arithmetik
  5. Einheitswurzel: Wikipedia - Einheitswurzel

Problem Özeti

Sabit bir \(B\) için yalnızca \(p \le B\) olan asallar önemlidir. Bu asallar, 3 moduna göre sıfır olmayan iki sınıfa ayrılır:

$$\mathcal P_1=\{p \le B : p \text{ asal},\ p \equiv 1 \pmod{3}\}, \qquad \mathcal P_2=\{p \le B : p \text{ asal},\ p \equiv 2 \pmod{3}\}.$$

Asal \(3\) özellikle dışarıda bırakılır; çünkü bu iki sınıftan hiçbirine ait değildir ve sayımı belirleyen kongruansı değiştirmez. Herhangi bir \(n\) tam sayısı için

$$\omega_1(n)=\#\{p\in \mathcal P_1 : p\mid n\},\qquad \omega_2(n)=\#\{p\in \mathcal P_2 : p\mid n\}$$

tanımını yapalım. Hesaplanan büyüklük

$$F(N,B)=\#\{1\le n\le N : \omega_1(n)+2\omega_2(n)\equiv0\pmod3\}$$

olur. Eşdeğer olarak, \(\omega_1(n)\equiv\omega_2(n)\pmod3\) olan tam sayılar sayılır. Burada yalnızca ilgili asal bölenlerin varlığı önemlidir; üsleri hiçbir rol oynamaz. Sorunun gerçek girdisinde \(N=10^{18}\) ve \(B=120\) olduğundan, bütün sayıları tek tek taramak mümkün değildir.

Matematiksel Yaklaşım

Çözüm, ilgili asalların altkümeleri üzerinde iki boyutlu bir dahil etme-hariç tutma hesabıdır. Amaç, \(\omega_1(n)\) ve \(\omega_2(n)\) üzerindeki mod 3 koşulunu, karesiz asal çarpımları üzerinden ağırlıklı bir toplama dönüştürmektir.

İlgili asal imzası

Şunu tanımlayalım:

$$\mathcal P_B=\mathcal P_1\cup\mathcal P_2.$$

Her \(S\subseteq \mathcal P_B\) altkümesi için

$$d(S)=\prod_{p\in S} p,\qquad a(S)=\#(S\cap\mathcal P_1),\qquad b(S)=\#(S\cap\mathcal P_2)$$

olsun. Eğer \(d(S)\mid n\) ise, \(S\) içindeki bütün asallar \(n\)'yi böler. Tanım yalnızca hangi ilgili asalların \(n\)'yi böldüğüne baktığı için her tam sayı, \((\omega_1(n),\omega_2(n))\) imzası ile özetlenir. \(3\)'ün kuvvetleri, zaten var olan asalların daha yüksek üsleri ve \(B\)'den büyük asal çarpanlar bu imzayı değiştirmez.

Kongruans koşulunu ayıklamak

Gösterge fonksiyonunu

$$G(x,y)=\begin{cases}1,& x+2y\equiv0\pmod3,\\0,& \text{aksi halde}.\end{cases}$$

şeklinde tanımlayalım. Bu koşulu ayıklamanın temiz bir yolu, kübik birim kök filtresidir. \(\zeta\neq1\) ve \(\zeta^3=1\) ise

$$G(x,y)=\frac{1+\zeta^{x+2y}+\zeta^{2x+y}}{3}$$

olur. Böylece yalnızca sayıların mod 3 kalanı ile ilgilendiğimiz açıkça görülür. Uygulamalar karmaşık sayılarla çalışmaz; fakat bu kimlik, önceden hesaplanan ağırlık tablosunun neden yalnızca iki sınıftaki sayılara bağlı olduğunu açıklar.

Binom terslemesi ile tam sayılı ağırlıklar üretmek

Diyelim ki bir \(n\) sayısının \(\mathcal P_1\) içinden tam \(x\), \(\mathcal P_2\) içinden de tam \(y\) tane ilgili asal böleni var. Şöyle bir \(w(a,b)\) tablosu istiyoruz:

$$G(x,y)=\sum_{a=0}^{x}\sum_{b=0}^{y}\binom{x}{a}\binom{y}{b}w(a,b).$$

Bunun kombinatoryal anlamı açıktır: \(x\) adet sınıf-\(1\) asal içinden herhangi \(a\) tanesini, \(y\) adet sınıf-\(2\) asal içinden de herhangi \(b\) tanesini seçiyoruz. Her seçim, \(n\)'nin ilgili asal bölenlerinden bir altkümeye karşılık gelir.

İki değişkenli binom terslemesi uygulanınca

$$w(a,b)=\sum_{i=0}^{a}\sum_{j=0}^{b}(-1)^{a-i+b-j}\binom{a}{i}\binom{b}{j}G(i,j)$$

elde edilir. \(G(i,j)\) yalnızca \(i+2j\equiv0\pmod3\) olduğunda sıfır olmadığı için bu formül

$$w(a,b)=(-1)^{a+b}\sum_{\substack{0\le i\le a\\0\le j\le b\\i+2j\equiv0\!\!\!\pmod3}}(-1)^{i+j}\binom{a}{i}\binom{b}{j}$$

biçimine iner. C++, Python ve Java uygulamalarının önceden hesapladığı tablo tam olarak budur.

Sayımı karesiz çarpımlara indirmek

Her \(n\) için, \(n\)'yi bölen ilgili asallar kümesini \(R_B(n)\subseteq\mathcal P_B\) ile gösterelim; \(\mathbf{1}_B(n)\) ise \(n\)'nin \(B\)-trivisible olup olmadığını gösteren fonksiyon olsun. Tersleme formülü

$$\mathbf{1}_B(n)=\sum_{S\subseteq R_B(n)} w\big(a(S),b(S)\big)$$

sonucunu verir. Şimdi bunu tüm \(n\le N\) için toplayıp toplamların sırasını değiştirirsek, sabit bir \(S\) altkümesi her \(d(S)\) katı için bir kez katkı yapar. Yani toplam katkısı

$$\left\lfloor\frac{N}{d(S)}\right\rfloor$$

olur. Böylece

$$\boxed{F(N,B)=\sum_{S\subseteq\mathcal P_B} w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor.}$$

temel formülüne ulaşırız. Başlangıçta \(1,\dots,N\) üzerinde görünen problem, artık yalnızca \(B\)'ye kadar olan asalların altkümeleri üzerinde bir toplamdır. Her asal ya seçilir ya seçilmez olduğu için yalnızca karesiz çarpımlar ortaya çıkar.

Çalışılmış örnek: \(B=10\)

\(B=10\) için ilgili asallar

$$\mathcal P_1=\{7\},\qquad \mathcal P_2=\{2,5\}$$

olur. Koşul \(\omega_1(n)+2\omega_2(n)\equiv0\pmod3\) şeklindedir. \(1\le n\le10\) aralığında \(1,3,9\) sayılarının imzası \((0,0)\) olduğu için bunlar sayılır. \(2,4,5,6,8\) sayılarının imzası \((0,1)\) olduğundan koşulu sağlamazlar. \(7\) için imza \((1,0)\), \(10\) için imza \((0,2)\) olur; bunlar da dışarıda kalır. Dolayısıyla

$$F(10,10)=3$$

elde edilir; bu da uygulamalardaki küçük kontrol değeriyle aynıdır. Aynı mantık \(F(10,4)=5\) sonucunu da açıklar: ilgili tek asal \(2\) ise, sayılan tam sayılar tam olarak tek sayılardır.

Kod Nasıl Çalışır

C++, Python ve Java uygulamaları bu türetmeyi neredeyse doğrudan uygular. Önce asal sınıfları hazırlanır, sonra tam sayılı ağırlık tablosu kurulur, en son da asal altkümeleri özyineli olarak dolaşılır.

İlgili asalları üretmek

İlk aşamada \(B\)'ye kadar bir asal eleği çalıştırılır, sonra \(3\) listeden çıkarılır. Kalan asalların kaç tanesinin \(\mathcal P_1\), kaç tanesinin \(\mathcal P_2\) içinde olduğu sayılır. Bu sayılar ağırlık tablosunun boyutlarını belirler.

Pascal üçgeni ve ağırlık tablosu

Daha sonra \(\binom{n}{k}\) katsayıları Pascal kuralı ile üretilir. Bu katsayılar kullanılarak \(w(a,b)\) tablosu yukarıdaki tersleme formülüyle doldurulur. Birim kök filtresi sadece matematiksel açıklama sağlar; gerçek hesap tamamen tam sayılarla yapılır.

Karesiz çarpımların özyineli taranması

Ana özyineleme, her ilgili asalı mevcut altkümeye katıp katmamaya karar verir. Bu sırada güncel çarpım \(d(S)\), \(\mathcal P_1\)'den seçilen asal sayısı ve \(\mathcal P_2\)'den seçilen asal sayısı tutulur. Bir yaprakta

$$w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor$$

ifadeleri sonuca eklenir. Sıradaki asalı eklemek çarpımı \(N\)'yi aşıracaksa, o dal hemen budanır; çünkü daha derin hiçbir uzantı katkı yapamaz. C++ ve Java sürümleri son değeri yazmadan önce birkaç küçük kontrol örneğini de doğrular.

Karmaşıklık Analizi

\(P=|\mathcal P_B|\), \(M_1=|\mathcal P_1|\), \(M_2=|\mathcal P_2|\) ve \(D=\max(M_1,M_2)\) olsun. Binom tablosunu kurmanın maliyeti zaman ve bellek açısından \(O(D^2)\)'dir. Doğrudan tersleme formülü ile ağırlık tablosunu doldurmak \(O(M_1^2M_2^2)\) zaman ve \(O(M_1M_2)\) bellek ister. \(B=120\) için bu tablolar çok küçüktür.

Altküme özyinelemesinin teorik üst sınırı \(O(2^P)\)'dir; çünkü her ilgili asal ya vardır ya yoktur. Pratikte dallar, çarpım \(N\)'yi aşar aşmaz budandığı için gerçek maliyet çok daha düşüktür. Bu problemde \(P=29\) olduğundan, algoritma \(N=10^{18}\) ölçeğinde sayıları dolaşmak yerine yönetilebilir büyüklükte bir arama ağacı üzerinde çalışır.

Dipnotlar ve Kaynaklar

  1. Problem sayfası: Project Euler 967 - B-Trivisible Numbers
  2. Dahil etme-hariç tutma ilkesi: Wikipedia - Dahil etme-hariç tutma prensibi
  3. Binom dönüşümü ve tersleme: Wikipedia - Binomial transform
  4. Modüler aritmetik: Wikipedia - Modüler aritmetik
  5. Birim kökler: Wikipedia - Birliğin kökleri

Resumen del Problema

Para un \(B\) fijo, solo importan los primos \(p \le B\). Esos primos se separan en las dos clases residuales no nulas módulo 3:

$$\mathcal P_1=\{p \le B : p \text{ primo},\ p \equiv 1 \pmod{3}\}, \qquad \mathcal P_2=\{p \le B : p \text{ primo},\ p \equiv 2 \pmod{3}\}.$$

El primo \(3\) se excluye a propósito, porque no pertenece a ninguna de esas dos clases y no altera la congruencia que gobierna el conteo. Para cualquier entero \(n\), definimos

$$\omega_1(n)=\#\{p\in \mathcal P_1 : p\mid n\},\qquad \omega_2(n)=\#\{p\in \mathcal P_2 : p\mid n\}.$$

Entonces la cantidad buscada es

$$F(N,B)=\#\{1\le n\le N : \omega_1(n)+2\omega_2(n)\equiv0\pmod3\}.$$

De forma equivalente, se cuentan los enteros para los que \(\omega_1(n)\equiv\omega_2(n)\pmod3\). Solo importa si un primo relevante divide o no a \(n\); las potencias no influyen. Con \(N=10^{18}\) y \(B=120\), una exploración directa de todos los enteros es imposible.

Enfoque Matemático

La solución es una inclusión-exclusión bidimensional sobre subconjuntos de los primos relevantes. La idea es convertir la condición modular sobre \(\omega_1(n)\) y \(\omega_2(n)\) en una suma ponderada sobre productos libres de cuadrados.

La firma prima relevante

Sea

$$\mathcal P_B=\mathcal P_1\cup\mathcal P_2.$$

Para un subconjunto \(S\subseteq \mathcal P_B\), escribimos

$$d(S)=\prod_{p\in S} p,\qquad a(S)=\#(S\cap\mathcal P_1),\qquad b(S)=\#(S\cap\mathcal P_2).$$

Si \(d(S)\mid n\), entonces todos los primos de \(S\) dividen a \(n\). Como la definición depende solo de qué primos relevantes dividen a \(n\), cada entero queda resumido por la firma \((\omega_1(n),\omega_2(n))\). Los factores \(3\), las potencias de primos ya presentes y los primos mayores que \(B\) son invisibles para esta firma.

Cómo aislar la clase de congruencia

Definimos el indicador

$$G(x,y)=\begin{cases}1,& x+2y\equiv0\pmod3,\\0,& \text{en otro caso}.\end{cases}$$

Una forma elegante de detectar esta condición es usar un filtro con raíces cúbicas de la unidad. Si \(\zeta\neq1\) y \(\zeta^3=1\), entonces

$$G(x,y)=\frac{1+\zeta^{x+2y}+\zeta^{2x+y}}{3}.$$

La identidad deja claro que solo importan los recuentos módulo 3. Las implementaciones no usan números complejos, pero este filtro explica la estructura de la tabla de pesos.

La inversión binomial produce pesos enteros

Supongamos que un entero \(n\) tiene exactamente \(x\) primos relevantes de \(\mathcal P_1\) y \(y\) de \(\mathcal P_2\). Queremos pesos \(w(a,b)\) tales que

$$G(x,y)=\sum_{a=0}^{x}\sum_{b=0}^{y}\binom{x}{a}\binom{y}{b}w(a,b).$$

La interpretación combinatoria es directa: elegimos cualesquiera \(a\) de los \(x\) primos de clase \(1\) y cualesquiera \(b\) de los \(y\) primos de clase \(2\). Cada elección corresponde a un subconjunto de los primos relevantes que dividen a \(n\).

Aplicando inversión binomial en ambas variables se obtiene

$$w(a,b)=\sum_{i=0}^{a}\sum_{j=0}^{b}(-1)^{a-i+b-j}\binom{a}{i}\binom{b}{j}G(i,j).$$

Como \(G(i,j)\) solo es no nula cuando \(i+2j\equiv0\pmod3\), resulta

$$w(a,b)=(-1)^{a+b}\sum_{\substack{0\le i\le a\\0\le j\le b\\i+2j\equiv0\!\!\!\pmod3}}(-1)^{i+j}\binom{a}{i}\binom{b}{j}.$$

Esa es exactamente la fórmula que las implementaciones precalculan.

Reducir el conteo a productos libres de cuadrados

Para cada entero \(n\), sea \(R_B(n)\subseteq\mathcal P_B\) el conjunto de primos relevantes que dividen a \(n\), y sea \(\mathbf{1}_B(n)\) el indicador de que \(n\) es \(B\)-trivisible. La inversión implica

$$\mathbf{1}_B(n)=\sum_{S\subseteq R_B(n)} w\big(a(S),b(S)\big).$$

Ahora sumamos sobre todos los \(n\le N\) e intercambiamos el orden de la suma. Un subconjunto fijo \(S\) aporta una vez por cada múltiplo de \(d(S)\), así que su contribución total es

$$\left\lfloor\frac{N}{d(S)}\right\rfloor.$$

Por tanto,

$$\boxed{F(N,B)=\sum_{S\subseteq\mathcal P_B} w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor.}$$

Esta es la reducción clave. El problema ya no es recorrer \(1,\dots,N\), sino sumar sobre subconjuntos de primos hasta \(B\). Solo aparecen productos libres de cuadrados porque cada primo se elige una sola vez o no se elige.

Ejemplo trabajado: \(B=10\)

Cuando \(B=10\), los primos relevantes son

$$\mathcal P_1=\{7\},\qquad \mathcal P_2=\{2,5\}.$$

La condición es \(\omega_1(n)+2\omega_2(n)\equiv0\pmod3\). Entre \(1\le n\le10\), los números \(1,3,9\) tienen firma \((0,0)\) y sí cuentan. Los números \(2,4,5,6,8\) tienen firma \((0,1)\), así que fallan la congruencia. El \(7\) tiene firma \((1,0)\) y el \(10\) tiene firma \((0,2)\); ambos quedan fuera. Por lo tanto,

$$F(10,10)=3,$$

exactamente el mismo valor que aparece en la comprobación pequeña de las implementaciones. La misma idea explica \(F(10,4)=5\): si el único primo relevante es \(2\), los enteros contados son simplemente los impares.

Cómo Funciona el Código

Las implementaciones en C++, Python y Java siguen esta derivación casi línea por línea. Primero preparan las clases de primos, luego construyen la tabla de pesos enteros y finalmente recorren de forma recursiva los subconjuntos de primos.

Generación de los primos relevantes

La implementación comienza con una criba hasta \(B\), elimina el primo \(3\) y cuenta cuántos primos restantes pertenecen a \(\mathcal P_1\) y cuántos a \(\mathcal P_2\). Esos conteos fijan el tamaño de la tabla de pesos.

Construcción del triángulo de Pascal y de la tabla de pesos

Después se calculan los coeficientes binomiales \(\binom{n}{k}\) con la regla de Pascal. Con ellos se rellena la tabla \(w(a,b)\) usando la fórmula de inversión anterior. La interpretación con raíces de la unidad sirve para justificar el método, pero el cálculo real se hace por completo con enteros.

Enumeración recursiva de productos libres de cuadrados

La recursión principal decide, para cada primo relevante, si se excluye del subconjunto actual o si se incluye. Durante el recorrido se mantienen el producto actual \(d(S)\), el número de primos elegidos de \(\mathcal P_1\) y el número de primos elegidos de \(\mathcal P_2\). En una hoja se suma

$$w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor$$

al acumulado. Si incluir el siguiente primo haría que el producto superara \(N\), esa rama se poda inmediatamente, porque cualquier extensión posterior aportaría cero. Las versiones en C++ y Java además verifican unos pocos valores pequeños antes de imprimir el resultado final.

Análisis de Complejidad

Sean \(P=|\mathcal P_B|\), \(M_1=|\mathcal P_1|\), \(M_2=|\mathcal P_2|\) y \(D=\max(M_1,M_2)\). Construir la tabla binomial cuesta \(O(D^2)\) en tiempo y memoria. Rellenar la tabla de pesos mediante la fórmula directa de inversión cuesta \(O(M_1^2M_2^2)\) en tiempo y \(O(M_1M_2)\) en memoria; para \(B=120\), esto es muy pequeño.

La exploración recursiva de subconjuntos tiene tamaño máximo \(O(2^P)\), porque cada primo relevante puede estar o no estar. En la práctica es bastante menor, ya que las ramas se podan cuando el producto supera \(N\). En el caso real, \(P=29\), así que el coste queda dominado por un árbol de búsqueda manejable y no por el tamaño gigantesco de \(N=10^{18}\).

Notas y Referencias

  1. Página del problema: Project Euler 967 - B-Trivisible Numbers
  2. Principio de inclusión-exclusión: Wikipedia - Principio de inclusión-exclusión
  3. Transformada binomial e inversión: Wikipedia - Binomial transform
  4. Aritmética modular: Wikipedia - Aritmética modular
  5. Raíces de la unidad: Wikipedia - Raíz de la unidad

问题概述

固定一个上界 \(B\) 之后,真正相关的只有所有满足 \(p \le B\) 的素数。实现把这些素数按模 3 的两个非零剩余类分开:

$$\mathcal P_1=\{p \le B : p \text{ 是素数},\ p \equiv 1 \pmod{3}\}, \qquad \mathcal P_2=\{p \le B : p \text{ 是素数},\ p \equiv 2 \pmod{3}\}.$$

素数 \(3\) 被刻意排除,因为它不属于这两个非零剩余类,也不会改变最终要检测的同余条件。对任意整数 \(n\),定义

$$\omega_1(n)=\#\{p\in \mathcal P_1 : p\mid n\},\qquad \omega_2(n)=\#\{p\in \mathcal P_2 : p\mid n\}.$$

程序实际计算的是

$$F(N,B)=\#\{1\le n\le N : \omega_1(n)+2\omega_2(n)\equiv0\pmod3\}.$$

这也等价于要求 \(\omega_1(n)\equiv\omega_2(n)\pmod3\)。这里关心的只是某个相关素数是否整除 \(n\),完全不关心它出现了几次,所以素因子的指数在整个公式中都不会出现。对题目的实际参数 \(N=10^{18}\)、\(B=120\) 来说,逐个检查 \(1\) 到 \(N\) 的整数显然不可行。

数学方法

整个解法可以看成是对相关素数子集做一个二维容斥。核心目标是把 \(\omega_1(n)\) 与 \(\omega_2(n)\) 上的模 3 条件,改写成对平方自由素数乘积的加权求和。

相关素因子签名

$$\mathcal P_B=\mathcal P_1\cup\mathcal P_2.$$

对任意子集 \(S\subseteq \mathcal P_B\),定义

$$d(S)=\prod_{p\in S} p,\qquad a(S)=\#(S\cap\mathcal P_1),\qquad b(S)=\#(S\cap\mathcal P_2).$$

如果 \(d(S)\mid n\),就说明 \(S\) 里的每个素数都整除 \(n\)。由于题目关心的性质只取决于哪些相关素数整除 \(n\),每个整数都可以用签名 \((\omega_1(n),\omega_2(n))\) 来概括。\(3\) 的因子、已经出现过的素数的更高次幂,以及所有大于 \(B\) 的素因子,对这个签名都没有影响。

用三次单位根筛出同余条件

定义指示函数

$$G(x,y)=\begin{cases}1,& x+2y\equiv0\pmod3,\\0,& \text{否则}.\end{cases}$$

一个很自然的写法是使用三次单位根筛。若 \(\zeta\neq1\) 且 \(\zeta^3=1\),那么

$$G(x,y)=\frac{1+\zeta^{x+2y}+\zeta^{2x+y}}{3}.$$

这条恒等式直接说明:真正重要的只是 \(x\) 与 \(y\) 对 3 取模后的结果。实现本身并不真的用复数去算,但这层观点解释了为什么后面的权重表只需要依赖两类素因子的个数。

用二项反演构造整数权重

设某个整数 \(n\) 恰好有 \(x\) 个来自 \(\mathcal P_1\) 的相关素因子,以及 \(y\) 个来自 \(\mathcal P_2\) 的相关素因子。我们希望找到一个权重表 \(w(a,b)\),使得

$$G(x,y)=\sum_{a=0}^{x}\sum_{b=0}^{y}\binom{x}{a}\binom{y}{b}w(a,b).$$

它的组合意义很直接:从这 \(x\) 个第一类素数里任选 \(a\) 个,再从这 \(y\) 个第二类素数里任选 \(b\) 个。每一种选择都对应于 \(n\) 的相关素因子集合的一个子集。

对两个变量同时做二项反演,就得到

$$w(a,b)=\sum_{i=0}^{a}\sum_{j=0}^{b}(-1)^{a-i+b-j}\binom{a}{i}\binom{b}{j}G(i,j).$$

由于只有在 \(i+2j\equiv0\pmod3\) 时 \(G(i,j)\) 才不为零,上式可化为

$$w(a,b)=(-1)^{a+b}\sum_{\substack{0\le i\le a\\0\le j\le b\\i+2j\equiv0\!\!\!\pmod3}}(-1)^{i+j}\binom{a}{i}\binom{b}{j}.$$

这正是三种实现预先计算出来的权重公式。所有复杂的模 3 条件,都被压缩进这张小表里。

把计数压缩为平方自由乘积求和

对每个整数 \(n\),记 \(R_B(n)\subseteq\mathcal P_B\) 为所有整除 \(n\) 的相关素数集合,再记 \(\mathbf{1}_B(n)\) 为 \(n\) 是否为 \(B\)-trivisible 的指示函数。由反演公式立刻得到

$$\mathbf{1}_B(n)=\sum_{S\subseteq R_B(n)} w\big(a(S),b(S)\big).$$

接着对所有 \(n\le N\) 求和,并交换求和顺序。固定一个子集 \(S\) 之后,它会对每个 \(d(S)\) 的倍数贡献一次,因此总贡献正好是

$$\left\lfloor\frac{N}{d(S)}\right\rfloor.$$

于是得到核心公式

$$\boxed{F(N,B)=\sum_{S\subseteq\mathcal P_B} w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor.}$$

原本看起来必须在 \(1,\dots,N\) 上逐个判断的任务,就这样被压缩成了对所有相关素数子集的求和。每个素数只会被选一次或者不选,因此出现的都是平方自由乘积。

例子:\(B=10\)

当 \(B=10\) 时,相关素数是

$$\mathcal P_1=\{7\},\qquad \mathcal P_2=\{2,5\}.$$

条件是 \(\omega_1(n)+2\omega_2(n)\equiv0\pmod3\)。在 \(1\le n\le10\) 中,\(1,3,9\) 的签名都是 \((0,0)\),因此会被计入。\(2,4,5,6,8\) 的签名为 \((0,1)\),不满足条件。\(7\) 的签名是 \((1,0)\),\(10\) 的签名是 \((0,2)\),也都不满足。于是

$$F(10,10)=3,$$

这恰好与实现里的小样例校验一致。同样地,\(F(10,4)=5\) 也很好理解:若唯一相关素数只有 \(2\),那么被计入的就是所有奇数。

代码如何工作

C++、Python 和 Java 实现与上面的推导几乎是一一对应的。它们先预处理相关素数,再建立整数权重表,最后递归枚举素数子集。

生成相关素数

第一步先用筛法找出不超过 \(B\) 的所有素数,然后去掉 \(3\)。接着统计剩余素数中有多少个属于 \(\mathcal P_1\),多少个属于 \(\mathcal P_2\)。这些数量决定了后面权重表的大小。

构造杨辉三角和权重表

然后,程序按 Pascal 递推构造二项式系数 \(\binom{n}{k}\)。借助这些系数,按前面的反演公式填出 \(w(a,b)\) 整数表。三次单位根的观点只出现在数学说明里,真正的程序始终只做整数运算。

递归枚举平方自由乘积

主递归对每个相关素数做一次“选或不选”的分支。递归状态包含当前乘积 \(d(S)\)、从 \(\mathcal P_1\) 中选了多少个素数、以及从 \(\mathcal P_2\) 中选了多少个素数。走到叶子时,就把

$$w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor$$

加到答案里。如果继续加入下一个素数会使乘积超过 \(N\),这条分支会立刻剪掉,因为后续任何扩展的贡献都只能是零。C++ 和 Java 版本还会在输出最终结果前先验证几个小规模检查值。

复杂度分析

记 \(P=|\mathcal P_B|\)、\(M_1=|\mathcal P_1|\)、\(M_2=|\mathcal P_2|\),以及 \(D=\max(M_1,M_2)\)。构造二项式系数表需要 \(O(D^2)\) 时间和内存。按直接反演公式填充权重表需要 \(O(M_1^2M_2^2)\) 时间和 \(O(M_1M_2)\) 内存;对于 \(B=120\) 来说,这部分开销非常小。

递归枚举子集的最坏规模是 \(O(2^P)\),因为每个相关素数都可以选或不选。但实际运行要快得多,因为一旦当前乘积超过 \(N\),整条分支就会被剪掉。对本题而言,\(P=29\),因此整体代价由一个可控的深度优先搜索树决定,而不是由 \(N=10^{18}\) 的数量级决定。

注释与参考资料

  1. 题目页面: Project Euler 967 - B-Trivisible Numbers
  2. 容斥原理: Wikipedia - 容斥原理
  3. 二项变换与反演: Wikipedia - Binomial transform
  4. 模运算: Wikipedia - 模算术
  5. 单位根: Wikipedia - 单位根

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

Для фиксированного \(B\) важны только простые числа \(p \le B\). Они разбиваются на два ненулевых класса вычетов по модулю 3:

$$\mathcal P_1=\{p \le B : p \text{ простое},\ p \equiv 1 \pmod{3}\}, \qquad \mathcal P_2=\{p \le B : p \text{ простое},\ p \equiv 2 \pmod{3}\}.$$

Простое \(3\) специально исключается, потому что оно не принадлежит ни одному из этих двух классов и не влияет на проверяемую конгруэнцию. Для любого целого \(n\) зададим

$$\omega_1(n)=\#\{p\in \mathcal P_1 : p\mid n\},\qquad \omega_2(n)=\#\{p\in \mathcal P_2 : p\mid n\}.$$

Тогда вычисляется величина

$$F(N,B)=\#\{1\le n\le N : \omega_1(n)+2\omega_2(n)\equiv0\pmod3\}.$$

То же самое можно записать как условие \(\omega_1(n)\equiv\omega_2(n)\pmod3\). Важно, что учитывается только сам факт делимости на релевантный простой, а не его степень. Поэтому задача для \(N=10^{18}\) и \(B=120\) должна решаться без перебора всех чисел до \(N\).

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

Вся идея состоит в двумерном включении-исключении по подмножествам релевантных простых. Нужно превратить модульное условие на \(\omega_1(n)\) и \(\omega_2(n)\) в взвешенную сумму по бесквадратным произведениям простых.

Релевантная простая сигнатура

Положим

$$\mathcal P_B=\mathcal P_1\cup\mathcal P_2.$$

Для подмножества \(S\subseteq \mathcal P_B\) обозначим

$$d(S)=\prod_{p\in S} p,\qquad a(S)=\#(S\cap\mathcal P_1),\qquad b(S)=\#(S\cap\mathcal P_2).$$

Если \(d(S)\mid n\), то каждое простое из \(S\) делит \(n\). Так как свойство зависит только от того, какие релевантные простые делят \(n\), любое число можно сжать до сигнатуры \((\omega_1(n),\omega_2(n))\). Множители \(3\), более высокие степени уже присутствующих простых и простые делители больше \(B\) на эту сигнатуру не влияют.

Как выделяется нужная конгруэнция

Введем индикатор

$$G(x,y)=\begin{cases}1,& x+2y\equiv0\pmod3,\\0,& \text{иначе}.\end{cases}$$

Удобно выразить его через фильтр на кубических корнях единицы. Если \(\zeta\neq1\) и \(\zeta^3=1\), то

$$G(x,y)=\frac{1+\zeta^{x+2y}+\zeta^{2x+y}}{3}.$$

Эта формула сразу показывает, что важны только остатки \(x\) и \(y\) по модулю 3. В программе комплексные числа не используются, но именно этот взгляд объясняет структуру таблицы весов.

Биномиальная инверсия и таблица весов

Пусть число \(n\) имеет ровно \(x\) релевантных простых делителей из \(\mathcal P_1\) и ровно \(y\) из \(\mathcal P_2\). Требуется построить веса \(w(a,b)\), удовлетворяющие равенству

$$G(x,y)=\sum_{a=0}^{x}\sum_{b=0}^{y}\binom{x}{a}\binom{y}{b}w(a,b).$$

Комбинаторный смысл здесь прямой: мы выбираем любые \(a\) из \(x\) простых первого класса и любые \(b\) из \(y\) простых второго класса. Каждому такому выбору соответствует некоторое подмножество релевантных простых делителей числа \(n\).

Двумерная биномиальная инверсия дает

$$w(a,b)=\sum_{i=0}^{a}\sum_{j=0}^{b}(-1)^{a-i+b-j}\binom{a}{i}\binom{b}{j}G(i,j).$$

Поскольку \(G(i,j)\neq0\) только при \(i+2j\equiv0\pmod3\), формула упрощается до

$$w(a,b)=(-1)^{a+b}\sum_{\substack{0\le i\le a\\0\le j\le b\\i+2j\equiv0\!\!\!\pmod3}}(-1)^{i+j}\binom{a}{i}\binom{b}{j}.$$

Именно такую таблицу весов заранее строят реализации на C++, Python и Java.

Сведение подсчета к бесквадратным произведениям

Для каждого \(n\) обозначим через \(R_B(n)\subseteq\mathcal P_B\) множество релевантных простых делителей числа \(n\), а через \(\mathbf{1}_B(n)\) индикатор того, что \(n\) является \(B\)-trivisible. Тогда из инверсии следует

$$\mathbf{1}_B(n)=\sum_{S\subseteq R_B(n)} w\big(a(S),b(S)\big).$$

Теперь суммируем по всем \(n\le N\) и меняем порядок суммирования. Фиксированное подмножество \(S\) вносит вклад один раз для каждого кратного \(d(S)\), то есть суммарно

$$\left\lfloor\frac{N}{d(S)}\right\rfloor.$$

Получаем основную формулу

$$\boxed{F(N,B)=\sum_{S\subseteq\mathcal P_B} w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor.}$$

Именно здесь задача радикально упрощается: вместо перебора чисел от \(1\) до \(N\) остается сумма по подмножествам простых не больше \(B\). Появляются только бесквадратные произведения, потому что каждый простой либо выбран, либо нет.

Разобранный пример: \(B=10\)

При \(B=10\) релевантные простые таковы:

$$\mathcal P_1=\{7\},\qquad \mathcal P_2=\{2,5\}.$$

Нужно проверить условие \(\omega_1(n)+2\omega_2(n)\equiv0\pmod3\). Среди чисел \(1\le n\le10\) числа \(1,3,9\) имеют сигнатуру \((0,0)\) и засчитываются. Числа \(2,4,5,6,8\) имеют сигнатуру \((0,1)\), поэтому не подходят. У числа \(7\) сигнатура \((1,0)\), у числа \(10\) сигнатура \((0,2)\); они тоже не засчитываются. Следовательно,

$$F(10,10)=3,$$

что в точности совпадает с контрольным значением в реализации. Тот же аргумент объясняет и \(F(10,4)=5\): если единственный релевантный простой равен \(2\), то подходят ровно нечетные числа.

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

Реализации на C++, Python и Java очень близко следуют описанной выше математике. Сначала они подготавливают классы простых, затем строят целочисленную таблицу весов и после этого рекурсивно перебирают подмножества простых.

Построение списка релевантных простых

Сначала запускается решето до \(B\), затем из списка удаляется простой \(3\). После этого считается, сколько оставшихся простых лежит в \(\mathcal P_1\) и сколько в \(\mathcal P_2\). Эти значения задают размеры таблицы весов.

Таблицы биномиальных коэффициентов и весов

Далее по правилу Паскаля строятся коэффициенты \(\binom{n}{k}\). На их основе по формуле инверсии заполняется таблица \(w(a,b)\). Подход с корнями единицы используется только как математическое объяснение; сама программа выполняет лишь целочисленные операции.

Рекурсивный обход бесквадратных произведений

Основная рекурсия для каждого релевантного простого выбирает один из двух вариантов: не брать его в текущее подмножество или взять. В состоянии хранятся текущее произведение \(d(S)\), число выбранных простых из \(\mathcal P_1\) и число выбранных простых из \(\mathcal P_2\). В листе к ответу прибавляется

$$w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor.$$

Если добавление следующего простого уже делает произведение больше \(N\), ветвь немедленно отсекается, потому что любой ее дальнейший потомок даст нулевой вклад. Версии на C++ и Java также проверяют несколько небольших контрольных значений перед выводом окончательного ответа.

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

Пусть \(P=|\mathcal P_B|\), \(M_1=|\mathcal P_1|\), \(M_2=|\mathcal P_2|\), а \(D=\max(M_1,M_2)\). Построение таблицы биномиальных коэффициентов требует \(O(D^2)\) времени и памяти. Заполнение таблицы весов по прямой формуле инверсии требует \(O(M_1^2M_2^2)\) времени и \(O(M_1M_2)\) памяти; при \(B=120\) это очень мало.

Рекурсивный перебор подмножеств имеет худшую оценку \(O(2^P)\), потому что каждый релевантный простой либо входит в подмножество, либо нет. На практике работа заметно быстрее, так как ветви обрезаются, как только их произведение превосходит \(N\). Для данной задачи \(P=29\), поэтому основная стоимость связана с вполне управляемым деревом поиска, а не с самим масштабом \(N=10^{18}\).

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

  1. Страница задачи: Project Euler 967 - B-Trivisible Numbers
  2. Принцип включения-исключения: Wikipedia - Принцип включения-исключения
  3. Биномиальное преобразование и инверсия: Wikipedia - Binomial transform
  4. Модульная арифметика: Wikipedia - Модульная арифметика
  5. Корни из единицы: Wikipedia - Корни из единицы

ملخص المسألة

عند تثبيت حد أعلى \(B\)، فإن الأعداد الأولية الوحيدة التي تهم هي \(p \le B\). تقسم هذه الأوليات إلى فئتين غير صفريتين بترديد 3:

$$\mathcal P_1=\{p \le B : p \text{ prime},\ p \equiv 1 \pmod{3}\}, \qquad \mathcal P_2=\{p \le B : p \text{ prime},\ p \equiv 2 \pmod{3}\}.$$

يستبعد العدد الأولي \(3\) عمدًا، لأنه لا ينتمي إلى أي من هاتين الفئتين ولا يغير شرط التطابق الذي يحكم العد. ولكل عدد صحيح \(n\) نعرّف

$$\omega_1(n)=\#\{p\in \mathcal P_1 : p\mid n\},\qquad \omega_2(n)=\#\{p\in \mathcal P_2 : p\mid n\}.$$

وعندئذ تكون الكمية المطلوبة هي

$$F(N,B)=\#\{1\le n\le N : \omega_1(n)+2\omega_2(n)\equiv0\pmod3\}.$$

وبصورة مكافئة، فإن الأعداد المحسوبة هي التي تحقق \(\omega_1(n)\equiv\omega_2(n)\pmod3\). المهم هنا هو مجرد وجود القاسم الأولي ذي الصلة أو عدم وجوده؛ أما الأسس فلا تؤثر إطلاقًا. لذلك لا يمكن التعامل مع الحالة الفعلية \(N=10^{18}\) و\(B=120\) بفحص الأعداد واحدًا واحدًا.

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

الفكرة كلها هي شمول واستبعاد ثنائي الأبعاد على مجموعات فرعية من الأوليات ذات الصلة. نريد تحويل الشرط المعياري على \(\omega_1(n)\) و\(\omega_2(n)\) إلى مجموع موزون على جداءات أولية خالية من المربعات.

التوقيع الأولي ذي الصلة

لنكتب

$$\mathcal P_B=\mathcal P_1\cup\mathcal P_2.$$

ولكل مجموعة فرعية \(S\subseteq \mathcal P_B\) نعرّف

$$d(S)=\prod_{p\in S} p,\qquad a(S)=\#(S\cap\mathcal P_1),\qquad b(S)=\#(S\cap\mathcal P_2).$$

إذا كان \(d(S)\mid n\)، فهذا يعني أن كل أولي في \(S\) يقسم \(n\). وبما أن الخاصية تعتمد فقط على أي الأوليات ذات الصلة تقسم \(n\)، فإن كل عدد يمكن تلخيصه بالتوقيع \((\omega_1(n),\omega_2(n))\). عوامل \(3\)، والقوى الأعلى للأوليات نفسها، والأوليات الأكبر من \(B\) كلها لا تؤثر في هذا التوقيع.

استخراج شرط التطابق

عرّف دالة المؤشر

$$G(x,y)=\begin{cases}1,& x+2y\equiv0\pmod3,\\0,& \text{otherwise}.\end{cases}$$

وهناك طريقة أنيقة لعزل هذا الشرط باستعمال مرشح جذور الوحدة التكعيبية. إذا كان \(\zeta\neq1\) و\(\zeta^3=1\)، فإن

$$G(x,y)=\frac{1+\zeta^{x+2y}+\zeta^{2x+y}}{3}.$$

وهذا يوضح مباشرة أن المهم هو العدّان بترديد 3 فقط. التنفيذ لا يستعمل أعدادًا عقدية فعلًا، لكن هذه الهوية تشرح لماذا يكفي جدول أوزان يعتمد على عدد الأوليات في كل فئة.

العكس الثنائي الحد يبني جدول الأوزان الصحيح

افترض أن عددًا \(n\) يملك بالضبط \(x\) من الأوليات ذات الصلة من \(\mathcal P_1\)، و\(y\) من \(\mathcal P_2\). نريد أوزانًا \(w(a,b)\) بحيث

$$G(x,y)=\sum_{a=0}^{x}\sum_{b=0}^{y}\binom{x}{a}\binom{y}{b}w(a,b).$$

المعنى التوافقي مباشر: نختار أي \(a\) من الأوليات \(x\) في الفئة الأولى، وأي \(b\) من الأوليات \(y\) في الفئة الثانية. وكل اختيار من هذا النوع يقابل مجموعة فرعية من الأوليات ذات الصلة التي تقسم \(n\).

بتطبيق العكس الثنائي الحد في المتغيرين نحصل على

$$w(a,b)=\sum_{i=0}^{a}\sum_{j=0}^{b}(-1)^{a-i+b-j}\binom{a}{i}\binom{b}{j}G(i,j).$$

وبما أن \(G(i,j)\) لا يكون غير صفري إلا عندما \(i+2j\equiv0\pmod3\)، فإن الصيغة تصبح

$$w(a,b)=(-1)^{a+b}\sum_{\substack{0\le i\le a\\0\le j\le b\\i+2j\equiv0\!\!\!\pmod3}}(-1)^{i+j}\binom{a}{i}\binom{b}{j}.$$

وهذه هي بالضبط الصيغة التي تسبقها التطبيقات في البناء المسبق لجدول الأوزان.

اختزال العد إلى جداءات خالية من المربعات

لكل عدد \(n\)، ليكن \(R_B(n)\subseteq\mathcal P_B\) مجموعة الأوليات ذات الصلة التي تقسمه، ولتكن \(\mathbf{1}_B(n)\) دالة المؤشر لكون \(n\) عددًا \(B\)-trivisible. من صيغة العكس ينتج أن

$$\mathbf{1}_B(n)=\sum_{S\subseteq R_B(n)} w\big(a(S),b(S)\big).$$

إذا جمعنا الآن على جميع \(n\le N\) وبدلنا ترتيب الجمع، فإن المجموعة الفرعية الثابتة \(S\) تسهم مرة واحدة لكل مضاعف لـ \(d(S)\)، ولذلك يكون إسهامها الكلي

$$\left\lfloor\frac{N}{d(S)}\right\rfloor.$$

ومن ثم نحصل على الصيغة المركزية

$$\boxed{F(N,B)=\sum_{S\subseteq\mathcal P_B} w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor.}$$

هنا تتحول المسألة من فحص جميع الأعداد من \(1\) إلى \(N\) إلى مجموع على مجموعات فرعية من الأوليات حتى \(B\). ولا تظهر إلا الجداءات الخالية من المربعات لأن كل أولي يختار مرة واحدة أو لا يختار.

مثال محلول: \(B=10\)

عندما يكون \(B=10\)، فإن الأوليات ذات الصلة هي

$$\mathcal P_1=\{7\},\qquad \mathcal P_2=\{2,5\}.$$

والشرط هو \(\omega_1(n)+2\omega_2(n)\equiv0\pmod3\). بين الأعداد \(1\le n\le10\)، تمتلك الأعداد \(1,3,9\) التوقيع \((0,0)\) ولذلك تُحسب. أما الأعداد \(2,4,5,6,8\) فتمتلك التوقيع \((0,1)\) فلا تحقق الشرط. والعدد \(7\) له التوقيع \((1,0)\)، والعدد \(10\) له التوقيع \((0,2)\)، ولذلك يستبعدان أيضًا. إذًا

$$F(10,10)=3,$$

وهو نفس مقدار التحقق الصغير الموجود في التطبيقات. وبالمنطق نفسه نفهم أيضًا أن \(F(10,4)=5\): إذا كان الأولي الوحيد ذي الصلة هو \(2\)، فإن الأعداد المحسوبة هي الأعداد الفردية تمامًا.

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

تتبع تطبيقات C++ وPython وJava هذا الاشتقاق بصورة شبه مباشرة. فهي تجهز أولًا فئات الأوليات، ثم تبني جدول الأوزان الصحيح، وبعد ذلك تعدد المجموعات الفرعية للأوليات بطريقة تكرارية قائمة على الاستدعاء الذاتي.

توليد الأوليات ذات الصلة

تبدأ الشيفرة بمنخل حتى \(B\)، ثم تحذف العدد الأولي \(3\)، ثم تحصي كم أولي متبقٍ ينتمي إلى \(\mathcal P_1\) وكم أولي ينتمي إلى \(\mathcal P_2\). وهذه الأعداد هي التي تحدد أبعاد جدول الأوزان.

بناء مثلث باسكال وجدول الأوزان

بعد ذلك تُبنى معاملات ثنائية الحد \(\binom{n}{k}\) باستعمال قاعدة باسكال. ومن خلال هذه المعاملات يملأ الجدول \(w(a,b)\) وفق صيغة العكس السابقة. منظور جذور الوحدة يبقى مجرد تفسير رياضي، أما الحساب الفعلي فيبقى كله حسابًا صحيحًا على الأعداد الصحيحة.

الاستعراض التكراري للجداءات الخالية من المربعات

الجزء الرئيسي من الخوارزمية يقرر لكل أولي ذي صلة هل يدخل في المجموعة الفرعية الحالية أم لا. وأثناء هذا الاستعراض يحتفظ التنفيذ بالجداء الحالي \(d(S)\)، وعدد الأوليات المختارة من \(\mathcal P_1\)، وعدد الأوليات المختارة من \(\mathcal P_2\). وعند الوصول إلى ورقة في شجرة البحث يضاف المقدار

$$w\big(a(S),b(S)\big)\left\lfloor\frac{N}{d(S)}\right\rfloor$$

إلى الجواب. وإذا كان إدخال الأولي التالي سيجعل الجداء يتجاوز \(N\)، فإن ذلك الفرع يقص فورًا لأن أي امتداد أعمق له سيعطي مساهمة صفرية. كما أن نسختي C++ وJava تتحققان من بعض القيم الصغيرة قبل طباعة النتيجة النهائية.

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

لنكتب \(P=|\mathcal P_B|\)، و\(M_1=|\mathcal P_1|\)، و\(M_2=|\mathcal P_2|\)، و\(D=\max(M_1,M_2)\). بناء جدول معاملات ثنائية الحد يكلف \(O(D^2)\) زمنًا وذاكرة. أما ملء جدول الأوزان وفق صيغة العكس المباشرة فيكلف \(O(M_1^2M_2^2)\) زمنًا و\(O(M_1M_2)\) ذاكرة، وهو صغير جدًا في حالة \(B=120\).

الاستعراض التكراري للمجموعات الفرعية له حد أسوأ يساوي \(O(2^P)\)، لأن كل أولي ذي صلة قد يكون موجودًا أو غير موجود. عمليًا يكون أسرع بكثير، لأن الفروع تقص عندما يتجاوز الجداء \(N\). في المسألة الفعلية لدينا \(P=29\)، ولذلك تهيمن شجرة بحث قابلة للإدارة على زمن التنفيذ بدلًا من الحجم الهائل لـ \(N=10^{18}\).

هوامش ومراجع

  1. صفحة المسألة: Project Euler 967 - B-Trivisible Numbers
  2. مبدأ الشمول والاستبعاد: Wikipedia - مبدأ الشمول والاستبعاد
  3. التحويل الثنائي الحد والعكس: Wikipedia - Binomial transform
  4. الحسابيات المعيارية: Wikipedia - حسابيات معيارية
  5. جذور الوحدة: Wikipedia - جذور الوحدة