Detailed mathematical approach

Problem Summary

The process starts with \(n\) stones at position \(0\). Positions are then processed in increasing order. If the current position \(i\) contains \(m\) stones, one singleton is left behind exactly when \(m\) is odd; in all cases \(\lfloor m/2\rfloor\) stones are sent to \(i+1\) and another \(\lfloor m/2\rfloor\) stones are sent to \(i+3\).

Only the positions where the current pile was odd matter for this problem. A singleton is lonely when no other singleton is within distance \(2\), so a sad integer is one whose singleton positions are pairwise separated by at least \(3\). The function \(S(k)\) asks for the sum of all sad starting values whose complete singleton set is contained in \(0\le i\lt k\). The target is \(S(80)\), and the problem gives the checks \(S(14)=159\) and \(S(30)=33438\).

The main challenge is that many starting values do not terminate. Even when the process keeps moving stones forever, it can still leave only finitely many singleton positions. So the solution cannot rely on simulating until all stones disappear; it must characterize the finite parity trace left by the process.

Mathematical Approach

Viewing the process as a carrying rule

At each position the process removes pairs of stones and turns each pair at \(i\) into one stone at \(i+1\) and one stone at \(i+3\). The possible leftover digit at \(i\) is therefore only \(0\) or \(1\). This is analogous to carrying in a positional numeral system, except the carry rule is unusual:

$$2\cdot \text{position }i \quad\longrightarrow\quad \text{position }(i+1)+\text{position }(i+3).$$

If \(s_i\in\{0,1\}\) records whether a singleton was left at \(i\), then the process deterministically normalizes the starting pile \(n\) into a digit sequence \(s_0,s_1,s_2,\dots\). The sadness condition is then a condition on the support of this digit sequence.

An algebraic invariant for the carry rule

Attach a formal weight \(x^i\) to position \(i\). For the carry rule to preserve weight, one pair at position \(i\) must have the same value as the two carried stones:

$$2x^i=x^{i+1}+x^{i+3}.$$

After division by \(x^i\), this becomes

$$x^3+x-2=0.$$

The polynomial factors as

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

The factor \(x-1\) only tracks total mass, which is not enough here because the process may continue forever while leaving finitely many singletons. The useful finite trace is captured by the quadratic factor, so the implementation works in the ring where

$$x^2+x+2=0,\qquad\text{equivalently}\qquad x^2=-x-2.$$

Reducing each position to two coefficients

Modulo \(x^2+x+2\), every power has the form

$$x^i=a_i+b_i x.$$

The first two values are \(x^0=1\) and \(x^1=x\). If \(x^{i-1}=a+b x\), then

$$x^i=x(a+b x)=a x+b x^2=a x+b(-x-2)=-2b+(a-b)x.$$

Thus the code stores only two integers per position and advances by

$$a_i=-2b_{i-1},\qquad b_i=a_{i-1}-b_{i-1}.$$

For a finite singleton set \(P\), define

$$V(P)=\sum_{i\in P}x^i=A(P)+B(P)x.$$

If \(P\) is the complete finite singleton set for some starting value \(n\), then its value in this quadratic quotient must equal the initial value \(n\), which has no \(x\)-part. Therefore the required algebraic test is

$$B(P)=0,\qquad n=A(P),\qquad A(P)\gt0.$$

Why this enumerates the right objects

The left-to-right process is deterministic: once \(n\) is fixed, the parity digit \(s_i\) and the outgoing carry from position \(i\) are fixed. Conversely, a finite digit set \(P\) whose quadratic value is a positive integer is a valid finite trace for that integer under this carrying rule. The program uses this as a bijection between sad starting values counted by \(S(k)\) and finite subsets \(P\subseteq\{0,\dots,k-1\}\) satisfying the algebraic test and the loneliness gap.

This is also why the final answer is a sum of values \(A(P)\), not a count of sets. Each accepted set \(P\) represents one starting value \(n=A(P)\), and \(S(k)\) asks for the sum of those starting values.

The loneliness condition as an independent-set condition

A finite set \(P\) of singleton positions is lonely exactly when

$$|i-j|\ge 3\qquad\text{for all distinct }i,j\in P.$$

Equivalently, in the path graph on positions \(0,1,\dots,k-1\), where vertices at distance \(1\) or \(2\) conflict, \(P\) is an independent set. The problem has now become a finite combinatorial enumeration: list gap-\(3\) subsets, reduce their algebraic value to \(A+Bx\), and retain only positive integer values with \(B=0\).

