Detailed mathematical approach

Problem Summary

The task is to evaluate

$B(x,y,n)=\sum_{\substack{p\ \text{prime}\\ x\le p\le x+y}} \left(a_n \bmod p\right),$

for the specific target \(B(10^9,10^7,10^{15})\), where the sequence is defined by

$a_1=1,\qquad a_{t+1}=6a_t^2+10a_t+3.$

The recurrence is explosive: every step squares the previous value, so computing \(a_n\) directly is hopeless when \(n=10^{15}\). The solution therefore avoids building the integer \(a_n\) itself and instead works modulo each prime in the interval.

Mathematical Approach

The central idea is to turn the quadratic recurrence into a Lucas-sequence problem whose enormous index can be reduced modulo a small group order.

Step 1: Affine Change of Variables

Introduce

$b_t=6a_t+5.$

Then the recurrence becomes

$b_{t+1}=6a_{t+1}+5=6(6a_t^2+10a_t+3)+5=(6a_t+5)^2-2=b_t^2-2,$

with initial value

$b_1=6\cdot 1+5=11.$

This is the key simplification. Instead of a quadratic polynomial with three terms, we now have repeated application of the much cleaner map \(u\mapsto u^2-2\).

Step 2: Recognize the Lucas Structure

Let

$\alpha,\beta=\frac{11\pm\sqrt{117}}{2},\qquad \alpha+\beta=11,\qquad \alpha\beta=1.$

Define the Lucas \(V\)-sequence with parameters \(P=11\) and \(Q=1\) by

$V_m=\alpha^m+\beta^m.$

Then

$V_0=2,\qquad V_1=11,\qquad V_{m+1}=11V_m-V_{m-1}.$

Because \(\alpha\beta=1\), we also have the doubling identity

$V_{2m}=\alpha^{2m}+\beta^{2m}=(\alpha^m+\beta^m)^2-2(\alpha\beta)^m=V_m^2-2.$

Since \(b_1=V_1\) and both sequences obey the same doubling rule, induction gives

$b_t=V_{2^{t-1}}.$

Therefore the original problem is equivalent to evaluating

$a_n \bmod p \quad \text{from} \quad V_{2^{n-1}} \bmod p.$

Step 3: Reduce the Giant Index

For the primes relevant to the target sum, \(p>13\), so the discriminant

$\Delta=11^2-4=117$

is nonzero modulo \(p\), and \(6\) is invertible modulo \(p\).

If \(117\) is a quadratic residue modulo \(p\), then \(\sqrt{117}\in\mathbb{F}_p\), so \(\alpha,\beta\in\mathbb{F}_p^\times\). Since \(\beta=\alpha^{-1}\), the value

$V_k=\alpha^k+\alpha^{-k}$

depends only on \(k\bmod(p-1)\).

If \(117\) is a quadratic nonresidue modulo \(p\), then \(\alpha\) and \(\beta\) live in \(\mathbb{F}_{p^2}\), and the Frobenius map swaps them:

$\alpha^p=\beta=\alpha^{-1}.$

Hence

$\alpha^{p+1}=1,$

so \(V_k\) depends only on \(k\bmod(p+1)\).

Thus the enormous index

$k=2^{n-1}$

is reduced by

$k\equiv 2^{n-1}\pmod{p-1}\quad\text{if }\left(\frac{117}{p}\right)=1,$

$k\equiv 2^{n-1}\pmod{p+1}\quad\text{if }\left(\frac{117}{p}\right)=-1.$

This turns an index of size roughly \(2^{10^{15}}\) into one smaller than \(p+1\).

Step 4: Fast Doubling for the Lucas Term

After the index has been reduced, the Lucas value is computed with standard doubling identities:

$V_{2m}=V_m^2-2,$

$V_{2m+1}=V_mV_{m+1}-11,$

$V_{2m+2}=V_{m+1}^2-2.$

Processing the bits of \(k\) from top to bottom keeps a pair of consecutive Lucas values and reaches \(V_k \bmod p\) in \(O(\log k)\) arithmetic steps.

Once \(b_n \equiv V_k \pmod p\) is known, we recover the original sequence via

$a_n\equiv (b_n-5)\cdot 6^{-1}\pmod p.$

Step 5: Sum Over the Prime Interval

The interval \([x,x+y]\) is handled with a segmented sieve. First, all base primes up to \(\sqrt{x+y}\) are generated. Then their multiples are marked inside the target interval, leaving exactly the primes that contribute to the sum.

For each prime \(p\) in the segment, the algorithm performs three modular tasks:

$\text{determine whether }117\text{ is a residue mod }p,$

$\text{compute }2^{n-1}\text{ modulo }p-1\text{ or }p+1,$

$\text{evaluate }V_k\text{ and convert it back to }a_n\bmod p.$

Adding those residues yields the required value of \(B(x,y,n)\).

Worked Example

A small example shows the reduction clearly. Take \(n=4\) and the interval \([5,7]\), so the contributing primes are \(5\) and \(7\).

For \(p=5\), \(117\equiv 2\pmod 5\), which is a nonresidue, so the index is reduced modulo \(p+1=6\):

$k\equiv 2^{4-1}=8\equiv 2\pmod 6.$

Then

$V_2=11^2-2=119\equiv 4\pmod 5,$

so

$a_4\equiv (4-5)\cdot 6^{-1}\equiv (-1)\cdot 1\equiv 4\pmod 5.$

For \(p=7\), \(117\equiv 5\pmod 7\), also a nonresidue, so the index is reduced modulo \(p+1=8\):

$k\equiv 8\equiv 0\pmod 8.$

Therefore

$V_0=2,\qquad a_4\equiv (2-5)\cdot 6^{-1}\equiv (-3)\cdot 6\equiv 3\pmod 7.$

Hence

$B(5,2,4)=4+3=7.$

How the Code Works

The C++, Python, and Java implementations follow the same pipeline. They first build a simple sieve up to \(\sqrt{x+y}\), then use those base primes to mark a segmented interval \([x,x+y]\). Every unmarked entry is a prime that contributes one term to the final sum.

For each such prime, the implementation computes a Legendre-style residue test for \(117\), chooses \(p-1\) or \(p+1\) as the period modulus, and evaluates \(2^{n-1}\) modulo that period with fast modular exponentiation. It then performs Lucas fast doubling on a pair of consecutive values until it reaches the reduced index \(k\), converts the result back from \(b_n\) to \(a_n\), and adds the residue to the running total.

No large integer value of \(a_n\) is ever constructed. Only modular values are propagated, which is why the method stays fast even when \(n=10^{15}\).

Complexity Analysis

Let \(H=x+y\), and let \(\pi(x,y)\) denote the number of primes in \([x,x+y]\). Building the base sieve up to \(\sqrt{H}\) costs \(O(\sqrt{H}\log\log H)\) time and \(O(\sqrt{H})\) memory. Marking the segment costs \(O(y\log\log H)\) time and \(O(y)\) memory.

For each prime in the segment, the modular residue test, the reduction of \(2^{n-1}\), and the Lucas doubling stage together require \(O(\log n+\log p)\) arithmetic operations. The overall running time is therefore

$O\!\left(\sqrt{H}\log\log H+y\log\log H+\pi(x,y)(\log n+\log H)\right),$

with memory usage

$O(\sqrt{H}+y).$

Footnotes and References

1. Project Euler Problem 492
2. Wikipedia - Lucas sequence
3. Wikipedia - Finite field
4. cp-algorithms - Fibonacci numbers and fast doubling
5. cp-algorithms - Sieve of Eratosthenes