Solving a big system of equations using certain conditions

Question

I have the following problem: I want to solve a big system of equations for a certain number of variables which I call $a_k$, the variables need to satisfy the following conditions: all $a_k\in\mathbb{N}^+$ so they have to be a natural number bigger than $0$, besides that all the variables can not be equal to each other so: $a_k\ne a_{k+1}$ for all $k$. The value of $n$ is chosen in the problem.

So, as an example I wrote the following code:

n = 15; Solve[{n == a1 + a2 + a3, n == a4 + a5 + a6, 
  n == a7 + a8 + a9, n == a1 + a4 + a7, n == a2 + a5 + a8, 
  n == a3 + a6 + a9, n == a1 + a5 + a9, n == a3 + a5 + a7, 
  1 <= a1 < n && 1 <= a2 < n && 1 <= a3 < n && 1 <= a4 < n && 
   1 <= a5 < n && 1 <= a6 < n && 1 <= a7 < n && 1 <= a8 < n && 
   1 <= a9 < n && a1 != a2 && a1 != a3 && a1 != a4 && a1 != a5 && 
   a1 != a6 && a1 != a7 && a1 != a8 && a1 != a9 && a2 != a3 && 
   a2 != a4 && a2 != a5 && a2 != a6 && a2 != a7 && a2 != a8 && 
   a2 != a9 && a3 != a4 && a3 != a5 && a3 != a6 && a3 != a7 && 
   a3 != a8 && a3 != a9 && a4 != a5 && a4 != a6 && a4 != a7 && 
   a4 != a8 && a4 != a9 && a5 != a6 && a5 != a7 && a5 != a8 && 
   a5 != a9 && a6 != a7 && a6 != a8 && a6 != a9 && a7 != a8 && 
   a7 != a9 && a8 != a9 && a1 \[Element] Integers && 
   a2 \[Element] Integers && a3 \[Element] Integers && 
   a4 \[Element] Integers && a5 \[Element] Integers && 
   a6 \[Element] Integers && a7 \[Element] Integers && 
   a8 \[Element] Integers && a9 \[Element] Integers}, {a1, a2, a3, a4,
   a5, a6, a7, a8, a9}]

But my question is: is there a way to program way more easily than all the variables have to take another value so no number shows up twice or more as a solution? In other words, how can I make my code more compact? The second question is: the output I get is basically the same but it shows me:

{{a1 -> 2, a2 -> 7, a3 -> 6, a4 -> 9, a5 -> 5, a6 -> 1, a7 -> 4, 
  a8 -> 3, a9 -> 8}, {a1 -> 2, a2 -> 9, a3 -> 4, a4 -> 7, a5 -> 5, 
  a6 -> 3, a7 -> 6, a8 -> 1, a9 -> 8}, {a1 -> 4, a2 -> 3, a3 -> 8, 
  a4 -> 9, a5 -> 5, a6 -> 1, a7 -> 2, a8 -> 7, a9 -> 6}, {a1 -> 4, 
  a2 -> 9, a3 -> 2, a4 -> 3, a5 -> 5, a6 -> 7, a7 -> 8, a8 -> 1, 
  a9 -> 6}, {a1 -> 6, a2 -> 1, a3 -> 8, a4 -> 7, a5 -> 5, a6 -> 3, 
  a7 -> 2, a8 -> 9, a9 -> 4}, {a1 -> 6, a2 -> 7, a3 -> 2, a4 -> 1, 
  a5 -> 5, a6 -> 9, a7 -> 8, a8 -> 3, a9 -> 4}, {a1 -> 8, a2 -> 1, 
  a3 -> 6, a4 -> 3, a5 -> 5, a6 -> 7, a7 -> 4, a8 -> 9, 
  a9 -> 2}, {a1 -> 8, a2 -> 3, a3 -> 4, a4 -> 1, a5 -> 5, a6 -> 9, 
  a7 -> 6, a8 -> 7, a9 -> 2}}

These are all the same solution, so how can I make my code such that I only see one of these same solutions?

Answer 1

One way I prefer to go about such problems is using 0-1 programming. So each of the 9 variables gets split into 15 component sub-variables. The components are all 0 except for one which is 1. The all-different constraint is handled by insisting that the sum over each common sub-component index is between 0 and 1.

Here is code for this.

n = 15;
nvars = 9;
vars = Array[a, {nvars, n - 1}];
fvars = Flatten[vars];
c1 = Map[0 <= # <= 1 &, fvars];
c2 = Map[Total[#] == 1 &, vars];
c3 = Map[0 <= Total[#] <= 1 &, Transpose@vars];
vals = vars.Range[n - 1];
c4 = {n == Total[vals[[1 ;; 3]]], n == Total[vals[[4 ;; 6]]], 
   n == Total[vals[[7 ;; 9]]], n == Total[vals[[1 ;; -1 ;; 3]]], 
   n == Total[vals[[2 ;; -1 ;; 3]]], 
   n == Total[vals[[3 ;; -1 ;; 3]]],
   n == Total[vals[[1 ;; -1 ;; 4]]], 
   n == Total[vals[[3 ;; 7 ;; 2]]]};
constraints = Join[c1, c2, c3, c4];

Now solve the system.

Timing[solns = Solve[constraints, fvars, Integers];]
Length[solns]

(* Out[198]= {0.492, Null}

Out[199]= 8 *)

As for dealing with symmetries, I do not have a good way to figure out what they are. If I did, however, I would impose ordering on all variables that are in an equivalence class by using a chain of <- inequalities on the corresponding components of vals.

Answer 2

You get a quick solution, if you apply a condition for quadratic deviation from mean value of allowed numbers.

Only numbers 1 to 9 are allowed, because the sum of all rows or collumns may not exceed 45. Quadratic deviation from mean value 5 is 60.

s1 = First@Select[Subsets[Range[15], {9}], Total[#] <= 45 &]

(*   {1, 2, 3, 4, 5, 6, 7, 8, 9}   *)

tot = Total[(# - 5)^2 & /@ Range[9]]

(*   60   *)

There are a lot of combinations with tot == 60, but this condition plus the row sums, column sums and diagonal sum are sufficient to get the desired solution.

ta = Total[(# - 5)^2 & /@ {a1, a2, a3, a4, a5, a6, a7, a8, a9}]

n = 15; sol = {{a1, a2, a3}, {a4, a5, a6}, {a7, a8, a9}} /. 
   Solve[{ta == 60, n == a1 + a2 + a3, n == a4 + a5 + a6, 
          n == a7 + a8 + a9, n == a1 + a4 + a7, n == a2 + a5 + a8, 
          n == a3 + a6 + a9, n == a1 + a5 + a9, n == a3 + a5 + a7}, 
   {a1, a2,a3, a4, a5, a6, a7, a8, a9}, Integers]

(*   {{{2, 7, 6}, {9, 5, 1}, {4, 3, 8}}, 
      {{2, 9, 4}, {7, 5, 3}, {6, 1, 8}}, 
      {{4, 3, 8}, {9, 5, 1}, {2, 7, 6}}, 
      {{4, 9, 2}, {3, 5, 7}, {8, 1, 6}}, 
      {{6, 1, 8}, {7, 5, 3}, {2, 9, 4}}, 
      {{6, 7, 2}, {1, 5, 9}, {8, 3, 4}}, 
      {{8, 1, 6}, {3, 5, 7}, {4, 9, 2}}, 
      {{8, 3, 4}, {1, 5, 9}, {6, 7, 2}}}   *)