/* 

  Consecutive integer partition problem in Picat.

  From Tadao Kitzawa "Arithmetical, Geometrical and Combinatorial Puzzles from Japan", 
  problem 2, page 1
  """
  Problem 2
  The consecutive integers from 1 to n are partitioned into several groups such
  that in each group, the largest number is equal to the sum of the remaining
  numbers.
  (a) How many groups are there if n = 12?
  (b) What are the numbers of groups that may work if n = 14?
  """

  Also see http://hakank.org/cpmpy/consecutive_integer_partition_problem.py

  Note: Here I don't assume anything about the sizes of the groups, except
  that a group must have at least three values, otherwise there can be no 
  maximum value and remainder numbers that adds to the maximum number. 
  There can, for example, be groups with different number of elements:
  for n = 12 all groups are of size 3 but for n=16 there are groups of 
  size 3 and of size 4.


  a) For n = 12, the number of groups are 4
     Without symmetry breaking, there are 192 different solutions.
     With symmetry breaking (see the model), there are only these 3 solutions:

     m = 4
     maxes = [11,7,9,12]
     [1,2,3,4,2,3,2,4,3,1,1,4]
     1 = [1,10,11]
     2 = [2,5,7]
     3 = [3,6,9]
     4 = [4,8,12]

     m = 4
     maxes = [6,11,10,12]
     [1,2,3,4,1,1,3,4,2,3,2,4]
     1 = [1,5,6]
     2 = [2,9,11]
     3 = [3,7,10]
     4 = [4,8,12]

     m = 4
     maxes = [10,6,11,12]
     [1,2,3,2,4,2,4,3,1,1,3,4]
     1 = [1,9,10]
     2 = [2,4,6]
     3 = [3,8,11]
     4 = [5,7,12]


  b) For n=14 there is no solution. (Why?)


  Here are some explorations.

  * go2/0 
    A little more investigation give the following number of solutions for N = 2..20
    (using symmetry breaking):
    2 = 0  
    3 = 1 (w/o symmetry breaking: 1)
    4 = 0
    5 = 0
    6 = 0
    7 = 0
    8 = 0
    9 = 0
    10 = 0
    11 = 0
    12 = 3 (w/o symmetry breaking: 192)
    13 = 0
    14 = 0
    15 = 8 (w/o symmetry breaking: 2520) 
    16 = 31 (w/o symmetry breaking: 13200) 
    17 = 0
    18 = 0
    19 = 215
    20 = 530

  * go3/0 
    Here are the number of groups for each of these successful cases, including
    some assignments. Here we test different values of the number of groups in case
    there might be more than one valid value of the number of groups, which I
    have not found.

    x = [1,1,1]
    maxes = [3]
    1 = [1,2,3]

    n = 12
    m = 4
    x = [1,2,3,4,2,3,2,4,3,1,1,4]
    maxes = [11,7,9,12]
    1 = [1,10,11]
    2 = [2,5,7]
    3 = [3,6,9]
    4 = [4,8,12]

    n = 15
    m = 5
    x = [1,2,3,4,5,5,1,1,3,4,5,3,2,4,2]
    maxes = [8,15,12,14,11]
    1 = [1,7,8]
    2 = [2,13,15]
    3 = [3,9,12]
    4 = [4,10,14]
    5 = [5,6,11]

    n = 16
    m = 5
    x = [1,2,3,4,5,2,4,2,3,5,4,3,1,1,5,2]
    maxes = [14,16,12,11,15]
    1 = [1,13,14]
    2 = [2,6,8,16]
    3 = [3,9,12]
    4 = [4,7,11]
    5 = [5,10,15]

    n = 19
    m = 6
    x = [1,1,2,4,3,6,5,4,6,5,1,4,3,1,6,2,5,3,2]
    maxes = [14,19,18,12,17,15]
    1 = [1,2,11,14]
    2 = [3,16,19]
    3 = [5,13,18]
    4 = [4,8,12]
    5 = [7,10,17]
    6 = [6,9,15]

    n = 20
    m = 6
    x = [1,1,1,2,5,4,6,3,3,4,2,6,5,1,2,4,3,5,6,1]
    maxes = [20,15,17,16,18,19]
    1 = [1,2,3,14,20]
    2 = [4,11,15]
    3 = [8,9,17]
    4 = [6,10,16]
    5 = [5,13,18]
    6 = [7,12,19]

    n = 23
    m = 7
    x = [1,1,3,3,2,5,4,6,3,4,7,7,5,6,2,3,4,1,5,2,1,6,7]
    maxes = [21,20,16,17,19,22,23]
    1 = [1,2,18,21]
    2 = [5,15,20]
    3 = [3,4,9,16]
    4 = [7,10,17]
    5 = [6,13,19]
    6 = [8,14,22]
    7 = [11,12,23]

    n = 24
    m = 8
    x = [1,2,3,4,1,1,7,8,6,5,8,6,5,3,7,2,3,2,8,4,6,7,5,4]
    maxes = [6,18,17,24,23,21,22,19]
    1 = [1,5,6]
    2 = [2,16,18]
    3 = [3,14,17]
    4 = [4,20,24]
    5 = [10,13,23]
    6 = [9,12,21]
    7 = [7,15,22]
    8 = [8,11,19]

    n = 27
    m = 9
    x = [1,2,3,4,4,5,6,7,4,6,8,9,8,2,9,2,6,7,5,3,1,1,3,8,5,7,9]
    maxes = [22,16,23,9,25,17,26,24,27]
    1 = [1,21,22]
    2 = [2,14,16]
    3 = [3,20,23]
    4 = [4,5,9]
    5 = [6,19,25]
    6 = [7,10,17]
    7 = [8,18,26]
    8 = [11,13,24]
    9 = [12,15,27]

    n = 28
    m = 9
    x = [1,2,3,4,5,6,6,7,7,6,8,2,9,2,9,8,7,1,1,4,5,3,6,4,3,5,8,9]
    maxes = [19,14,25,24,26,23,17,27,28]
    1 = [1,18,19]
    2 = [2,12,14]
    3 = [3,22,25]
    4 = [4,20,24]
    5 = [5,21,26]
    6 = [6,7,10,23]
    7 = [8,9,17]
    8 = [11,16,27]
    9 = [13,15,28]

    n = 31
    m = 10
    x = [1,2,3,4,5,6,7,1,8,7,9,2,8,2,10,10,7,9,1,4,6,8,3,4,5,3,6,1,9,5,10]
    maxes = [28,14,26,24,30,27,17,22,29,31]
    1 = [1,8,19,28]
    2 = [2,12,14]
    3 = [3,23,26]
    4 = [4,20,24]
    5 = [5,25,30]
    6 = [6,21,27]
    7 = [7,10,17]
    8 = [9,13,22]
    9 = [11,18,29]
    10 = [15,16,31]


  Some findings:

  - It seems that the number of groups is max([1,n div 3]) for those n's that has a solution.
    But which n has a solution?

  - By running the model, I've found the following instances to be solvable:

    3,12,15,16,19,20,23,24,27,28,31,32

    What is the closed formula for this?

    Using Mathematica's FindSequenceFunction:
     ff = FindSequenceFunction[{3, 12, 15, 16, 19, 20, 23, 24, 27, 28, 31, 32}, x]

     This returns a DifferenceRoot  (-19-4 N+Y[N]+Y[1+N]==0, Y[1]==3, Y[2]==12)
     Rewritten as
        Y[1+N]==4*N+Y[N]-19, Y[1]==3, Y[2]==12

     ff  /. x -> Range[40]

     give these values:
     3, 12, 15, 16, 19, 20, 23, 24, 27, 28, 31, 32, 35, 36, 39, 40, 43, 44, 47, 48,
     51, 52, 55, 56, 59, 60, 63, 64, 67, 68, 71, 72, 75, 76, 79, 80, 83, 84, 87, 88

     See go4/0 for running these cases and with m = max(1,n div 3)   

     If we skip the first two instances (3,12) there is a simpler pattern:
       s = {15, 16, 19, 20, 23, 24, 27, 28, 31, 32}
       tt = FindSequenceFunction[s,n]
         1/2 (-1)^n (-1 + 25 (-1)^n + 4 (-1)^n n)

       tt /. n -> Range[4]

    X = [15,16,19,20,23,24,27,28,31,32,35,36,39,40,43,44,47,48,51,52]

    Here is a simpler version: checking modulo 4 == 0 or 3, i.e. for n=12.. 

    Picat> X= [I : I in 1..40, member(I mod 4,[0,3])]
    X = [3,4,7,8,11,12,15,16,19,20,23,24,27,28,31,32,35,36,39,40]

    Above 12 this seems to work, and below 12 there are no solutions for 
    4, 7, 8, 11. Why? 
    There must be at least three elements in each group but we cannot create
    groups for N = 4 (and 1+2+3 != 4)

    But why is there no solution of N=7,8,and 11?
    For 7 there _could_ be one group of 3 and one group of 4. But there isn't.
    For 8 there _could_ be one group of 4 and one group of 4. But there isn't.
    For 11 there _could_ be one group of 3, and two groups of 4. But there isn't.

    What 3 and 4 numbers satisfies the max + sum constraint for these N's?
    See go5/0 for some explorations on this.
    Well, it seems to be that one simply cannot find any disjoint max-sum groups 
    of size >= 3.
    
  * There are solutions with groups of size 5, for example for N=20

      n = 20
      m = 6
      maxes = [15,17,18,16,19,20]
      [1,2,2,2,3,4,5,2,6,4,6,5,3,1,1,4,2,3,5,6]
      1 = [1,14,15]
      2 = [2,3,4,8,17]
      3 = [5,13,18]
      4 = [6,10,16]
      5 = [7,12,19]
      6 = [9,11,20]

    Is it possible with groups of size 6? If not, why?


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

*/

