Problem 507: Shortest Lattice Vector

Project Euler Problem 507 Solution

EulerSolve provides an optimized solution for Project Euler Problem 507, Shortest Lattice Vector, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Let the Tribonacci sequence be defined by $$t_1=0,\qquad t_2=t_3=1,\qquad t_{k+3}\equiv t_{k+2}+t_{k+1}+t_k \pmod{10^7}.$$ For each $n\ge 1$, take the 12 consecutive values $$r_k=t_{12n-12+k}\qquad (1\le k\le 12),$$ and form two vectors in $\mathbb{Z}^3$: $$v_n=(r_1-r_2,\ r_3+r_4,\ r_5r_6),\qquad w_n=(r_7-r_8,\ r_9+r_{10},\ r_{11}r_{12}).$$ We must compute $$s(n)=\min_{(a,b)\in \mathbb{Z}^2\setminus\{(0,0)\}} \left\|a v_n+b w_n\right\|_1,$$ where $$\left\|(x,y,z)\right\|_1=|x|+|y|+|z|,$$ and then sum these values: $$S(N)=\sum_{n=1}^{N} s(n).$$ The target computation uses $N=20{,}000{,}000$, so the solution has to generate the sequence and solve each rank-2 lattice problem without storing huge tables. Mathematical Approach The lattice generated by $v_n$ and $w_n$ is $$L_n=\{a v_n+b w_n:\ a,b\in\mathbb{Z}\}.$$ Because the lattice has rank 2, the implementations use an $L_1$-analogue of two-vector reduction: repeatedly shorten the longer generator by subtracting an integer multiple of the shorter one. The critical point is that the best multiple can be found by checking only a handful of integer candidates. Step 1: Generate the Tribonacci data in blocks of 12 Each lattice instance uses exactly 12 consecutive Tribonacci values....

Detailed mathematical approach

Problem Summary

Let the Tribonacci sequence be defined by

$$t_1=0,\qquad t_2=t_3=1,\qquad t_{k+3}\equiv t_{k+2}+t_{k+1}+t_k \pmod{10^7}.$$

For each $n\ge 1$, take the 12 consecutive values

$$r_k=t_{12n-12+k}\qquad (1\le k\le 12),$$

and form two vectors in $\mathbb{Z}^3$:

$$v_n=(r_1-r_2,\ r_3+r_4,\ r_5r_6),\qquad w_n=(r_7-r_8,\ r_9+r_{10},\ r_{11}r_{12}).$$

We must compute

$$s(n)=\min_{(a,b)\in \mathbb{Z}^2\setminus\{(0,0)\}} \left\|a v_n+b w_n\right\|_1,$$

where

$$\left\|(x,y,z)\right\|_1=|x|+|y|+|z|,$$

and then sum these values:

$$S(N)=\sum_{n=1}^{N} s(n).$$

The target computation uses $N=20{,}000{,}000$, so the solution has to generate the sequence and solve each rank-2 lattice problem without storing huge tables.

Mathematical Approach

The lattice generated by $v_n$ and $w_n$ is

$$L_n=\{a v_n+b w_n:\ a,b\in\mathbb{Z}\}.$$

Because the lattice has rank 2, the implementations use an $L_1$-analogue of two-vector reduction: repeatedly shorten the longer generator by subtracting an integer multiple of the shorter one. The critical point is that the best multiple can be found by checking only a handful of integer candidates.

Step 1: Generate the Tribonacci data in blocks of 12

Each lattice instance uses exactly 12 consecutive Tribonacci values. There is no need to materialize the whole sequence up to index $12N$; we only need a rolling state for three consecutive terms.

Once the current state $(t_k,t_{k+1},t_{k+2})$ is known, the next term is

$$t_{k+3}\equiv t_{k+2}+t_{k+1}+t_k \pmod{10^7}.$$

After every 12 streamed values, we obtain one pair $(v_n,w_n)$. This turns the outer problem into a sequence of independent shortest-vector queries in lattices of the form $\mathbb{Z}v+\mathbb{Z}w\subseteq \mathbb{Z}^3$.

Step 2: Rewrite the shortest-vector search as a basis reduction problem

For fixed vectors $v,w$, we want

$$m(v,w)=\min_{(a,b)\in \mathbb{Z}^2\setminus\{(0,0)\}} \left\|av+bw\right\|_1.$$

