Problem 818: SET

Project Euler Problem 818 Solution

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

Problem Summary The 81 cards of the SET deck can be modeled as the 81 points of $D=\mathbb F_3^4$. For an $n$-card subset $A\subseteq D$, let $X(A)$ be the number of SETs contained entirely in $A$. The problem defines $$F(n)=\sum_{\substack{A\subseteq D\\|A|=n}} X(A)^4$$ and asks for $F(12)$, with the checkpoints $F(3)=1080$ and $F(6)=159690960$. The implemented solution does not enumerate the $\binom{81}{12}$ subsets. Instead, it rewrites the fourth moment in terms of ordered quadruples of affine lines and counts those quadruples by the size of their union. Mathematical Approach The four attributes of a SET card each have three possible values, so $\mathbb F_3^4$ is the natural ambient space. In that model, every valid SET is exactly a 3-point affine line. Step 1: Encode the deck as affine geometry over $\mathbb F_3$ Encode each attribute value by $0$, $1$, or $2$. Then a card becomes a vector in $\mathbb F_3^4$. Three cards $u,v,w$ form a SET exactly when, coordinate by coordinate, the values are either all equal or all distinct....

Detailed mathematical approach

Problem Summary

The 81 cards of the SET deck can be modeled as the 81 points of $D=\mathbb F_3^4$. For an $n$-card subset $A\subseteq D$, let $X(A)$ be the number of SETs contained entirely in $A$. The problem defines

$$F(n)=\sum_{\substack{A\subseteq D\\|A|=n}} X(A)^4$$

and asks for $F(12)$, with the checkpoints $F(3)=1080$ and $F(6)=159690960$. The implemented solution does not enumerate the $\binom{81}{12}$ subsets. Instead, it rewrites the fourth moment in terms of ordered quadruples of affine lines and counts those quadruples by the size of their union.

Mathematical Approach

The four attributes of a SET card each have three possible values, so $\mathbb F_3^4$ is the natural ambient space. In that model, every valid SET is exactly a 3-point affine line.

Step 1: Encode the deck as affine geometry over $\mathbb F_3$

Encode each attribute value by $0$, $1$, or $2$. Then a card becomes a vector in $\mathbb F_3^4$. Three cards $u,v,w$ form a SET exactly when, coordinate by coordinate, the values are either all equal or all distinct. Over $\mathbb F_3$, both cases are equivalent to

$$u+v+w=0.$$

So once two cards are known, the third one is forced:

$$w=-(u+v).$$

Geometrically, every SET is an affine line of the form

$$\{a,\ a+d,\ a+2d\},\qquad d\ne 0.$$

The number of possible directions is the number of 1-dimensional subspaces of $\mathbb F_3^4$:

$$\frac{3^4-1}{3-1}=40.$$

Each direction partitions the 81 points into $81/3=27$ parallel lines, hence the total number of SETs is

$$L=40\cdot 27=1080.$$

Equivalently, every point lies on one line in each direction, so every card belongs to exactly $40$ SETs.

Step 2: Rewrite the fourth moment as a sum over ordered line quadruples

For any $n$-card subset $A\subseteq D$, $X(A)$ counts how many affine lines are fully contained in $A$. Therefore

$$X(A)^4=\sum_{(\ell_1,\ell_2,\ell_3,\ell_4)} \mathbf{1}_{\ell_1\cup\ell_2\cup\ell_3\cup\ell_4\subseteq A},$$

where the sum runs over ordered quadruples of SET-lines. Summing over all $A$ with $|A|=n$ yields

$$F(n)=\sum_{(\ell_1,\ell_2,\ell_3,\ell_4)} \#\{A\subseteq D:\ |A|=n,\ \ell_1\cup\ell_2\cup\ell_3\cup\ell_4\subseteq A\}.$$

This is the key transformation: instead of choosing subsets first and counting SETs afterward, we choose ordered quadruples of SETs first and then ask how many subsets contain them.

Step 3: Group everything by the union size

For one ordered quadruple $(\ell_1,\ell_2,\ell_3,\ell_4)$, define