import sat. % faster
% import cp.
% import smt.

main => go.

go ?=>
  nolog,
  N = 12, % OK
  % N = 14, % No solutions
  println(n=N),
  % Here we check all possible group sizes (M)
  % In other explorations we assume that the group size is N div 3
  % M = N div 3,
  member(M,1..N),
  println(m=M),
  SymmetryBreaking = true,
  consecutive_integers(N,SymmetryBreaking, M,X,Maxes),
  println(m=M),
  println(maxes=Maxes),
  println(X),
  foreach(G in 1..M)
    println(G=[J : J in 1..N, X[J]==G])
  end,
  nl,
  flush(stdout),
  fail,
  
  nl.
  
go => true.

%
% Number of solutions (with symmetry breaking)
%
go2 ?=>
  nolog,
  SymmetryBreaking = true,
  member(N,2..50),
  % foreach(N in 2..50)
  M = max([1,N div 3]),
  Count = count_all(consecutive_integers(N,SymmetryBreaking,M,_X,_Maxes)),
  println(N=Count),
  fail,
  % end,
  nl.
go2 => true.

%
% First solution for N = 2..50
%
go3 ?=>
  nolog,
  SymmetryBreaking = true,  
  % Here we also test all possible number of groups (M)
  foreach(N in 2..50)
    foreach(M in 1..1+(N div 3)) 
      if consecutive_integers(N,SymmetryBreaking,M,X,Maxes) then
        println(n=N),
        println(m=M),
        println(x=X),
        println(maxes=Maxes),
        foreach(G in 1..M)
          println(G=[J : J in 1..N, X[J]==G])
        end,
        nl
      end
    end
  end,
  nl.
