2:- use_module(library(mcintyre)).    3
    4:- if(current_predicate(use_rendering/1)).    5:- use_rendering(c3).    6:- use_rendering(graphviz).    7:- endif.    8
    9
   10
   11:- mc.   12
   13:- begin_lpad.   14
   15a(X):user(X,binomial_user(20,0.5)).
   16
   17:- end_lpad.   18
   19
   20binomial_user(N,1.0,N):-!.
   21
   22binomial_user(N,P,X):-
   23  Pr0 is (1-P)^N,
   24  random(U),
   25  binomial_cycle(0,X,N,P,Pr0,Pr0,U).
   26
   27binomial_cycle(X,X,_N,_P,_Pr,CPr,U):-
   28  U=<CPr,!.
   29
   30binomial_cycle(X0,X,N,P,Pr0,CPr0,U):-
   31  X1 is X0+1,
   32  Pr is Pr0*P*(N-X0)/(X1*(1-P)),
   33  CPr is CPr0+Pr,
   34  binomial_cycle(X1,X,N,P,Pr,CPr,U).
   35
   36disc(binomial_user).

?- mc_expectation(a(X),1000,X,E).