Both this and the Mathematica solution linked below [1] are brute force, which strikes me as far less interesting than using linear programming. Of course, at this size of problem brute force works fine. You can write a brute force version easily in any language, including, of course, Julia:
julia> function brute_force(P,t)
for x1=P[1,:], x2=P[2,:], x3=P[3,:], x4=P[4,:], x5=P[5,:], x6=P[6,:]
x1 + x2 + x3 + x4 + x5 + x6 == t || continue
println("solution: $x1 $x2 $x3 $x4 $x5 $x6")
end
end
brute_force (generic function with 1 method)
julia> brute_force(P,419)
solution: 90 68 91 66 41 63
solution: 90 48 91 66 41 83
solution: 90 64 67 66 49 83
solution: 88 68 91 40 49 83
The interesting thing about Iain's linear programming approach is that you can write the problem down the way it is defined and get a solution efficiently.
I also find it fascinating to read the algebraic solutions of others. This algorithm seems to rely on an open source integer problem solver called "CBC" [1]. But is it necessary for a general 'searching solver' of that kind to be used for this particular problem?
The way in which the problem is framed suggests there is a unique solution. The model is presented as a square matrix subject to linear constraints. There are 36 unknowns, being x[1][1] through to x[36][36], each of which is 0 or 1. There is the constraint that the sum of the unknowns is 419. And there are two constraints to ensure that only one number is chosen from each row and each column, which together give rise to 36 constraint equations.
If any one of the latter constraints is replaced by the sum constraint, there is a linear system of 36 equations in 36 unknowns. Is it possible for that system to be constructed and solved directly through matrix inversion?
(For those trying the exercise at home, the number at P[1,2] appears as 90 but from the full problem statement I think it should be 6.)
Unfortunately not. Matrix inversion (or any type of linear algebra, really) deals with real numbers. This problem is nearly a linear program, which can be solved very efficiently, except that some variables are constrained to { 0, 1 }. As a result, it's an integer program, which is NP-hard. Intuitively, this is hard because there are no derivatives (e.g. gradient, Hessian) to exploit in discrete optimization problems.
One heuristic is to solve an IP is to relax the integer constraint to inequality constraints, solve the LP, and round the results. However, this can do arbitrarily poorly on most problems.
