Problem 688: Piles of Plates

Project Euler Problem 688 Solution

EulerSolve provides an optimized solution for Project Euler Problem 688, Piles of Plates, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Consider plate layouts whose pile heights are consecutive positive integers. If a layout uses exactly $k$ piles and the smallest pile has height $a \ge 1$, then the heights are $$a,\ a+1,\ a+2,\ \dots,\ a+k-1.$$ The total number of plates in that layout is therefore $$ka+\frac{k(k-1)}{2}.$$ For a limit $n$, let $f(n,k)$ count the valid layouts with exactly $k$ piles and total at most $n$. Then $$f(n,k)=\max\left(0,\left\lfloor\frac{n-\frac{k(k-1)}{2}}{k}\right\rfloor\right).$$ Now define $$F(n)=\sum_{k\ge 1} f(n,k),\qquad S(N)=\sum_{n=1}^{N}F(n).$$ The goal is to compute $S(10^{16}) \bmod (10^9+7)$. A direct enumeration over $n$, $k$, and the starting height $a$ would be far too slow, so the solution rewrites the count as a closed floor-sum over $k$. Mathematical Approach Write $$\Delta_k=\frac{k(k-1)}{2}.$$ This is the extra number of plates forced by the consecutive increase $0,1,\dots,k-1$ once the smallest pile height is fixed. Step 1: Count layouts for fixed $n$ and fixed $k$ If the smallest pile has height $a$, the total is $ka+\Delta_k$....

Detailed mathematical approach

Problem Summary

Consider plate layouts whose pile heights are consecutive positive integers. If a layout uses exactly $k$ piles and the smallest pile has height $a \ge 1$, then the heights are

$$a,\ a+1,\ a+2,\ \dots,\ a+k-1.$$

The total number of plates in that layout is therefore

$$ka+\frac{k(k-1)}{2}.$$

For a limit $n$, let $f(n,k)$ count the valid layouts with exactly $k$ piles and total at most $n$. Then

$$f(n,k)=\max\left(0,\left\lfloor\frac{n-\frac{k(k-1)}{2}}{k}\right\rfloor\right).$$

Now define

$$F(n)=\sum_{k\ge 1} f(n,k),\qquad S(N)=\sum_{n=1}^{N}F(n).$$

The goal is to compute $S(10^{16}) \bmod (10^9+7)$. A direct enumeration over $n$, $k$, and the starting height $a$ would be far too slow, so the solution rewrites the count as a closed floor-sum over $k$.

Mathematical Approach

Write

$$\Delta_k=\frac{k(k-1)}{2}.$$

This is the extra number of plates forced by the consecutive increase $0,1,\dots,k-1$ once the smallest pile height is fixed.

Step 1: Count layouts for fixed $n$ and fixed $k$

If the smallest pile has height $a$, the total is $ka+\Delta_k$. Therefore the admissible values of $a$ satisfy

$$ka+\Delta_k\le n.$$

Equivalently,

$$a\le \frac{n-\Delta_k}{k}.$$

Since $a$ must be a positive integer, the number of choices is exactly

$$f(n,k)=\max\left(0,\left\lfloor\frac{n-\Delta_k}{k}\right\rfloor\right).$$

This already explains why large $k$ values disappear quickly: once $\Delta_k$ is close to $n$, there is no room left for even the smallest possible first pile.

Step 2: Swap the order of summation

Instead of computing $F(n)$ for each $n$ and then adding, sum over the number of piles first:

$$S(N)=\sum_{n=1}^{N}\sum_{k\ge 1} f(n,k)=\sum_{k\ge 1}\sum_{n=1}^{N} f(n,k).$$

For a fixed $k$, substitute the formula from Step 1:

$$\sum_{n=1}^{N} f(n,k)=\sum_{n=1}^{N}\max\left(0,\left\lfloor\frac{n-\Delta_k}{k}\right\rfloor\right).$$

Now shift the index by setting $m=n-\Delta_k$. Terms with $m<k$ contribute $0$, so the whole block becomes

$$\sum_{m=0}^{N-\Delta_k}\left\lfloor\frac{m}{k}\right\rfloor,$$

