Problem 251: Cardano Triplets

Project Euler Problem 251 Solution

EulerSolve provides an optimized solution for Project Euler Problem 251, Cardano Triplets, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary A Cardano triplet is a triple of positive integers $(a,b,c)$ satisfying $$\sqrt[3]{a+b\sqrt{c}}+\sqrt[3]{a-b\sqrt{c}}=1.$$ The task is to count all such triplets with $a+b+c\le L$. The implementation does not search directly in $(a,b,c)$-space; it first converts the radical identity into a Diophantine equation and then enumerates the resulting factorization patterns. Mathematical Approach Let $$x=\sqrt[3]{a+b\sqrt{c}},\qquad y=\sqrt[3]{a-b\sqrt{c}}.$$ Then $x+y=1$. The whole solution comes from extracting arithmetic consequences of this identity. From the Cube-Root Identity to an Integer Equation Using $$x^3+y^3=(x+y)^3-3xy(x+y),$$ we get $$2a=x^3+y^3=1-3xy,$$ so $$xy=\frac{1-2a}{3}.$$ Multiplying the two cube-root arguments also gives $$x^3y^3=(a+b\sqrt{c})(a-b\sqrt{c})=a^2-b^2c,$$ hence $$a^2-b^2c=\left(\frac{1-2a}{3}\right)^3.$$ After rearranging, $$27b^2c=27a^2+(2a-1)^3=(a+1)^2(8a-1).$$ This factorization is the key algebraic simplification behind the code. Why $a=3u-1$ The left-hand side is divisible by $27$, so $(a+1)^2(8a-1)$ must also be divisible by $27$....

Detailed mathematical approach

Problem Summary

A Cardano triplet is a triple of positive integers $(a,b,c)$ satisfying

$$\sqrt[3]{a+b\sqrt{c}}+\sqrt[3]{a-b\sqrt{c}}=1.$$

The task is to count all such triplets with $a+b+c\le L$. The implementation does not search directly in $(a,b,c)$-space; it first converts the radical identity into a Diophantine equation and then enumerates the resulting factorization patterns.

Mathematical Approach

Let

$$x=\sqrt[3]{a+b\sqrt{c}},\qquad y=\sqrt[3]{a-b\sqrt{c}}.$$

Then $x+y=1$. The whole solution comes from extracting arithmetic consequences of this identity.

From the Cube-Root Identity to an Integer Equation

Using

$$x^3+y^3=(x+y)^3-3xy(x+y),$$

we get

$$2a=x^3+y^3=1-3xy,$$

$$xy=\frac{1-2a}{3}.$$

Multiplying the two cube-root arguments also gives

$$x^3y^3=(a+b\sqrt{c})(a-b\sqrt{c})=a^2-b^2c,$$

hence

$$a^2-b^2c=\left(\frac{1-2a}{3}\right)^3.$$

After rearranging,

$$27b^2c=27a^2+(2a-1)^3=(a+1)^2(8a-1).$$

This factorization is the key algebraic simplification behind the code.

Why $a=3u-1$

The left-hand side is divisible by $27$, so $(a+1)^2(8a-1)$ must also be divisible by $27$. Modulo $3$, this forces

$$a\equiv 2\pmod 3.$$

Therefore we write

$$a=3u-1,\qquad u\ge 1.$$

Then $a+1=3u$ and $8a-1=3(8u-3)$, so the previous equation becomes

$$b^2c=u^2(8u-3).$$

Thus every Cardano triplet corresponds to a positive integer $u$ together with a decomposition of $u^2(8u-3)$ into a square part and a residual part.

Complete Parameterization Counted by the Implementation

Set

$$m=8u-3.$$

We need positive integers $b,c$ such that

$$b^2c=u^2m.$$

Let

$$g=\gcd(b,u),\qquad u=gh,\qquad b=gb_1,$$

so that $\gcd(h,b_1)=1$. Substituting into $b^2c=u^2m$ gives

$$g^2b_1^2c=g^2h^2m\qquad\Longrightarrow\qquad b_1^2c=h^2m.$$

Because $\gcd(h,b_1)=1$, every prime dividing $b_1$ must come from $m$, so $b_1^2\mid m$. Write

