%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% AGOSTINO DOVIER
%%% TESTED IN B PROLOG
%%% AUGUST 31 2012
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%


%%% Simplifying assumptions:
%%% Assume Gs all of the same length
%%% Assume all of them have at least a "2"
%%% Assume L > 1

haplo_decision(Gs, L):-
        statistics(runtime,_),
		size(Gs,Ng,Lg),
        Nh is 2*Lg,
        generate_haplo(Hs,Nh,Ng),
        constraint(Gs,Hs),
        size_fd(Hs,SizeHs),
        SizeHs #=< L,
        term_variables(Hs, Vars), 
		labeling([ffc], Vars), !,
        statistics(runtime,[_,T]),
		format("Sol in ~d ms ~n",[T]),
        print_sol(Hs).
        
haplo_decision(_,_) :-         
        statistics(runtime,[_,T]),
		format("Unsat, stated in ~d ms ~n",[T]).

%
size([G|Gs],Ng,Lg) :-
     length(G,Ng),
     length([G|Gs],Lg).

% Count the different haplotypes employed

size_fd([], 0).
size_fd([H|Hs], N) :-
   size_fd(Hs, N1),
   memberflag(H,Hs,F),
   F1 #<=> (F#=0),
   N #= F1 + N1.
   
memberflag(_,[],0).
memberflag(A,[B|C],F) :-
   uguali(A,B,F1),
   FAB #= min(F1),
   memberflag(A,C,F2), 
   F #= FAB + F2.

uguali([],[],[]).
uguali([A|R],[B|S],[F|Fs]) :-
   F #<=> (A #= B),
   uguali(R,S,Fs).

%
generate_haplo(Hs,L,N) :-
     length(Hs,L),
     generate_haplo_aux(Hs,N).      

generate_haplo_aux([],_).        
generate_haplo_aux([H | Hs],N) :-
     length(H,N),domain(H,0,1),
     generate_haplo_aux(Hs,N).

constraint(Gs,Hs) :-
	 slightly_sorted(Hs),
     constraint_aux(Gs, Hs).  
           
    
slightly_sorted([]).
slightly_sorted([A,B|R]) :-
     %lex_chain([A,B],[op(#<)]),
     strict_order_aux(A,B,1),
     slightly_sorted(R).
     
%%% Not built-in in B Prolog

strict_order([]).
strict_order([_]).
strict_order([A,B|R]) :-
     strict_order_aux(A,B,1),
     strict_order([B|R]).

strict_order_aux([],[],0).
strict_order_aux([A|B],[C|D],Fo) :-
     strict_order_aux(B,D,Fi),
     ordtable(T),
     [(A,C,Fi,Fo)] in T.
     %A #= C #=> Fo #= Fi,
     %A #< C #=> Fo #= 1,
     %A #> C #=> Fo #= 0.
 
ordtable( [
% A=C
(0,0, 0,0),     
(0,0, 1,1),     
(1,1, 0,0),     
(1,1, 1,1),     
% A < C
(0,1, 0,1),     
(0,1, 1,1),
% A > C    
(1,0, 0,0),     
(1,0, 1,0)
]).
    
%% Apply constraints to all genotypes %%
		
constraint_aux([], _ ).
constraint_aux([G|Gs], [H1,H2|Hs] ):-
		solve(G, H1, H2, FLAGS),
        1 #= min(FLAGS), 
        constraint_aux(Gs, Hs ).
		
%% solve a single genotype
		
solve([], [], [], []).
solve([G|Gs], [H|Hs], [K|Ks], [F|Flags]):- 
		solve(Gs, Hs, Ks, Flags),
%%%% 1 DIRECT
        (G=2 -> ( F #<=>  H #\= K);
		 G<2 -> ( F #<=> (H#=G #/\ K#=G))).
%%%% 2 Model sum
%        F #<=> (G #= H+K).


%%%%% PRINTING AND EXAMPLES
				
print_sol(Hs) :-
        sort(Hs,Hsorted), 
        myprint(Hsorted).
        
myprint(Hs) :-  
      %foreach(H,Hs) do (print(H),nl).
      forall(member(H,Hs),(print(H),nl)).

%% 1: length 5 (no), 6 (yes)

input(1,Gs) :-
Gs = [ [2,1,1,0,2,0,2],
       [2,0,0,2,1,1,0],
       [2,0,1,2,2,0,2],
       [2,2,2,2,1,2,2]].
   
input(2,Gs) :-
Gs = [ [2,2,2,2],
       [2,2,0,0],
       [0,2,2,0],
       [0,0,2,2] ].

input(3,Gs) :-
Gs = [ [2,2,2,2,0],
       [2,2,0,0,0],
       [0,2,2,0,1],
       [0,1,2,0,1],
       [0,0,2,2,1] ].

%%% 10 no, 11 yes

input(4,Gs) :-
Gs = [ [2,1,1,0,2,0,2],
       [2,0,0,2,1,1,0],
       [2,0,1,2,2,0,2],
       [2,2,2,2,1,2,2],
       [1,2,0,1,2,0,1],
       [1,2,0,2,2,0,1],
       [1,2,2,2,1,2,2],
       [1,1,2,2,1,2,2],
       [1,1,2,1,1,2,2],
       [1,0,2,2,1,0,1] ].

