% Algorithm to find the strongly connected components of a directed graph
% Taken from "Graph Algorithms" by Shimon Even, Pitman (1979), page 65.
%
% The algorithm has complexity k1.E + k2.V + k3 where E is number
% of edges and V is number of vertices (see R. Tarjan:  SIAM J. Computing;
% Vol 1. No 2. p. 146-160).
%
% (Note: the complexity of this version may be a bit worse because it uses
% some less efficient data structures).
%
% The graph is given by the dependency relation of a logic program.
% I.e. the vertices are predicate names and there is an edge p-q if 
% q occurs in the body of a clause with head p.
%
% JPG: 16.9.93   please report bugs etc to john@compsci.bristol.ac.uk
%
% strong_connected_components(Ps,Cs):   
% strong_connected_components(+,-):
%	Ps:  a list of vertices (the edges are generated by procedure step1/2).
%	Cs:  a list of components, each labelled as recursive or non-recursive
%
% Description of data structures (see above references for more details):
%
% vertex(V,Es,K,L,F):  	V a vertex
%			Es is list of vertices U such that V-U is an edge
%			K is the number of V
%			L is the lowpoint of V (see algorithm description)
%			F is either "undef" or the node by which V was reached.
%
% state(G,V,S,I):	G is the current graph 
%			V the current vertex
%			S the stack of visited vertices
%			I the counter
%			
%
% Note:  should run in Sicstus or ProLog by BIM.
% 	 if running in Sicstus, include the last definitions (see comments
%		at end of this program).
%	 if running in BIM, include the BIM definitions at the end
%		of the program, and run the query ?- bim_options.
%
% To use the procedure, the clauses of the object program should be dynamic.
% Read in the file using sp_consult(Filename).
% This asserts the program clauses and a fact program_preds(P),
% where P is a list of the predicates in Filename.
% (Note, at the moment only the head predicates are included).
%
% example:
%
%  ?- sp_consult('../Exs/press.pl'),
%     program_preds(P),
%     strong_connected_components(P,C).
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

strong_connected_components([],[]) :-
	!.
strong_connected_components(Ps,Cs) :-
	step1(Ps,S0),
	step2(S0,Cs),
	recursive_classify(Cs).

step1(Ps,state([V|Vs],V,[],0)) :-
	init_depends(Ps,[V|Vs]).

step2(S,Cs) :-
	next_vertex(S,S1),
	step3(S1,Cs).

step3(S1,Cs) :-
	no_unused_edges(S1),
	!,
	step7(S1,Cs).
step3(S1,Cs) :-
	next_edge(S1,U,S2),
	step4(S2,U,Cs).

step4(S1,U,Cs) :-
	new_vertex(U,S1,S2),
	!,
	step2(S2,Cs).
step4(S1,U,Cs) :-
	step5(S1,U,Cs).

step5(S1,U,Cs) :-
	forward_edge(U,S1),
	!,
	step3(S1,Cs).
step5(S1,U,Cs) :-
	different_component(S1,U),
	!,
	step3(S1,Cs).
step5(S1,U,Cs) :-
	step6(S1,U,Cs).

step6(S1,U,Cs) :-
	newLV(S1,U,S2),
	step3(S2,Cs).

step7(S1,[(_,C)|Cs]) :-
	lowpoint(S1,V),
	!,
	component(S1,V,S2,C),
	step8(S2,Cs).
step7(S1,Cs) :-
	step8(S1,Cs).

step8(S1,Cs) :-
	definedFV(S1),
	!,
	newLV1(S1,S2),
	step3(S2,Cs).
step8(S1,Cs) :-
	step9(S1,Cs).

step9(S1,Cs) :-
	unused_vertex(S1,S2),
	!,
	step2(S2,Cs).
step9(S1,Cs) :-
	step10(S1,Cs).

step10(_,[]).

%--------------------------------------------------
next_vertex(state(G,vertex(V,Es,_,_,F),S,I),
		state(G1,vertex(V,Es,I1,I1,F),[V|S],I1)) :-
	I1 is I+1,
	update_vertex(vertex(V,Es,I1,I1,F),G,G1).
	
no_unused_edges(state(_,vertex(_,[],_,_,_),_,_)).
	
