We must count integers in the half-open range \(1 \le x < N\) with \(N=10^{16}\) such that \(x\) is divisible by at least four distinct primes below \(100\). Multiplicity does not matter: a prime like \(2\) contributes only once, whether \(x\) is divisible by \(2\) or by \(2^{20}\).
The C++ program exposes the bound as --limit=<N>, so the same method works for any \(N\ge 1\).
Let
$$\mathcal P=\{p \text{ prime}: p<100\}.$$
There are \(25\) such primes. For a subset \(S\subseteq\mathcal P\), define its squarefree product
$$q_S=\prod_{p\in S}p.$$
An integer \(x\) is divisible by every prime in \(S\) exactly when \(q_S\mid x\). Therefore the number of integers \(x<N\) divisible by all primes of \(S\) is
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
The code uses \(N-1\) because the searched range is \(1\le x<N\), not \(1\le x\le N\).
For each \(k\ge4\), define
$$T_k=\sum_{|S|=k}\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
So \(T_k\) counts, with multiplicity, how many numbers are divisible by every prime in some \(k\)-element subset of \(\mathcal P\).
If a number \(x\) has exactly \(r\) distinct prime divisors from \(\mathcal P\), then:
1. it appears in \(T_k\) exactly \(\binom{r}{k}\) times, because we may choose any \(k\) of those \(r\) primes,
2. it contributes nothing to \(T_k\) when \(k>r\).
So the whole problem is to find coefficients \(a_k\) such that
$$\sum_{k=4}^{r} a_k\binom{r}{k}= \begin{cases} 0,&r<4,\\ 1,&r\ge4. \end{cases}$$
The program uses
$$a_k=(-1)^{k-4}\binom{k-1}{3}.$$
Thus the final answer is
$$\sum_{k=4}^{25} a_kT_k =\sum_{k=4}^{25}(-1)^{k-4}\binom{k-1}{3}T_k.$$
The first few coefficients are
$$a_4=1,\qquad a_5=-4,\qquad a_6=10,\qquad a_7=-20,\dots$$
So the formula starts as
$$T_4-4T_5+10T_6-20T_7+\cdots.$$
Take a number with exactly \(r\) eligible prime factors. Its total weight becomes
$$w(r)=\sum_{k=4}^{r}(-1)^{k-4}\binom{k-1}{3}\binom{r}{k}.$$
For \(r<4\), this sum is empty, so \(w(r)=0\), exactly as required.
For \(r\ge4\), there is a neat closed-form evaluation. First use
$$\binom{r}{k}\binom{k-1}{3}=\binom{r}{4}\binom{r-4}{k-4}\frac{4}{k}.$$
Now write \(j=k-4\). Then
$$w(r)=4\binom{r}{4}\sum_{j=0}^{r-4}(-1)^j\binom{r-4}{j}\frac{1}{j+4}.$$
Using \(\frac{1}{j+4}=\int_0^1 x^{j+3}\,dx\), we get
$$w(r)=4\binom{r}{4}\int_0^1 x^3(1-x)^{r-4}\,dx.$$
This is a beta integral:
$$\int_0^1 x^3(1-x)^{r-4}\,dx=\frac{3!(r-4)!}{r!}.$$
Therefore
$$w(r)=4\cdot\frac{r!}{4!(r-4)!}\cdot\frac{3!(r-4)!}{r!}=1.$$
So every number with at least four distinct eligible primes contributes exactly once, and every number with fewer than four contributes zero.
The number \(2310=2\cdot3\cdot5\cdot7\cdot11\) has \(r=5\) eligible prime factors.
It appears in \(T_4\) exactly \(\binom54=5\) times, once for each choice of four of its five primes. It appears in \(T_5\) exactly \(\binom55=1\) time. It appears in no \(T_k\) with \(k\ge6\).
Its final weighted contribution is therefore
$$1\cdot5+(-4)\cdot1=1.$$
This is the simplest way to see why the coefficients correct the heavy overcount from \(T_4\).
The program never enumerates all \(2^{25}\) subsets explicitly. Instead it performs a depth-first search over increasing prime indices.
At a DFS state, the code knows:
1. the next prime index it may use,
2. how many primes have already been taken,
3. the current squarefree product \(q_S\).
If taken >= 4, the code immediately adds
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor$$
to subset_sums[taken].
Before multiplying by a new prime \(p\), it checks
$$q_S \le \frac{N-1}{p},$$
implemented as product > max_value / p. This avoids overflow and prunes any branch whose product already exceeds \(N-1\).
The first checkpoint is
$$\texttt{solve}(1000)=23.$$
Indeed, the numbers below \(1000\) with at least four distinct primes below \(100\) are exactly \(23\) in number; examples include \(210,330,390,462,\dots,990\).
The second checkpoint compares
$$\texttt{solve}(100000)$$
against a brute-force routine that explicitly counts distinct prime divisors for every \(x<100000\). This confirms that the weighted subset formula and the DFS implementation match the direct definition.
primes_below_100() generates the 25 primes less than \(100\) with a small linear sieve.
solve(limit) sets max_value = limit - 1, runs the DFS, stores all \(T_k\)-type totals in subset_sums, then forms the weighted sum with coefficients \((-1)^{k-4}\binom{k-1}{3}\).
binom is only used for those small coefficients, so 64-bit signed integers are sufficient.
The command-line interface accepts --limit=<N> and --skip-checkpoints.
The theoretical subset space has size \(2^{25}\), but in practice most branches are cut immediately once the squarefree product becomes too large. The main cost is therefore the number of valid subset products \(q_S\le N-1\), not the full power set.
The final weighted sum over \(k=4,\dots,25\) is tiny compared with the DFS itself, and memory usage is \(O(25)\) beyond the recursion stack because only the prime list and the 26-element accumulator array are stored.
Gesucht ist die Anzahl der Zahlen im halboffenen Bereich \(1\le x<N\) mit \(N=10^{16}\), die durch mindestens vier verschiedene Primzahlen kleiner als \(100\) teilbar sind. Vielfachheiten spielen keine Rolle: Ob \(2\) einmal oder zwanzigmal vorkommt, zählt nur als ein verschiedener Primfaktor.
Das C++-Programm parametrisiert die Grenze über --limit=<N>.
Sei
$$\mathcal P=\{p \text{ prim}: p<100\}.$$
Diese Menge hat \(25\) Elemente. Für eine Teilmenge \(S\subseteq\mathcal P\) definieren wir das quadratfreie Produkt
$$q_S=\prod_{p\in S}p.$$
Eine Zahl \(x\) ist genau dann durch alle Primzahlen aus \(S\) teilbar, wenn \(q_S\mid x\). Also ist die Anzahl solcher \(x<N\)
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
Das \(N-1\) erscheint, weil der Suchbereich \(1\le x<N\) ist.
Für \(k\ge4\) setzen wir
$$T_k=\sum_{|S|=k}\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
Hat eine Zahl \(x\) genau \(r\) verschiedene Primfaktoren aus \(\mathcal P\), dann erscheint sie in \(T_k\) genau \(\binom{r}{k}\)-mal.
Wir brauchen also Koeffizienten \(a_k\) mit
$$\sum_{k=4}^{r}a_k\binom{r}{k}= \begin{cases} 0,&r<4,\\ 1,&r\ge4. \end{cases}$$
Der Code verwendet
$$a_k=(-1)^{k-4}\binom{k-1}{3}.$$
Damit lautet die Endformel
$$\sum_{k=4}^{25}a_kT_k =\sum_{k=4}^{25}(-1)^{k-4}\binom{k-1}{3}T_k.$$
Die ersten Gewichte sind
$$a_4=1,\qquad a_5=-4,\qquad a_6=10,\qquad a_7=-20,\dots$$
also beginnt der Ausdruck als
$$T_4-4T_5+10T_6-20T_7+\cdots.$$
Für eine Zahl mit genau \(r\) zulässigen Primfaktoren ist das Gesamtgewicht
$$w(r)=\sum_{k=4}^{r}(-1)^{k-4}\binom{k-1}{3}\binom{r}{k}.$$
Für \(r<4\) ist die Summe leer, also \(w(r)=0\).
Für \(r\ge4\) benutzen wir die Identität
$$\binom{r}{k}\binom{k-1}{3}=\binom{r}{4}\binom{r-4}{k-4}\frac{4}{k}.$$
Mit \(j=k-4\) folgt
$$w(r)=4\binom{r}{4}\sum_{j=0}^{r-4}(-1)^j\binom{r-4}{j}\frac{1}{j+4}.$$
Nun \(\frac{1}{j+4}=\int_0^1 x^{j+3}\,dx\), also
$$w(r)=4\binom{r}{4}\int_0^1 x^3(1-x)^{r-4}\,dx.$$
Das Integral ist eine Beta-Funktion:
$$\int_0^1 x^3(1-x)^{r-4}\,dx=\frac{3!(r-4)!}{r!}.$$
Daraus ergibt sich
$$w(r)=4\cdot\frac{r!}{4!(r-4)!}\cdot\frac{3!(r-4)!}{r!}=1.$$
Also zählt jede Zahl mit mindestens vier zulässigen Primfaktoren genau einmal.
Die Zahl \(2310=2\cdot3\cdot5\cdot7\cdot11\) hat \(r=5\) zulässige Primfaktoren.
Sie erscheint \(\binom54=5\)-mal in \(T_4\) und \(\binom55=1\)-mal in \(T_5\). Höhere \(T_k\) tragen nichts bei. Das Gesamtgewicht ist also
$$1\cdot5+(-4)\cdot1=1.$$
Genau so korrigieren die Gewichte die starke Überzählung in \(T_4\).
Das Programm läuft nicht naiv über alle \(2^{25}\) Teilmengen. Stattdessen benutzt es eine Tiefensuche über aufsteigende Primindizes.
Ein DFS-Zustand enthält
1. den nächsten erlaubten Primindex,
2. die bereits gewählte Anzahl von Primzahlen,
3. das aktuelle Produkt \(q_S\).
Sobald taken >= 4 gilt, wird
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor$$
zu subset_sums[taken] addiert.
Vor der Multiplikation mit einer neuen Primzahl \(p\) prüft der Code
$$q_S\le\frac{N-1}{p},$$
implementiert als product > max_value / p. Das verhindert Überlauf und schneidet unmögliche Äste sofort ab.
Der erste Kontrollpunkt ist
$$\texttt{solve}(1000)=23.$$
Tatsächlich gibt es unter \(1000\) genau \(23\) Zahlen mit mindestens vier verschiedenen Primfaktoren kleiner als \(100\); Beispiele sind \(210,330,390,462,\dots,990\).
Der zweite Kontrollpunkt vergleicht
$$\texttt{solve}(100000)$$
mit einer Brute-Force-Funktion, die für jedes \(x<100000\) die verschiedenen Primfaktoren explizit zählt. So wird bestätigt, dass gewichtete Teilmengenformel und direkte Definition übereinstimmen.
primes_below_100() erzeugt mit einem kleinen linearen Sieb die \(25\) Primzahlen kleiner als \(100\).
solve(limit) setzt max_value = limit - 1, führt die DFS aus, sammelt alle \(T_k\)-Werte in subset_sums und bildet anschließend die gewichtete Summe mit \((-1)^{k-4}\binom{k-1}{3}\).
binom wird nur für diese kleinen Koeffizienten benötigt, daher reichen 64-Bit-Integer.
Die Kommandozeile akzeptiert --limit=<N> und --skip-checkpoints.
Theoretisch gibt es \(2^{25}\) Teilmengen, praktisch werden aber die meisten Zweige schon sehr früh wegen der Produktgrenze abgeschnitten. Entscheidend ist also die Anzahl gültiger Produkte \(q_S\le N-1\), nicht die volle Potenzmenge.
Die abschließende Summe über \(k=4,\dots,25\) ist gegenüber der DFS vernachlässigbar, und der Speicherverbrauch ist bis auf den Rekursionsstack \(O(25)\), da nur die Primliste und ein Akkumulatorarray der Länge \(26\) gehalten werden.
\(N=10^{16}\) için \(1\le x<N\) aralığındaki, \(100\)'den küçük en az dört farklı asal ile bölünebilen sayıların sayısı istenir. Asal üsleri önemli değildir; örneğin \(2\) sayısı \(x\)'i bir kez de bölse, on kez de bölse yalnızca tek bir farklı asal bölen sayılır.
C++ programı üst sınırı --limit=<N> ile parametreleştirir.
$$\mathcal P=\{p \text{ asal}: p<100\}$$
olsun. Bu kümede \(25\) asal vardır. Her altküme \(S\subseteq\mathcal P\) için karesiz çarpımı
$$q_S=\prod_{p\in S}p$$
olarak tanımlarız. Bir \(x\) sayısı, \(S\) içindeki bütün asallara tam bölünüyorsa ve ancak o zaman \(q_S\mid x\) olur. Dolayısıyla \(x<N\) aralığında bu koşulu sağlayan sayı sayısı
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor$$
şeklindedir. Burada \(N-1\) kullanılmasının sebebi aralığın \(1\le x<N\) olmasıdır.
\(k\ge4\) için
$$T_k=\sum_{|S|=k}\left\lfloor\frac{N-1}{q_S}\right\rfloor$$
tanımını yapalım. Eğer bir sayı \(\mathcal P\) içinden tam \(r\) farklı asal bölen içeriyorsa, \(T_k\) içinde tam \(\binom{r}{k}\) kez görünür; çünkü bu \(r\) asaldan herhangi \(k\) tanesi seçilebilir.
O halde amacımız, şu koşulu sağlayan ağırlıklar bulmaktır:
$$\sum_{k=4}^{r}a_k\binom{r}{k}= \begin{cases} 0,&r<4,\\ 1,&r\ge4. \end{cases}$$
Kodun kullandığı ağırlıklar
$$a_k=(-1)^{k-4}\binom{k-1}{3}$$
şeklindedir. Böylece nihai formül
$$\sum_{k=4}^{25}a_kT_k =\sum_{k=4}^{25}(-1)^{k-4}\binom{k-1}{3}T_k$$
olur. İlk katsayılar
$$a_4=1,\qquad a_5=-4,\qquad a_6=10,\qquad a_7=-20,\dots$$
yani ifade
$$T_4-4T_5+10T_6-20T_7+\cdots$$
şeklinde başlar.
Tam \(r\) adet uygun asal bölen içeren bir sayının toplam ağırlığı
$$w(r)=\sum_{k=4}^{r}(-1)^{k-4}\binom{k-1}{3}\binom{r}{k}$$
olur. \(r<4\) için toplam boştur, dolayısıyla \(w(r)=0\).
\(r\ge4\) için şu özdeşliği kullanırız:
$$\binom{r}{k}\binom{k-1}{3}=\binom{r}{4}\binom{r-4}{k-4}\frac{4}{k}.$$
\(j=k-4\) yazarsak
$$w(r)=4\binom{r}{4}\sum_{j=0}^{r-4}(-1)^j\binom{r-4}{j}\frac{1}{j+4}.$$
Burada \(\frac{1}{j+4}=\int_0^1 x^{j+3}\,dx\) olduğundan
$$w(r)=4\binom{r}{4}\int_0^1 x^3(1-x)^{r-4}\,dx.$$
Bu integral bir beta integrali olup
$$\int_0^1 x^3(1-x)^{r-4}\,dx=\frac{3!(r-4)!}{r!}$$
değerine eşittir. Sonuç olarak
$$w(r)=4\cdot\frac{r!}{4!(r-4)!}\cdot\frac{3!(r-4)!}{r!}=1.$$
Yani en az dört uygun asal böleni olan her sayı tam bir kez, daha azı ise sıfır kez sayılır.
\(2310=2\cdot3\cdot5\cdot7\cdot11\) sayısının \(r=5\) uygun asal böleni vardır.
Bu sayı, \(T_4\) içinde \(\binom54=5\) kez, \(T_5\) içinde \(\binom55=1\) kez görünür; daha büyük \(T_k\) toplamlarına katkı yapmaz. Toplam ağırlık bu yüzden
$$1\cdot5+(-4)\cdot1=1$$
olur. Böylece \(T_4\)'teki fazla sayım tam olarak düzeltilmiş olur.
Program bütün \(2^{25}\) altkümeyi körlemesine gezmez. Bunun yerine artan asal indeksleri üzerinde derinlik öncelikli arama yapar.
Bir DFS durumunda kod şu bilgileri taşır:
1. sıradaki kullanılabilir asal indeksi,
2. şimdiye kadar seçilen asal sayısı,
3. mevcut karesiz çarpım \(q_S\).
Eğer taken >= 4 ise
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor$$
değeri doğrudan subset_sums[taken] içine eklenir.
Yeni bir asal \(p\) ile çarpmadan önce
$$q_S\le\frac{N-1}{p}$$
kontrol edilir; kodda bu, product > max_value / p testiyle yapılır. Böylece taşma önlenir ve artık ürün sınırını aşan dallar erken kesilir.
İlk checkpoint
$$\texttt{solve}(1000)=23$$
eşitliğidir. Gerçekten de \(1000\)'den küçük olup en az dört farklı uygun asal bölen içeren tam \(23\) sayı vardır; örnekler arasında \(210,330,390,462,\dots,990\) bulunur.
İkinci checkpoint ise
$$\texttt{solve}(100000)$$
değerini, her \(x<100000\) için farklı asal bölen sayısını açıkça hesaplayan brute-force yordamıyla karşılaştırır. Böylece hem ağırlıklı altküme formülü hem de DFS uygulaması doğrudan tanımla çapraz doğrulanmış olur.
primes_below_100(), küçük bir lineer elekle \(100\)'den küçük 25 asalı üretir.
solve(limit), max_value = limit - 1 değerini kurar, DFS'i çalıştırır, bütün \(T_k\)-tipi toplamları subset_sums dizisinde biriktirir ve sonunda \((-1)^{k-4}\binom{k-1}{3}\) katsayılarıyla ağırlıklı toplamı alır.
binom fonksiyonu yalnızca bu küçük katsayılar için kullanıldığı için 64-bit signed tamsayı yeterlidir.
Komut satırında --limit=<N> ve --skip-checkpoints seçenekleri vardır.
Kuramsal altküme uzayı \(2^{25}\) olsa da, karesiz çarpım sınırı yüzünden dalların büyük kısmı çok erken budanır. Bu yüzden gerçek maliyeti belirleyen şey tam kuvvet kümesi değil, \(q_S\le N-1\) koşulunu sağlayan geçerli altküme çarpımlarının sayısıdır.
Son aşamadaki \(k=4,\dots,25\) ağırlıklı toplaması DFS yanında ihmal edilebilir. Bellek kullanımı da özyineleme yığını dışında \(O(25)\) düzeyindedir; çünkü yalnızca asal listesi ile 26 elemanlı bir sayaç dizisi tutulur.
Debemos contar los enteros del intervalo semiabierto \(1\le x<N\), con \(N=10^{16}\), que sean divisibles por al menos cuatro primos distintos menores que \(100\). Las multiplicidades no importan: que \(2\) divida una vez o muchas veces solo cuenta como un primo distinto.
El programa de C++ expone el límite mediante --limit=<N>.
Sea
$$\mathcal P=\{p \text{ primo}: p<100\}.$$
Hay \(25\) primos en ese conjunto. Para cada subconjunto \(S\subseteq\mathcal P\), definimos el producto libre de cuadrados
$$q_S=\prod_{p\in S}p.$$
Un entero \(x\) es divisible por todos los primos de \(S\) si y solo si \(q_S\mid x\). Por tanto, el número de enteros \(x<N\) con esa propiedad es
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
Para \(k\ge4\), definimos
$$T_k=\sum_{|S|=k}\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
Si un número \(x\) tiene exactamente \(r\) primos distintos de \(\mathcal P\), entonces aparece en \(T_k\) exactamente \(\binom{r}{k}\) veces.
Necesitamos coeficientes \(a_k\) tales que
$$\sum_{k=4}^{r}a_k\binom{r}{k}= \begin{cases} 0,&r<4,\\ 1,&r\ge4. \end{cases}$$
El código usa
$$a_k=(-1)^{k-4}\binom{k-1}{3}.$$
Así, la cuenta final es
$$\sum_{k=4}^{25}a_kT_k =\sum_{k=4}^{25}(-1)^{k-4}\binom{k-1}{3}T_k.$$
Los primeros coeficientes son
$$a_4=1,\qquad a_5=-4,\qquad a_6=10,\qquad a_7=-20,\dots$$
de modo que la fórmula empieza como
$$T_4-4T_5+10T_6-20T_7+\cdots.$$
Si un número tiene exactamente \(r\) primos válidos, su peso total es
$$w(r)=\sum_{k=4}^{r}(-1)^{k-4}\binom{k-1}{3}\binom{r}{k}.$$
Para \(r<4\), la suma es vacía, así que \(w(r)=0\).
Para \(r\ge4\), usamos
$$\binom{r}{k}\binom{k-1}{3}=\binom{r}{4}\binom{r-4}{k-4}\frac{4}{k}.$$
Escribiendo \(j=k-4\), obtenemos
$$w(r)=4\binom{r}{4}\sum_{j=0}^{r-4}(-1)^j\binom{r-4}{j}\frac{1}{j+4}.$$
Ahora \(\frac{1}{j+4}=\int_0^1 x^{j+3}\,dx\), así que
$$w(r)=4\binom{r}{4}\int_0^1 x^3(1-x)^{r-4}\,dx.$$
El integral beta vale
$$\int_0^1 x^3(1-x)^{r-4}\,dx=\frac{3!(r-4)!}{r!},$$
y por tanto
$$w(r)=4\cdot\frac{r!}{4!(r-4)!}\cdot\frac{3!(r-4)!}{r!}=1.$$
Cada número con al menos cuatro primos válidos cuenta exactamente una vez.
El número \(2310=2\cdot3\cdot5\cdot7\cdot11\) tiene \(r=5\) primos válidos.
Aparece \(\binom54=5\) veces en \(T_4\) y \(\binom55=1\) vez en \(T_5\). No contribuye a \(T_k\) con \(k\ge6\). Su peso final es
$$1\cdot5+(-4)\cdot1=1.$$
Eso muestra cómo los coeficientes corrigen el sobreconteo bruto de \(T_4\).
El programa no recorre ingenuamente los \(2^{25}\) subconjuntos. En su lugar usa una búsqueda en profundidad sobre índices de primos crecientes.
En cada estado DFS, el código conoce
1. el siguiente índice disponible,
2. cuántos primos se han tomado,
3. el producto actual \(q_S\).
Si taken >= 4, añade de inmediato
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor$$
a subset_sums[taken].
Antes de multiplicar por un nuevo primo \(p\), verifica
$$q_S\le\frac{N-1}{p},$$
implementado como product > max_value / p. Eso evita overflow y corta cualquier rama cuyo producto ya sea demasiado grande.
El primer checkpoint es
$$\texttt{solve}(1000)=23.$$
En efecto, por debajo de \(1000\) hay exactamente \(23\) números con al menos cuatro primos distintos menores que \(100\); por ejemplo \(210,330,390,462,\dots,990\).
El segundo checkpoint compara
$$\texttt{solve}(100000)$$
con una rutina brute-force que cuenta explícitamente los factores primos distintos de cada \(x<100000\). Así se valida que la fórmula ponderada y la implementación DFS coinciden con la definición directa.
primes_below_100() genera los 25 primos menores que \(100\) mediante una criba lineal pequeña.
solve(limit) fija max_value = limit - 1, ejecuta el DFS, acumula todos los totales de tipo \(T_k\) en subset_sums y luego forma la suma ponderada con coeficientes \((-1)^{k-4}\binom{k-1}{3}\).
binom solo se usa para esos coeficientes pequeños, así que los enteros con signo de 64 bits son suficientes.
La interfaz acepta --limit=<N> y --skip-checkpoints.
El espacio teórico de subconjuntos es \(2^{25}\), pero en la práctica la mayoría de ramas se eliminan muy pronto por el límite sobre el producto. Por eso el coste real depende del número de productos válidos \(q_S\le N-1\), no del conjunto potencia completo.
La suma final sobre \(k=4,\dots,25\) es despreciable frente al DFS, y la memoria adicional es \(O(25)\) aparte de la pila recursiva, porque solo se guardan la lista de primos y un arreglo acumulador de longitud \(26\).
我们要统计半开区间 \(1\le x<N\) 中满足条件的整数个数，其中 \(N=10^{16}\)。条件是：\(x\) 至少能被四个 不同的、且小于 \(100\) 的素数整除。素因子的重数不重要，例如 \(2\mid x\) 和 \(2^{20}\mid x\) 都只算一个不同素因子。
C++ 程序通过 --limit=<N> 传入上界，因此同样的方法适用于任意 \(N\ge1\)。
设
$$\mathcal P=\{p \text{ 为素数}: p<100\}.$$
这个集合共有 \(25\) 个素数。对任意子集 \(S\subseteq\mathcal P\)，定义其无平方部分乘积
$$q_S=\prod_{p\in S}p.$$
一个整数 \(x\) 同时被 \(S\) 中所有素数整除，当且仅当 \(q_S\mid x\)。因此满足 \(x<N\) 且被 \(S\) 全部素数整除的整数个数为
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
这里出现 \(N-1\)，是因为程序统计的是 \(1\le x<N\)，不是 \(1\le x\le N\)。
对每个 \(k\ge4\)，定义
$$T_k=\sum_{|S|=k}\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
如果某个数 \(x\) 在 \(\mathcal P\) 中恰好有 \(r\) 个不同素因子，那么它会在 \(T_k\) 中出现恰好 \(\binom{r}{k}\) 次。
因此我们需要找到系数 \(a_k\)，使得
$$\sum_{k=4}^{r}a_k\binom{r}{k}= \begin{cases} 0,&r<4,\\ 1,&r\ge4. \end{cases}$$
程序使用的系数是
$$a_k=(-1)^{k-4}\binom{k-1}{3}.$$
于是最终答案写成
$$\sum_{k=4}^{25}a_kT_k =\sum_{k=4}^{25}(-1)^{k-4}\binom{k-1}{3}T_k.$$
前几个权重为
$$a_4=1,\qquad a_5=-4,\qquad a_6=10,\qquad a_7=-20,\dots$$
所以表达式开头就是
$$T_4-4T_5+10T_6-20T_7+\cdots.$$
若某个整数恰好有 \(r\) 个符合条件的素因子，它的总权重为
$$w(r)=\sum_{k=4}^{r}(-1)^{k-4}\binom{k-1}{3}\binom{r}{k}.$$
当 \(r<4\) 时，这个和为空，因此 \(w(r)=0\)。
对 \(r\ge4\)，先用恒等式
$$\binom{r}{k}\binom{k-1}{3}=\binom{r}{4}\binom{r-4}{k-4}\frac{4}{k}.$$
令 \(j=k-4\)，得到
$$w(r)=4\binom{r}{4}\sum_{j=0}^{r-4}(-1)^j\binom{r-4}{j}\frac{1}{j+4}.$$
再用 \(\frac{1}{j+4}=\int_0^1 x^{j+3}\,dx\)，可写成
$$w(r)=4\binom{r}{4}\int_0^1 x^3(1-x)^{r-4}\,dx.$$
这个 beta 积分等于
$$\int_0^1 x^3(1-x)^{r-4}\,dx=\frac{3!(r-4)!}{r!},$$
因此
$$w(r)=4\cdot\frac{r!}{4!(r-4)!}\cdot\frac{3!(r-4)!}{r!}=1.$$
所以只要一个数至少拥有四个符合条件的不同素因子，它最终就恰好贡献 \(1\)；否则贡献 \(0\)。
数 \(2310=2\cdot3\cdot5\cdot7\cdot11\) 有 \(r=5\) 个符合条件的素因子。
它会在 \(T_4\) 中出现 \(\binom54=5\) 次，在 \(T_5\) 中出现 \(\binom55=1\) 次，在更高的 \(T_k\) 中不再出现。因此它的最终权重是
$$1\cdot5+(-4)\cdot1=1.$$
这正好展示了权重如何修正 \(T_4\) 带来的严重重复计数。
程序不会显式枚举全部 \(2^{25}\) 个子集，而是按素数下标递增的顺序做深度优先搜索。
一个 DFS 状态里保存着：
1. 下一个可选素数的位置，
2. 已经选了多少个素数，
3. 当前的无平方乘积 \(q_S\)。
一旦 taken >= 4，程序就立刻把
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor$$
加到 subset_sums[taken] 中。
在尝试乘上新素数 \(p\) 之前，程序会检查
$$q_S\le\frac{N-1}{p},$$
代码里写成 product > max_value / p。这既避免了整数溢出，也把乘积已经超界的分支立刻剪掉。
第一个 checkpoint 是
$$\texttt{solve}(1000)=23.$$
事实上，小于 \(1000\) 且至少有四个不同合格素因子的数恰好有 \(23\) 个，例如 \(210,330,390,462,\dots,990\)。
第二个 checkpoint 把
$$\texttt{solve}(100000)$$
与一个 brute-force 函数比较。后者对每个 \(x<100000\) 显式统计不同素因子个数，从而验证加权子集公式与 DFS 实现确实和原始定义一致。
primes_below_100() 通过一个小型线性筛生成 100 以下的 25 个素数。
solve(limit) 先设定 max_value = limit - 1，运行 DFS，把所有 \(T_k\) 型总和存入 subset_sums，最后按 \((-1)^{k-4}\binom{k-1}{3}\) 的权重求和。
binom 只用于这些很小的组合数权重，因此 64 位有符号整数足够。
命令行接口支持 --limit=<N> 和 --skip-checkpoints。
理论上的子集空间大小是 \(2^{25}\)，但实际运行时绝大多数分支会因为乘积过大而很早被剪掉。因此真正决定运行时间的不是完整幂集，而是满足 \(q_S\le N-1\) 的有效子集乘积个数。
最后对 \(k=4,\dots,25\) 的加权求和与 DFS 相比几乎可以忽略。除递归栈外，额外内存仅为 \(O(25)\)，因为程序只保存素数表和一个长度为 \(26\) 的累加数组。
Нужно посчитать числа в полуинтервале \(1\le x<N\), где \(N=10^{16}\), которые делятся как минимум на четыре различных простых числа меньше \(100\). Кратности не важны: простое \(2\) считается один раз, даже если входит в разложение с большой степенью.
Программа на C++ принимает верхнюю границу через --limit=<N>.
Пусть
$$\mathcal P=\{p \text{ простое}: p<100\}.$$
В этом множестве \(25\) простых. Для подмножества \(S\subseteq\mathcal P\) определим квадратсвободное произведение
$$q_S=\prod_{p\in S}p.$$
Число \(x\) делится на все простые из \(S\) тогда и только тогда, когда \(q_S\mid x\). Значит, количество таких \(x<N\) равно
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
Для \(k\ge4\) зададим
$$T_k=\sum_{|S|=k}\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
Если число \(x\) имеет ровно \(r\) различных простых делителей из \(\mathcal P\), то в \(T_k\) оно появляется ровно \(\binom{r}{k}\) раз.
Нужны коэффициенты \(a_k\), для которых
$$\sum_{k=4}^{r}a_k\binom{r}{k}= \begin{cases} 0,&r<4,\\ 1,&r\ge4. \end{cases}$$
Код использует
$$a_k=(-1)^{k-4}\binom{k-1}{3}.$$
Итоговая формула имеет вид
$$\sum_{k=4}^{25}a_kT_k =\sum_{k=4}^{25}(-1)^{k-4}\binom{k-1}{3}T_k.$$
Первые коэффициенты:
$$a_4=1,\qquad a_5=-4,\qquad a_6=10,\qquad a_7=-20,\dots$$
То есть начало суммы выглядит так:
$$T_4-4T_5+10T_6-20T_7+\cdots.$$
Для числа с ровно \(r\) допустимыми простыми делителями общий вес равен
$$w(r)=\sum_{k=4}^{r}(-1)^{k-4}\binom{k-1}{3}\binom{r}{k}.$$
Если \(r<4\), сумма пуста и \(w(r)=0\).
Если \(r\ge4\), используем тождество
$$\binom{r}{k}\binom{k-1}{3}=\binom{r}{4}\binom{r-4}{k-4}\frac{4}{k}.$$
Положив \(j=k-4\), получаем
$$w(r)=4\binom{r}{4}\sum_{j=0}^{r-4}(-1)^j\binom{r-4}{j}\frac{1}{j+4}.$$
Так как \(\frac{1}{j+4}=\int_0^1 x^{j+3}\,dx\), имеем
$$w(r)=4\binom{r}{4}\int_0^1 x^3(1-x)^{r-4}\,dx.$$
Это бета-интеграл:
$$\int_0^1 x^3(1-x)^{r-4}\,dx=\frac{3!(r-4)!}{r!}.$$
Следовательно,
$$w(r)=4\cdot\frac{r!}{4!(r-4)!}\cdot\frac{3!(r-4)!}{r!}=1.$$
Значит, каждое число с хотя бы четырьмя допустимыми простыми делителями учитывается ровно один раз.
Число \(2310=2\cdot3\cdot5\cdot7\cdot11\) имеет \(r=5\) допустимых простых делителей.
Оно входит в \(T_4\) ровно \(\binom54=5\) раз и в \(T_5\) ровно \(\binom55=1\) раз. В более высоких \(T_k\) его уже нет. Поэтому итоговый вес равен
$$1\cdot5+(-4)\cdot1=1.$$
Так коэффициенты исправляют сильный переучёт в \(T_4\).
Программа не перебирает все \(2^{25}\) подмножеств напрямую. Вместо этого используется поиск в глубину по возрастающим индексам простых.
Состояние DFS хранит:
1. индекс следующего допустимого простого,
2. сколько простых уже выбрано,
3. текущее квадратсвободное произведение \(q_S\).
Как только taken >= 4, в subset_sums[taken] сразу добавляется
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
Перед умножением на новый простой \(p\) код проверяет
$$q_S\le\frac{N-1}{p},$$
что реализовано как product > max_value / p. Это одновременно защищает от переполнения и отсеивает бесполезные ветви.
Первый checkpoint:
$$\texttt{solve}(1000)=23.$$
Действительно, ниже \(1000\) существует ровно \(23\) числа хотя бы с четырьмя различными простыми меньше \(100\); примеры: \(210,330,390,462,\dots,990\).
Второй checkpoint сравнивает
$$\texttt{solve}(100000)$$
с brute-force процедурой, которая для каждого \(x<100000\) напрямую подсчитывает количество различных простых делителей. Это подтверждает совпадение формулы с прямым определением.
primes_below_100() генерирует 25 простых меньше \(100\) с помощью небольшого линейного решета.
solve(limit) задаёт max_value = limit - 1, запускает DFS, накапливает все суммы типа \(T_k\) в массиве subset_sums, а затем берёт взвешенную сумму с коэффициентами \((-1)^{k-4}\binom{k-1}{3}\).
binom нужен только для этих маленьких коэффициентов, поэтому 64-битных целых со знаком достаточно.
Поддерживаются параметры --limit=<N> и --skip-checkpoints.
Теоретический размер пространства подмножеств равен \(2^{25}\), но на практике большая часть ветвей быстро отсекается ограничением на произведение. Поэтому реальная стоимость определяется числом допустимых произведений \(q_S\le N-1\), а не полной мощностью множества.
Финальная сумма по \(k=4,\dots,25\) ничтожна по сравнению с DFS, а дополнительная память равна \(O(25)\) помимо рекурсивного стека, поскольку хранятся только список простых и массив-накопитель длины \(26\).
المطلوب هو عدّ الأعداد في المجال شبه المفتوح \(1\le x<N\) حيث \(N=10^{16}\)، والتي تقبل القسمة على أربعة أعداد أولية مختلفة على الأقل من بين الأوليات الأصغر من \(100\). تكرار العامل الأولي لا يهم؛ فوجود \(2\) مرة أو مرات كثيرة يبقى عاملًا أوليًا مميزًا واحدًا فقط.
برنامج C++ يمرر الحد الأعلى عبر --limit=<N>.
لنضع
$$\mathcal P=\{p \text{ أولي}: p<100\}.$$
وفيها \(25\) عددًا أوليًا. لكل مجموعة جزئية \(S\subseteq\mathcal P\) نعرّف حاصل الضرب الخالي من التكرار
$$q_S=\prod_{p\in S}p.$$
يكون العدد \(x\) قابلًا للقسمة على جميع أوليات \(S\) إذا وفقط إذا كان \(q_S\mid x\). لذلك فإن عدد الأعداد \(x<N\) التي تحقق ذلك هو
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
ويظهر \(N-1\) لأن المجال الذي يعدّه الكود هو \(1\le x<N\).
لكل \(k\ge4\) نعرّف
$$T_k=\sum_{|S|=k}\left\lfloor\frac{N-1}{q_S}\right\rfloor.$$
إذا كان لعدد ما \(r\) عوامل أولية مميزة من \(\mathcal P\)، فإنه يظهر في \(T_k\) بعدد مرات يساوي \(\binom{r}{k}\).
لذلك نحتاج إلى معاملات \(a_k\) تحقق
$$\sum_{k=4}^{r}a_k\binom{r}{k}= \begin{cases} 0,&r<4,\\ 1,&r\ge4. \end{cases}$$
يستخدم الكود
$$a_k=(-1)^{k-4}\binom{k-1}{3}.$$
ومن ثم يكون الجواب النهائي
$$\sum_{k=4}^{25}a_kT_k =\sum_{k=4}^{25}(-1)^{k-4}\binom{k-1}{3}T_k.$$
وأول هذه المعاملات هو
$$a_4=1,\qquad a_5=-4,\qquad a_6=10,\qquad a_7=-20,\dots$$
أي أن الصيغة تبدأ بـ
$$T_4-4T_5+10T_6-20T_7+\cdots.$$
إذا كان لعدد ما بالضبط \(r\) عوامل أولية صالحة، فإن وزنه الكلي يساوي
$$w(r)=\sum_{k=4}^{r}(-1)^{k-4}\binom{k-1}{3}\binom{r}{k}.$$
عندما \(r<4\) يكون المجموع فارغًا، وبالتالي \(w(r)=0\).
أما عندما \(r\ge4\)، فنستعمل الهوية
$$\binom{r}{k}\binom{k-1}{3}=\binom{r}{4}\binom{r-4}{k-4}\frac{4}{k}.$$
وبوضع \(j=k-4\) نحصل على
$$w(r)=4\binom{r}{4}\sum_{j=0}^{r-4}(-1)^j\binom{r-4}{j}\frac{1}{j+4}.$$
ثم نستخدم \(\frac{1}{j+4}=\int_0^1 x^{j+3}\,dx\)، فنصل إلى
$$w(r)=4\binom{r}{4}\int_0^1 x^3(1-x)^{r-4}\,dx.$$
وهذا تكامل بيتا يساوي
$$\int_0^1 x^3(1-x)^{r-4}\,dx=\frac{3!(r-4)!}{r!}.$$
ومن ثم
$$w(r)=4\cdot\frac{r!}{4!(r-4)!}\cdot\frac{3!(r-4)!}{r!}=1.$$
إذن كل عدد يملك أربعة عوامل أولية صالحة أو أكثر يُحسب مرة واحدة تمامًا، وما دون ذلك لا يُحسب مطلقًا.
العدد \(2310=2\cdot3\cdot5\cdot7\cdot11\) له \(r=5\) عوامل أولية صالحة.
لذلك يظهر في \(T_4\) بعدد \(\binom54=5\) مرات، وفي \(T_5\) مرة واحدة، ولا يظهر في أي \(T_k\) أكبر. ومن ثم يكون وزنه النهائي
$$1\cdot5+(-4)\cdot1=1.$$
وهكذا نرى بوضوح كيف تصحح الأوزان العد الزائد الناتج عن \(T_4\).
لا يمر البرنامج على جميع المجموعات الجزئية \(2^{25}\) بشكل ساذج. بدلًا من ذلك، يستخدم بحثًا بالعمق على فهارس الأوليات بترتيب متزايد.
تحتوي حالة DFS على:
1. فهرس العدد الأولي التالي المسموح به،
2. عدد الأوليات المختارة حتى الآن،
3. حاصل الضرب الحالي \(q_S\).
متى ما تحقق taken >= 4، يضيف الكود مباشرة
$$\left\lfloor\frac{N-1}{q_S}\right\rfloor$$
إلى subset_sums[taken].
وقبل الضرب في أولي جديد \(p\)، يتحقق من الشرط
$$q_S\le\frac{N-1}{p},$$
والمكتوب في الكود على صورة product > max_value / p. وهذا يمنع تجاوز السعة ويقطع الفروع التي تجاوز حاصل ضربها الحد مبكرًا.
نقطة التحقق الأولى هي
$$\texttt{solve}(1000)=23.$$
وفعلًا يوجد تحت \(1000\) عددها \(23\) فقط من الأعداد التي تملك أربعة عوامل أولية مميزة صالحة على الأقل؛ من أمثلتها \(210,330,390,462,\dots,990\).
أما نقطة التحقق الثانية فتقارن
$$\texttt{solve}(100000)$$
مع دالة brute-force تعدّ مباشرة عدد العوامل الأولية المميزة لكل \(x<100000\). وهذا يؤكد أن صيغة الأوزان وتنفيذ DFS متطابقان مع التعريف الأصلي للمسألة.
primes_below_100() يولد الأعداد الأولية الخمسة والعشرين الأصغر من \(100\) باستخدام منخل خطي صغير.
solve(limit) يحدد max_value = limit - 1، ثم يشغل DFS، ويجمع جميع القيم من نوع \(T_k\) في subset_sums، ثم يأخذ المجموع الموزون بمعاملات \((-1)^{k-4}\binom{k-1}{3}\).
تُستخدم الدالة binom فقط لهذه المعاملات الصغيرة، لذا فإن الأعداد الصحيحة ذات 64 بت تكفي تمامًا.
تقبل الواجهة الوسيطين --limit=<N> و--skip-checkpoints.
الحجم النظري لفضاء المجموعات الجزئية هو \(2^{25}\)، لكن عمليًا تُحذف أغلب الفروع مبكرًا بسبب حد حاصل الضرب. لذلك تعتمد الكلفة الفعلية على عدد الحواصل \(q_S\le N-1\) الصالحة، لا على كامل مجموعة القوى.
أما الجمع النهائي على \(k=4,\dots,25\) فتكلفته صغيرة جدًا مقارنة بـ DFS، واستهلاك الذاكرة الإضافي هو \(O(25)\) عدا مكدس الاستدعاء، لأن البرنامج يخزن فقط قائمة الأوليات ومصفوفة تراكم بطول \(26\).