$$m=b_1^2r.$$

Then automatically

$$c=h^2r,\qquad b=\frac{u}{h}b_1.$$

Hence every solution enumerated by the code has the form

$$a=3u-1,\qquad m=8u-3=b_1^2r,\qquad h\mid u,\qquad \gcd(h,b_1)=1,$$

$$b=\frac{u}{h}b_1,\qquad c=h^2r.$$

Conversely, any such choice satisfies $b^2c=u^2m$, so this parameterization is complete.

Local Bounds from the Budget Constraint

For fixed $u$, the value of $a=3u-1$ is already known. Define the remaining budget

$$R=L-a.$$

Any candidate must satisfy

$$b+c=\frac{u}{h}b_1+h^2r\le R.$$

The code derives two strong bounds for $h$. Since $u/h\ge 1$, we have $b\ge b_1$, hence

$$c=h^2r\le R-b_1,$$

which yields

$$h\le\left\lfloor\sqrt{\frac{R-b_1}{r}}\right\rfloor.$$

Also $h\ge 1$ implies $c=h^2r\ge r$, so $b\le R-r$. Therefore

$$\frac{u}{h}b_1\le R-r\qquad\Longrightarrow\qquad h\ge\left\lceil\frac{ub_1}{R-r}\right\rceil$$

whenever $R-r\gt 0$. The implementation checks only divisors $h\mid u$ inside this interval.

A Global Upper Bound for $u$

Before the main enumeration starts, the code computes a safe maximal $u$ by lower-bounding $b+c$. Let

$$x=\frac{u}{h}b_1,\qquad y=h^2r.$$

Then $x+y=b+c$ and

$$x^2y=u^2b_1^2r=u^2(8u-3).$$

For fixed $x^2y=K$, the minimum of $x+y$ is $3\sqrt[3]{K/4}$, by calculus or the weighted AM-GM inequality. Therefore every valid triplet must satisfy

$$a+b+c\ge 3u-1+3\sqrt[3]{\frac{u^2(8u-3)}{4}}.$$

If this lower bound already exceeds $L$, then that $u$ cannot contribute. Since the bound grows with $u$, binary search yields the cutoff used by max_u_bound.

Worked Example: $(2,1,5)$

Take $u=1$. Then

$$a=3\cdot 1-1=2,\qquad m=8\cdot 1-3=5.$$

The only square divisor of $m=5$ is $b_1=1$, and the only divisor of $u=1$ is $h=1$. Thus

$$b=\frac{1}{1}\cdot 1=1,\qquad c=1^2\cdot 5=5.$$

This reproduces the example $(2,1,5)$, and the parameterization counts it exactly once.

How the Code Works

The implementation first computes max_u_bound(limit), then builds an SPF table for fast factorization of every $u$, and a prime table up to $\sqrt{8u_{\max}-3}$ for factoring $m=8u-3$. For each $u$, it enumerates all divisors $h\mid u$, all square-divisor choices $b_1$ of $m$, applies the interval bounds for $h$, checks $\gcd(h,b_1)=1$, and finally verifies $b+c\le R$. The outer loop over $u$ is split across threads because different $u$-values are independent.

Complexity Analysis

Precomputation costs $O(u_{\max})$ time and memory for the SPF table, plus a smaller prime sieve up to $\sqrt{8u_{\max}-3}$. The counting phase is dominated by divisor generation: for each $u$, the algorithm enumerates divisors of $u$ and square divisors of $8u-3$. A simple closed form is awkward because divisor counts fluctuate, but the effective runtime is far below naive search thanks to the global $u$-bound, the $h$-interval pruning, and the coprimality filter.

Problem 251 source code

C++

#include <algorithm>
#include <cmath>
#include <cstdint>
#include <iostream>
#include <numeric>
#include <string>
#include <thread>
#include <utility>
#include <vector>