next_edge(state(G,vertex(V,[U|Us],K,L,F),S,I),Urec,
		state(G1,vertex(V,Us,K,L,F),S,I)) :-
	get_vertex(U,G,Urec),
	update_vertex(vertex(V,Us,K,L,F),G,G1).

new_vertex(vertex(U,Es,0,L,_),state(G,vertex(V,_,_,_,_),S,I),
		state(G1,vertex(U,Es,0,L,V),S,I)) :-
	update_vertex(vertex(U,Es,0,L,V),G,G1).

get_vertex(U,[vertex(U,Es,K,L,F)|_],vertex(U,Es,K,L,F)) :-
	!.
get_vertex(U,[_|G],Urec) :-
	get_vertex(U,G,Urec).


forward_edge(vertex(_,_,KU,_,_),state(_,vertex(_,_,KV,_,_),_,_)) :-
	KU > KV.

different_component(state(_,_,S,_),vertex(U,_,_,_,_)) :-
	\+ memb1(U,S).

newLV(state(G,vertex(V,Es,KV,LV,FV),S,I),vertex(_,_,KU,_,_),
		state(G1, vertex(V,Es,KV,LV1,FV),S,I)) :-
	minimum(LV,KU,LV1),
	update_vertex(vertex(V,Es,KV,LV1,FV),G,G1).

lowpoint(state(_,vertex(V,_,LV,LV,_),_,_),V).

component(state(G,V,S,I),V1,state(G,V,S1,I),C) :-
	pop(S,V1,S1,C).

pop([V|S],V,S,[V]) :-
	!.
pop([U|S],V,S1,[U|C]) :-
	pop(S,V,S1,C).

definedFV(state(_,vertex(_,_,_,_,FV),_,_)) :-
	\+ FV = undef.

newLV1(state(G,vertex(_,_,_,L,F),S,I),state(G1,vertex(F,FEs,KF,L1,FF),S,I)) :-
	get_vertex(F,G,vertex(F,FEs,KF,LF,FF)),
	minimum(LF,L,L1),
	update_vertex(vertex(F,FEs,KF,L1,FF),G,G1).
	
unused_vertex(state(G,_,S,I),state(G,U,S,I)) :-
	newstart(G,U).

newstart([vertex(V,Es,0,L,F)|_],vertex(V,Es,0,L,F)) :-
	!.
newstart([_|G],V) :-
	newstart(G,V).
	
update_vertex(vertex(V,Es,K,L,F),[vertex(V,_,_,_,_)|G],
			[vertex(V,Es,K,L,F)|G]) :-
	!.
update_vertex(Y,[V|G],[V|G1]) :-
	update_vertex(Y,G,G1).

minimum(X,Y,X) :-
	X =< Y.
minimum(X,Y,Y) :-
	X > Y.

%-----------------------------------------------

init_depends([P|Ps],[vertex(P,Es,0,0,undef)|Vs]) :-
	immed_depends(1,P,[],Es),
	init_depends(Ps,Vs).
init_depends([],[]).
	

immed_depends(0,_,S,S).
immed_depends(K,P/N,Rs,Rs1) :-
	K > 0,
	user_clauses(P/N,Cls),
	body_preds(K,Cls,Rs,Rs1).


body_preds(K,[(_ :- B)|Cs],S,S2) :-	
	bodylits(K,B,S,S1),
	body_preds(K,Cs,S1,S2).
body_preds(_,[],S,S).


bodylits(K,(B,Bs),S,S3) :-
        \+ sp_builtin(B),
	!,
        functor(B,T,N),
	insertp(T/N,S,S1),
	K1 is K-1,
	immed_depends(K1,T/N,S1,S2),
	bodylits(K,Bs,S2,S3).
bodylits(K,(_,Bs),S,S1) :-
	!,
	bodylits(K,Bs,S,S1).
bodylits(K,B,S,S2) :-
        \+ sp_builtin(B),
	!,
        functor(B,T,N),
        insertp(T/N,S,S1),
	K1 is K-1,
	immed_depends(K1,T/N,S1,S2).
bodylits(_,_,S,S).

insertp(X,L,L) :-
	memb1(X,L),
	!.
