Detailed mathematical approach

Problem Summary

Problem 984, Knights and Horses, ultimately asks for the value of a specific integer sequence at the enormous index \(n=10^{18}\), reduced modulo \(10^9+7\). The implementations make that sequence concrete through its first terms \(1,4,9,92,903,4411,14959,41083,98200,212418,425756,803074,1441065,2479669\).

The key point is that the computation is not done by any direct combinatorial search at size \(10^{18}\). Instead, the sequence has already been compressed into two exact algebraic descriptions: a parity-sensitive closed form for the tail of the sequence, and a linear-recurrence evaluator that can reproduce the remaining cases modulo the prime modulus.

Mathematical Approach

The three implementations expose the same structure: after the first few exceptional indices, the Knights and Horses sequence is a quasi-polynomial in \(n\) with period 2. The recurrence engine and the explicit even formula are two faces of the same underlying mathematics.

The Sequence Encoded by the Implementations

From the stored initial values and the signed recurrence used in the exact checkpoint, the tail of the sequence satisfies

$ f_n=7f_{n-1}-19f_{n-2}+21f_{n-3}+6f_{n-4}-42f_{n-5}+42f_{n-6}-6f_{n-7}-21f_{n-8}+19f_{n-9}-7f_{n-10}+f_{n-11}, \qquad n\ge 15. $

This is the actual algebraic law extracted from the implementations. The final Project Euler answer is simply the value of this sequence at a huge even index.

The Transition Polynomial and Its Factorization

The generic evaluator stores fourteen basis values, so its transition polynomial is

$ \chi(x)=x^{14}-7x^{13}+19x^{12}-21x^{11}-6x^{10}+42x^9-42x^8+6x^7+21x^6-19x^5+7x^4-x^3. $

Factoring it gives

$ \chi(x)=x^3(x-1)^9(x+1)^2. $

The factor \(x^3\) explains why the smallest indices are treated separately. After that transient prefix, the behavior is governed by the roots \(1\) and \(-1\): multiplicity 9 at \(1\), and multiplicity 2 at \(-1\).

Why the Tail Becomes a Parity-Dependent Polynomial

A root \(1\) of multiplicity 9 contributes a polynomial of degree at most 8, while a root \(-1\) of multiplicity 2 contributes \((-1)^n\) times a polynomial of degree at most 1. Therefore, for \(n\ge 4\), the sequence must have the form

$ f(n)=A_8(n)+(-1)^nB_1(n), $

with \(\deg A_8\le 8\) and \(\deg B_1\le 1\). Matching this shape against the values encoded in the implementations yields

$ f(n)=\frac{31}{40320}n^8+\frac{31}{3360}n^7+\frac{67}{1440}n^6+\frac{41}{320}n^5+\frac{313}{1440}n^4-\frac{5699}{240}n^3+\frac{16049}{420}n^2+\frac{941251}{4480}n-\frac{107261}{256}-(-1)^n\frac{2n+3}{256}, \qquad n\ge 4. $

This single formula explains why the code can switch between a polynomial shortcut and a generic recurrence solver without changing the underlying sequence.

Recovering the Even Closed Form

The required query is \(n=10^{18}\), which is even. Setting \((-1)^n=1\) collapses the quasi-polynomial to the degree-8 polynomial used directly by the implementations:

$ P(n)=\frac{31}{40320}n^8+\frac{31}{3360}n^7+\frac{67}{1440}n^6+\frac{41}{320}n^5+\frac{313}{1440}n^4-\frac{5699}{240}n^3+\frac{16049}{420}n^2+\frac{29413}{140}n-419, \qquad n>3,\ n\equiv 0\pmod 2. $

Because the modulus \(10^9+7\) is prime and larger than every denominator that appears here, each division is implemented as multiplication by a modular inverse.

Worked Example: Evaluating \(f(16)\)

The parity reduction is already enough to show the method on a concrete value. Since \(16\) is even and greater than 3, we use the closed polynomial:

$ f(16)=P(16)=6623481. $

The recurrence branch produces the same value from earlier terms, so the polynomial and the recurrence are not independent tricks; they are two exact descriptions of the same tail sequence.

How the Code Works

Closed-Form Evaluation for the Target Index

The C++, Python, and Java implementations first notice that the requested index \(10^{18}\) is even. They compute \(n,n^2,\dots,n^8\) modulo \(10^9+7\), replace each rational denominator by its modular inverse, and evaluate the polynomial \(P(n)\) term by term. That is enough to obtain the final answer for the problem instance.

Generic Linear-Recurrence Evaluator

The implementations also keep a general recurrence solver. It represents \(x^{n-1}\) modulo the transition polynomial \(\chi(x)\), using repeated squaring on coefficient vectors. After each coefficient-vector multiplication, higher powers are reduced with the recurrence coefficients, which is the standard quotient-ring or Kitamasa viewpoint for linear recurrences. Although fourteen basis slots are stored, the last three recurrence coefficients are zero, so the genuine tail relation is the 11-term recurrence written above.

Validation Checkpoints

One implementation performs explicit sanity checks before the final evaluation. It verifies small values such as \(f(3)=9\) and \(f(5)=903\), checks the exact integer value \(f(100)=8658918531876\) using signed 128-bit recurrence arithmetic, and confirms the modular checkpoint \(f(10000)\equiv 377956308\pmod{10^9+7}\). These checks protect against indexing mistakes, sign errors, and incorrect modular reductions.

Complexity Analysis

For the actual Project Euler query, the runtime is \(O(1)\): the target index is even, so only a fixed number of modular multiplications and inverses are needed to evaluate the degree-8 polynomial. Memory usage is also \(O(1)\).

The retained generic evaluator is more general. With a fixed state size \(K=14\), repeated squaring of reduced coefficient vectors costs \(O(K^2\log n)\) time and \(O(K)\) memory. Since only 11 recurrence coefficients are nonzero, the constant factor is small.

Footnotes and References

1. Project Euler problem page: https://projecteuler.net/problem=984
2. Linear recurrences with constant coefficients: Wikipedia - Linear recurrence with constant coefficients
3. Characteristic polynomial: Wikipedia - Characteristic polynomial
4. Finite differences and polynomial sequences: Wikipedia - Finite difference
5. Modular multiplicative inverse: Wikipedia - Modular multiplicative inverse
6. Exponentiation by squaring: Wikipedia - Exponentiation by squaring