/* 

  6 digit numbers in Picat.

  From Gunnar Blom, Lars Holst, Dennis Sandell:
  "Problems and Snapshots from the World of Probability"
  Page 19f, Problem 2.5 Problems concerning random numbers

  Given the 6 digits numbers:

  a) Problem 1
     find the probability that at least one of the digits 0..9
     appears exactly twice.
     Answer: 2943/4000 ~ 0.7358

  b) Problem 2
     find the probability that at least two of the digits 0..9
     appears exactly once.
     Answer: 1179/1250 ~ 0.9432


  Cf my Gamble model gamble_6_digit_numbers.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.

/*
  The p variable is the answer of the problem (k_occ is shown just 
  for fun :-)).

  Part 1:
  At least one (m) of the digits 0..9 appears exactly twice (k). (6,2,1)
  var : k occ
  Probabilities:
  1: 0.5020500000000000
  0: 0.2634000000000000
  2: 0.2250000000000000
  3: 0.0095500000000000
  mean = 0.9807

  var : p
  Probabilities:
  true: 0.7366000000000000
  false: 0.2634000000000000
  mean = [true = 0.7366,false = 0.2634]


  Part 2:
  At least two (m) of the digits 0..9 appears exactly one (k). (6,1,2)
  var : k occ
  Probabilities:
  4: 0.4500500000000000
  2: 0.2402500000000000
  6: 0.1522500000000000
  3: 0.1016000000000000
  1: 0.0433000000000000
  0: 0.0125500000000000
  mean = 3.5423

  var : p
  Probabilities:
  true: 0.9441500000000000
  false: 0.0558500000000000
  mean = [true = 0.94415,false = 0.05585]

*/
go ?=>
  println("Part 1:"),
  println("At least one (m) of the digits 0..9 appears exactly twice (k). (6,2,1)"),
  reset_store,
  run_model(20_000,$model(6,2,1),[show_probs_trunc,mean]),
  nl,
  println("Part 2:"),
  println("At least two (m) of the digits 0..9 appears exactly one (k). (6,1,2)"),
  reset_store,
  run_model(20_000,$model(6,1,2),[show_probs_trunc,mean]),
  nl.
go => true.

  
model(N,K,M) =>
  Digits = random_integer_n(10,N),

  % Number of occurrences of each digit in digits
  Freq = get_freq(Digits),
  
  % Number of digits that occurs exactly k times
  KOcc = get_freq(Freq.values).get(K,0),

  % Probabily that this occurs at least m times
  P = check(KOcc >= M),
  
  add("k occ",KOcc),
  add("p",P).