/* 

  Coupon collector problem in Picat.

  There are N different collecter's cards hidden in a package, but we don't
  know which card there is in the package we buy.
  We want to collect all of them, how many packages must one buy to collect
  all the different cards?

  See https://en.wikipedia.org/wiki/Coupon_collector%27s_problem
  """
  In probability theory, the coupon collector's problem describes 'collect all coupons and win' 
  contests. It asks the following question: If each box of a brand of cereals contains a 
  coupon, and there are n different types of coupons, what is the probability that more 
  than t boxes need to be bought to collect all n coupons? 

  An alternative statement is: Given n coupons, how many coupons do you expect you need 
  to draw with replacement before having drawn each coupon at least once? The mathematical 
  analysis of the problem reveals that the expected number of trials needed grows as 
  Θ(n log(n).
  For example, when n = 50 it takes about 225[b] trials on average to collect all 50 coupons. 

  ...

  [b]: E(50) = 50(1 + 1/2 + 1/3 + ... + 1/50) = 224.9603, the expected number of trials to 
  collect all 50 coupons. 

   The approximation n*log(n) + γ*n + 1/2 for this expected number gives in this case 
   ≈ 195.6011 + 28.8608 + 0.5 ≈ 224.9619.  [log is the natural logarithm]
  """ 

  Cf 
  - ppl_coupon_collectors3.pi
  - ppl_coupon_collectors_problem.pi
  - my Gamble model gamble_coupon_collector2.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.

main => go.

/*

  [n = 10,approx = 29.298,exact = 29.2897]
  var : len
  Probabilities (truncated):
  25: 0.0462000000000000
  21: 0.0459000000000000
  24: 0.0454000000000000
  23: 0.0450000000000000
  .........
  92: 0.0001000000000000
  87: 0.0001000000000000
  83: 0.0001000000000000
  81: 0.0001000000000000
  mean = 29.3405


  [n = 40,approx = 171.144,exact = 171.142]
  var : len
  Probabilities (truncated):
  139: 0.0111000000000000
  151: 0.0110000000000000
  132: 0.0108000000000000
  147: 0.0107000000000000
  .........
  81: 0.0001000000000000
  78: 0.0001000000000000
  76: 0.0001000000000000
  64: 0.0001000000000000
  mean = 170.063


  [n = 100,approx = 518.739,exact = 518.738]
  var : len
  Probabilities (truncated):
  484: 0.0050000000000000
  463: 0.0049000000000000
  483: 0.0046000000000000
  468: 0.0046000000000000
  .........
  267: 0.0001000000000000
  266: 0.0001000000000000
  257: 0.0001000000000000
  238: 0.0001000000000000
  mean = 521.346

*/
go ?=>
  member(N,[10,40,100]),
  println([n=N,approx=coupon_collectors_problem_approx(N),exact=coupon_collectors_problem_theoretical(N)]),
  reset_store,
  run_model(10_000,$model(N),[show_probs_trunc,mean]),
  nl,
  % show_store_lengths,nl,
  fail,
  nl.
go => true.
  

collect_coupons(L,N) = Res =>
  if N == L.remove_dups.len then
    Res = L
  else
    Res = collect_coupons(L++[random_integer1(N)],N)
  end.
  
model(N) =>
  
  Len = collect_coupons([],N).len,
  add("len",Len).