with the understanding that this contribution is zero if $N\le \Delta_k$.

Step 3: Evaluate the floor-sum in closed form

Let

$$M=N-\Delta_k.$$

If $M>0$, divide it by $k$:

$$M=uk+v,\qquad 0\le v<k.$$

Then the values of $\left\lfloor m/k\right\rfloor$ appear in blocks:

$$\underbrace{0,\dots,0}_{k\text{ terms}},\ \underbrace{1,\dots,1}_{k\text{ terms}},\ \dots,\ \underbrace{u-1,\dots,u-1}_{k\text{ terms}},\ \underbrace{u,\dots,u}_{v+1\text{ terms}}.$$

Therefore

$$\sum_{m=0}^{M}\left\lfloor\frac{m}{k}\right\rfloor =k\sum_{j=0}^{u-1} j + u(v+1) =k\frac{u(u-1)}{2}+u(v+1).$$

This removes the inner sum completely.

Step 4: Assemble the final formula

For every $k$ with $\Delta_k<N$, define $u_k$ and $v_k$ by

$$N-\Delta_k=u_k k+v_k,\qquad 0\le v_k<k.$$

Then

$$\boxed{S(N)=\sum_{\substack{k\ge 1\\ \Delta_k<N}}\left(k\frac{u_k(u_k-1)}{2}+u_k(v_k+1)\right)\pmod{10^9+7}.}$$

Once $\Delta_k\ge N$, all later terms vanish. Since $\Delta_k\sim k^2/2$, only $O(\sqrt{N})$ values of $k$ need to be examined.

Worked Example

For a local example, take $n=14$ and $k=3$. The possible layouts have heights

$$1,2,3;\qquad 2,3,4;\qquad 3,4,5,$$

with totals $6,9,12$. The next starting height would give $4+5+6=15>14$, so

$$f(14,3)=3.$$

The formula agrees:

$$\Delta_3=\frac{3\cdot 2}{2}=3,\qquad \left\lfloor\frac{14-3}{3}\right\rfloor=3.$$

For a cumulative example, take $N=10$. The nonzero contributions are

$$\begin{aligned} k=1&:&&\Delta_1=0,\ M=10=10\cdot 1+0,\ &&55,\\ k=2&:&&\Delta_2=1,\ M=9=4\cdot 2+1,\ &&20,\\ k=3&:&&\Delta_3=3,\ M=7=2\cdot 3+1,\ &&7,\\ k=4&:&&\Delta_4=6,\ M=4=1\cdot 4+0,\ &&1. \end{aligned}$$

For $k\ge 5$, the contribution is zero. Hence

$$S(10)=55+20+7+1=83.$$

How the Code Works

The C++, Python, and Java implementations all follow the same formula. They iterate over $k=1,2,3,\dots$, compute the triangular offset $\Delta_k=k(k-1)/2$, and stop as soon as $\Delta_k\ge N$.

For each surviving $k$, the implementation forms $M=N-\Delta_k$, performs one integer division by $k$, and obtains the quotient and remainder needed for

$$k\frac{u(u-1)}{2}+u(v+1).$$

Because the final answer is required modulo $10^9+7$, the division by $2$ is handled by multiplying with the modular inverse of $2$. Every contribution is reduced modulo $10^9+7$ before being added to the running total.

The C++ implementation also uses a wider integer type for intermediate products so that multiplying values near the $10^{16}$ scale remains safe before the modular reduction. In addition, one implementation checks the closed formula against a tiny brute-force computation at $N=100$, where both methods produce $12656$.

Complexity Analysis

The loop stops when $k(k-1)/2\ge N$, so the number of iterations is $O(\sqrt{N})$. Each iteration performs only constant-time arithmetic: one triangular-number computation, one subtraction, one division with remainder, a handful of modular multiplications, and one modular addition. Therefore the total running time is $O(\sqrt{N})$, and the extra memory usage is $O(1)$.

Footnotes and References

Problem page: https://projecteuler.net/problem=688
Triangular number: Wikipedia - Triangular number
Arithmetic progression: Wikipedia - Arithmetic progression
Floor and ceiling functions: Wikipedia - Floor and ceiling functions
Modular arithmetic: Wikipedia - Modular arithmetic

