#### Steiner Triple Systems

This is based on a description given in "A Short Course in Combinatorial Designs" by Ian Anderson and Iiro Honkala.

```/*
* Steiner Triple Systems
*
* A Steiner triple system of order v, denoted STS(v), is a pair (X,B) where
* X is a set of cardinality v, and B is a set of 3-subsets of X (called
* blocks), with the property that any 2-subset of X is contained in a unique
* block. STS(v) exists if and only if v=1 or 3 (mod 6).
* There are v(v-1)/6 blocks, and each point occurs (v-1)/2 times.
*/

% e.g. sts(7, Blocks) gives
%        Blocks=[[1,2,7],[1,3,5],[1,4,6],[2,3,4],[2,5,6],[3,6,7],[4,5,7]]
sts(V, Blocks):-
% The Bose construction...
3 =:= V mod 6, !,
T is V // 6,
findall(Element, element_3(T ,Element), BaseSet),
findall(Block, block_3(T, BaseSet, Block), Blocks0),
sort_elements(Blocks0, Blocks1),
sort(Blocks1, Blocks).
sts(V, Blocks):-
% The Skolem construction...
1 =:= V mod 6, !,
T is V // 6,
findall(Element, element_1(T ,Element), BaseSet),
findall(Block, block_1(T, BaseSet, Block), Blocks0),
sort_elements(Blocks0, Blocks1),
sort(Blocks1, Blocks).

% If v=6t+3, the base set is (i,k) for 0<=i<=2t, 0<=k<=2.
element_3(T, e(I,K)):-
TwoT is 2 * T,
for(0, TwoT, I),
for(0, 2, K).

% If v=6t+3, the blocks are:
%   {(i,0), (i,1), (i,2)} for all 0<=i<=2t and
%   {(i,k), (j,k), ((t+1)(i+j),k+1)} for all 0<=i< j <=2t, 0<=k<=2.
block_3(T, BaseSet, [P,Q,R]):-
TwoT is 2 * T,
for(0, TwoT, I),
nth_member(BaseSet, P, e(I,0)),
nth_member(BaseSet, Q, e(I,1)),
nth_member(BaseSet, R, e(I,2)).
block_3(T, BaseSet, [P,Q,R]):-
TwoT is 2 * T,
for(0, TwoT, I),
J0 is I + 1,
for(J0, TwoT, J),
for(0, 2, K),
L is ((T + 1) * (I + J)) mod (TwoT + 1),
M is (K + 1) mod 3,
nth_member(BaseSet, P, e(I,K)),
nth_member(BaseSet, Q, e(J,K)),
nth_member(BaseSet, R, e(L,M)).

% If v=6t+1, the base set is (0,0) and (i,k) for 1<=i<=2t, 0<=k<=2.
element_1(_, e(0,0)).
element_1(T, e(I,K)):-
TwoT is 2 * T,
for(1, TwoT, I),
for(0, 2, K).

% If v=6t+1, the blocks are:
%   {(i,0), (i,1), (i,2)} for 1<=i<=t,
%   {(i,k), (i-t,k+1), (0,0))} for t+1<=i<=2t, 0<=k<=2, and
%   {(i,k), (j,k), (L(i,j),k+1)} for 1<=i< j <=2t, 0<=k<=2
block_1(T, BaseSet, [P,Q,R]):-
for(1, T, I),
nth_member(BaseSet, P, e(I,0)),
nth_member(BaseSet, Q, e(I,1)),
nth_member(BaseSet, R, e(I,2)).
block_1(T, BaseSet, [P,Q,R]):-
T_1 is T + 1,
TwoT is 2 * T,
for(T_1, TwoT, I),
for(0, 2, K),
nth_member(BaseSet, P, e(I,K)),
I_T is I - T,
K_1 is (K + 1) mod 3,
nth_member(BaseSet, Q, e(I_T,K_1)),
nth_member(BaseSet, R, e(0,0)).
block_1(T, BaseSet, [P,Q,R]):-
TwoT is 2 * T,
for(1, TwoT, I),
J0 is I + 1,
for(J0, TwoT, J),
for(0, 2, K),
nth_member(BaseSet, P, e(I,K)),
nth_member(BaseSet, Q, e(J,K)),
K_1 is (K + 1) mod 3,
l(I, J, TwoT, L),
nth_member(BaseSet, R, e(L,K_1)).

% L(i,j)=1+p((i+j-2) mod 2t)
l(I, J, TwoT, L):-
X is (I + J - 2) mod TwoT,
p(X, TwoT, L).

% p(x)=x/2 if x is even; p(x)=(x+2t-1)/2 if x is odd
p(X, _, L):-
0 =:= X mod 2, !,
L is X // 2 + 1.
p(X, TwoT, L):-
L is (X + TwoT - 1) // 2 + 1.

/* for(I, J, K) is true if K is an integer between I and J inclusive.        */
for(I, J, I):-I =< J.
for(I, J, K):-I < J, I1 is I + 1, for(I1, J, K).

/* nth_member(+Xs, ?N, ?X) is true if X is the N-th (base 1) element of      */
/*   the list Xs.                                                            */
nth_member(Xs, N, X):-nth_member_1(Xs, X, 1, N).

nth_member_1([X|_], X, I, I).
nth_member_1([_|Xs], X, I0, I):-
I1 is I0 + 1,
nth_member_1(Xs, X, I1, I).

sort_elements([], []).
sort_elements([X|Xs], [Y|Ys]):-
sort(X, Y),
sort_elements(Xs, Ys).

/*
If c(v-1,k,t,t) is a Steiner system then c(v,k,t,t) = L(v,k,t) where L is
the Schonheim lower bound.
c(v,k,t,m) <= c(v-1,k,t,m)+c(v-1,k-1,t-1,m-1)
c(v,k,2,2) <= c(v-1,k,2,2)+c(v-1,k-1,1,1)
c(v,k,2,2) <= c(v-1,k,2,2)+ceil((v-1)/(k-1))
c(v,k,2,2) = c(v-1,k,2,2)+ceil((v-1)/(k-1)) if c(v-1,k,2,2) is an STS(v-1)
*/

% Constructs a c(V,3,2,2) cover equal to the Schonheim lower bound, where
%   c(V-1,3,2,2) is a Steiner Triple System.
% e.g. v322(8, Blocks) gives
%   Blocks=[[1,2,7],[1,2,8],[1,3,5],[1,4,6],[2,3,4],[2,5,6],
%           [3,4,8],[3,6,7],[4,5,7],[5,6,8],[6,7,8]]
v322(V, Blocks):-
V_1 is V - 1,
sts(V_1, STS),
setof(Block, v322_1(V, 3, STS, Block), Blocks).

% This works for any K.
v322_1(_, _, Cover, Block):-
member(Block, Cover).             % Block extracted from c(V-1,K,2,2)
v322_1(V, K, _, Block):- % V=11, K=4
V_1 is V - 1,
for(1, V_1, I),        % 1..10
I mod (K - 1) =:= 0,   % 3,6,9
J is I - K + 2,        % 1,4,7
findall(X, for(J, I, X), Block0), % Block derived from c(V-1,K-1,1,1)
append(Block0, [V], Block).
v322_1(V, K, _, Block):- % V=11, K=4
I is V - 1,
I mod (K - 1) > 0,     % 10
J is I - K + 2,        % 8
findall(X, for(J, I, X), Block0), % Block derived from c(V-1,K-1,1,1)
append(Block0, [V], Block).

/* append(Xs, Ys, Zs) is true if Zs is the result of appending the list Xs   */
/*   to the list Ys.                                                         */
%append([], Ys, Ys).
%append([X|Xs], Ys, [X|Zs]):-
%  append(Xs, Ys, Zs).

/* member(X, Xs) is true if the element X is contained in the list Xs.       */
%member(X, [X|_]).
%member(X, [_|Xs]):-member(X, Xs).
```