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)
.