/* 

  Prouhe-Tarry-Escott problem in Picat.

  https://en.wikipedia.org/wiki/Prouhet%E2%80%93Tarry%E2%80%93Escott_problem
  """
  In mathematics, the Prouhet–Tarry–Escott problem asks for two disjoint sets 
  A and B of n integers each, such that:

    \sum_{a\in A} a^i = \sum_{b\in B} b^i

  for each integer i from 1 to a given k.

  This problem was named after Eugène Prouhet, who studied it in the early 1850s, and 
  Gaston Tarry and Escott, who studied it in the early 1910s. The problem originates 
  from letters of Christian Goldbach and Leonhard Euler (1750/1751).

  ...

  Examples:
  It has been shown that n must be strictly greater than k. The largest value of k for 
  which a solution with n = k+1 is known is given by A = {±22, ±61, ±86, ±127, ±140, ±151}, 
  B = {±35, ±47, ±94, ±121, ±146, ±148} for which k = 11.

  An example for n = 6 and k = 5 is given by the two sets { 0, 5, 6, 16, 17, 22 } and 
  { 1, 2, 10, 12, 20, 21 }, because:

    0^1 + 5^1 + 6^1 + 16^1 + 17^1 + 22^1 = 1^1 + 2^1 + 10^1 + 12^1 + 20^1 + 21^1

    0^2 + 5^2 + 6^2 + 16^2 + 17^2 + 22^2 = 1^2 + 2^2 + 10^2 + 12^2 + 20^2 + 21^2

    0^3 + 5^3 + 6^3 + 16^3 + 17^3 + 22^3 = 1^3 + 2^3 + 10^3 + 12^3 + 20^3 + 21^3

    0^4 + 5^4 + 6^4 + 16^4 + 17^4 + 22^4 = 1^4 + 2^4 + 10^4 + 12^4 + 20^4 + 21^4

    0^5 + 5^5 + 6^5 + 16^5 + 17^5 + 22^5 = 1^5 + 2^5 + 10^5 + 12^5 + 20^5 + 21^5. 
  """


  [n = 5,k = 4,limit = 18]: 2.7s
  a = [0,4,8,16,17]
  b = [1,2,10,14,18]


  N = 6, K=4:  (0.6s)
  a = [0,5,6,7,13,14]
  b = [1,2,8,9,10,15]


  [n = 6,k = 5,limit = 22]: 85.78s
  a = [0,5,6,16,17,22]
  b = [1,2,10,12,20,21]


  This Picat model was created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Picat page: http://www.hakank.org/picat/

*/


import cp.

main => go.

go =>

  time(test(5,4)),
  nl,
  time(test(6,4)),
  nl,
  time(test(6,5)),
  nl,
  nl.


test(N,K) => 
  println($test(N,K)),
  member(Limit,1..22),
  % println(limit=Limit),
  % println([n=N,k=K,limit=Limit]),

  p(N,K,Limit, A, B),
  println(limit=Limit),
  println(a=A),
  println(b=B),
  nl.


p(N,K,Limit, A, B) => 
  A = new_list(N),
  A :: 0..Limit,
  B = new_list(N),
  B :: 0..Limit,

  % A = [0,5,6,16,17,22],
  % B = [1,2,10,12,20,21],

  Vars = A ++ B,
  all_different(Vars),
  foreach(I in 1..K)
    sum([A[J]**I : J in 1..N]) #= sum([B[J]**I : J in 1..N])
  end,

  lex_lt(A,B),
  increasing(A),
  increasing(B),

  solve($[ffc,split], Vars).