If we replace $w$ by

$$w'=w-tv\qquad (t\in\mathbb{Z}),$$

then the lattice does not change, because

$$\mathbb{Z}v+\mathbb{Z}w=\mathbb{Z}v+\mathbb{Z}(w-tv).$$

This is the same unimodular basis change used in Gauss reduction. Therefore we are free to shear one generator by an integer multiple of the other whenever that makes the basis shorter in the $L_1$ norm.

Step 3: Minimize one shear step by studying a convex piecewise-linear function

Suppose $v$ is the shorter current vector and we want the best replacement for $w$. Define

$$f(t)=\left\|w-tv\right\|_1=\sum_{i=1}^{3}\left|w_i-t v_i\right|,\qquad t\in\mathbb{Z}.$$

Each summand $\left|w_i-t v_i\right|$ is a convex piecewise-linear function of $t$. Its slope changes only when $t=w_i/v_i$ for a coordinate with $v_i\ne 0$. Therefore $f(t)$ is also convex and piecewise linear, with all breakpoints coming from those coordinate ratios.

Over the integers, a minimizer of a convex piecewise-linear function must occur at an integer adjacent to one of the breakpoints, so it is enough to test

$$\left\lfloor\frac{w_i}{v_i}\right\rfloor,\qquad \left\lfloor\frac{w_i}{v_i}\right\rfloor+1$$

for every coordinate with $v_i\ne 0$, together with $t=0$. In three dimensions this gives at most seven candidates.

That is why the implementations can find the best integer shear in constant time per reduction step instead of searching over all integers.

Step 4: Repeat the reduction until no further improvement is possible

The reduction loop always keeps the shorter vector first. It computes the best integer $t$, forms $w-tv$, and accepts the change only if the $L_1$ norm strictly decreases:

$$\left\|w-tv\right\|_1<\left\|w\right\|_1.$$

Every accepted step preserves the lattice and strictly lowers the norm of one basis vector, so the process must terminate. When no tested integer multiple makes the longer vector shorter, the pair is locally reduced for this $L_1$ shear operation.

The implementations then return the shorter nonzero vector of the reduced pair. For the rank-2 lattices arising here, this is the quantity used throughout the computation of $s(n)$.

Step 5: Jump to any starting index with a companion matrix

The C++ and Java implementations do not restart the Tribonacci sequence from the beginning whenever they need a later index. Instead they use the companion matrix

$$A=\begin{pmatrix} 1 & 1 & 1\\ 1 & 0 & 0\\ 0 & 1 & 0 \end{pmatrix}.$$

This matrix advances the state vector by one step:

$$A\begin{pmatrix}t_{k+2}\\ t_{k+1}\\ t_k\end{pmatrix}=\begin{pmatrix}t_{k+3}\\ t_{k+2}\\ t_{k+1}\end{pmatrix}.$$

So exponentiation by squaring gives any required starting state in $O(\log k)$ time. After that jump, the next terms are streamed linearly.

Worked Example: The first lattice pair

The first 12 Tribonacci values are

$$0,1,1,2,4,7,13,24,44,81,149,274.$$

Hence

$$v=(-1,3,28),\qquad w=(-11,125,40826).$$

The initial norms are

$$\|v\|_1=|-1|+|3|+|28|=32,$$

$$\|w\|_1=|-11|+|125|+|40826|=40962.$$

To reduce $w$, examine the coordinate ratios

$$\frac{-11}{-1}=11,\qquad \frac{125}{3}\approx 41.67,\qquad \frac{40826}{28}\approx 1458.07.$$

The best tested integer is $t=1458$, giving

$$w-1458v=(1447,-4249,2),$$

with

$$\|w-1458v\|_1=1447+4249+2=5698.$$

This is a major improvement over $40962$, but it is still larger than $\|v\|_1=32$. No further integer shear shortens the second vector, so the shortest nonzero lattice vector has $L_1$ length

$$s(1)=32.$$

This matches the first checkpoint used by the implementation. A second checkpoint verifies

$$S(10)=130762273722.$$

How the Code Works

The C++ and Java implementations share the same numerical method. They compute the Tribonacci values on the fly, buffer 12 consecutive terms, build the two 3-dimensional generators, run the reduction loop described above, and add the resulting shortest $L_1$ length to a running total.

