```% Card generation

% Classic setup
setup([[red, green, blue],
[1, 2, 3],
[circle, triangle, square],
[empty, dashed, filled]]).

% Alternative setup
% setup([[red, green, blue, black],
%        [1, 2, 3, 4],
%        [circle, triangle, square, pentagon]]).

card(C) :- setup(L), maplist(member, C, L).

cards(L) :- findall(C, card(C), L).

% Set property

set([[]|_]).
set(L) :-
( all_same(H) ; all_diff(H) ),
set(T).

all_same([_]) :- !.
all_same([X,X|Xs]) :- all_same([X|Xs]).

all_diff([_]) :- !.
all_diff([X|Xs]) :- forall(member(Y, Xs), X \= Y), all_diff(Xs).

% Probability experiments

% What is the probability P of having a set in K cards?
%   (using M Monte-Carlo iterations)
probability(K, M, P) :-
cards(C), length(C, N),
monte_carlo(M, C, N, K, S),
P is S / M.

monte_carlo(0, _, _, _, 0) :- !.
monte_carlo(M, Cards, NumCards, K, SetCount) :-
M > 0, M1 is M - 1,
monte_carlo(M1, Cards, NumCards, K, S),
randset(K, NumCards, Indices),
findall(X, (member(I, Indices), nth1(I, Cards, X)), L),
( contains_set(L) -> SetCount is S + 1 ; SetCount = S ).

contains_set(L) :-
setup([X|_]), length(X, N),
select_n(N, L, S),
set(S).

select_n(0, _, []).
select_n(N, [X|Xs], [X|L]) :- N > 0, N1 is N - 1, select_n(N1, Xs, L).
select_n(N, [_|Xs], L) :- N > 0, select_n(N, Xs, L).

% Test

% Classic setup
% ?- probability(12, 10000, P).
% P = 0.9665.
%
% ?- probability(15, 10000, P).
% P = 0.9998.
%
% Alternative setup
% ?- probability(12, 10000, P).
% P = 0.5199.
%
% ?- probability(16, 10000, P).
% P = 0.936.
%
% ?- probability(18, 10000, P).
% P = 0.9882.
%
% ?- probability(20, 10000, P).
% P = 0.9987.

% So we should use 16 of the alternative cards,
% adding another 4 when there is no set.
```