go3 => true.


%
% Solutions for the instances that Mathematica's FindSequenceFunction generated
% (see above) and only for m = max(1,n div 3).
%
go4 ?=>
  nolog,
  SymmetryBreaking = true,
  L = [3, 12, 15, 16, 19, 20, 23, 24, 27, 28, 31, 32, 35, 36, 39, 40, 43, 44, 47, 48,
      51, 52, 55, 56, 59, 60, 63, 64, 67, 68, 71, 72, 75, 76, 79, 80, 83, 84, 87, 88],
  % L = [closed_formula(I) : I in 1..20],
  println(L),
  % NoSolutions = [],
  % foreach(N in L)
  member(N,L),
  garbage_collect,
  M = max([1,N div 3]),
  println([n=N,m=M]),
  flush(stdout),
  if time(consecutive_integers(N,SymmetryBreaking,M,X,Maxes)) then
    println(n=N),
    println(m=M),
    println(x=X),
    println(maxes=Maxes),
    foreach(G in 1..M)
      println(G=[J : J in 1..N, X[J]==G])
    end,
    nl
  else
    println("No solution!")
  end,
  nl,
  flush(stdout),
  fail,
  nl.
go4 => true.


%
% Which 3 and 4 number combinations satisfies the
% main constraint that:
%    Max = max(X)
%    Max = sum(X) - Max
% For N = 7 these are the only possible max-sum groups >= 3
% [3,4,7]
% [1,6,7]
% [1,4,5]
% [1,5,6]
% [1,2,3]
% [1,3,4]
% [2,4,6]
% [2,5,7]
% [2,3,5]
% [1,2,4,7]
% [1,2,3,6]
% and we cannot create two disjoint groups with these.
%
go5 ?=>
  nolog,
  N=7,
  member(M,3..7),
  X = new_list(M),
  X :: 1..N,

  increasing_strict(X),
  Max #= X[M],
  sum(X[1..M-1]) #= Max,

  solve(X),
  println(X),
  fail,
 
  nl.