Because the third coordinate is a product of two Tribonacci values, intermediate vector norms can become large. The C++ implementation therefore uses 128-bit integer arithmetic for norm calculations and final accumulation, while the Java implementation uses arbitrary-precision integers for the same reason.

The C++ implementation also splits the range $1\le n\le N$ into contiguous chunks when multiple hardware threads are available. Each chunk computes its own starting Tribonacci state via matrix exponentiation, then streams its local terms independently and returns a partial sum. Those partial sums are added at the end.

The Python implementation is intentionally thin: it ensures that the C++ solver is available, runs that compiled program, and returns the printed decimal answer. In other words, the mathematical work is the same, while Python serves as a bridge layer.

Before the full target run, the implementation checks the first generated vector pair, confirms $s(1)=32$, and confirms $S(10)=130762273722$. Those checkpoints protect against indexing mistakes in the sequence stream and against incorrect reduction behavior.

Complexity Analysis

Generating the 12 Tribonacci values for one lattice instance costs $O(1)$ time and $O(1)$ memory. Each reduction iteration evaluates at most seven candidate integers, so one iteration is also $O(1)$. In practice the number of accepted reductions is small for these 3-dimensional rank-2 bases, which makes the average cost per $n$ effectively constant.

As a result, the total sequential running time is linear in the number of instances:

$$O(N).$$

If the work is split across $P$ independent chunks, each chunk pays an additional $O(\log N)$ initialization cost for the matrix jump, so the total work stays linear and the wall-clock time is close to

$$O\!\left(\frac{N}{P}+\log N\right)$$

per chunk, with $O(1)$ memory per chunk.

Footnotes and References

Problem 507 source code

C++

#include <algorithm>
#include <array>
#include <cstdint>
#include <cstdlib>
#include <iostream>
#include <string>
#include <thread>
#include <vector>

using namespace std;

