# Solving a system of equations using a list of possible solutions

Well, I have for example the following input and output of a code:

In[1]:=Clear["Global*"];
a = 3;
b = 2;
n = 21609;
Select[PowersRepresentations[n, a, b],
DuplicateFreeQ[#] && ! MemberQ[#, 0] &]

Out[1]={{1, 38, 142}, {1, 82, 122}, {2, 17, 146}, {2, 34, 143}, {2, 74,
127}, {2, 94, 113}, {10, 22, 145}, {10, 97, 110}, {12, 27,
144}, {12, 99, 108}, {17, 34, 142}, {17, 58, 134}, {17, 86,
118}, {22, 26, 143}, {22, 31, 142}, {22, 65, 130}, {22, 79,
122}, {22, 95, 110}, {27, 96, 108}, {28, 35, 140}, {28, 56,
133}, {28, 91, 112}, {31, 82, 118}, {36, 72, 123}, {36, 93,
108}, {38, 47, 134}, {38, 79, 118}, {38, 86, 113}, {42, 63,
126}, {45, 72, 120}, {46, 58, 127}, {46, 82, 113}, {47, 50,
130}, {47, 74, 118}, {50, 65, 122}, {56, 77, 112}, {58, 74,
113}, {58, 94, 97}, {69, 72, 108}, {74, 82, 97}}


Now, how can I find solutions to the following system of equations:

Solve[{a1^2 + a2^2 + a3^2 == n, a4^2 + a5^2 + a6^2 == n,
a7^2 + a8^2 + a9^2 == n, a1^2 + a4^2 + a7^2 == n,
a2^2 + a5^2 + a8^2 == n, a3^2 + a6^2 + a9^2 == n}, {a1, a2, a3, a4,
a5, a6, a7, a8, a9}]


Using the possible solutions given in the output of the code from above? I think that I need to check every set of possible solutions but I do not see how I can code that?

Thanks for any help and advice.

• n = 21609; Solve[{a1^2 + a2^2 + a3^2 == n, a4^2 + a5^2 + a6^2 == n, a7^2 + a8^2 + a9^2 == n, a1^2 + a4^2 + a7^2 == n, a2^2 + a5^2 + a8^2 == n, a3^2 + a6^2 + a9^2 == n}, {a1, a2, a3, a4, a5, a6, a7, a8, a9}, PositiveIntegers] produces 1278 solutions in short time. Dec 8, 2021 at 17:16
• I want to use the possible ones that are produced by the first code. Dec 8, 2021 at 17:18
• Sorry, don't understand you. Could you elaborate your comment? Dec 8, 2021 at 17:23
• @user64494 Sorry, English is not my native language so sorry for the misunderstanding. Maybe you can check out the answer given below and the discussion below the answer. I think that my question becomes clear after reading that. Thanks for the time and effort. Dec 8, 2021 at 17:25
• Sorry, I don't see any answer below. I repeat could you elaborate your comment? Dec 8, 2021 at 17:33

Since exploring all the combinations of triples and their inner ordering is very slow, this method uses randomization and ResourceFunction["BacktrackSearch"] instead to come up with solutions from your set, where the inner order of the triples is randomized. Change the SeedRandom to get other solutions.

Remove["Global*"];

a = 3;
b = 2;
n = 21609;
triples =
Select[PowersRepresentations[n, a, b],
DuplicateFreeQ[#] && ! MemberQ[#, 0] &];
system = {a1^2 + a2^2 + a3^2 == n, a4^2 + a5^2 + a6^2 == n,
a7^2 + a8^2 + a9^2 == n, a1^2 + a4^2 + a7^2 == n,
a2^2 + a5^2 + a8^2 == n, a3^2 + a6^2 + a9^2 == n};

vars = {{a1, a2, a3}, {a4, a5, a6}, {a7, a8, a9}};
(** assigns all the trips encountered so far -
this will either reduce to False,
or an expression with outstanding variables,
or True if all variables are assigned **)
test[trips_] := With[{expr = And @@ system},
expr /.
]

SeedRandom[123];
maxsols = 5;
ResourceFunction["BacktrackSearch"][
(* randomize the order of all the triples and make three choice sets for each block of assignments *)
Map[RandomSample, Table[triples, {3}], {2}],
(* backtrack if test is False i.e some assignment definitely doesn't work *)
test[#] =!= False &,
(* solution is correct if test is exactly True *)
test[#] === True &, maxsols
]

(** results:
{
{{127, 74, 2}, {58, 97, 94}, {46, 82, 113}},
{{2, 113, 94}, {74, 82, 97}, {127, 46, 58}},
{{12, 99, 108}, {27, 108, 96}, {144, 12, 27}},
{{22, 95, 110}, {95, 110, 22}, {110, 22, 95}},
{{118, 82, 31}, {82, 31, 118}, {31, 118, 82}}}
**)


• I verified that results are all in the triples by doing Map[MemberQ[triples, Sort[#]] &, results, {2}] Dec 8, 2021 at 17:47