Project Euler Problem 515
Let d(p, n, 0) be the multiplicative inverse of n modulo prime p, defined as n times d(p, n, 0) = 1 bmod p.
Solution
Answer: 2422639000800
Let
$$d(p,n,0)\equiv n^{-1}\pmod p$$
for prime $p$, and recursively
$$d(p,n,k)=\sum_{i=1}^{n} d(p,i,k-1).$$
We want
$$D(a,b,k)=\sum_{a\le p<a+b\atop p\text{ prime}} \left(d(p,p-1,k)\bmod p\right)$$
for $a=10^9,;b=10^5,;k=10^5$.
1. Mathematical analysis
Define
$$S_k(n)=d(p,n,k).$$
By repeatedly expanding the recurrence,
$$S_k(n) = \sum_{i=1}^{n} \binom{n-i+k-1}{k-1}, i^{-1} \pmod p.$$
We need $n=p-1$, so:
$$S_k(p-1) = \sum_{i=1}^{p-1} \binom{p+k-i-2}{k-1} i^{-1} \pmod p.$$
Let $m=k-1$. Since in this problem $p>k$, we have $m<p$.
Then
$$\binom{p+m-i-1}{m} \equiv \binom{m-i-1}{m} \pmod p$$
because the $p$-part disappears modulo $p$.
Using the negative binomial identity
$$\binom{m-i-1}{m} = (-1)^m \binom{i}{m},$$
we get
$$S_k(p-1) \equiv (-1)^m \sum_{i=1}^{p-1} \binom{i}{m}\frac1i.$$
Now use
$$\frac{1}{i}\binom{i}{m} = \frac{1}{m}\binom{i-1}{m-1}.$$
Hence
$$S_k(p-1) \equiv \frac{(-1)^m}{m} \sum_{i=1}^{p-1} \binom{i-1}{m-1}.$$
Apply the hockey-stick identity:
$$\sum_{i=m}^{p-1} \binom{i-1}{m-1} = \binom{p-1}{m}.$$
So
$$S_k(p-1) \equiv \frac{(-1)^m}{m} \binom{p-1}{m}.$$
For prime $p$,
$$\binom{p-1}{m} \equiv (-1)^m \pmod p.$$
Therefore
$$S_k(p-1) \equiv \frac{1}{m} = (k-1)^{-1} \pmod p.$$
This is the key simplification:
$$d(p,p-1,k)\bmod p = (k-1)^{-1}\pmod p.$$
So the whole problem becomes
$$D(a,b,k) = \sum_{a\le p<a+b} (k-1)^{-1}\pmod p.$$
For $k=10^5$,
$$k-1=99999.$$
Thus we sum the modular inverse of $99999$ modulo each prime $p\in[10^9,10^9+10^5)$.
2. Python implementation
from sympy import primerange
def D(a, b, k):
m = k - 1
total = 0
for p in primerange(a, a + b):
total += pow(m, -1, p)
return total
print(D(101, 1, 10)) # 45
print(D(10**3, 10**2, 10**2)) # 8334
print(D(10**6, 10**3, 10**3)) # 38162302
print(D(10**9, 10**5, 10**5))
3. Code walkthrough
m = k - 1
From the proof,
$$d(p,p-1,k)\equiv (k-1)^{-1}\pmod p.$$
primerange(a, a+b)
Iterates over all primes $p$ in the interval.
pow(m, -1, p)
Python computes the modular inverse:
$$m^{-1}\pmod p.$$
- Add all inverses into
total.
Verification against examples
The code reproduces:
- $D(101,1,10)=45$
- $D(10^3,10^2,10^2)=8334$
- $D(10^6,10^3,10^3)=38162302$
matching the problem statement exactly.
4. Final verification
There are about $\frac{10^5}{\log(10^9)}\approx 4800$ primes in the interval, and each inverse is roughly size $p/2\sim 5\times 10^8$, so a total around $2\times 10^{12}$ is plausible.
The computed value is:
$$2422639000800$$
which is consistent in magnitude.
Answer: 2422639000800