namespace {

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

struct Options {
    i64 limit = 110000000;
    int threads = 0;
    bool run_checkpoints = true;
};

bool parse_i64_after_prefix(const std::string& arg, const std::string& prefix, i64& value) {
    if (arg.rfind(prefix, 0U) != 0U) {
        return false;
    }
    const std::string tail = arg.substr(prefix.size());
    if (tail.empty()) {
        return false;
    }

    i64 parsed = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10 + static_cast<i64>(c - '0');
    }
    value = parsed;
    return true;
}

bool parse_int_after_prefix(const std::string& arg, const std::string& prefix, int& value) {
    if (arg.rfind(prefix, 0U) != 0U) {
        return false;
    }
    const std::string tail = arg.substr(prefix.size());
    if (tail.empty()) {
        return false;
    }

    int parsed = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10 + static_cast<int>(c - '0');
    }
    value = parsed;
    return true;
}

bool parse_arguments(int argc, char** argv, Options& options) {
    for (int i = 1; i < argc; ++i) {
        const std::string arg(argv[i]);
        if (arg == "--skip-checkpoints") {
            options.run_checkpoints = false;
            continue;
        }
        if (parse_i64_after_prefix(arg, "--limit=", options.limit) ||
            parse_int_after_prefix(arg, "--threads=", options.threads)) {
            continue;
        }

        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }
    return options.limit >= 1 && options.threads >= 0;
}

std::vector<int> sieve_primes(const int limit) {
    std::vector<bool> is_prime(static_cast<std::size_t>(limit + 1), true);
    is_prime[0] = false;
    if (limit >= 1) {
        is_prime[1] = false;
    }

    for (int i = 2; static_cast<i64>(i) * i <= limit; ++i) {
        if (!is_prime[static_cast<std::size_t>(i)]) {
            continue;
        }
        for (int j = i * i; j <= limit; j += i) {
            is_prime[static_cast<std::size_t>(j)] = false;
        }
    }

    std::vector<int> primes;
    for (int i = 2; i <= limit; ++i) {
        if (is_prime[static_cast<std::size_t>(i)]) {
            primes.push_back(i);
        }
    }
    return primes;
}

std::vector<int> build_spf(const int n) {
    std::vector<int> spf(static_cast<std::size_t>(n + 1), 0);
    std::vector<int> primes;
    spf[0] = 0;
    if (n >= 1) {
        spf[1] = 1;
    }

    for (int i = 2; i <= n; ++i) {
        if (spf[static_cast<std::size_t>(i)] == 0) {
            spf[static_cast<std::size_t>(i)] = i;
            primes.push_back(i);
        }
        for (int p : primes) {
            const i64 v = static_cast<i64>(p) * static_cast<i64>(i);
            if (v > n || p > spf[static_cast<std::size_t>(i)]) {
                break;
            }
            spf[static_cast<std::size_t>(v)] = p;
        }
    }

    return spf;
}

void factorize_with_spf(int x, const std::vector<int>& spf, std::vector<std::pair<int, int>>& factors) {
    factors.clear();
    while (x > 1) {
        const int p = spf[static_cast<std::size_t>(x)];
        int e = 0;
        do {
            x /= p;
            ++e;
        } while (x > 1 && spf[static_cast<std::size_t>(x)] == p);
        factors.push_back({p, e});
    }
}

void factorize_square_part(i64 m, const std::vector<int>& primes, std::vector<std::pair<int, int>>& square_factors) {
    square_factors.clear();
    i64 x = m;

    for (int p : primes) {
        if (static_cast<i64>(p) * static_cast<i64>(p) > x) {
            break;
        }
        if ((x % p) != 0) {
            continue;
        }

        int e = 0;
        do {
            x /= p;
            ++e;
        } while ((x % p) == 0);

        if (e >= 2) {
            square_factors.push_back({p, e / 2});
        }
    }
}

void generate_divisors(const std::vector<std::pair<int, int>>& factors, std::vector<i64>& divisors) {
    divisors.clear();
    divisors.push_back(1);

    for (const auto& [p, e] : factors) {
        const std::size_t base_size = divisors.size();
        i64 pe = 1;
        for (int k = 1; k <= e; ++k) {
            pe *= static_cast<i64>(p);
            for (std::size_t i = 0; i < base_size; ++i) {
                divisors.push_back(divisors[i] * pe);
            }
        }
    }
}

