% Balanced Incomplete Block Designed
% prob028 CSPLib (http://4c.ucc.ie/~tw/csplib/prob/prob028/)
% Implementation : guillaume.arnaud@inria.fr
% Date : 22/10/2003

bibd(Solution,V,B,R,K,L):-
	% Solution is a list of V objects, each object is represented 
	% by a list of B blocks. Each block has value 1 if it contains
	% the object, O otherwise.
	length(Solution,V),
	init(Solution,B,R),
	flatten(Solution,AllPairs),
	% constraints: each object occurs in R blocks
	occurrenceObjects(Solution,K,B),
	% constraints: scalar product
	scalarProductAll(Solution,L),
	fd_labeling(AllPairs),
	write_sq(Solution)
	.

% Initializes variables
init([],_,_).
init([Object|S2],B,R) :-
	length(Object,B),
	% init domain
	fd_domain(Object,[0,1]),
	% constraints: the current object is contains in exactly 
	% R blocks
	fd_exactly(R,Object,1),
	init(S2,B,R).
	
% Constraints number of occurences of each object
occurrenceObjects(S,K,B) :-
	occurenceObjectsRec(S,K,1,B).
	
occurenceObjectsRec(_,_,I,B):-
	I>B,!.
occurenceObjectsRec(S,K,I,B):-
	takeCol(S,I,Block),
	fd_exactly(K,Block,1),
	Ip1 is I+1,
	occurenceObjectsRec(S,K,Ip1,B).	

% Constraint onto product scalar between each object
scalarProductAll([],_).
scalarProductAll([Object1|Objects],L):-
	fd_domain(PS,[L]),
	scalarProductAllRec(Object1,Objects,PS,L),
	scalarProductAll(Objects,L).

scalarProductAllRec(_,[],_,_).
scalarProductAllRec(Object1,[Object2|Objects],PS,L):-
%   Comment one of this two lines for use "product scalar" constraint
%	or fd_cardinality/2.
%	scalarProductConst(Object1,Object2,PS,L),
	productListConst(Object1,Object2,L),
	scalarProductAllRec(Object1,Objects,PS,L).

% Scalar product
scalarProductConst([],[],PS,_):- PS#=0.
scalarProductConst([X1|L1],[X2|L2],PS,L):-
	fd_domain(P,[0,1]),
	fd_domain(PSNext,0,L),
	X1*X2#=P, % X1*X2=P
	scalarProductConst(L1,L2,PSNext,L),
	PS#=P+PSNext.

% Alternative of scalar product 
productListConst(Object1,Object2,PS):-
	productList(Object1,Object2,LProduct),
	fd_cardinality(LProduct,PS).

% Product lists
productList([],[],[]).
productList([X|LX],[Y|LY],[X*Y#=1|L]):-
	productList(LX,LY,L).

% Displays a square
write_sq([]).
write_sq([Lg|R]):-
	write(Lg),nl,
	write_sq(R).

% Find instance
findInstance(L):-
	L=[V,B,R,K,Ld],
	fd_domain(L,1,30),
	V#>2,
	Ld#<5,
	R*V#=B*K,
	Ld*(V-1)#=R*(K-1),
	B#>=V,
	randomize,
	fd_labeling(L,[value_method(random)]).

% Functions for manipulate lists

flatten([],[]):-!.
flatten([L1|Rest],L):-
	append(L1,L2,L),
	flatten(Rest,L2).

% Extract a column from a matrix with the feature [[...],[...],...]
takeCol([],_,[]).
takeCol([Ligne|Matrice],I,[X|Col]) :-
	nth(I,Ligne,X),
	takeCol(Matrice,I,Col).