$$m=\left|\ell_1\cup\ell_2\cup\ell_3\cup\ell_4\right|.$$

If that union has size $m$, then any $n$-card subset containing it is formed by choosing the remaining $n-m$ cards from the other $81-m$ cards. Thus each such quadruple contributes

$$\binom{81-m}{n-m}$$

to $F(n)$. Now define the histogram

$$N_m=\#\left\{(\ell_1,\ell_2,\ell_3,\ell_4):\left|\ell_1\cup\ell_2\cup\ell_3\cup\ell_4\right|=m\right\}.$$

Because each line has exactly three points, the union size can only range from $3$ to $12$. Therefore

$$\boxed{F(n)=\sum_{m=3}^{12} N_m\binom{81-m}{n-m}.}$$

The original problem has now become a finite-geometry counting problem: determine the ten integers $N_3,N_4,\dots,N_{12}$.

Step 4: Use symmetry to count $N_m$ in $O(L^3)$ time

A brute-force count over all ordered quadruples would require $L^4=1080^4$ cases. The implementations avoid that by representing each line as an 81-bit mask with exactly three set bits. Then the union size of several lines is just the popcount of the bitwise OR of those masks.

The crucial simplification is symmetry: the affine automorphism group of $\mathbb F_3^4$ sends any line to any other line, so the union-size distribution does not depend on which first line is chosen. One may therefore fix a single line $\ell_1$, enumerate all ordered triples $(\ell_2,\ell_3,\ell_4)$, record the resulting union sizes, and multiply the histogram by $L=1080$ at the end.

This turns the dominant work into $O(L^3)$ union operations instead of $O(L^4)$, while preserving the exact values of all $N_m$.

Worked Example: Why the checkpoint $F(3)=1080$ is immediate

If $n=3$, then a 3-card subset $A$ either is a SET or it is not. In the first case $X(A)=1$; in the second case $X(A)=0$. Hence $X(A)^4=X(A)$, and so

$$F(3)=\sum_{\substack{A\subseteq D\\|A|=3}} X(A)=1080.$$

This is exactly the first validation checkpoint used by the implementation.

The same union-size formula also explains the range of $m$. If all four ordered lines are identical, then $m=3$ and each quadruple contributes $\binom{78}{n-3}$. If the four lines are pairwise disjoint, then $m=12$ and each quadruple contributes $\binom{69}{n-12}$. All intermediate overlap patterns fall into the bins $m=4,5,\dots,11$.

A second checkpoint, obtained from the same histogram formula after summing all bins, is

$$F(6)=159690960.$$

Matching both checkpoints is strong evidence that the union-size counts are correct before evaluating the target case $n=12$.

How the Code Works

The C++, Python, and Java implementations all follow the same mathematics. They first list the 81 points of $\mathbb F_3^4$, then enumerate SETs by choosing pairs of points, computing the unique third point $w=-(u+v)$, sorting each triple, and discarding duplicates. That produces exactly 1080 distinct lines.

Each line is then stored as an 81-bit mask split across machine words. The implementation fixes one reference line, scans all ordered triples of additional lines, computes the popcount of the OR-union, and increments the histogram entry for that union size. After multiplying by 1080, it has all values $N_m$.

Finally, the implementation evaluates

$$F(n)=\sum_{m=3}^{12} N_m\binom{81-m}{n-m}$$

with exact integer arithmetic. The C++ and Java implementations parallelize the outer part of the triple scan across threads. The Python implementation delegates to the same compiled counting strategy. All three implementations also verify structural identities such as the total number of lines, the 40 lines through each point, and the checkpoints $F(3)$ and $F(6)$.

Complexity Analysis

Generating the line list is inexpensive: it comes from scanning unordered point pairs, so that stage is bounded by $O(81^2)$ work plus deduplication. The dominant phase is the fixed-line histogram count, which performs $O(L^3)$ iterations with $L=1080$. Each iteration uses only a constant number of bitwise OR and popcount operations, so the method is fast enough in practice despite the large constant $1080^3$.

Memory usage is $O(L)$: the program stores one bit mask per line and a small histogram for the possible union sizes. No large table over subsets is needed.

Footnotes and References

