utils_me.pi

/*

  Some general utilities in Picat.

  Note: Some of these are experimental.

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

*/

module utils_me.

import util.
import cp.

% main => go.

test => 
   list_test,
   test3.

get_n_solutions(Goal, N) =>
  printf("Get %d solutions:\n", N),
  get_global_map().put(solcount,1), 
  time2(Goal),
  C = get_global_map().get(solcount),
  if C < N then get_global_map().put(solcount,C+1), fail end.


%
% Calls get_n_solutions/3 and fetch the
% solutions (in get_global_map().get(sols))
%
get_n_solutions(Goal,N,Arg,Sols) =>
  M = get_global_map(), 
  M.put(sols,[]), % reset the solution list
  get_n_solutions(Goal,N,Arg),
  Sols = M.get(sols).

%
% get_n_solutions(Goal, N, Arg)
% get max N solutions for goal Goal and put the Arg'th element
% in the global map list (sols).
% 
% Note: the call will fail after N solutions so be aware that
% it might affect other code.
%
get_n_solutions(Goal, N, Arg) =>
  printf("Get %d solutions:\n", N),
  get_global_map().put(solcount,1),
  % time2(Goal), % If we want the time for the calls
  call(Goal),
  M = get_global_map(),
  M.put(sols, M.get(sols)++[Goal[Arg]]),
  C = M.get(solcount),
  if C < N then M.put(solcount,C+1), fail end.



% Flatten a list/array and remove duplicates
term_variables(L) = Flatten =>
   Flatten1 = [],
   foreach(LL in L) 
      Flatten1 := Flatten1 ++ LL
   end,
   Flatten2 = remove_dups(Flatten1),
   Flatten = Flatten2.


next_permutation(P) = Perm =>
   Perm1 = P,
   N = Perm1.length,
   K = N - 1,
   while (Perm1[K] > Perm1[K+1], K >= 0) 
      K := K - 1
   end,
   if K > 0 then
      J = N,
      while (Perm1[K] > Perm1[J])  J := J - 1 end,      
      Tmp := Perm1[K],
      Perm1[K] := Perm1[J],
      Perm1[J] := Tmp,
      R = N, 
      S = K + 1,
      while (R > S) 
         Tmp := Perm1[R],
         Perm1[R] := Perm1[S],
         Perm1[S] := Tmp,
         R := R - 1, 
         S := S + 1
      end
   end,
   Perm = Perm1.

% From T. Van Le, "Techniques of Prolog Programming", page 100f
next_higher_permutation(L,L1) => 
   reverse(L,[],L2),
   append(A,[X,Y|B],L2), X > Y,
   append(A,[X],C),
   append(A1,[U|B1],C), U > Y,
   append(A1,[Y|B1], B2),
   reverse([U|B], B2,L1).

% From T. Van Le, "Techniques of Prolog Programming", page 99
reverse([],R1,R2) => R1=R2.
reverse([H|T],R,L1) => reverse(T,[H|R],L1).

random_perm(L,N) = Perm => 
  Perm = L,
  Len = L.length,
  _ = random2(),
  foreach(_I in 1..N) 
     R1 = 1+(random() mod Len),
     R2 = 1+(random() mod Len),
     T = Perm[R1],
     Perm[R1] := Perm[R2],
     Perm[R2] := T
  end.