insertp(X,L,[X|L]).

recursive_classify([]).
recursive_classify([(recursive,[_,_|_])|Cs]) :-
	!,
	recursive_classify(Cs).
recursive_classify([(recursive,[P])|Cs]) :-
	direct_recursive(P),
	!,
	recursive_classify(Cs).
recursive_classify([(non_recursive,_)|Cs]) :-
	recursive_classify(Cs).

direct_recursive(P/N) :-
	functor(H,P,N),
	clause(H,B),
	rec_body(P/N,B).

rec_body(P/N,(B,_)) :-
	functor(B,P,N),
	!.
rec_body(P/N,(_,Bs)) :-
	!,
	rec_body(P/N,Bs).
rec_body(P/N,B) :-
	functor(B,P,N).
	
user_clauses(P/N,Cls) :-
        functor(H,P,N),
        findall((H:-B),clause(H,B),Cls),
        !.
user_clauses(_,[]).

memb1(X,[X|_]) :-
	!.
memb1(X,[_|Y]) :-
	memb1(X,Y).



/* read in a program and assert clauses, keeping a list of the predicates */

:- dynamic(program_preds/1).

sp_consult(F) :-
	nl,
	see(F),
	abolish_all_existing,
	read(C),
	read_clauses(C,Ps),
	assert(program_preds(Ps)),
	seen,
	nl,
	write('Finished reading '),write(F),
	nl,
	!.
sp_consult(F) :-
	seen,
	write('Cannot find file '), write(F),
	write(', or problems while reading.'), nl,nl.

read_clauses(C,Ps) :- 
	end_of_file(C),
	!,
	close_list(Ps).
read_clauses((:- Dir),Ps) :- 
	!,
	call(Dir),
	read(C1),
	!,
	read_clauses(C1,Ps).
read_clauses(C,Ps) :-
        get_pred_name(C,Pred,Bodypreds),
        each_memb1([Pred|Bodypreds],Ps),
        assert(C),
        read(C1),
        !,
        read_clauses(C1,Ps).

get_pred_name((H :- B),P/N,BPs) :-
        !,
        functor(H,P,N),
        body_preds(B,BPs).
get_pred_name(H ,P/N,[]) :-
        functor(H,P,N).

body_preds(true,[]) :-
        !.
body_preds((\+ B,Bs),Ps) :-
        !,
        body_preds((B,Bs),Ps).
body_preds((B,Bs),Ps) :-
        sp_builtin(B),
        !,
        body_preds(Bs,Ps).
body_preds((B,Bs),[P/N|Ps]) :-
        !,
        functor(B,P,N),
        body_preds(Bs,Ps).
body_preds(\+ B,Ps) :-
        !,
        body_preds(B,Ps).
body_preds(B,[]) :-
        sp_builtin(B),
        !.
body_preds(B,[P/N]) :-
        functor(B,P,N).

each_memb1([],_).
each_memb1([P|Ps],S) :-
	memb1(P,S),
	each_memb1(Ps,S).


close_list([]) :-
	!.
close_list([_|X]) :-
	close_list(X).

abolish_all_existing :-
	program_preds(Ps),
	!,
	write('Clearing workspace....'),
	abolish_each(Ps),
	abolish_pred_list.
abolish_all_existing.

abolish_each([]) :-
	nl.
abolish_each([P/N|Ps]) :-
	abolish_pred(P/N),
	!,
	abolish_each(Ps).

% BIM-specific - remove these when using Sicstus

%bim_options :-
%	please(readeoffail,off),
%	please(readeofatom,end_of_file).
%
end_of_file(end_of_file).	% for compatibility with Sicstus
%
%abolish_pred(P/N) :-
%	functor(H,P,N),
%	abolish(H).
%
%abolish_pred_list :-
%	abolish(program_preds(_)).
%
% end of BIM-specific

% the following to be used with Sicstus

abolish_pred(P) :-
	abolish(P).

abolish_pred_list :-
	abolish(program_preds/1).

sp_builtin(G) :-
        predicate_property(G,built_in).

test_scc(F) :-
	time(sp_consult(F)),
	program_preds(P),
	time(strong_connected_components(P,_)).

