/* 

  Half time in Picat.

  https://brainstellar.com/puzzles/probability/13
  """
  The probability of having accidents on a road in one hour is 3/4. 
  What is the probability of accidents in half an hour?

  Answer: 1/2
  Solution: Probability of no accident in 1 hour = 
    P(no accidents in the first half an hour ) * P(no accidents in the next half an hour)
   --> (1 - 3 / 4) = p^2 -> p = 1/2
  
  """

  Here's a simple approach (though not as simple as the reasoning in the Solution part):
  - Split into 60 minutes and generate the number of accidents. Sum the total number of accidents.
  - Take the first 30 minutes of accidents and sum these.
  - Divide the number of 30 minute accidents with the number of 60 minutes accidents


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

/*

  var : sum t1
  Probabilities (truncated):
  45: 0.0620000000000000
  44: 0.0588000000000000
  43: 0.0560000000000000
  42: 0.0556000000000000
  .........
  24: 0.0002000000000000
  75: 0.0001000000000000
  71: 0.0001000000000000
  22: 0.0001000000000000
  mean = 45.0647

  var : sum t2
  Probabilities (truncated):
  21: 0.0848000000000000
  23: 0.0822000000000000
  22: 0.0798000000000000
  20: 0.0790000000000000
  .........
  40: 0.0003000000000000
  8: 0.0003000000000000
  41: 0.0002000000000000
  7: 0.0001000000000000
  mean = 22.5127

  var : t3
  Probabilities (truncated):
  0.5: 0.0624000000000000
  0.428571428571429: 0.0130000000000000
  0.444444444444444: 0.0104000000000000
  0.4: 0.0102000000000000
  .........
  0.25: 0.0001000000000000
  0.236842105263158: 0.0001000000000000
  0.235294117647059: 0.0001000000000000
  0.212121212121212: 0.0001000000000000
  mean = 0.499376

*/
go ?=>
  reset_store,
  run_model(10_000,$model,[show_probs_trunc,mean]),
  nl,
  % show_store_lengths,nl,
  % fail,
  nl.
go => true.
  
  
model() =>
  Lambda = 3/4,
  T1 = poisson_dist_n(Lambda,60),
  T2 = T1[1..30],
  T3 = T2.sum / T1.sum,
  
  add("sum t1",T1.sum),
  add("sum t2",T2.sum),
  add("t3",T3).