go5=>true.


%
% Solve the consecutive_integers partition problem
%
consecutive_integers(N,SymmetryBreaking, M,X,Maxes) =>
  if var(M) then
    member(M,1..1+(N div 3)) % number of groups
  end,

  X = new_list(N),
  X :: 1..M,

  Maxes = new_list(M),
  Maxes :: M..N,

  foreach(G in 1..M)
    count(G,X,#>=,3), % The group size must be at least 3 elements
    Maxes[G] #= max([J*(X[J] #= G) : J in 1..N]),
    Maxes[G] #= (sum([J*(X[J] #= G) : J in 1..N]) - Maxes[G])
  end,

  if SymmetryBreaking then
    % Symmetry breaking: assign new group ids in order
    % foreach(G in 1..M)
    %   X[G] #<= G
    % end,
    value_precede_chain(1..M,X), % Faster

    % this is slower than value_precede_chain/2 and don't
    % remove as many symmetries.
    % increasing(Maxes), 

    % N must be one of the maximum numbers, and we break
    % symmetry by decide it to be the last (M'th) group
    Maxes[M] #= N
  end,

  Vars = X ++ Maxes,
  % Vars = Maxes ++ X,
  solve($[degree,updown],Vars).

%
%  The global constraint
%      value_precede_chain(C, X)
%  ensures that the value C[I-1] precedes the value C[I] is the array X
%  if both C[I-1] and C[I] are in X.
%
value_precede_chain(C, X) =>
  foreach(I in 2..C.length)
    value_precede(C[I-1], C[I], X)
  end.

% This definition is inspired by 
% MiniZinc definition value_precede_int.mzn
value_precede(S,T,X) =>
   XLen = X.length,
   B = new_list(XLen+1),
   B :: 0..1,
   foreach(I in 1..XLen)
     % Xis :: 0..1,
     Xis #= (X[I] #= S),
     (Xis #=> (B[I+1] #= 1))
     #/\ ((#~ Xis #= 1) #=> (B[I] #= B[I+1]))
     #/\ ((#~ B[I] #= 1) #=> (X[I] #!= T))
   end,
   B[1] #= 0.

% 
% This is the closed formula function given by Mathematica's FindSequenceFunction
% when excluding n = 3 and 12, i.e. 15, 16, 19, 20, 23, 24, 27, 28, 31, 32.
% See more discussions above.
% 
% Picat> X = [closed_formula(I) : I in 1..20]
% X = [15,16,19,20,23,24,27,28,31,32,35,36,39,40,43,44,47,48,51,52]
%
% This happens to be the same as checking for modulo 4:
% Picat> X = [I : I in 12..50, member(I mod 4,[0,3])]
% X = [12,15,16,19,20,23,24,27,28,31,32,35,36,39,40,43,44,47,48]
%
closed_formula(N) = (-1)**N * (-1 + 25*(-1)**N + 4*(-1)**N*N) div 2.