namespace {

constexpr int64_t MOD = 10000000LL;

struct Vec3 {
    int64_t x = 0;
    int64_t y = 0;
    int64_t z = 0;
};

inline __int128 abs128(int64_t v) {
    return v < 0 ? -static_cast<__int128>(v) : static_cast<__int128>(v);
}

inline __int128 l1(const Vec3& v) {
    return abs128(v.x) + abs128(v.y) + abs128(v.z);
}

int64_t floor_div(int64_t a, int64_t b) {
    int64_t q = a / b;
    int64_t r = a % b;
    if (r != 0 && ((r > 0) != (b > 0))) --q;
    return q;
}

__int128 l1_diff(const Vec3& w, const Vec3& v, int64_t t) {
    __int128 sum = 0;
    __int128 val = static_cast<__int128>(w.x) - static_cast<__int128>(t) * v.x;
    sum += val < 0 ? -val : val;
    val = static_cast<__int128>(w.y) - static_cast<__int128>(t) * v.y;
    sum += val < 0 ? -val : val;
    val = static_cast<__int128>(w.z) - static_cast<__int128>(t) * v.z;
    sum += val < 0 ? -val : val;
    return sum;
}

int64_t best_t(const Vec3& v, const Vec3& w) {
    int64_t cand[8];
    int count = 0;
    cand[count++] = 0;

    auto add = [&](int64_t t) { cand[count++] = t; };
    if (v.x != 0) {
        int64_t q = floor_div(w.x, v.x);
        add(q);
        add(q + 1);
    }
    if (v.y != 0) {
        int64_t q = floor_div(w.y, v.y);
        add(q);
        add(q + 1);
    }
    if (v.z != 0) {
        int64_t q = floor_div(w.z, v.z);
        add(q);
        add(q + 1);
    }

    int64_t best = cand[0];
    __int128 best_val = l1_diff(w, v, best);
    for (int i = 1; i < count; ++i) {
        __int128 val = l1_diff(w, v, cand[i]);
        if (val < best_val) {
            best_val = val;
            best = cand[i];
        }
    }
    return best;
}

__int128 shortest_l1(Vec3 a, Vec3 b) {
    while (true) {
        if (l1(b) < l1(a)) swap(a, b);
        int64_t t = best_t(a, b);
        Vec3 nb{
            static_cast<int64_t>(static_cast<__int128>(b.x) - static_cast<__int128>(t) * a.x),
            static_cast<int64_t>(static_cast<__int128>(b.y) - static_cast<__int128>(t) * a.y),
            static_cast<int64_t>(static_cast<__int128>(b.z) - static_cast<__int128>(t) * a.z),
        };
        if (l1(nb) < l1(b)) {
            b = nb;
        } else {
            break;
        }
    }

    __int128 la = l1(a);
    __int128 lb = l1(b);
    if (la == 0) return lb;
    if (lb == 0) return la;
    return la < lb ? la : lb;
}

struct Mat3 {
    int64_t m[3][3];
};

Mat3 mul(const Mat3& a, const Mat3& b) {
    Mat3 r{};
    for (int i = 0; i < 3; ++i) {
        for (int j = 0; j < 3; ++j) {
            __int128 sum = 0;
            for (int k = 0; k < 3; ++k) {
                sum += static_cast<__int128>(a.m[i][k]) * b.m[k][j];
            }
            r.m[i][j] = static_cast<int64_t>(sum % MOD);
        }
    }
    return r;
}

Mat3 mat_pow(Mat3 base, int64_t exp) {
    Mat3 res{};
    for (int i = 0; i < 3; ++i) {
        for (int j = 0; j < 3; ++j) res.m[i][j] = (i == j) ? 1 : 0;
    }
    while (exp > 0) {
        if (exp & 1LL) res = mul(res, base);
        base = mul(base, base);
        exp >>= 1;
    }
    return res;
}

struct Triple {
    int64_t tn = 0;
    int64_t tn1 = 0;
    int64_t tn2 = 0;
};

Triple trib_state(int64_t n) {
    if (n == 1) return {0, 0, 1};
    if (n == 2) return {1, 0, 0};
    Mat3 base{{
        {1, 1, 1},
        {1, 0, 0},
        {0, 1, 0},
    }};
    Mat3 p = mat_pow(base, n - 2);
    return {p.m[0][0], p.m[1][0], p.m[2][0]};
}

__int128 sum_range(int64_t n_start, int64_t n_end) {
    if (n_start > n_end) return 0;
    int64_t idx_start = 12 * n_start - 11;
    int64_t idx_end = 12 * n_end;

    Triple state = trib_state(idx_start);
    int64_t a = state.tn2;
    int64_t b = state.tn1;
    int64_t c = state.tn;

    int64_t total = idx_end - idx_start + 1;
    int64_t rvals[12];
    int pos = 0;
    __int128 sum = 0;

    for (int64_t i = 0; i < total; ++i) {
        rvals[pos++] = c;
        if (pos == 12) {
            Vec3 v{
                rvals[0] - rvals[1],
                rvals[2] + rvals[3],
                rvals[4] * rvals[5],
            };
            Vec3 w{
                rvals[6] - rvals[7],
                rvals[8] + rvals[9],
                rvals[10] * rvals[11],
            };
            sum += shortest_l1(v, w);
            pos = 0;
        }

        int64_t next = a + b + c;
        if (next >= MOD) {
            next -= MOD;
            if (next >= MOD) next -= MOD;
        }
        a = b;
        b = c;
        c = next;
    }

    return sum;
}

__int128 compute_sum(int64_t n, unsigned threads) {
    if (n <= 0) return 0;
    if (threads == 0) threads = 1;
    if (n < 1000) threads = 1;
    if (threads > static_cast<unsigned>(n)) threads = static_cast<unsigned>(n);

    vector<__int128> partial(threads, 0);
    vector<thread> pool;
    pool.reserve(threads);

    int64_t base = n / threads;
    int64_t rem = n % threads;
    int64_t start = 1;

    for (unsigned i = 0; i < threads; ++i) {
        int64_t count = base + (i < static_cast<unsigned>(rem) ? 1 : 0);
        int64_t n_start = start;
        int64_t n_end = start + count - 1;
        start = n_end + 1;

        pool.emplace_back([&, i, n_start, n_end, count]() {
            if (count == 0) {
                partial[i] = 0;
                return;
            }
            partial[i] = sum_range(n_start, n_end);
        });
    }

    for (auto& th : pool) th.join();

    __int128 total = 0;
    for (const auto& v : partial) total += v;
    return total;
}

string to_string128(__int128 value) {
    if (value == 0) return "0";
    bool neg = value < 0;
    if (neg) value = -value;
    string out;
    while (value > 0) {
        int digit = static_cast<int>(value % 10);
        out.push_back(static_cast<char>('0' + digit));
        value /= 10;
    }
    if (neg) out.push_back('-');
    reverse(out.begin(), out.end());
    return out;
}

bool validate() {
    int64_t r[13];
    r[0] = 0;
    r[1] = 0;
    r[2] = 1;
    for (int i = 3; i <= 12; ++i) {
        int64_t next = r[i - 1] + r[i - 2] + r[i - 3];
        if (next >= MOD) {
            next -= MOD;
            if (next >= MOD) next -= MOD;
        }
        r[i] = next;
    }

    Vec3 v{r[1] - r[2], r[3] + r[4], r[5] * r[6]};
    Vec3 w{r[7] - r[8], r[9] + r[10], r[11] * r[12]};
    if (!(v.x == -1 && v.y == 3 && v.z == 28 && w.x == -11 && w.y == 125 && w.z == 40826)) {
        cerr << "Validation failed: first vector pair mismatch." << '\n';
        return false;
    }

    __int128 s1 = compute_sum(1, 1);
    if (s1 != 32) {
        cerr << "Validation failed: S(1)=" << to_string128(s1) << ", expected 32." << '\n';
        return false;
    }

    __int128 s10 = compute_sum(10, 1);
    if (s10 != 130762273722LL) {
        cerr << "Validation failed: sum S(1..10)=" << to_string128(s10)
             << ", expected 130762273722." << '\n';
        return false;
    }
    return true;
}

} // namespace