%
% Do a weighted random choice based on the user's previous choices.
% Input is a map with the choice=counts, e.g.
%   PlayerMap = new_map([paper=1,rock=2,scissors=300]),
%
% (See http://hakank.org/picat/rock_paper_scissors.pi for an example)
%
weighted_choice(Map) = Choice => 
  Map2 = [(V+1)=K : K=V in Map].sort, % ensure that all choices can be made
  % Prepare the probability matrix M
  Total = sum([P : P=_ in Map2]),
  Len = Map2.len,
  M = new_array(Len,2),
  T = new_list(Len),
  foreach({I,P=C} in zip(1..Len,Map2))
    if I == 1 then
      M[I,1] := 1,
      M[I,2] := P
    else
      M[I,1] := M[I-1,2]+1,
      M[I,2] := M[I,1]+P-1
    end,
    T[I] := C
  end,
  M[Len,2] := Total,

  % Pick a random number in 1..Total
  R = 1 + random() mod Total,
  Choice = _,
  % Check were R match
  foreach(I in 1..Len, var(Choice))
    if M[I,1] <= R, M[I,2] >= R then
      Choice := T[I]
    end
  end.


%sublist(String1, String2):
%   succeeds if S2 is a sublist of S
% (nondet.)
sublist(S1,S2) =>
    append(_, S, S1), append(S2, _, S).



% all the substrings of string List
all_substrings(List) = L, list(List) =>
     L = findall(S2, $sublist(List,S2)).remove_dups().


%
% Some list utils etc.
%

% drop(Xs,0) = Xs.
% drop([],_N) = [].
% drop([_X|Xs],N) = drop(Xs,N-1).

drop_while([],_P) = [].
drop_while([X|Xs],P) = cond(call(P,X), drop_while(Xs,P), [X] ++ Xs).


% L=take(1..10, 4)
% take(_Xs,0) = [].
% take([],_N) = [].
% take([X|Xs],N) = [X] ++ take(Xs,N-1).

% take_while
% L=take_while([2,4,6,7,8,9,16],even) -> [2,4,6]
take_while([],_P) = [].
take_while([X|Xs],P) = cond(call(P,X),
                           [X] ++ take_while(Xs,P), 
                           []).


butfirst(L) = drop(L,1).
butlast(L) = take(L,L.length-1).

butfirst2(S) = [S[I] : I in 2..S.length].
butlast2(S) = [S[I] : I in 1..S.length-1].

% sort/3 using B-Prolog's sort/3
sort_order(Order, L1) = L2 => bp.sort(Order,L1,L2).

%
% sort_map(Map,keys)  : sort Map on keys
% sort_map(Map,values): sort Map on values
%
sort_map(Map,values) = [K=V:_=(K=V) in sort([V=(K=V): K=V in Map])].
sort_map(Map,keys)   = sort([K=V:K=V in Map]).
sort_map(Map) = sort_map(Map,keys).


swap(A,B) = [B,A].

% sorted(L) =>
%   foreach(I in 2..L.length) L[I-1] =< L[I] end.


% sort list L using position Pos in the list of lists
sort_by(L,Pos) = [L[I] : _=I in sort(Perm)] =>
  Perm = [L[I,Pos]=I : I in 1..L.length].

sort_down_by(L,Pos) = [L[I] : _=I in sort_down(Perm)] =>
  Perm = [L[I,Pos]=I : I in 1..L.length].

% get the max value in position Pos
% max/1 only works for numbers
maximal_by(L,Pos) = L[argmax([L[I,Pos] : I in 1..L.length]).first()], integer(L[1,Pos]) => true.
maximal_by(L,Pos) = sort_by(L,Pos).last().

% get the min value in position Pos
% min/1 only works for numbers
minimal_by(L,Pos) = L[argmin([L[I,Pos] : I in 1..L.length]).first()], integer(L[1,Pos]) => true.
minimal_by(L,Pos) = sort_by(L,Pos).first().


%
% The scan family.
% Note that the list is the first argument so we can chain, e.g. 
%    Picat> L=1..4, L2=L.scanl1(+).scanl1(+).scanl1(+).reduce(+) 
%    L = [1,2,3,4]
%    L2 = 56
%   
%
% scanl(1..10,+,0)   [0,1,3,6,10,15,21,28,36,45,55]
% scanl(1..10,*,1)   [1,1,2,6,24,120,720,5040,40320,362880,3628800]
% scanl(1..4,/,1)    [1,1.0,0.5,0.166666666666667,0.041666666666667]
% scanl([[1],[2],[3]],++,[])    [[],[1],[1,2],[1,2,3]]
% scanl1(1..10,+)    [1,3,6,10,15,21,28,36,45,55]
% scanl1(1..10,*)    [1,2,6,24,120,720,5040,40320,362880,3628800]
% scanl1(2..5,**)    [2,8,4096,1152921504606846976]   
% scanr(1..10,+,0)   [55,54,52,49,45,40,34,27,19,10,0]
% scanr(1..10,*,1)   [3628800,3628800,1814400,604800,151200,30240,5040,720,90,10,1]
% scanr([[1],[2],[3]],++,[])   [[1,2,3],[2,3],[3],[]]
% scanr1(1..10,+)    [55,54,52,49,45,40,34,27,19,10]
% scanr1(1..10,*)    [3628800,3628800,1814400,604800,151200,30240,5040,720,90,10]
%

% scanl
scanl([],_F, E) = [E].
scanl(Xs,F,E) = [E] ++ scanlT(Xs,F,E).
scanlT([],_F,_A) = [].
scanlT([Y|Ys],F,A) = scanl(Ys, F,apply(F,A,Y)).

% scanl1
scanl1([],_F) = _ => throw $error(empty_list,scanl1,[]).
scanl1([X|Xs],F) = scanl(Xs,F,X).

% scanr
scanr([],_F,E) = [E].
scanr([X|Xs],F,E) = [apply(F,X,head(Ys))] ++ Ys =>
                     Ys = scanr(Xs,F,E).

% scanr1
scanr1([],_F) = _ => throw $error(empty_list,scanr1,[]).
scanr1([X],_F) = [X].
scanr1([X,Y|Xs],F) = [apply(F,X,head(Zs))] ++ Zs =>
                     Zs = scanr1([Y|Xs],F).


% L=span([2,3,5,6,7,8,9],prime) -> [[2,3,5],[6,7,8,9]]
span([],_P) = [[],[]].
span([X|Xs],P) = cond(call(P,X),
                      [[X|Ys],Zs],
                      [[],[X|Xs]]) =>
                 [Ys,Zs] = span(Xs,P).

%
% L=split_at(1..10,3) -> [[1,2],[3,4,5,6,7,8,9,10]]
% Note: 1-based
split_at(Xs,1) = [[],Xs].
split_at([],_N) = [[],[]].
split_at([X|Xs],N) = [[X|Ys],Zs] => [Ys,Zs] = split_at(Xs,N-1).


%
% group(List, F) 
% groups together elements in List  according to function F.
% Note: F must be a defined function.
%
% Example:
%   Picat> L=1..10, G=L.group(mod3)      
%   L = [1,2,3,4,5,6,7,8,9,10]
%   G = (map)[0=[3,6,9],1=[1,4,7,10],2=[2,5,8]]
%
group(List, F) = P, list(List) =>
   P = new_map(),
   foreach(E in List) 
      V = apply(F,E),
      P.put(V, P.get(V,[]) ++ [E])
   end.

%
% group/2 with function at first argument.
%
% Picat>  L=1..10, F=mod3,  G=F.group(L)
% L = [1,2,3,4,5,6,7,8,9,10]
% F = mod3
% G = (map)[0=[3,6,9],1=[1,4,7,10],2=[2,5,8]]
% 
group(F,List) = P, list(List) =>
   P = new_map(),
   foreach(E in List) 
      V = apply(F,E),
      P.put(V, P.get(V,[]) ++ [E])
   end.


% example function for group/2:
mod3(I) = I mod 3.   

% partition/3
% Example:
%   partition(1..10, prime)   -> [[2,3,5,7],[1,4,6,8,9,10]]
partition([],_P) = [[],[]].
partition([X|Xs],P) = cond(call(P,X),
                            [[X|Ys],Zs],
                            [Ys,[X|Zs]]) =>
                       [Ys,Zs] = partition(Xs,P).

% partition/4
% Example:
%   partition(1..10, prime,Primes, NoPrimes)
%    Primes = [2,3,5,7]
%    NoPrimes = [1,4,6,8,9,10] ?;
% (nondet)
%
partition([],_, LHS,RHS) => 
  LHS = [],
  RHS = [].

partition([Head|Tail], Pred, HeadLHS, RHS) ?=>
  HeadLHS = [Head|LHS],
  call(Pred,Head),
  partition(Tail, Pred,LHS, RHS).

partition([Head|Tail], Pred, LHS, HeadRHS) =>
  HeadRHS = [Head|RHS],
  not call(Pred,Head),
  partition(Tail, Pred,LHS, RHS).





% More utils



% Inspired by http://www.sci.brooklyn.cuny.edu/~zhou/comp_lang_prolog.pl
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Cartesian product of lists Options 
% Posted July 1, 2011
cartprod(L,R) =>
    L = [L1,L2],    
    R = [[X1,X2] : X1 in L1, X2 in L2].
cartprod(L1,R) =>
    L1 = [L|Ls],
    cartprod(Ls,R1), 
    R = [[X|P] : X in L, P in R1].

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Inspired by http://www.sci.brooklyn.cuny.edu/~zhou/comp_lang_prolog.pl
% All possible pairs
% Posted Nov. 25, 2010
all_pairs(List,Pairs) =>
    Pairs = [[I=J] : I in List, J in List, I < J].

%
% From http://www.probp.com/publib/listut.html
% 
%   correspond(X, Xlist, Ylist, Y)
%   is true when Xlist and Ylist are lists, X is an element of Xlist, Y is
%   an element of Ylist, and X and Y are in similar places in their lists.
%
correspond(X, X1, Y1, Y) ?=> X1=[X|_], Y1=[Y|_].
correspond(X, X1, Y1, Y) ?=>
        X1 = [_|T],
        Y1 = [_|U],
	correspond(X, T, U, Y).


% Inspired from SWI-Prolog's lib/lists.pl
%%	prefix(?Part, ?Whole)
%
%	True iff Part is a leading substring of Whole.  This is the same
%	as append(Part, _, Whole).
prefix([],X) ?=> X = _.
prefix([E|T0], [E|T]) =>
	prefix(T0, T).


%
% Return a random element from List
%
random_list(List) = List[random(1,List.length)].

% Get one element of list L
oneof(L) = L[random(1,L.length)].


%
% Shuffles the list List.
% 
shuffle(List) = List2 => 
  List2 = List,
  Len = List.length,
  _ = random2(),
  foreach(I in 1..Len) 
    R2 = random(1,Len),
    List2 := swap(List2,I,R2)
  end.

% pick an element from the list List
choice(List) = List[choice(List.length)], list(List) => true.

% pick a number from 1..N
choice(N) = random2(1,N), integer(N) => true.

% pick a number from 0..N-1
choice0(N) = random(0,N-1).



statistics_all_map() = new_map([ K=V : [K,V] in statistics_all()]).

%
% General functions for the simple integer operations:
% 
%  plus(A,B)
%  minus(A,B)
%  multiply(A,B)
%  divide(A,B)
%  power(A,B)
% 
% where the parameters can be
%   integer op integer
%   integer op list
%   list op integer
%   list op list
%
% In the latter case, if the two lists are of unequal length we pad the shorter
% to left with a default value, see pad/3 and pad_default/2.
% 
% Also, we define self operators
% 
%  plus(A) => plus(A,A)
%  minus(A) => minus(A,A)
%  divide(A) => divide(A,A)
%  plus(A) => plus(A,A)
%  power(A) => power(A,A)
%
arithm(Op,Int1,Int2) = apply(Op, Int1,Int2), integer(Int1), integer(Int2) => true.
arithm(Op,Int1,List2) = [apply(Op,Int1,E2) : E2 in List2], integer(Int1), list(List2) => true.
arithm(Op,List1,Int2) = [apply(Op,E1,Int2): E1 in List1], list(List1), integer(Int2)  => true.
arithm(Op,List1,List2) = [apply(Op,E1,E2): {E1,E2} in zip(List1,List2)], list(List1), list(List2), length(List1) == length(List2) => true.
% lists of unequal length
arithm(Op,List1,List2) = Res, list(List1), list(List2), length(List1) != length(List2) =>
   MaxLen = max(List1.length,List2.length),
   List1b = pad(List1,MaxLen,Op),
   List2b = pad(List2,MaxLen,Op),
   Res = [apply(Op, E1,E2): {E1,E2} in zip(List1b,List2b)].

% operate on self
arithm(Op,List) = arithm(Op,List,List), list(List) => true.

% 
% Picat> L=plus(2,3).plus(4).multiply(4).power(4)
% L = 1679616
% 
plus(A,B) = arithm(+,A,B).
minus(A,B) = arithm(-,A,B).
multiply(A,B) = arithm(*,A,B).
divide(A,B) = arithm(/,A,B).
power(A,B) = arithm(**,A,B).
modulus(A,B) = arithm(mod,A,B).

plus(A) = arithm(+,A,A).
minus(A) = arithm(-,A,A). % quite silly...
multiply(A) = arithm(*,A,A).
divide(A) = arithm(/,A,A). % quite silly...
power(A) = arithm(**,A,A). % beware of this!
modulus(A) = arithm(modulus,A,A).

% We always pad to left
pad(List, N, Op) = Res => 
   Res = List,
   Len = List.length,
   if N > Len then
      (pad_default(Op,Default); Default = 1),
      Res := [Default : _I in 1..N-Len] ++ List
   end.

% the padding values
index(+,+)
pad_default(+,0).
pad_default(-,0).
pad_default(*,1).
pad_default(/,1).
pad_default(**,1).
pad_default(mod,1).


%
% factors of N
%
factors(N) = Factors =>
     M = N,
     Factors1 = [],
     while (M mod 2 == 0) 
         Factors1 := Factors1 ++ [2],  
         M := M div 2 
     end,
     T = 3,
     
     while (M > 1, T < 1+(sqrt(M)))
        if M mod T == 0 then
           [Divisors, NewM] = alldivisorsM(M, T),
           Factors1 := Factors1 ++ Divisors,
           M := NewM
        end,
        T := T + 2
        % next_prime(T, T2),
        % T := T2
     end,
     if M > 1 then Factors1 := Factors1 ++ [M] end,
     Factors = Factors1.

alldivisorsM(N,Div) = [Divisors,NewN] =>
   M = N,
   Divisors1 = [],
   while (M mod Div == 0) 
      Divisors1 := Divisors1 ++ [Div],
      M := M div Div
   end,
   NewN := M,
   Divisors = Divisors1.


table
prime_cached(N) => prime(N).


inc(X) = X+1.
inc(X,N) = X+N.
dec(X) = X-1.
dec(X,N) = X-N.

inc2(N1,N2) =>
  between(N1,99999999,N2).

% with step Step
inc2(N1,Step,N2) =>
  between(N1,99999999,Tmp),
  N2 = (Tmp-N1)*Step+N1.

between(From,Step,To,N) =>
  between(From,To,Tmp),
  N = (Tmp-From)*Step+From.


to_string2(Expr) = [E.to_string() : E in Expr].flatten().

% E=["4","+","3","-","2"].flatten().eval()
% E=["4","-","(","3","-","2",")"].eval()
eval(Expr) = parse_term(Expr.flatten()).apply().

% Replacing all occurrences of the string/list Old in List with New.
% Both Old and New must be lists (or strings).
replace2(List,Old,New) = Res =>
  Res = List,
  while (find(Res,Old,_,_))
    once(append(Before,Old,After,Res)),
    Res := Before ++ New ++ After
  end.

%
% Get first Len digits of Pi
%
pi(Len) = Res => 
   Ndigits = 0,
   Q = 1,
   R = 0,
   T = Q,
   K = Q,
   N = 3,
   L = N,

   Res = "",
   C = 1, % counter for presentation
   while (C <= Len)
       if 4 * Q + R - T < N * T then
           Res := Res ++ [N],
           C := C + 1,
           Ndigits := Ndigits + 1,
           NR := 10 * (R - N * T),
           N := ((10 * (3 * Q + R)) // T) - 10 * N,
           Q := Q* 10,
           R := NR
       else
           NR := (2    * Q + R) * L,
           NN := (Q * (7 * K + 2) + R * L) // (T * L),
           Q := Q * K,
           T := T* L,
           L := L + 2,
           K := K+1,
           N := NN,
           R := NR
       end
   end,
   nl.


%
% Swap position I <=> J in list L
%
swap(L,I,J) = L2, list(L) =>
  L2 = L,
  T = L2[I],
  L2[I] := L2[J],
  L2[J] := T.

num_fixpoints(List) = sum([1 :  I in 1..List.length, List[I] = I]).

%
% occurrences(List): 
% returns a map with all the keys and the 
% number of occurrences of each elements.
%
occurrences(List) = Map =>
  Map = new_map(),
  foreach(E in List)
     Map.put(E, cond(Map.has_key(E),Map.get(E)+1,1))
  end.

push(L,E) = L ++ [E].

%
% This is a bit strange, but might work...
%
% Picat> L=1..10, L:=L.pop(E1).pop(E2).pop(E3)
% L = [4,5,6,7,8,9,10]
% E1 = 1
% E2 = 2
% E3 = 3
%
% Note: This is actually a quite useful idiom since Picat don't
%       have global variables.
%
pop(L,First) = L.tail() => First = L.first().


%
% scalar produt of L1 and L2
%
scalar(L1,L2) = [L1I*L2I : {L1I,L2I} in zip(L1,L2)].sum().


table
binomial(N,K) = Res =>
   R = 1,
   foreach(I in 1..K) 
     R := floor(R * ((N-I+1)/I))
   end,
   Res = R.


% From 
% Mihaela Malita: "Library for Logic Puzzles"
% http://www.anselm.edu/internet/compsci/faculty_staff/mmalita/HOMEPAGE/logic/bibmm.txt
% 
% Note: I have changed some of the predicate names.
%

/* mem(Lr,L). Elements from Lr are all members in L. Same as subset(Lr,L).
	?- mem([X,Y],[a,b,c]).
	X = Y = a ;
	X = a ,Y = b ;
	X = a ,Y = c ;
	X = b ,Y = a ;	
	X = Y = b ;
	X = b ,Y = c ;
	X = c ,Y = a ;
	X = c ,Y = b ;
	X = Y = c ;
*/
mem([],_) => true.
mem([H|T],Y) => member(H,Y),mem(T,Y).

subset(X,Y) => mem(X,Y).

/* all(N,L,X). All possible pairs with elements from [a,b,c,d]:
	| ?- all(2,[a,b,c],I).
	I = [a,a] ;
	I = [a,b] ;
	I = [a,c] ;
	I = [b,a] ;
	I = [b,b] ;
	I = [b,c] ;
	I = [c,a] ;
	I = [c,b] ;
	I = [c,c] ;
*/
subsets(N,L,X) => X = new_list(N),mem(X,L).

/* list_all(N,L,R).  
	?- list_all(2,[a,b,c],R).
	R = [[a,a],[a,b],[a,c],[b,a],[b,b],[b,c],[c,a],[c,b],[c,c]]
*/
all_subsets(N,L) = findall(X,subsets(N,L,X)).


/* 
 set(L). Tests if a list has no duplicates. 
*/
set([]) => true.
set([H|T]) => not(member(H,T)),set(T).

no_duplicates(M) => set(M).
distinct(M) => set(M).
% aldifferent(M) =>  set(M).

/*
 Picat> N=3, R=all_subsets(N,1..N).filter2(distinct) 
 R = [[1,2,3],[1,3,2],[2,1,3],[2,3,1],[3,1,2],[3,2,1]]
*/

filter1(_P,[]) = [].
filter1(P,[X|Xs]) = cond(call(P,X), [X] ++ filter1(P,Xs), filter1(P,Xs)).

filter2([],_P) = [].
filter2([X|Xs],P) = cond(call(P,X), [X] ++ filter2(Xs,P), filter2(Xs,P)).

/* arange(N,L,R). Lists have to be sets! That is: elements do not repeat.
	?- arange(2,[a,b,c],R).    
	R = [a,b] ;
	R = [a,c] ;
	R = [b,a] ;
	R = [b,c] ;	
	R = [c,a] ;
	R = [c,b] ;
*/
arange(N,L,X) => X = new_list(N), mem(X,L),no_duplicates(X).
list_arange(N,L) = findall(X,arange(N,L,X)).



/* mem1(Lr,L). For comb/3. Same as mem/2 but does not generate [a,b] and [b,a]. 	
	?- mem1([X,Y],[a,b,c]),write([X,Y]),fail.
	[a,b][a,c][b,a][b,c][c,a][c,b]no
*/
mem1([],_Y) ?=> true.
mem1([H|T],Y) => member(H,Y),rest(H,Y,New),mem1(T,New).

/* rest(A,L,R). For comb/3. Return the rest of the list after the first occurrence of A. 
	| ?- rest(a,[a,b,c,d],I).
	I = [b,c,d]
	| ?- rest(a,[b,c,a,d],I).
	I = [d]
	| ?- rest(a,[b,c,d],I).
	I = []
*/
rest(_X,[],T) => T = [].
rest(X,[X|T],TT) => TT=T.
rest(X,XT,R) => XT = [_|T], rest(X,T,R).

/* comb(N,L,Res). Combinations. Arangements without " order".	
	| ?- comb(2,[a,b,c],I).
	I = [a,b] ;
	I = [a,c] ;
	I = [b,c] ;
*/
comb(N,L,X) => X = new_list(N),mem1(X,L).

/* list_comb(N,L,Res).	
	?-  list_comb(2,[a,b,c,d],L).
	L = [[a,b],[a,c],[a,d],[b,c],[b,d],[c,d]]
*/
% Picat> L = [C : C in all_comb(4,1..10), sum(C) = 30]
% L = [[3,8,9,10],[4,7,9,10],[5,6,9,10],[5,7,8,10],[6,7,8,9]]
%
all_comb(N,L) = findall(X,comb(N,L,X)).


/* Order in a list: 
	List= ....X,Y,... means: Y is in the right of X
	?- right(a,b,[c,a,b,m,n]).
	yes
	?- right(a,X,[c,a,b,m,n]).
	X = b ;
	?- right(a,m,[c,a,b,m,n]).
	no
*/
right(X,Y,L) => append(_,[X,Y|_],L).


/* next(X,Y,L). If X and Y are next to each other in the list.	
	?- next(X,Y,[a,b,c]).
	X = a ,Y = b ;
	X = b ,Y = c ;
	X = b ,Y = a ;
	X = c ,Y = b ;
	?- next(a,b,[m,a,b,c]).
	yes
	?- next(a,c,[m,a,b,c,d]).
	yes
	?- next(a,c,[m,a,b,c,d]).
	no
	?- next(a,X,[m,a,b,c]).
	X = b ;
	X = m ;
*/
next(X,Y,L) => right(X,Y,L) ; right(Y,X,L).



/* neighbor(+X,+Y,X1,Y1,+S). Two cells are neighbors in an array size S
- starts from 0. Assume X and Y are in the range 0 - S.
integer_bound(0,I,5) means 0 <= I <= 5
	?-neighbor(1,1,2,2,3). 
	yes
	?-neighbor(1,2,I,J,2).  
	I = J = 2 ;
	I = 0 ,J = 2 ;
	I = J = 1 ;
	I = 2 ,	J = 1 ;
	I = 0 ,	J = 1 ;
*/
neighbor(X,Y,X1,Y1,S) => ((X1 is X+1, Y=Y1);
                      (X1 is X-1, Y=Y1);
		      (X1=X, Y1 is Y+1);
                      (X1=X, Y1 is Y-1);
		      (X1 is X+1, Y1 is Y-1);
		      (X1 is X-1, Y1 is Y-1);
		      (X1 is X+1, Y1 is Y+1);
		      (X1 is X-1, Y1 is Y+1)), 
		      between(0,S,X1),between(0,S,Y1).


/* before (X,Y,List). Checks if X is before Y in the List.
	Starts from Left to right (normal order..).
	?-before(a,c,[m,a,v,c,d]).
	yes
*/
before(X,Y,L) => append(_,[X|T],L), member(Y,T).

%
% length/2
% as Prolog standard "two-way" length/2
%
length(X,Len), var(X), integer(Len) => X = new_list(Len).
length(X,Len), var(Len), list(X) => Len = X.length.

%% This don't work from Len -> List
% lengthX([],Len) ?=> Len = 0.
% lengthX([_|T], L) =>
%   lengthX(T,L1),
%   L = L1 + 1.

% This works (with accumulator)
length2(List,Len) => 
  length2(List,0,Len).
length2(List,Len0,Len) ?=>  List = [], Len = Len0.
length2(List, Len0,Len) =>
  List = [_|T],
  Len1 = Len0 + 1,
  length2(T,Len1,Len).

% This works as well, but is not tail recursive.
length3(List,Len) ?=> 
  List = [],
  Len = 0.
length3(List,Len) => 
  List = [_|T],
  length3(T,Len1),
  Len = Len1+1.



% lcm/2
lcm(X,Y)= abs(X*Y)//gcd(X,Y).
% lcm/3
lcm(X,Y,LCM) => LCM = abs(X*Y)//gcd(X,Y).

lcm(List) = fold(lcm,1,List).

gcd(List) = fold(gcd,List.first(),List.tail()).


%
% argmax:
% find the index/indices for the max value(s) of L
%
argmax(L) = MaxIxs =>
  Max = max(L),
  MaxIxs = [I : I in 1..L.length, L[I] == Max].

%
% argmin:
% find the index/indices for the min value(s) of L
%
argmin(L) = MinIxs =>
  Min = min(L),
  MinIxs = [I : I in 1..L.length, L[I] == Min].


%
% argsort(L) = P:
% P contains the indices to permute L to a sorted list.
% L[P[1]] represents the smallest value in L,
% L[P[2]] represents the second smallest value in L, etc.
% 
% Picat> L=[2,3,1,5,4], L.argsort=T
% L = [2,3,1,5,4]
% T = [3,1,2,5,4]
%
% T[1] = 3 -> A[3] = 1
% T[2] = 1 -> A[1] = 2
% T[3] = 2 -> A[2] = 3
% T[4] = 5 -> A[5] = 4
% T[5] = 4 -> A[4] = 5
%
% Restore the order (in A)
% Picat>  L=[2,3,1,5,4],P=L.argsort, A=[L[I] : I in P]        
% L = [2,3,1,5,4]
% P = [3,1,2,5,4]
% A = [1,2,3,4,5]
%
% Duplicates is handled as expected (here there are two 1s):
% Picat> L = [2,1,3,4,5,1].argsort
% L = [2,6,1,3,4,5]
%
% Get the 3 smallest values:
% Picat> L=[2,1,3,4,5,1,4,7,5,3,89,5,3,4], B=[L[I] : I in L.argsort.take(3)]
% L = [2,1,3,4,5,1,4,7,5,3,89,5,3,4]
% B = [1,1,2]
% 3 largest values:
% Picat> L=[2,1,3,4,5,1,4,7,5,3,89,5,3,4], B=[L[I] : I in L.argsort.reverse.take(3)]
% L = [2,1,3,4,5,1,4,7,5,3,89,5,3,4]
% B = [89,7,5]
%
% Cf permutation3/3 in cp_utils.pi
% What permutation is needed to convert [2,3,1,5,4] to [1,2,3,4,5]:
% Picat> P=new_list(5),P::1..5,permutation3([2,3,1,5,4],P,[1,2,3,4,5])
% P = [3,1,2,5,4]
%
argsort(L) = zip(L,1..L.len).sort.map(second).

% Fast fibonacci function.
table
fibt(0)=1.
fibt(1)=1.
fibt(N)=F,N>1 => F=fibt(N-1)+fibt(N-2).

%
% rotate_left/1 and rotate_right/2:
% rotate list L to the left N times
%
rotate_left(L) = rotate_left(L,1).
rotate_left(L,N) = slice(L,N+1,L.length) ++ slice(L,1,N).

rotate_left_all(L) = [ rotate_left(L,I-1): I in 1..L.len].
 
%
% rotate_right/1 and rotate_right/2:
% rotate list L to the right N times 
%
rotate_right(L) = rotate_right(L,1).
rotate_right(L,N) = Rot =>
  Len = L.length,
  Rot=slice(L,Len-N+1,Len) ++ slice(L,1,Len-N).

rotate_right_all(L) = [ rotate_right(L,I-1): I in 1..L.len].


%
% Rotate the list L N steps (either positive or negative N)
%   rotate(1..10,3) -> [4,5,6,7,8,9,10,1,2,3]
%   rotate(1..10,-3) -> [8,9,10,1,2,3,4,5,6,7]
%
rotate_n(L,N) = Rot => 
  Len = L.length,
  R = cond(N < 0, Len + N, N),
  Rot = [L[I] : I in (R+1..Len) ++ 1..R].


%
% nest/3:
% apply F to Expr N times, return a single value
%
nest(F,Expr,N) = Nest =>
   Nest = apply(F,Expr),
   foreach(_I in 1..N-1)
      Nest := apply(F,Nest)
   end.

% recursive version of nest/3
nest2(_F,Expr0,0,Expr1) => 
   Expr1=Expr0.
nest2(F,Expr0,N,Expr1) => 
   ExprNew = apply(F,Expr0),
   nest2(F,ExprNew,N-1,Expr1).


% 
% next_list/3:
% apply F to Expr N times, return the complete history
%
nest_list(F,Expr,N) = Nest =>
   L = [Expr] ++ [apply(F,Expr)],
   foreach(_I in 1..N-1)
      L := L ++ [apply(F,L.last())]
   end,
   Nest = L.

% recursive version
nest_list2(F,Expr,N,Nest) =>
   nest_list2(F,Expr,N,[Expr],Nest).

nest_list2(_F,_Expr,0,Nest0,Nest) =>
   Nest = Nest0.reverse().
nest_list2(F,Expr,N,Nest0,Nest) =>
   Expr1 = apply(F,Expr),
   nest_list2(F,Expr1,N-1,[Expr1|Nest0],Nest).


%
% apply F to Expr until no difference
%
fixpoint(F,Expr) = Fix =>
   E1 = Expr, 
   E2 = apply(F,Expr),
   L = [E1,E2],
   while (E1 != E2) 
      E1 := E2,
      E2 := apply(F,E1),
      if (E1 != E2) then
         L := L ++ [E2]
      end
   end,
   Fix = L.

dec_to_base(N, Base) = [Alpha[D+1] : D in reverse(Res)] =>
  Alpha = "0123456789abcdefghijklmnopqrstuvwxyz",
  Res = [],
  while (N > 0) 
    R := N mod Base,
    N := N div Base,
    Res := Res ++ [R]
  end.


base_to_dec(N, Base) = base_to_dec(N.to_string(), Base), integer(N) => true.
base_to_dec(N, Base) = Res =>
  println($base_to_dec(N, Base)),
  Alpha = "0123456789abcdefghijklmnopqrstuvwxyz",
  Map = new_map([A=I : {A,I} in zip(Alpha,0..length(Alpha)-1)]),
  Len = N.length,
  Res = sum([Map.get(D)*Base**(Len-I) : {D,I} in zip(N,1..N.length)]).


%
% Day of week, Sakamoto's method
% http://en.wikipedia.org/wiki/Weekday_determination#Sakamoto.27s_Method
%
dow(Y, M, D) = R =>
  T = [0, 3, 2, 5, 0, 3, 5, 1, 4, 6, 2, 4],
  if M < 3 then
     Y := Y - 1
  end,
  R = (Y + Y // 4 - Y // 100 + Y // 400 + T[M] + D) mod 7.

% Maximum days in a month
max_days_in_month(Year,Month) = Days => 
  if member(Month, [1,3,5,7,8,10,12]) then 
    Days = 31
  elseif member(Month,[4,6,9,11]) then
    Days = 30
  else
    if leap_year(Year) then
     Days = 29
   else
     Days = 28
   end
  end.

% is Year a leap year?
leap_year(Year) => 
  (Year mod 4 == 0, Year mod 100 != 0) ; Year mod 400 == 0.


%
% http://en.wikipedia.org/wiki/Julian_day
% Gregorian date -> Julian Day
gregorian2julian(Year,Month,Day) = JD =>
  A = floor((14-Month) / 12), % 1 for Jan or Feb, 0 for other months
  Y = Year + 4800 - A,
  M = Month + 12*A - 3, % 0 for Mars, 11 for Feb
  JD = Day + floor( (153*M + 2) / 5) + 365*Y + floor(Y/4) -
       floor(Y / 100) + floor(Y / 400) - 32045.

% Julian Day -> Gregorian date
julian2gregorian(JD) = Date =>
  Y=4716,
  V=3,
  J=1401,
  U=5,
  M=2,
  S=153,
  N=12,
  W=2,
  R=4,
  B=274277,
  P=1461,
  C= -38,
  F = JD + J + (((4 * JD + B) div 146097) * 3) div 4 + C,
  E = R * F + V,
  G = mod(E, P) div R,
  H = U * G + W,
  Day = (mod(H, S)) div U + 1,
  Month = mod(H div S + M, N) + 1,
  Year = (E div P) - Y + (N + M - Month) div N,
  Date = [Year,Month,Day].



%
% select_many(List, N NewL, Selected)
%
% Selects N elements from the list L:
%  - Selected: the selected elements
%  - NewL: The new L after Selected has been removed
% 
% (used in bridge_and_torch_problem.pi)
%
select_many(L,N, NewL, Selected) => 
  select_many(L,N, 0, NewL, [], Selected).

select_many(L,N, I, NewL,Selected0, Selected), I = N =>
  Selected = Selected0,
  NewL = L.

select_many(L,N, I, NewL,Selected0, Selected) =>
  select(E,L,Rest),
  select_many(Rest,N,I+1,NewL,Selected0++[E], Selected).



% Get the first N numbers that satisfies function F, starting with S
take_n(F,N,S) = L =>
  I = S,
  C = 0,
  L = [],
  while(C < N)
    if call(F,I) then
       L := L ++ [I],
       C := C + 1
    end,
    I := I + 1
  end.


%
% time2 + time_out as a function.
%
time2f(Goal,Timeout) = [End,Backtracks,Status] =>
    statistics(runtime,_),
    statistics(backtracks, Backtracks1),
    time_out(Goal,Timeout,Status),
    statistics(backtracks, Backtracks2),
    statistics(runtime, [_,End]),
    Backtracks = Backtracks2 - Backtracks1.


% 
% Generate a random matrix of size N x N with values of 1..MaxVal,
% default 0 and with a probability of a matrix of Pct %.
%
generate_matrix(N,MaxVal,Pct,Randomize) = Matrix =>
  println($generate_matrix(N,MaxVal,Pct,Randomize)),
  Matrix1 = new_array(N,N),
  bind_vars(Matrix1,0),
  Gen = N*N div Pct,
  % println(gen=Gen),
  if Randomize then _ = random2() end,
  foreach(I in 1..N) 
    Rand = [[1 + random() mod N, 1 + random() mod MaxVal] :  _ in 1..Gen],
    % println(rand=Rand),
    foreach([J,Num] in Rand)
      Matrix1[I,J] := Num
    end
  end,
  Matrix = Matrix1.array_matrix_to_list_matrix().


% From Bob Welham, Lawrence Byrd, R.A.O'Keefe
% (listut.pl)
nth0(0, [Head|_], Head1) => Head1 = Head.

nth0(N, [_|Tail], Elem) =>
    nonvar(N),
    M is N-1,
    nth0(M, Tail, Elem).

nth0(N,[_|T],Item) =>       % Clause added KJ 4-5-87 to allow mode
    var(N),         % nth0(-,+,+)
    nth0(M,T,Item),
    N is M + 1.


pad_string(S,Len,Fill) = Padded =>
   if S.length < Len then
      Padded = [Fill : _I in 1..Len-S.length].flatten().to_string() ++ S
   else
      Padded = S.to_string()
   end.

%
% new_matrix(L,Rows,Cols,M)
%
% M is a Rows x Cols matrix representing the list L.
%
new_matrix(Rows,Cols,L) = M =>
  M = new_array(Rows,Cols),
  foreach(I in 0..Rows-1, J in 0..Cols-1)
    M[I+1,J+1] := L[I*Cols+J+1]
  end.



%
% Get all diagonals on a square matrix of size NxN.
%
all_diagonals(X) = Diagonals =>
  N = X.len,
  % There are in total 2 * NumDiagonals (from left and from right)  
  NumDiagonals = (N*2-1),
  Diagonals = [],
  foreach(K in 1..NumDiagonals)
    Diagonals := Diagonals ++ [[X[I,J] : I in 1..N, J in 1..N, I+J == K+1]],
    Diagonals := Diagonals ++ [[X[J,I] : I in 1..N, J in 1..N, I+(N-J+1) == K+1]]
  end.


%
% Convert a list representing an Rows x Cols matrix (e.g. from MiniZinc) to a proper array matrix.
%
% to a list matrix
convert_board(R,C,Board0) = Board.array_matrix_to_list_matrix =>
  Board = new_array(R,C),
  foreach(I in 0..R-1, J in 0..C-1)
    Board[I+1,J+1] := Board0[C*I+J+1]
  end.
% to an array matrix
convert_board_array(R,C,Board0) = Board =>
  Board = new_array(R,C),
  foreach(I in 0..R-1, J in 0..C-1)
    Board[I+1,J+1] := Board0[C*I+J+1]
  end.

%
% clumped/2 (from SWI-Prolog)
% 
clumped(Items, Counts) :-
    clump(Items, Counts).
clump([], []).
clump([H|T0], [H-C|T]) :-
    ccount(T0, H, T1, 1, C),
    clump(T1, T).
ccount([H|T0], E, T, C0, C) :-
    E==H,
    !,
    C1 is C0+1,
    ccount(T0, E, T, C1, C).
ccount(List, _, List, C, C).

clumped(Items) = Counts =>
  clumped(Items,Counts).

%
% make/0
%
% Reloads all loaded modules.
%
make =>
  BaseModules=[base,sys,io,math],
  foreach(M in loaded_modules(), not member(M,BaseModules))
    % println(M),
    import(M)
  end.


list_test =>
    writeln(map=map(to_integer,to_string(1234567890))),
    writeln(map=1234567890.to_string().map(to_integer)),

    L = 1..10,
    writeln(butfirst=butfirst(L)),
    writeln(butlast=butlast(L)),

    writeln(sort_order=sort_order(>=,L)),
    writeln(swap=swap(1,9999)),
    writeln(swap=swap([1,2,3],[10,11])),

    writeln(group=(1..40).group(mod3)),
    writeln(rotate_left_all=rotate_left_all(1..5)),
    nl.


test3 =>
   cartprod([[a,b,d],[1,2,3,4]],CartProd),
   writeln(cartprod=CartProd),

   all_pairs([1,2,3], Pairs),
   writeln(all_pairs=Pairs),

   if prefix("euler","euler10.pi") then
     writeln($prefix("euler","euler10.pi"))
   end,

   correspond(CorrespondX,[6,2,3,4,5], [1,2,7,5],CorrespondY),
   writeln(correspond=[CorrespondX,CorrespondY]),
   CorrespondAll = findall([CX,CY], $correspond(CX,[6,2,3,4,5], [1,2,3,7,5],CY)),
   writeln(correspond=CorrespondAll),

   writeln(factors=factors(600851475143)),
   % foreach(I in 61..100)         
   %     time(F=factors(2**I-1)),         
   %     writeln([I,2**I-1, F]) 
   % end,

   println("arithmetic (integer) functions"),
   writeln(plus2=plus(1,2)),
   writeln(multiply2=multiply(1..10,2)),
   writeln(divide2=divide(1..10,1..3)),
   writeln(power2=power(2,1..10)),

   writeln(power1x2=power(1..3).power()),
   writeln(mult_divide=multiply(1..20).divide(2)),
   writeln(pi100=pi(100)),

   nl,

   writeln(shuffle=(1..20).shuffle()),
   writeln(num_fixpoints=num_fixpoints((1..20).shuffle())),

   nl,

   P = 1..10,
   writeln(push=P.push(11)),
   writeln(P.shuffle().pop(E1).pop(E2).pop(E3)),
   writeln([e1=E1,e2=E2,e3=E3]),


   nl.