Problem 782: Distinct Rows and Columns
View on Project EulerProject Euler Problem 782 Solution
EulerSolve provides an optimized solution for Project Euler Problem 782, Distinct Rows and Columns, with C++, Python, Java, and a step-by-step mathematical explanation.
Problem Summary For an \(n\times n\) binary matrix \(A\), let \(r(A)\) be the number of distinct rows and \(c(A)\) the number of distinct columns. The matrix complexity is $$\kappa(A)=\max(r(A),c(A)).$$ For each integer \(k\in\{0,1,\dots,n^2\}\), define $$m_n(k)=\min\{\kappa(A): A \text{ is an } n\times n \text{ binary matrix with exactly } k \text{ ones}\}.$$ The required quantity is $$C(n)=\sum_{k=0}^{n^2} m_n(k).$$ The implementation does not search over matrices directly. Instead, it marks which values of \(k\) are achievable with complexity \(1\), \(2\), or \(3\), and then every unmarked value contributes \(4\). Mathematical Approach The key idea is to treat the problem as a reachability problem on the interval \(0,1,\dots,n^2\). Each construction marks certain counts of ones, and the final sum only depends on how many counts fall into each complexity class. Step 1: Complexity 1 Complexity \(1\) means every row is identical and every column is identical. For a binary square matrix that happens only for the all-zero matrix and the all-one matrix, so $$m_n(0)=m_n(n^2)=1.$$ These two endpoints are the only values with complexity \(1\), and they account for the constant term in the final counting identity....
Detailed mathematical approach
Problem Summary
For an \(n\times n\) binary matrix \(A\), let \(r(A)\) be the number of distinct rows and \(c(A)\) the number of distinct columns. The matrix complexity is
$$\kappa(A)=\max(r(A),c(A)).$$
For each integer \(k\in\{0,1,\dots,n^2\}\), define
$$m_n(k)=\min\{\kappa(A): A \text{ is an } n\times n \text{ binary matrix with exactly } k \text{ ones}\}.$$
The required quantity is
$$C(n)=\sum_{k=0}^{n^2} m_n(k).$$
The implementation does not search over matrices directly. Instead, it marks which values of \(k\) are achievable with complexity \(1\), \(2\), or \(3\), and then every unmarked value contributes \(4\).
Mathematical Approach
The key idea is to treat the problem as a reachability problem on the interval \(0,1,\dots,n^2\). Each construction marks certain counts of ones, and the final sum only depends on how many counts fall into each complexity class.
Step 1: Complexity 1
Complexity \(1\) means every row is identical and every column is identical. For a binary square matrix that happens only for the all-zero matrix and the all-one matrix, so
$$m_n(0)=m_n(n^2)=1.$$
These two endpoints are the only values with complexity \(1\), and they account for the constant term in the final counting identity.
Step 2: Complexity 2
For the two-type constructions used by the implementations, the \(n\) indices are split into two groups of sizes \(a\) and \(b\) with
$$a+b=n.$$
The admissible block patterns reduce to three independent weights:
$$a^2,\qquad b^2,\qquad 2ab.$$
Hence, for each split \(a+b=n\), the implementation enumerates the seven nonempty subset sums
$$a^2,\ b^2,\ 2ab,\ a^2+b^2,\ a^2+2ab,\ b^2+2ab,\ a^2+b^2+2ab=n^2,$$
and then discards the endpoints \(0\) and \(n^2\), which were already handled in the complexity-1 case. The union over all \(a\) gives the full set of counts with complexity \(2\).
Step 3: A Simple Family with Complexity at Most 3
Before the full three-type analysis, the implementations mark the explicit product family
$$k=ab,\qquad k=n^2-ab,\qquad 0\le a,b\le n.$$
The second value follows from complementation: replacing every \(0\) by \(1\) preserves the numbers of distinct rows and columns, so achievable counts come in complementary pairs \(k\leftrightarrow n^2-k\). This inexpensive family already covers many values with complexity at most \(3\).
Step 4: Complexity at Most 3 from Three Types
The main construction splits \(n\) into three parts
$$a+b+c=n,\qquad a\le b\le c.$$
A three-type pattern is encoded by a binary \(3\times3\) template. The implementations enumerate all \(2^9=512\) templates, keep only those whose three rows are pairwise distinct, whose three columns are pairwise distinct, and whose row-pattern multiset matches the column-pattern multiset. After deduplication, this leaves exactly \(46\) distinct coefficient tuples.
For every surviving template, the number of ones has the quadratic form
$$k=d_1a^2+d_2b^2+d_3c^2+e_{ab}ab+e_{ac}ac+e_{bc}bc,$$
where each diagonal coefficient satisfies \(d_i\in\{0,1\}\) and each mixed coefficient satisfies \(e_{ab},e_{ac},e_{bc}\in\{0,1,2\}\). The diagonal terms count within-type blocks, while the mixed terms count cross-type blocks. Evaluating all \(46\) forms over all ordered partitions \(a\le b\le c\) marks every value reached by this three-type construction.
Step 5: Turn Reachability into \(C(n)\)
Let
$$T=n^2+1,$$
$$N_2=\#\{k: m_n(k)=2\},\qquad N_{\le 3}=\#\{k: m_n(k)\le 3\},\qquad N_4=T-N_{\le 3}.$$
Since exactly two values have complexity \(1\), we obtain
$$C(n)=2+2N_2+3(N_{\le 3}-N_2-2)+4N_4.$$
Simplifying gives the identity used by the implementations:
$$C(n)=3T-4-N_2+N_4.$$
So the task is reduced to marking the complexity-2 set and the complexity-\(\le 3\) set, then counting how many integers lie in each.
Worked Example: \(n=5\)
For \(n=5\), we have \(T=5^2+1=26\). The complexity-2 construction marks exactly
$$\{1,4,8,9,12,13,16,17,21,24\},$$
so \(N_2=10\). The combined complexity-\(\le 3\) construction covers all values \(0,1,\dots,25\), hence \(N_4=0\).
Therefore
$$C(5)=3\cdot 26-4-10+0=64.$$
The same result can be read directly from the complexity classes:
$$1+1+10\cdot 2+14\cdot 3=64.$$
So only \(0\) and \(25\) have complexity \(1\), ten values have complexity \(2\), and the remaining fourteen values have complexity \(3\).
How the Code Works
The C++, Python, and Java implementations use the same pipeline. First they generate the finite catalogue of \(46\) three-type coefficient tuples from the \(512\) binary \(3\times3\) templates. Next they allocate two reachability containers over the interval \(0\) through \(n^2\): one container for values with complexity at most \(3\), and one container for values with complexity \(2\). The C++ and Java implementations pack these states into bitsets, while the Python implementation stores them in byte arrays.
After that, the implementation marks the simple family \(ab\) and \(n^2-ab\), precomputes the squares \(0^2,1^2,\dots,n^2\), and loops over all ordered partitions \(a\le b\le c\) with \(a+b+c=n\). For each partition it evaluates the \(46\) quadratic forms and marks the resulting counts. A separate loop over all \(a+b=n\) marks the seven nonempty subset sums built from \(a^2\), \(b^2\), and \(2ab\), ignoring \(0\) and \(n^2\). Finally it counts the marked positions and applies
$$C(n)=3(n^2+1)-4-N_2+\bigl((n^2+1)-N_{\le 3}\bigr).$$
The C++ implementation also checks the known values \(C(5)=64\), \(C(10)=274\), and \(C(20)=1150\) before evaluating the full target input.
Complexity Analysis
The template generation is constant work because the search space has only \(512\) patterns. Marking the product family costs \(O(n^2)\) time. The main partition loop has \(O(n^2)\) states, since \(a\le b\le c\) and \(a+b+c=n\), and each state evaluates only \(46\) constant-size quadratic forms. The complexity-2 loop is \(O(n)\). Therefore the overall running time is \(O(n^2)\).
Memory usage is dominated by the two reachability containers on \(0,\dots,n^2\). The packed C++ and Java versions use \(\Theta(n^2)\) bits, while the Python version uses \(\Theta(n^2)\) bytes for the same logical state space.
Footnotes and References
- Problem page: https://projecteuler.net/problem=782
- Binary matrix: Wikipedia — Binary matrix
- Bit array: Wikipedia — Bit array
- Quadratic form: Wikipedia — Quadratic form
- Integer partition: Wikipedia — Integer partition
Problem 782 source code
C++
#include <algorithm>
#include <array>
#include <cstdint>
#include <iostream>
#include <set>
#include <thread>
#include <vector>
using std::int64_t;
namespace {
struct Coeff {
int c00;
int c11;
int c22;
int cab;
int cac;
int cbc;
};
class Bitset {
public:
explicit Bitset(size_t bits) : data_((bits + 63) / 64, 0) {}
inline void set(size_t idx) {
data_[idx >> 6] |= 1ULL << (idx & 63);
}
void or_with(const Bitset& other) {
const size_t n = data_.size();
for (size_t i = 0; i < n; ++i) {
data_[i] |= other.data_[i];
}
}
uint64_t count() const {
uint64_t total = 0;
for (uint64_t w : data_) {
total += static_cast<uint64_t>(__builtin_popcountll(w));
}
return total;
}
private:
std::vector<uint64_t> data_;
};
std::vector<Coeff> build_coeffs() {
std::set<std::array<int, 6>> uniq;
for (int bits = 0; bits < (1 << 9); ++bits) {
int B[3][3];
for (int i = 0; i < 3; ++i) {
for (int j = 0; j < 3; ++j) {
B[i][j] = (bits >> (i * 3 + j)) & 1;
}
}
int rows[3] = {
B[0][0] | (B[0][1] << 1) | (B[0][2] << 2),
B[1][0] | (B[1][1] << 1) | (B[1][2] << 2),
B[2][0] | (B[2][1] << 1) | (B[2][2] << 2),
};
int cols[3] = {
B[0][0] | (B[1][0] << 1) | (B[2][0] << 2),
B[0][1] | (B[1][1] << 1) | (B[2][1] << 2),
B[0][2] | (B[1][2] << 1) | (B[2][2] << 2),
};
if (!(rows[0] != rows[1] && rows[0] != rows[2] && rows[1] != rows[2])) {
continue;
}
if (!(cols[0] != cols[1] && cols[0] != cols[2] && cols[1] != cols[2])) {
continue;
}
std::array<int, 3> rsorted = {rows[0], rows[1], rows[2]};
std::array<int, 3> csorted = {cols[0], cols[1], cols[2]};
std::sort(rsorted.begin(), rsorted.end());
std::sort(csorted.begin(), csorted.end());
if (rsorted != csorted) {
continue;
}
const int c00 = B[0][0];
const int c11 = B[1][1];
const int c22 = B[2][2];
const int cab = B[0][1] + B[1][0];
const int cac = B[0][2] + B[2][0];
const int cbc = B[1][2] + B[2][1];
uniq.insert({c00, c11, c22, cab, cac, cbc});
}
std::vector<Coeff> coeffs;
coeffs.reserve(uniq.size());
for (const auto& v : uniq) {
coeffs.push_back({v[0], v[1], v[2], v[3], v[4], v[5]});
}
return coeffs;
}
int64_t compute_C(int n, unsigned threads, const std::vector<Coeff>& coeffs) {
const int64_t n2 = static_cast<int64_t>(n) * n;
const size_t bits = static_cast<size_t>(n2) + 1;
Bitset le3(bits);
for (int a = 0; a <= n; ++a) {
for (int b = 0; b <= n; ++b) {
const int64_t k = static_cast<int64_t>(a) * b;
le3.set(static_cast<size_t>(k));
le3.set(static_cast<size_t>(n2 - k));
}
}
std::vector<int64_t> sq(n + 1, 0);
for (int i = 0; i <= n; ++i) {
sq[i] = static_cast<int64_t>(i) * i;
}
const unsigned tcount = std::max(1u, threads);
if (tcount == 1) {
Bitset local(bits);
for (int a = 0; a <= n; ++a) {
const int bmax = (n - a) / 2;
if (bmax < a) continue;
const int64_t aa = sq[a];
for (int b = a; b <= bmax; ++b) {
const int c = n - a - b;
const int64_t bb = sq[b];
const int64_t cc = sq[c];
const int64_t ab = static_cast<int64_t>(a) * b;
const int64_t ac = static_cast<int64_t>(a) * c;
const int64_t bc = static_cast<int64_t>(b) * c;
for (const auto& cf : coeffs) {
const int64_t k = cf.c00 * aa + cf.c11 * bb + cf.c22 * cc
+ cf.cab * ab + cf.cac * ac + cf.cbc * bc;
local.set(static_cast<size_t>(k));
}
}
}
le3.or_with(local);
} else {
const unsigned cap_threads = std::min(threads, 8u);
std::vector<Bitset> locals;
locals.reserve(cap_threads);
for (unsigned t = 0; t < cap_threads; ++t) {
locals.emplace_back(bits);
}
std::vector<std::thread> workers;
workers.reserve(cap_threads);
for (unsigned t = 0; t < cap_threads; ++t) {
workers.emplace_back([&, t]() {
Bitset& local = locals[t];
for (int a = static_cast<int>(t); a <= n; a += static_cast<int>(cap_threads)) {
const int bmax = (n - a) / 2;
if (bmax < a) continue;
const int64_t aa = sq[a];
for (int b = a; b <= bmax; ++b) {
const int c = n - a - b;
const int64_t bb = sq[b];
const int64_t cc = sq[c];
const int64_t ab = static_cast<int64_t>(a) * b;
const int64_t ac = static_cast<int64_t>(a) * c;
const int64_t bc = static_cast<int64_t>(b) * c;
for (const auto& cf : coeffs) {
const int64_t k = cf.c00 * aa + cf.c11 * bb + cf.c22 * cc
+ cf.cab * ab + cf.cac * ac + cf.cbc * bc;
local.set(static_cast<size_t>(k));
}
}
}
});
}
for (auto& th : workers) th.join();
for (const auto& local : locals) {
le3.or_with(local);
}
}
Bitset s2(bits);
for (int a = 0; a <= n; ++a) {
const int b = n - a;
const int64_t w1 = sq[a];
const int64_t w2 = sq[b];
const int64_t w3 = 2LL * a * b;
const int64_t vals[7] = {
w1, w2, w3, w1 + w2, w1 + w3, w2 + w3, w1 + w2 + w3
};
for (int i = 0; i < 7; ++i) {
const int64_t k = vals[i];
if (k == 0 || k == n2) continue;
s2.set(static_cast<size_t>(k));
}
}
const uint64_t count_le3 = le3.count();
const uint64_t count_s2 = s2.count();
const uint64_t total = static_cast<uint64_t>(n2) + 1;
const uint64_t n4 = total - count_le3;
const uint64_t C = 3 * total - 4 - count_s2 + n4;
return static_cast<int64_t>(C);
}
} // namespace
int main() {
const auto coeffs = build_coeffs();
if (compute_C(5, 1, coeffs) != 64) {
std::cerr << "Validation failed for C(5)." << '\n';
return 1;
}
if (compute_C(10, 1, coeffs) != 274) {
std::cerr << "Validation failed for C(10)." << '\n';
return 1;
}
if (compute_C(20, 1, coeffs) != 1150) {
std::cerr << "Validation failed for C(20)." << '\n';
return 1;
}
const int n = 10000;
unsigned threads = std::thread::hardware_concurrency();
if (threads == 0) threads = 1;
threads = std::min(threads, 8u);
std::cout << compute_C(n, threads, coeffs) << '\n';
return 0;
}
Python
def solve():
n = 10000
# Build valid 3x3 coefficient templates
coeffs = set()
for bits in range(512):
B = [[(bits >> (i*3+j)) & 1 for j in range(3)] for i in range(3)]
rows = [B[i][0] | B[i][1]<<1 | B[i][2]<<2 for i in range(3)]
cols = [B[0][j] | B[1][j]<<1 | B[2][j]<<2 for j in range(3)]
if len(set(rows)) != 3 or len(set(cols)) != 3: continue
if sorted(rows) != sorted(cols): continue
c00, c11, c22 = B[0][0], B[1][1], B[2][2]
cab = B[0][1]+B[1][0]; cac = B[0][2]+B[2][0]; cbc = B[1][2]+B[2][1]
coeffs.add((c00, c11, c22, cab, cac, cbc))
coeffs = list(coeffs)
n2 = n*n; nbits = n2+1
# Use bytearray as bitset
le3 = bytearray(nbits); s2 = bytearray(nbits)
for a in range(n+1):
for b in range(n+1):
k = a*b
le3[k] = 1; le3[n2-k] = 1
sq = [i*i for i in range(n+1)]
for a in range(n+1):
bmax = (n-a)//2
if bmax < a: continue
for b in range(a, bmax+1):
c = n-a-b; aa = sq[a]; bb = sq[b]; cc = sq[c]
ab = a*b; ac = a*c; bc = b*c
for c00, c11, c22, cab, cac, cbc in coeffs:
k = c00*aa + c11*bb + c22*cc + cab*ab + cac*ac + cbc*bc
le3[k] = 1
for a in range(n+1):
b = n-a; w1 = sq[a]; w2 = sq[b]; w3 = 2*a*b
for v in (w1, w2, w3, w1+w2, w1+w3, w2+w3, w1+w2+w3):
if v != 0 and v != n2: s2[v] = 1
cnt_le3 = sum(le3); cnt_s2 = sum(s2); total = n2+1
n4 = total - cnt_le3
return str(3*total - 4 - cnt_s2 + n4)
if __name__ == '__main__':
print(solve())
Java
import java.util.ArrayList;
import java.util.Arrays;
import java.util.HashSet;
import java.util.List;
public class Euler782 {
static class Coeff {
int c00, c11, c22, cab, cac, cbc;
Coeff(int c00, int c11, int c22, int cab, int cac, int cbc) {
this.c00 = c00;
this.c11 = c11;
this.c22 = c22;
this.cab = cab;
this.cac = cac;
this.cbc = cbc;
}
@Override
public boolean equals(Object o) {
if (this == o)
return true;
if (!(o instanceof Coeff))
return false;
Coeff coeff = (Coeff) o;
return c00 == coeff.c00 && c11 == coeff.c11 && c22 == coeff.c22 &&
cab == coeff.cab && cac == coeff.cac && cbc == coeff.cbc;
}
@Override
public int hashCode() {
int h = Integer.hashCode(c00);
h = 31 * h + c11;
h = 31 * h + c22;
h = 31 * h + cab;
h = 31 * h + cac;
h = 31 * h + cbc;
return h;
}
}
static class Bitset {
long[] data;
Bitset(long bits) {
data = new long[(int) ((bits + 63) / 64)];
}
void set(long idx) {
data[(int) (idx >> 6)] |= (1L << (idx & 63));
}
void orWith(Bitset other) {
for (int i = 0; i < data.length; ++i) {
data[i] |= other.data[i];
}
}
long count() {
long total = 0;
for (long w : data) {
total += Long.bitCount(w);
}
return total;
}
}
static List<Coeff> buildCoeffs() {
HashSet<Coeff> uniq = new HashSet<>();
for (int bits = 0; bits < (1 << 9); ++bits) {
int[][] B = new int[3][3];
for (int i = 0; i < 3; ++i) {
for (int j = 0; j < 3; ++j) {
B[i][j] = (bits >> (i * 3 + j)) & 1;
}
}
int[] rows = {
B[0][0] | (B[0][1] << 1) | (B[0][2] << 2),
B[1][0] | (B[1][1] << 1) | (B[1][2] << 2),
B[2][0] | (B[2][1] << 1) | (B[2][2] << 2),
};
int[] cols = {
B[0][0] | (B[1][0] << 1) | (B[2][0] << 2),
B[0][1] | (B[1][1] << 1) | (B[2][1] << 2),
B[0][2] | (B[1][2] << 1) | (B[2][2] << 2),
};
if (!(rows[0] != rows[1] && rows[0] != rows[2] && rows[1] != rows[2]))
continue;
if (!(cols[0] != cols[1] && cols[0] != cols[2] && cols[1] != cols[2]))
continue;
int[] rsorted = rows.clone();
int[] csorted = cols.clone();
Arrays.sort(rsorted);
Arrays.sort(csorted);
if (!Arrays.equals(rsorted, csorted))
continue;
uniq.add(new Coeff(B[0][0], B[1][1], B[2][2],
B[0][1] + B[1][0], B[0][2] + B[2][0], B[1][2] + B[2][1]));
}
return new ArrayList<>(uniq);
}
static long computeC(int n, List<Coeff> coeffs) {
long n2 = (long) n * n;
long bits = n2 + 1;
Bitset le3 = new Bitset(bits);
for (int a = 0; a <= n; ++a) {
for (int b = 0; b <= n; ++b) {
long k = (long) a * b;
le3.set(k);
le3.set(n2 - k);
}
}
long[] sq = new long[n + 1];
for (int i = 0; i <= n; ++i) {
sq[i] = (long) i * i;
}
// Multithreading in Java
int numThreads = Runtime.getRuntime().availableProcessors();
numThreads = Math.min(numThreads, 8);
Bitset[] locals = new Bitset[numThreads];
Thread[] threadsArr = new Thread[numThreads];
for (int t = 0; t < numThreads; t++) {
final int tId = t;
final int tCount = numThreads;
locals[tId] = new Bitset(bits);
threadsArr[tId] = new Thread(() -> {
Bitset local = locals[tId];
for (int a = tId; a <= n; a += tCount) {
int bmax = (n - a) / 2;
if (bmax < a)
continue;
long aa = sq[a];
for (int b = a; b <= bmax; ++b) {
int c = n - a - b;
long bb = sq[b];
long cc = sq[c];
long ab = (long) a * b;
long ac = (long) a * c;
long bc = (long) b * c;
for (Coeff cf : coeffs) {
long k = cf.c00 * aa + cf.c11 * bb + cf.c22 * cc
+ cf.cab * ab + cf.cac * ac + cf.cbc * bc;
local.set(k);
}
}
}
});
threadsArr[tId].start();
}
for (int t = 0; t < numThreads; t++) {
try {
threadsArr[t].join();
le3.orWith(locals[t]);
} catch (InterruptedException e) {
e.printStackTrace();
}
}
Bitset s2 = new Bitset(bits);
for (int a = 0; a <= n; ++a) {
int b = n - a;
long w1 = sq[a];
long w2 = sq[b];
long w3 = 2L * a * b;
long[] vals = { w1, w2, w3, w1 + w2, w1 + w3, w2 + w3, w1 + w2 + w3 };
for (long k : vals) {
if (k == 0 || k == n2)
continue;
s2.set(k);
}
}
long countLe3 = le3.count();
long countS2 = s2.count();
long total = n2 + 1;
long n4 = total - countLe3;
return 3L * total - 4L - countS2 + n4;
}
public static String solve() {
List<Coeff> coeffs = buildCoeffs();
return Long.toString(computeC(10000, coeffs));
}
public static void main(String[] args) {
System.out.println(solve());
}
}