Project Euler Problem 234

For an integer n ge 4, we define the lower prime square root of n, denoted by operatorname{lps}(n), as the largest prime

Project Euler Problem 234

Solution

Answer: 1259187438574927161

Let $p<q$ be consecutive primes.

For every integer $n$ with

$$p^2 < n < q^2,$$

we have

$$\operatorname{lps}(n)=p,\qquad \operatorname{ups}(n)=q,$$

because there are no primes strictly between $p$ and $q$.

So inside the interval $(p^2,q^2)$, a number is semidivisible exactly when:

  • it is divisible by $p$ xor divisible by $q$,
  • i.e. divisible by exactly one of them.

That immediately turns the problem into a sum over consecutive prime pairs.


1. Mathematical analysis

For consecutive primes $p<q$, define

$$I_{p,q}={n : p^2<n<q^2}.$$

Inside this interval:

  • semidivisible numbers divisible by $p$:

$$p\cdot k,\qquad k=p+1,\dots,\left\lfloor\frac{q^2-1}{p}\right\rfloor$$

  • semidivisible numbers divisible by $q$:

$$q\cdot k,\qquad k=\left\lfloor\frac{p^2}{q}\right\rfloor+1,\dots,q-1$$

But numbers divisible by both $p$ and $q$ are multiples of $pq$, and these are not semidivisible, since both prime square roots divide them.

Therefore for each prime pair we use inclusion–exclusion:

$$S(p,q) = \sum_{\substack{p^2<n<q^2\ p\mid n}} n + \sum_{\substack{p^2<n<q^2\ q\mid n}} n - 2\sum_{\substack{p^2<n<q^2\ pq\mid n}} n.$$

The factor $2$ removes the doubly-counted multiples of $pq$ entirely.


Summing multiples efficiently

The sum of multiples of $m$ in an interval $[A,B]$ is:

$$m\sum_{k=\lceil A/m\rceil}^{\lfloor B/m\rfloor} k.$$

Using the arithmetic-series formula,

$$\sum_{k=L}^{R} k = \frac{(L+R)(R-L+1)}{2}.$$

So:

$$\operatorname{sumMultiples}(m,A,B) = m\cdot \frac{(L+R)(R-L+1)}{2},$$

where

$$L=\left\lceil\frac{A}{m}\right\rceil,\qquad R=\left\lfloor\frac{B}{m}\right\rfloor.$$


2. Python implementation

from math import isqrt

LIMIT = 999966663333

# ------------------------------------------------------------
# Generate all primes up to sqrt(LIMIT)
# ------------------------------------------------------------

max_prime = isqrt(LIMIT) + 1000

sieve = bytearray(b'\x01') * (max_prime + 1)
sieve[0] = sieve[1] = 0

for i in range(2, isqrt(max_prime) + 1):
    if sieve[i]:
        sieve[i * i : max_prime + 1 : i] = (
            b'\x00' * (((max_prime - i * i) // i) + 1)
        )

primes = [i for i, v in enumerate(sieve) if v]

# ------------------------------------------------------------
# Sum of multiples of m inside [A, B]
# ------------------------------------------------------------

def sum_multiples(m, A, B):
    L = (A + m - 1) // m
    R = B // m

    if L > R:
        return 0

    count = R - L + 1
    return m * (L + R) * count // 2

# ------------------------------------------------------------
# Main computation
# ------------------------------------------------------------

total = 0

for p, q in zip(primes, primes[1:]):

    # numbers with lps = p and ups = q
    A = p * p + 1
    B = min(q * q - 1, LIMIT)

    if A > B:
        break

    total += (
        sum_multiples(p, A, B)
        + sum_multiples(q, A, B)
        - 2 * sum_multiples(p * q, A, B)
    )

print(total)

3. Code walkthrough

Prime generation

max_prime = isqrt(LIMIT) + 1000

We only need primes up to slightly beyond $\sqrt{\text{LIMIT}}$, because the intervals are determined by consecutive primes around square roots.

The sieve computes all such primes efficiently.


Summing multiples

def sum_multiples(m, A, B):

This computes

$$\sum_{\substack{A \le n \le B \ m \mid n}} n$$

without iterating through every multiple individually.


Main loop

for p, q in zip(primes, primes[1:]):

We iterate through consecutive prime pairs.

For each pair:

A = p * p + 1
B = min(q * q - 1, LIMIT)

This is exactly the interval where

$$\operatorname{lps}(n)=p,\quad \operatorname{ups}(n)=q.$$

Then inclusion–exclusion gives the semidivisible contribution.


4. Verification on the examples

For $N=15$:

The semidivisible numbers are:

$$8,\ 10,\ 12$$

with sum

$$8+10+12=30.$$

The program reproduces this.

For $N=1000$:

The problem statement says:

  • there are $92$ semidivisible numbers,
  • their sum is $34825$.

The program also reproduces this exactly.


5. Final verification

The algorithm is mathematically exact:

  • every $n\ge4$ belongs to exactly one interval $(p^2,q^2)$,
  • semidivisibility reduces to divisibility by exactly one of $p,q$,
  • inclusion–exclusion removes the forbidden multiples of $pq$.

The computed value is:

$$1259187438574927161$$

which is of the expected magnitude ($\sim 10^{18}$) for a limit near $10^{12}$.

Answer: 1259187438574927161