1/*
    2A three-sided die is repeatedly thrown until the outcome is three.
    3on(T,F) means that on the Tth throw the face F came out.
    4From
    5J. Vennekens, S. Verbaeten, and M. Bruynooghe. Logic programs with annotated 
    6disjunctions. In International Conference on Logic Programming, 
    7volume 3131 of LNCS, pages 195-209. Springer, 2004.
    8*/
    9:- use_module(library(mcintyre)).   10
   11:- if(current_predicate(use_rendering/1)).   12:- use_rendering(c3).   13:- endif.   14
   15:- mc.   16
   17:- begin_lpad.   18
   19% on(T,F) means that the die landed on face F at time T
   20
   21on(0,1):1/3;on(0,2):1/3;on(0,3):1/3.
   22% at time 0 the dice lands on one of its faces with equal probability
   23
   24on(X,1):1/3;on(X,2):1/3;on(X,3):1/3:-
   25  X1 is X-1,X1>=0,
   26  on(X1,_),
   27  \+ on(X1,3).
   28% at time T the die lands on one of its faces with equal probability if
   29% at the previous time point it was thrown and it did not land on face 6
   30
   31evidence:-
   32  on(0,1),
   33  on(1,1).
   34
   35int(0).
   36
   37int(X1):-int(X),X1 is X+1.
   38
   39at_least_once_1:- int(X),((on(X,3),!,fail);on(X,1)).
   40
   41never_1:- \+ at_least_once_1.
   42
   43:- end_lpad.

?- mc_sample(at_least_once_1,1000,Prob,[successes(S),failures(F)]). % what is the probability that the die lands on face 1 at least once? % expected result 0.5 ?- mc_sample(never_1,1000,Prob,[successes(S),failures(F)]). % what is the probability that the die never lands on face 1? % expected result 0.5 ?- mc_prob(on(0,1),Prob),bar(Prob,C). % what is the probability that the die lands on face 1 at time 0? % expected result 0.333333333333333 ?- mc_prob(on(1,1),Prob),bar(Prob,C). % what is the probability that the die lands on face 1 at time 1? % expected result 0.222222222222222 ?- mc_prob(on(2,1),Prob),bar(Prob,C). % what is the probability that the die lands on face 1 at time 2? % expected result 0.148148147703704 ?- mc_mh_sample(on(2,1),on(1,1),1000,P,[mix(1000),successes(T),failures(F)]). % expected result 0.333333333333333 ?- mc_mh_sample(on(2,1),on(0,1),1000,P,[mix(1000),successes(T),failures(F)]). % expected result 0.222222222222222 ?- mc_gibbs_sample(on(0,1),1000,P). % expected result 0.333333333333333 ?- mc_gibbs_sample(on(1,1),1000,P). % expected result 0.222222222222222 ?- mc_gibbs_sample(on(2,1),1000,P). % expected result 0.148148147703704 ?- mc_gibbs_sample(on(2,1),on(1,1),1000,P,[mix(1000),successes(T),failures(F)]). % expected result 0.333333333333333 ?- mc_gibbs_sample(on(2,1),on(0,1),1000,P,[mix(1000),successes(T),failures(F)]). % expected result 0.222222222222222

*/