The examples fall out directly. \(\{0\}\) gives \(1\). The set \(\{2,5,8,13\}\) gives \(68\), and \(\{1,13\}\) gives \(90\). These are the sad contributions below \(k=14\), so \(S(14)=1+68+90=159\).

Why a direct subset search is too large

Let \(F(t)\) be the number of subsets of \(\{0,\dots,t-1\}\) with mutual distance at least \(3\). Looking at the last position gives

$$F(t)=F(t-1)+F(t-3).$$

The term \(F(t-1)\) covers subsets not using the last position; the term \(F(t-3)\) covers subsets using it, because the previous two positions are then forbidden. For \(t=80\), this recurrence already gives \(25058735850088\) candidate subsets. That is far beyond direct enumeration.

Splitting at \(m=\lfloor k/2\rfloor\) changes the scale. For \(k=80\), each half has length \(40\), and the number of valid half-subsets is only \(F(40)=5736961\). Millions of states are manageable; tens of trillions are not.

Meet-in-the-middle matching

Write the left-half value as \(A_L+B_Lx\) and the right-half value as \(A_R+B_Rx\). The union represents a positive integer exactly when

$$B_L+B_R=0,\qquad A_L+A_R\gt0.$$

So the right half only needs to be searchable by coefficient \(B_R\). For each left subset, the required right coefficient is \(-B_L\), and every compatible right subset with that coefficient contributes the starting value \(A_L+A_R\).

The C++ implementation stores right entries sorted by coefficient and keeps prefix sums of their constants. That allows it to support even the general case where many right subsets have the same coefficient: binary search locates the matching block, and a second binary search discards entries whose total constant is not positive. The Python and Java ports use a lighter coefficient-to-constant table because, for the \(k\le80\) instances used here, the right coefficients are unique inside each boundary class.

Handling the split boundary

The only missing condition after splitting is the distance rule across the cut. A left subset ending near \(m\) can conflict with a right subset beginning near \(m\). The code records exactly the information needed for this decision.

For the left half, the category is

$$0:\text{ no selected point at }m-2\text{ or }m-1,\qquad 1:\text{ last}=m-2,\qquad 2:\text{ last}=m-1.$$

For the right half, the category is

$$0:\text{ no selected point at }m\text{ or }m+1,\qquad 1:\text{ first}=m+1,\qquad 2:\text{ first}=m.$$

If the left category is \(0\), there is no left-side boundary conflict. If it is \(1\), only right category \(2\) is too close. If it is \(2\), both \(m\) and \(m+1\) are too close, so only right category \(0\) is allowed. This three-by-three compatibility rule is exactly the cross-boundary form of \(|i-j|\ge3\).

How the Code Works

build_powers constructs the pairs \((a_i,b_i)\) for \(x^i=a_i+b_i x\). evaluate_positions sums those pairs for a proposed singleton set. singleton_positions directly simulates the stone process for small checks, which verifies both the algebraic representation and the interpretation of the singleton trace.

enumerate_right recursively enumerates all right-half gap-\(3\) subsets and places them into the three boundary buckets. Each stored entry contains the reduced coefficient and constant. build_left_tasks creates shallow left prefixes; the C++ version uses them to distribute the deeper left search across worker threads.

combine_left_dfs finishes each left subset, computes its boundary category, and queries only compatible right buckets. A match is accepted when the coefficients cancel and the constant sum is positive. That positive constant is added to the answer because it is the represented starting value.

The checkpoint suite verifies the known traces for \(n=1\), \(n=68\), and \(n=90\), then verifies \(S(14)=159\) and \(S(30)=33438\). Only after those checks does the program evaluate \(S(80)\).

Complexity Analysis

The unsplit search would examine \(F(k)\) gap-valid subsets, which is not feasible for \(k=80\). The meet-in-the-middle method stores \(O(F(k/2))\) right-half states and performs \(O(F(k/2))\) left-half searches.

In the sorted-table C++ version, right-table construction costs \(O(F(k/2)\log F(k/2))\), and each completed left subset performs a constant number of binary-search lookups. The map-based Python and Java ports have expected \(O(F(k/2))\) table construction and expected constant-time lookups for the actual \(k\le80\) range. Memory usage is \(O(F(k/2))\).

Footnotes and References

1. Problem page: Project Euler 1003 - Lonely Singles
2. Algebraic number: Wikipedia - Algebraic number
3. Meet-in-the-middle algorithm: Wikipedia - Meet-in-the-middle
4. Independent set: Wikipedia - Independent set
5. Binary search: Wikipedia - Binary search