%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%% AGOSTINO DOVIER
%%% TESTED IN BPROLOG
%%% AUGUST 30, 2012
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

pf(Primary ) :- % Primary = [h,p,p,h,p,p,h,...] p/0 h/1
       statistics(runtime,_),
       constrain(Primary,Tertiary,Energy),
       labeling([ff,minimize(Energy)],Tertiary),
       statistics(runtime,[_,T]),
       format(" Running Time: ~d ms ~n Energy: ~d ~n Folding:~n",[T,Energy]),
       print_solution(Primary,Tertiary).
pf(_) :-
       statistics(runtime,[_,T]),
       format("No solutions: running time ~d ms~n",[T]).

% Tertiary = [X1,Y1,X2,Y2,X3,Y3,...]

constrain(Primary,Tertiary,Energy) :- 
       length(Primary,N),
       M is 2*N,
       length(Tertiary,M),
       Min is 0, %  N - integer(sqrt(N)), % 0
       Max is M, %  N + integer(sqrt(N)),   % 2*N
       domain(Tertiary,Min,Max),
       starting_point(Tertiary,N),
       avoid_self_loops(Tertiary),
       next_constraints(Tertiary),
       energy_constraint(Primary,Tertiary,Energy).
  
starting_point([N,N,N,N1|T],N) :-
      N1 is N + 1,
      % AVOID SYMMETRIES
      first_bend_right([N,N,N,N1|T],0).

first_bend_right([],_).
first_bend_right([_,_],_).
first_bend_right([X1,Y1,X2,Y2|R],F) :-
      F #= 1 #=> F1 #= 1, 
      (F #= 0 #/\ X1 #= X2) #=> F1 #= 0,
      (F #= 0 #/\ Y1 #= Y2) #=> (X1 #< X2 #/\ F1 #= 1),
      first_bend_right([X2,Y2|R],F1).


avoid_self_loops(Tertiary) :-
       positions_to_integers(Tertiary, ListaInteri),
       all_distinct(ListaInteri).

positions_to_integers([X,Y|R], [I|S]):-
       I  #= X*100+Y,
       positions_to_integers(R,S).
positions_to_integers([],[]).

next_constraints([_,_]) .
next_constraints([X1,Y1,X2,Y2|C]) :-
       next(X1,Y1,X2,Y2),
       next_constraints([X2,Y2|C]).

next(X1,Y1,X2,Y2):-
       abs(X1-X2)+ abs(Y1-Y2) #= 1.

energy_constraint([_],_,0).
energy_constraint([A,B|Primary],[XA,YA,XB,YB|Tertiary],E) :-
       energy_contribution_of_A(0,A,XA,YA,Primary,Tertiary,E1),
       energy_constraint([B|Primary],[XB,YB|Tertiary],E2),
       E #= E1 + E2.

energy_contribution_of_A(_,_,_,_,[],[],0).
energy_contribution_of_A(0,A,XA,YA,[_|Primary],[_,_|Tertiary],E):-
       energy_contribution_of_A(1,A,XA,YA,Primary,Tertiary,E).
energy_contribution_of_A(1,A,XA,YA,[B|Primary],[XB,YB|Tertiary],E):-
       energy(A,XA,YA,B,XB,YB,C),
       energy_contribution_of_A(0,A,XA,YA,Primary,Tertiary,E1),
       E #= E1 + C.

energy(h,XA,YA,h,XB,YB,C) :-
       C in -1..0,
       abs(XA - XB) + abs(YA - YB) #= 1  #<=> C #= -1.

energy(h,_,_,p,_,_,0).
energy(p,_,_,h,_,_,0).
energy(p,_,_,p,_,_,0).

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    
print_solution(Primary,Tertiary) :-
     length(Tertiary,N2),
     M is N2*N2,
     length(Matrix,M),
     scan(0,0,N2,Primary,Tertiary,Matrix),
     write_matrix(Matrix,N2).

scan(N,_,N,_,_,_) :- !.
scan(I,N,N,Primary,Tertiary,Matrix) :-
    !,
    I1 is I+1,
    scan(I1,0,N,Primary,Tertiary,Matrix).
scan(I,J,N,Primary,Tertiary,Matrix) :-
    N is 2*I,
    N is 2*J,
    !,  %%% caso particolare primo punto
    nth0(0,Primary,AMINO),
    (AMINO == h, POL = 'H';
     AMINO == p, POL = 'P'),
    K is (J*N+I), %%% 
    nth0(K,Matrix,POL),
    J1 is J+1,
    scan(I,J1,N,Primary,Tertiary,Matrix).
scan(I,J,N,Primary,Tertiary,Matrix) :-
    (member_pair(I,J,Tertiary,Val),
     nth0(Val,Primary,POL),
     !;
     POL = ' '),
    K is (J*N+I), %%% 
    nth0(K,Matrix,POL),
    J1 is J+1,
    scan(I,J1,N,Primary,Tertiary,Matrix).
        
%%%%
member_pair(I,J,[I,J|_],0) :- !.
member_pair(I,J,[_,_|R],K) :-
    member_pair(I,J,R,K1),
    K is K1 + 1.
      
%%%% 

write_matrix([],_) :-!.
write_matrix(Matrix,N) :-
   write_row(Matrix,Rest,N),
   write_matrix(Rest,N).

write_row(Matrix,Rest,N) :-
   length(First,N),
   append(First,Rest,Matrix),
   (null_row(First),!;
    write_vector(First),nl).   
   
write_vector([]) :-!.
write_vector([A|R]) :-
     write(A),
     write_vector(R).

null_row([]).
null_row([' '|R]) :-
     null_row(R).     
     