i64 floor_sqrt_i64(const i64 x) {
    if (x <= 0) {
        return 0;
    }
    i64 r = static_cast<i64>(std::sqrt(static_cast<long double>(x)));
    while (static_cast<i128>(r + 1) * static_cast<i128>(r + 1) <= x) {
        ++r;
    }
    while (static_cast<i128>(r) * static_cast<i128>(r) > x) {
        --r;
    }
    return r;
}

bool can_have_solution(const i64 u, const i64 limit) {
    const long double uu = static_cast<long double>(u);
    const long double mm = 8.0L * uu - 3.0L;
    const long double lower_bound = 3.0L * uu - 1.0L + 3.0L * std::cbrt((uu * uu * mm) / 4.0L);
    return lower_bound <= static_cast<long double>(limit) + 1e-9L;
}

i64 max_u_bound(const i64 limit) {
    i64 lo = 0;
    i64 hi = (limit + 1) / 3;
    while (lo < hi) {
        const i64 mid = lo + (hi - lo + 1) / 2;
        if (can_have_solution(mid, limit)) {
            lo = mid;
        } else {
            hi = mid - 1;
        }
    }
    return lo;
}

i64 solve(const i64 limit, int threads) {
    if (limit < 8) {
        return 0;
    }

    const i64 max_u = max_u_bound(limit);
    if (max_u <= 0) {
        return 0;
    }

    const int sqrt_max_m = static_cast<int>(std::sqrt(static_cast<long double>(8 * max_u - 3))) + 1;
    const std::vector<int> primes = sieve_primes(sqrt_max_m);
    const std::vector<int> spf = build_spf(static_cast<int>(max_u));

    if (threads <= 0) {
        threads = static_cast<int>(std::thread::hardware_concurrency());
        if (threads <= 0) {
            threads = 1;
        }
    }
    threads = std::max(1, std::min(threads, static_cast<int>(max_u)));

    std::vector<i64> local_counts(static_cast<std::size_t>(threads), 0);
    std::vector<std::thread> workers;
    workers.reserve(static_cast<std::size_t>(threads));

    for (int t = 0; t < threads; ++t) {
        workers.emplace_back([&, t]() {
            std::vector<std::pair<int, int>> factors_u;
            std::vector<std::pair<int, int>> square_factors_m;
            std::vector<i64> divisors_h;
            std::vector<i64> divisors_b1;
            factors_u.reserve(16);
            square_factors_m.reserve(8);
            divisors_h.reserve(256);
            divisors_b1.reserve(64);

            i64 subtotal = 0;
            for (i64 u = static_cast<i64>(t) + 1; u <= max_u; u += threads) {
                const i64 a = 3 * u - 1;
                const i64 remaining = limit - a;
                if (remaining <= 1) {
                    continue;
                }

                factorize_with_spf(static_cast<int>(u), spf, factors_u);
                generate_divisors(factors_u, divisors_h);

                const i64 m = 8 * u - 3;
                factorize_square_part(m, primes, square_factors_m);
                generate_divisors(square_factors_m, divisors_b1);

                for (const i64 b1 : divisors_b1) {
                    const i64 b1_sq = b1 * b1;
                    const i64 m_over = m / b1_sq;

                    if (b1 >= remaining || b1 + m_over > remaining) {
                        continue;
                    }

                    const i64 h_max = floor_sqrt_i64((remaining - b1) / m_over);
                    if (h_max <= 0) {
                        continue;
                    }

                    i64 h_min = 1;
                    const i64 denom = remaining - m_over;
                    if (denom > 0) {
                        h_min = (u * b1 + denom - 1) / denom;
                    }

                    for (const i64 h : divisors_h) {
                        if (h < h_min || h > h_max) {
                            continue;
                        }
                        if (std::gcd(h, b1) != 1) {
                            continue;
                        }

                        const i64 b = (u / h) * b1;
                        const i64 c = h * h * m_over;
                        if (static_cast<i128>(b) + static_cast<i128>(c) <= remaining) {
                            ++subtotal;
                        }
                    }
                }
            }
            local_counts[static_cast<std::size_t>(t)] = subtotal;
        });
    }

    for (std::thread& worker : workers) {
        worker.join();
    }

    i64 total = 0;
    for (const i64 v : local_counts) {
        total += v;
    }
    return total;
}