Problem page: https://projecteuler.net/problem=818
SET (card game): Wikipedia - SET (card game)
Affine space: Wikipedia - Affine space
Finite field: Wikipedia - Finite field
Binomial coefficient: Wikipedia - Binomial coefficient

Problem 818 source code

C++

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

#include <boost/multiprecision/cpp_int.hpp>

using namespace std;
using boost::multiprecision::cpp_int;

struct Mask {
    uint64_t lo;  // bits 0..63
    uint64_t hi;  // bits 64..80 stored in bits 0..16
};

static inline int popcnt_u64(uint64_t x) {
    return __builtin_popcountll(x);
}

static inline int mask_size(const Mask& m) {
    return popcnt_u64(m.lo) + popcnt_u64(m.hi);
}

static cpp_int binom_int(int n, int k) {
    if (k < 0 || k > n) return 0;
    k = min(k, n - k);
    cpp_int res = 1;
    for (int i = 1; i <= k; ++i) {
        res *= (n - k + i);
        res /= i;
    }
    return res;
}

int main(int argc, char** argv) {
    ios::sync_with_stdio(false);
    cin.tie(nullptr);

    int threads = static_cast<int>(thread::hardware_concurrency());
    if (threads <= 0) threads = 4;
    bool validate = true;

    for (int i = 1; i < argc; ++i) {
        string arg = argv[i];
        if (arg == "--no-validate") {
            validate = false;
        } else if (arg == "--validate") {
            validate = true;
        } else {
            threads = max(1, atoi(arg.c_str()));
        }
    }

    // Points are GF(3)^4, indexed by base-3 expansion.
    array<array<int, 4>, 81> coord{};
    for (int id = 0; id < 81; ++id) {
        int t = id;
        for (int d = 0; d < 4; ++d) {
            coord[id][d] = t % 3;
            t /= 3;
        }
    }
    const int pow3[4] = {1, 3, 9, 27};
    auto idx_from_coord = [&](const array<int, 4>& c) -> int {
        int id = 0;
        for (int d = 0; d < 4; ++d) id += c[d] * pow3[d];
        return id;
    };

    // Generate all affine lines (SETs) via point pairs.
    unordered_set<uint32_t> seen;
    seen.reserve(2000);
    vector<array<int, 3>> lines;
    lines.reserve(1080);

    for (int i = 0; i < 81; ++i) {
        for (int j = i + 1; j < 81; ++j) {
            array<int, 4> c{};
            for (int d = 0; d < 4; ++d) {
                int s = coord[i][d] + coord[j][d];
                c[d] = (3 - (s % 3)) % 3;
            }
            int k = idx_from_coord(c);

            array<int, 3> tri = {i, j, k};
            sort(tri.begin(), tri.end());
            uint32_t key = static_cast<uint32_t>(tri[0] | (tri[1] << 7) | (tri[2] << 14));
            if (seen.insert(key).second) {
                lines.push_back(tri);
            }
        }
    }

    if (validate && static_cast<int>(lines.size()) != 1080) {
        cerr << "ERROR: expected 1080 lines, got " << lines.size() << "\n";
        return 1;
    }

    vector<Mask> masks;
    masks.reserve(lines.size());
    for (const auto& tri : lines) {
        Mask m{0, 0};
        for (int p : tri) {
            if (p < 64) m.lo |= (1ULL << p);
            else m.hi |= (1ULL << (p - 64));
        }
        masks.push_back(m);

        if (validate && mask_size(m) != 3) {
            cerr << "ERROR: line mask size != 3\n";
            return 1;
        }
    }

    if (validate) {
        array<int, 81> point_deg{};
        for (const auto& tri : lines) {
            for (int p : tri) point_deg[p]++;
        }
        for (int p = 0; p < 81; ++p) {
            if (point_deg[p] != 40) {
                cerr << "ERROR: point " << p << " has degree " << point_deg[p]
                     << " (expected 40)\n";
                return 1;
            }
        }
    }

    const int L = static_cast<int>(masks.size());
    const Mask fixed = masks[0];

    vector<array<uint64_t, 13>> thread_counts(threads);
    auto worker = [&](int tid, int l2_begin, int l2_end) {
        array<uint64_t, 13> local{};
        for (int l2 = l2_begin; l2 < l2_end; ++l2) {
            const Mask& m2 = masks[l2];
            const uint64_t u12_lo = fixed.lo | m2.lo;
            const uint64_t u12_hi = fixed.hi | m2.hi;

            for (int l3 = 0; l3 < L; ++l3) {
                const Mask& m3 = masks[l3];
                const uint64_t u123_lo = u12_lo | m3.lo;
                const uint64_t u123_hi = u12_hi | m3.hi;

                for (int l4 = 0; l4 < L; ++l4) {
                    const Mask& m4 = masks[l4];
                    uint64_t u_lo = u123_lo | m4.lo;
                    uint64_t u_hi = u123_hi | m4.hi;
                    int sz = popcnt_u64(u_lo) + popcnt_u64(u_hi);
                    local[sz]++;
                }
            }
        }
        thread_counts[tid] = local;
    };

    vector<thread> pool;
    pool.reserve(threads);
    int chunk = (L + threads - 1) / threads;
    for (int t = 0; t < threads; ++t) {
        int b = t * chunk;
        int e = min(L, b + chunk);
        if (b >= e) {
            thread_counts[t].fill(0);
            continue;
        }
        pool.emplace_back(worker, t, b, e);
    }
    for (auto& th : pool) th.join();

    array<uint64_t, 13> counts_fixed{};
    for (int t = 0; t < threads; ++t) {
        for (int sz = 0; sz <= 12; ++sz) counts_fixed[sz] += thread_counts[t][sz];
    }

    if (validate) {
        const unsigned long long expected_triples = 1ULL * L * L * L;
        unsigned long long got_triples = 0;
        for (int sz = 0; sz <= 12; ++sz) got_triples += counts_fixed[sz];
        if (got_triples != expected_triples) {
            cerr << "ERROR: triple count mismatch: got " << got_triples
                 << ", expected " << expected_triples << "\n";
            return 1;
        }
    }

    array<uint64_t, 13> N{};
    for (int sz = 0; sz <= 12; ++sz) {
        N[sz] = counts_fixed[sz] * static_cast<uint64_t>(L);
    }

    if (validate) {
        const unsigned long long expected_quads = 1ULL * L * L * L * L;
        unsigned long long got_quads = 0;
        for (int sz = 0; sz <= 12; ++sz) got_quads += N[sz];
        if (got_quads != expected_quads) {
            cerr << "ERROR: quad count mismatch\n";
            return 1;
        }
    }

    if (validate) {
        cpp_int lhs = 0;
        for (int sz = 0; sz <= 12; ++sz) lhs += cpp_int(sz) * cpp_int(N[sz]);

        const long long diff = 74465;     // 27^4 - 26^4
        const long long forty4 = 2560000; // 40^4
        cpp_int rhs = cpp_int(81) * cpp_int(diff) * cpp_int(forty4);
        if (lhs != rhs) {
            cerr << "ERROR: expected union-size total mismatch\n";
            return 1;
        }
    }

    auto compute_F = [&](int n) -> cpp_int {
        cpp_int sum = 0;
        for (int m = 0; m <= 12; ++m) {
            if (N[m] == 0 || n < m) continue;
            sum += cpp_int(N[m]) * binom_int(81 - m, n - m);
        }
        return sum;
    };

    if (validate) {
        cpp_int F3 = compute_F(3);
        cpp_int F6 = compute_F(6);
        if (F3 != 1080) {
            cerr << "ERROR: F(3) mismatch: got " << F3 << ", expected 1080\n";
            return 1;
        }
        if (F6 != 159690960) {
            cerr << "ERROR: F(6) mismatch: got " << F6 << ", expected 159690960\n";
            return 1;
        }
    }

    cpp_int F12 = compute_F(12);
    cout << F12 << "\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.util.ArrayList;
import java.util.Arrays;
import java.util.HashSet;
import java.math.BigInteger;

public class Euler818 {

    static class Mask {
        long lo; // bits 0..63
        long hi; // bits 64..80
    }

    static BigInteger binomInt(int n, int k) {
        if (k < 0 || k > n)
            return BigInteger.ZERO;
        k = Math.min(k, n - k);
        BigInteger res = BigInteger.ONE;
        for (int i = 1; i <= k; ++i) {
            res = res.multiply(BigInteger.valueOf(n - k + i)).divide(BigInteger.valueOf(i));
        }
        return res;
    }

    public static String solve() throws InterruptedException {
        int[][] coord = new int[81][4];
        for (int id = 0; id < 81; ++id) {
            int t = id;
            for (int d = 0; d < 4; ++d) {
                coord[id][d] = t % 3;
                t /= 3;
            }
        }
        int[] pow3 = { 1, 3, 9, 27 };

        HashSet<Integer> seen = new HashSet<>();
        ArrayList<int[]> lines = new ArrayList<>();

        for (int i = 0; i < 81; ++i) {
            for (int j = i + 1; j < 81; ++j) {
                int[] c = new int[4];
                for (int d = 0; d < 4; ++d) {
                    int s = coord[i][d] + coord[j][d];
                    c[d] = (3 - (s % 3)) % 3;
                }
                int k = 0;
                for (int d = 0; d < 4; ++d)
                    k += c[d] * pow3[d];

                int[] tri = { i, j, k };
                Arrays.sort(tri);
                int key = tri[0] | (tri[1] << 7) | (tri[2] << 14);
                if (seen.add(key)) {
                    lines.add(tri);
                }
            }
        }

        ArrayList<Mask> masks = new ArrayList<>();
        for (int[] tri : lines) {
            Mask m = new Mask();
            for (int p : tri) {
                if (p < 64)
                    m.lo |= (1L << p);
                else
                    m.hi |= (1L << (p - 64));
            }
            masks.add(m);
        }

        int L = masks.size();
        Mask fixed = masks.get(0);

        int threads = Math.max(1, Runtime.getRuntime().availableProcessors());
        long[][] threadCounts = new long[threads][13];
        Thread[] pool = new Thread[threads];
        int chunk = (L + threads - 1) / threads;

        for (int t = 0; t < threads; ++t) {
            final int tid = t;
            final int b = t * chunk;
            final int e = Math.min(L, b + chunk);

            pool[t] = new Thread(() -> {
                long[] local = new long[13];
                for (int l2 = b; l2 < e; ++l2) {
                    Mask m2 = masks.get(l2);
                    long u12Lo = fixed.lo | m2.lo;
                    long u12Hi = fixed.hi | m2.hi;

                    for (int l3 = 0; l3 < L; ++l3) {
                        Mask m3 = masks.get(l3);
                        long u123Lo = u12Lo | m3.lo;
                        long u123Hi = u12Hi | m3.hi;

                        for (int l4 = 0; l4 < L; ++l4) {
                            Mask m4 = masks.get(l4);
                            long uLo = u123Lo | m4.lo;
                            long uHi = u123Hi | m4.hi;
                            int sz = Long.bitCount(uLo) + Long.bitCount(uHi);
                            local[sz]++;
                        }
                    }
                }
                threadCounts[tid] = local;
            });
            pool[t].start();
        }

        for (int t = 0; t < threads; ++t) {
            pool[t].join();
        }

        long[] countsFixed = new long[13];
        for (int t = 0; t < threads; ++t) {
            for (int sz = 0; sz <= 12; ++sz) {
                countsFixed[sz] += threadCounts[t][sz];
            }
        }

        BigInteger[] N = new BigInteger[13];
        for (int sz = 0; sz <= 12; ++sz) {
            N[sz] = BigInteger.valueOf(countsFixed[sz]).multiply(BigInteger.valueOf(L));
        }

        BigInteger F12 = BigInteger.ZERO;
        int nTarget = 12;
        for (int m = 0; m <= 12; ++m) {
            if (N[m].signum() == 0 || nTarget < m)
                continue;
            BigInteger add = N[m].multiply(binomInt(81 - m, nTarget - m));
            F12 = F12.add(add);
        }

        return F12.toString();
    }

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