/* 

  Birthday in Picat.

  From BLOG Swift example birthday.blog
  https://github.com/lileicc/swift/blob/master/example/birthday.blog
  """
  Distribution of values for 
   exists Person x (exists Person y ((!(x = y) & (Birthday(x) = Birthday(y)))))
	true	0.5446927374301685
	false	0.4553072625698328
  """

  Closed world = #t means that we observe that there are 23 people in the room.

  Cf my Gamble model gamble_birthday.rkt

  This program was created by Hakan Kjellerstrand, hakank@gmail.com
  See also my Picat page: http://www.hakank.org/picat/

*/

import ppl_distributions, ppl_utils.
import util.
% import ordset.

main => go.

/*

  closed_world = false
  var : num people
  Probabilities:
  12: 0.1122000000000000
  13: 0.1115000000000000
  14: 0.1077000000000000
  11: 0.0996000000000000
  15: 0.0842000000000000
  10: 0.0797000000000000
  16: 0.0722000000000000
  9: 0.0672000000000000
  17: 0.0578000000000000
  8: 0.0445000000000000
  18: 0.0402000000000000
  7: 0.0277000000000000
  19: 0.0274000000000000
  20: 0.0164000000000000
  6: 0.0152000000000000
  21: 0.0117000000000000
  22: 0.0083000000000000
  5: 0.0062000000000000
  23: 0.0044000000000000
  24: 0.0022000000000000
  4: 0.0019000000000000
  25: 0.0006000000000000
  26: 0.0004000000000000
  3: 0.0003000000000000
  29: 0.0002000000000000
  28: 0.0002000000000000
  31: 0.0001000000000000
  mean = 13.0769

  var : b
  Probabilities:
  0: 0.7953000000000000
  1: 0.1794000000000000
  2: 0.0230000000000000
  3: 0.0019000000000000
  4: 0.0004000000000000
  mean = 0.2327

  var : p
  Probabilities:
  false: 0.7953000000000000
  true: 0.2047000000000000
  mean = 0.2047

  closed_world = true
  min_accepted_samples = 1000
  var : num people
  Probabilities:
  23: 1.0000000000000000
  mean = 23.0

  var : b
  Probabilities:
  0: 0.4810000000000000
  1: 0.3680000000000000
  2: 0.1250000000000000
  3: 0.0240000000000000
  4: 0.0020000000000000
  mean = 0.698

  var : p
  Probabilities:
  true: 0.5190000000000000
  false: 0.4810000000000000
  mean = 0.519

*/
go ?=>
  member(ClosedWorld,[false,true]),
  println(closed_world=ClosedWorld),
  reset_store,
  Opts = [show_probs,mean],
  Options = cond(ClosedWorld==false,Opts,Opts ++ [min_accepted_samples=1000,show_accepted_samples=false]),
  run_model(10_000,$model(ClosedWorld),Options),
  nl,
  % show_store_lengths,nl,
  fail,
  nl.
go => true.
  
  
model(ClosedWorld) =>
  Day = 0..364,

  NumPeople = poisson_dist(13),

  Birthday = uniform_draw_n(Day,NumPeople),

  % Is there a duplicate of birthdays?
  B = NumPeople - Birthday.remove_dups.len,
  P = check(B > 0),

  if ClosedWorld == true then
    observe(NumPeople == 23)
  end,
  
  if ClosedWorld == false ; observed_ok then
    add("num people",NumPeople),
    add("b",B),
    add("p",P),
  end.