bool run_checkpoints() {
    if (solve(8, 1) != 1) {
        std::cerr << "Checkpoint failed for limit=8" << '\n';
        return false;
    }
    if (solve(100, 1) != 11) {
        std::cerr << "Checkpoint failed for limit=100" << '\n';
        return false;
    }
    if (solve(1000, 1) != 149) {
        std::cerr << "Checkpoint failed for limit=1000" << '\n';
        return false;
    }
    return true;
}

}  // namespace

int main(int argc, char** argv) {
    Options options;
    if (!parse_arguments(argc, argv, options)) {
        return 1;
    }
    if (options.run_checkpoints && !run_checkpoints()) {
        return 2;
    }

    std::cout << solve(options.limit, options.threads) << '\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 ""

    answer_candidates = []
    equal_candidates = []
    for line in lines:
        m1 = ANSWER_RE.search(line)
        if m1:
            answer_candidates.append(m1.group(1).strip())
        m2 = EQUAL_RE.search(line)
        if m2:
            equal_candidates.append(m2.group(1).strip())

    if answer_candidates:
        return answer_candidates[-1]
    if equal_candidates:
        return equal_candidates[-1]
    return lines[-1]


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 = subprocess.check_output([str(binary)], text=True)
    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.util.ArrayList;
import java.util.List;
import java.util.concurrent.*;

public class Euler251 {
    static final long LIMIT = 110000000L;

    static List<Integer> sievePrimes(int limit) {
        boolean[] isPrime = new boolean[limit + 1];
        for (int i = 2; i <= limit; i++)
            isPrime[i] = true;

        for (int i = 2; (long) i * i <= limit; ++i) {
            if (!isPrime[i])
                continue;
            for (int j = i * i; j <= limit; j += i) {
                isPrime[j] = false;
            }
        }

        List<Integer> primes = new ArrayList<>();
        for (int i = 2; i <= limit; ++i) {
            if (isPrime[i])
                primes.add(i);
        }
        return primes;
    }

    static int[] buildSPF(int n) {
        int[] spf = new int[n + 1];
        if (n >= 1)
            spf[1] = 1;
        List<Integer> primes = new ArrayList<>();

        for (int i = 2; i <= n; ++i) {
            if (spf[i] == 0) {
                spf[i] = i;
                primes.add(i);
            }
            for (int p : primes) {
                long v = (long) p * i;
                if (v > n || p > spf[i])
                    break;
                spf[(int) v] = p;
            }
        }
        return spf;
    }

    static class Factor {
        int p, e;

        Factor(int p, int e) {
            this.p = p;
            this.e = e;
        }
    }

    static void factorizeWithSPF(int x, int[] spf, List<Factor> factors) {
        factors.clear();
        while (x > 1) {
            int p = spf[x];
            int e = 0;
            do {
                x /= p;
                e++;
            } while (x > 1 && spf[x] == p);
            factors.add(new Factor(p, e));
        }
    }

    static void factorizeSquarePart(long m, List<Integer> primes, List<Factor> squareFactors) {
        squareFactors.clear();
        long x = m;
        for (int p : primes) {
            if ((long) p * p > x)
                break;
            if (x % p != 0)
                continue;
            int e = 0;
            do {
                x /= p;
                e++;
            } while (x % p == 0);
            if (e >= 2) {
                squareFactors.add(new Factor(p, e / 2));
            }
        }
    }

    static void generateDivisors(List<Factor> factors, List<Long> divisors) {
        divisors.clear();
        divisors.add(1L);
        for (Factor f : factors) {
            int baseSize = divisors.size();
            long pe = 1;
            for (int k = 1; k <= f.e; ++k) {
                pe *= f.p;
                for (int i = 0; i < baseSize; ++i) {
                    divisors.add(divisors.get(i) * pe);
                }
            }
        }
    }

    static long floorSqrt(long x) {
        if (x <= 0)
            return 0;
        long r = (long) Math.sqrt(x);
        while (BigIntegerUtil.multiply(r + 1, r + 1).compareTo(BigIntegerUtil.valueOf(x)) <= 0)
            r++;
        while (BigIntegerUtil.multiply(r, r).compareTo(BigIntegerUtil.valueOf(x)) > 0)
            r--;
        return r;
    }

    static boolean canHaveSolution(long u, long limit) {
        double uu = (double) u;
        double mm = 8.0 * uu - 3.0;
        double lowerBound = 3.0 * uu - 1.0 + 3.0 * Math.cbrt((uu * uu * mm) / 4.0);
        return lowerBound <= (double) limit + 1e-9;
    }

    static long maxUBound(long limit) {
        long lo = 0;
        long hi = (limit + 1) / 3;
        while (lo < hi) {
            long mid = lo + (hi - lo + 1) / 2;
            if (canHaveSolution(mid, limit)) {
                lo = mid;
            } else {
                hi = mid - 1;
            }
        }
        return lo;
    }

    static long gcd(long a, long b) {
        while (b != 0) {
            long t = b;
            b = a % b;
            a = t;
        }
        return a;
    }

    public static String solve() {
        if (LIMIT < 8)
            return "0";

        long maxU = maxUBound(LIMIT);
        if (maxU <= 0)
            return "0";

        int sqrtMaxM = (int) Math.sqrt(8.0 * maxU - 3.0) + 1;
        List<Integer> primes = sievePrimes(sqrtMaxM);
        int[] spf = buildSPF((int) maxU);

        int threads = Math.max(1, Runtime.getRuntime().availableProcessors());
        threads = (int) Math.min(threads, Math.max(1, maxU));

        ExecutorService executor = Executors.newFixedThreadPool(threads);
        List<Future<Long>> futures = new ArrayList<>();

        for (int t = 0; t < threads; ++t) {
            final int tId = t;
            final int tCount = threads;
            futures.add(executor.submit(() -> {
                List<Factor> factorsU = new ArrayList<>(16);
                List<Factor> squareFactorsM = new ArrayList<>(8);
                List<Long> divisorsH = new ArrayList<>(256);
                List<Long> divisorsB1 = new ArrayList<>(64);

                long subtotal = 0;
                for (long u = tId + 1; u <= maxU; u += tCount) {
                    long a = 3 * u - 1;
                    long remaining = LIMIT - a;
                    if (remaining <= 1)
                        continue;

                    factorizeWithSPF((int) u, spf, factorsU);
                    generateDivisors(factorsU, divisorsH);

                    long m = 8 * u - 3;
                    factorizeSquarePart(m, primes, squareFactorsM);
                    generateDivisors(squareFactorsM, divisorsB1);

                    for (long b1 : divisorsB1) {
                        long b1Sq = b1 * b1;
                        long mOver = m / b1Sq;

                        if (b1 >= remaining || b1 + mOver > remaining)
                            continue;

                        long hMax = floorSqrt((remaining - b1) / mOver);
                        if (hMax <= 0)
                            continue;

                        long hMin = 1;
                        long denom = remaining - mOver;
                        if (denom > 0) {
                            hMin = (u * b1 + denom - 1) / denom;
                        }

                        for (long h : divisorsH) {
                            if (h < hMin || h > hMax)
                                continue;
                            if (gcd(h, b1) != 1)
                                continue;

                            long b = (u / h) * b1;
                            long c = h * h * mOver;
                            if (BigIntegerUtil.valueOf(b).add(BigIntegerUtil.valueOf(c))
                                    .compareTo(BigIntegerUtil.valueOf(remaining)) <= 0) {
                                subtotal++;
                            }
                        }
                    }
                }
                return subtotal;
            }));
        }

        long total = 0;
        for (Future<Long> f : futures) {
            try {
                total += f.get();
            } catch (Exception e) {
            }
        }
        executor.shutdown();

        return String.valueOf(total);
    }

    static class BigIntegerUtil {
        static java.math.BigInteger valueOf(long val) {
            return java.math.BigInteger.valueOf(val);
        }

        static java.math.BigInteger multiply(long a, long b) {
            return valueOf(a).multiply(valueOf(b));
        }
    }

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