int main(int argc, char** argv) {
    if (!validate()) return 1;

    int64_t n = 20000000;
    if (argc > 1) n = max<int64_t>(1, atoll(argv[1]));
    unsigned threads = thread::hardware_concurrency();
    if (argc > 2) threads = max(1, atoi(argv[2]));

    __int128 result = compute_sum(n, threads);
    cout << to_string128(result) << '\n';
    return 0;
}

Python

from __future__ import annotations

import re
import shutil
import subprocess
from pathlib import Path

ANSWER_RE = re.compile(r"answer\s*:\s*(.+)$", re.IGNORECASE)
EQUAL_RE = re.compile(r"=\s*(.+)$")


def parse_output(stdout: str) -> str:
    lines = [line.strip() for line in stdout.splitlines() if line.strip()]
    if not lines:
        return ""
    answers = []
    equals = []
    for line in lines:
        m1 = ANSWER_RE.search(line)
        if m1:
            answers.append(m1.group(1).strip())
        m2 = EQUAL_RE.search(line)
        if m2:
            equals.append(m2.group(1).strip())
    if answers:
        return answers[-1]
    if equals:
        return equals[-1]
    return lines[-1]


def should_skip_cpp_checkpoints(src: Path) -> bool:
    try:
        text = src.read_text(encoding="utf-8", errors="ignore")
    except OSError:
        return False
    return "--skip-checkpoints" in text


def run_cpp(binary: Path, src: Path, root: Path) -> str:
    cmd = [str(binary)]
    if should_skip_cpp_checkpoints(src):
        cmd.append("--skip-checkpoints")

    try:
        return subprocess.check_output(cmd, text=True, cwd=root)
    except subprocess.CalledProcessError:
        return subprocess.check_output(cmd, text=True, cwd=src.parent)


def solve() -> str:
    problem_id = __file__.split("Euler")[-1].split(".")[0]
    root = Path(__file__).resolve().parent.parent
    src = root / "solutionsCpp" / f"Euler{problem_id}.cpp"
    binary = root / "solutionsCpp" / f".euler{problem_id}_py_bridge"

    if not binary.exists() or src.stat().st_mtime > binary.stat().st_mtime:
        compiler = shutil.which("clang++") or shutil.which("g++")
        if not compiler:
            raise RuntimeError("No C++ compiler found (clang++/g++).")
        subprocess.check_call([compiler, "-std=c++17", "-O2", str(src), "-o", str(binary)])

    output = run_cpp(binary=binary, src=src, root=root)
    parsed = parse_output(output)
    if not parsed:
        raise RuntimeError(f"Euler{problem_id} bridge produced empty output.")
    return parsed


if __name__ == "__main__":
    print(solve())

Java

import java.math.BigInteger;

public class Euler507 {
    static final long MOD = 10000000L;