Problem 688 source code

C++

#include <cassert>
#include <cstdint>
#include <iostream>

namespace {

using u64 = std::uint64_t;
using i64 = std::int64_t;
using i128 = __int128_t;

constexpr i64 MOD = 1'000'000'007LL;
constexpr i64 INV2 = (MOD + 1LL) / 2LL;

i64 norm(i64 x) {
    x %= MOD;
    if (x < 0) {
        x += MOD;
    }
    return x;
}

i64 mul_mod(i64 a, i64 b) {
    return static_cast<i64>((static_cast<i128>(a) * static_cast<i128>(b)) % MOD);
}

i64 S_mod(u64 N) {
    i64 ans = 0;

    for (u64 k = 1;; ++k) {
        const u64 c = k * (k - 1ULL) / 2ULL;
        if (c >= N) {
            break;
        }

        const u64 M = N - c;
        const u64 q = M / k;
        const u64 r = M % k;

        i64 term1 = mul_mod(static_cast<i64>(k % static_cast<u64>(MOD)),
                            mul_mod(static_cast<i64>(q % static_cast<u64>(MOD)),
                                    static_cast<i64>((q - 1ULL) % static_cast<u64>(MOD))));
        term1 = mul_mod(term1, INV2);

        i64 term2 = mul_mod(static_cast<i64>(q % static_cast<u64>(MOD)),
                            static_cast<i64>((r + 1ULL) % static_cast<u64>(MOD)));

        ans = norm(ans + term1 + term2);
    }

    return ans;
}

u64 f_small(u64 n, u64 k) {
    if (k == 0ULL) {
        return 0ULL;
    }
    if (n < k * (k + 1ULL) / 2ULL) {
        return 0ULL;
    }
    return (n - k * (k - 1ULL) / 2ULL) / k;
}

u64 F_small(u64 n) {
    u64 out = 0ULL;
    for (u64 k = 1ULL;; ++k) {
        const u64 v = f_small(n, k);
        if (v == 0ULL && k * (k + 1ULL) / 2ULL > n) {
            break;
        }
        out += v;
    }
    return out;
}

u64 S_small(u64 N) {
    u64 out = 0ULL;
    for (u64 n = 1ULL; n <= N; ++n) {
        out += F_small(n);
    }
    return out;
}

}  // namespace

int main() {
    assert(S_small(100ULL) == 12'656ULL);
    assert(S_mod(100ULL) == 12'656LL);

    std::cout << S_mod(10'000'000'000'000'000ULL) << "\n";
    return 0;
}

Python

def solve():
    MOD = 1000000007
    INV2 = (MOD + 1) // 2
    N = 10000000000000000

    ans = 0
    k = 1
    while True:
        c = k * (k - 1) // 2
        if c >= N: break
        M = N - c; q = M // k; r = M % k
        t1 = k % MOD * (q % MOD) % MOD * ((q - 1) % MOD) % MOD * INV2 % MOD
        t2 = q % MOD * ((r + 1) % MOD) % MOD
        ans = (ans + t1 + t2) % MOD
        k += 1

    return str(ans)

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

Java

public class Euler688 {
    public static String solve() {
        long N = 10000000000000000L;
        long MOD = 1000000007L;
        long INV2 = (MOD + 1L) / 2L;

        long ans = 0;

        for (long k = 1;; ++k) {
            long c = k * (k - 1) / 2;
            if (c >= N) {
                break;
            }

            long M = N - c;
            long q = M / k;
            long r = M % k;

            long qMod = q % MOD;

            long term1 = ((k % MOD) * qMod) % MOD;
            term1 = (term1 * ((q - 1 + MOD) % MOD)) % MOD;
            term1 = (term1 * INV2) % MOD;

            long term2 = (qMod * ((r + 1) % MOD)) % MOD;

            ans = (ans + term1 + term2) % MOD;
        }

        return Long.toString(ans);
    }

    public static void main(String[] args) {
        System.out.println(solve());
    }
}