    static long floorDiv(long a, long b) {
        long q = a / b, r = a % b;
        if (r != 0 && ((r > 0) != (b > 0)))
            q--;
        return q;
    }

    static long abs64(long v) {
        return v < 0 ? -v : v;
    }

    static BigInteger l1(long[] v) {
        return BigInteger.valueOf(abs64(v[0])).add(BigInteger.valueOf(abs64(v[1])))
                .add(BigInteger.valueOf(abs64(v[2])));
    }

    static BigInteger l1diff(long[] w, long[] v, long t) {
        BigInteger s = BigInteger.ZERO;
        for (int i = 0; i < 3; i++) {
            long val = w[i] - t * v[i];
            s = s.add(BigInteger.valueOf(abs64(val)));
        }
        return s;
    }

    static long bestT(long[] v, long[] w) {
        long[] cand = new long[8];
        int c = 0;
        cand[c++] = 0;
        for (int i = 0; i < 3; i++)
            if (v[i] != 0) {
                long q = floorDiv(w[i], v[i]);
                cand[c++] = q;
                cand[c++] = q + 1;
            }
        long best = cand[0];
        BigInteger bv = l1diff(w, v, best);
        for (int i = 1; i < c; i++) {
            BigInteger val = l1diff(w, v, cand[i]);
            if (val.compareTo(bv) < 0) {
                bv = val;
                best = cand[i];
            }
        }
        return best;
    }

    static BigInteger shortestL1(long[] a, long[] b) {
        while (true) {
            if (l1(b).compareTo(l1(a)) < 0) {
                long[] tmp = a;
                a = b;
                b = tmp;
            }
            long t = bestT(a, b);
            long[] nb = { b[0] - t * a[0], b[1] - t * a[1], b[2] - t * a[2] };
            if (l1(nb).compareTo(l1(b)) < 0)
                b = nb;
            else
                break;
        }
        BigInteger la = l1(a), lb = l1(b);
        if (la.signum() == 0)
            return lb;
        if (lb.signum() == 0)
            return la;
        return la.min(lb);
    }

    static long[][] matMul(long[][] a, long[][] b) {
        long[][] r = new long[3][3];
        for (int i = 0; i < 3; i++)
            for (int j = 0; j < 3; j++) {
                long s = 0;
                for (int k = 0; k < 3; k++)
                    s = (s + a[i][k] * b[k][j]) % MOD;
                r[i][j] = s;
            }
        return r;
    }

    static long[][] matPow(long[][] base, long exp) {
        long[][] res = { { 1, 0, 0 }, { 0, 1, 0 }, { 0, 0, 1 } };
        while (exp > 0) {
            if ((exp & 1) == 1)
                res = matMul(res, base);
            base = matMul(base, base);
            exp >>= 1;
        }
        return res;
    }

    static long[] tribState(long n) {
        if (n == 1)
            return new long[] { 0, 0, 1 };
        if (n == 2)
            return new long[] { 1, 0, 0 };
        long[][] base = { { 1, 1, 1 }, { 1, 0, 0 }, { 0, 1, 0 } };
        long[][] p = matPow(base, n - 2);
        return new long[] { p[0][0], p[1][0], p[2][0] };
    }

    static BigInteger sumRange(long nStart, long nEnd) {
        long idxStart = 12 * nStart - 11, idxEnd = 12 * nEnd;
        long[] st = tribState(idxStart);
        long a = st[2], b = st[1], c = st[0];
        long total = idxEnd - idxStart + 1;
        long[] rvals = new long[12];
        int pos = 0;
        BigInteger sum = BigInteger.ZERO;
        for (long i = 0; i < total; i++) {
            rvals[pos++] = c;
            if (pos == 12) {
                long[] v = { rvals[0] - rvals[1], rvals[2] + rvals[3], rvals[4] * rvals[5] };
                long[] w = { rvals[6] - rvals[7], rvals[8] + rvals[9], rvals[10] * rvals[11] };
                sum = sum.add(shortestL1(v, w));
                pos = 0;
            }
            long next = a + b + c;
            if (next >= MOD) {
                next -= MOD;
                if (next >= MOD)
                    next -= MOD;
            }
            a = b;
            b = c;
            c = next;
        }
        return sum;
    }

    public static void main(String[] args) {
        System.out.println(sumRange(1, 20000000));
    }
}