/* 

  Device life time (Pareto distribution) in Picat.

  From Mathematica ParetoDistribution
  """
  The lifetime of a device follows ParetoDistribution:
  D = ParetoDistribution(Quantity(1, "Years"), 1.23);

  Find the average lifetime of this device:
  Mean(D)
  -> Quantity(5.34783, "Years")

  Find the probability that the device will be operational for more than 6 years:
  Probability(x > Quantity(6, "Years"),  x -> D)
  -> 0.110376
  """

  Picat> X=pareto1_dist_mean(1,1.23)
  X = 5.347826086956522

  Picat> X=1-pareto1_dist_cdf(1,1.23,6)
  X = 0.110375824825248


  Cf my Gamble model gamble_pareto_device_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.
% import ordset.

main => go.

/*

  var : g
  Probabilities (truncated):
  2134.648891308180282: 0.0001000000000000
  1440.127349935320581: 0.0001000000000000
  1191.175460863742956: 0.0001000000000000
  1101.266700412873661: 0.0001000000000000
  .........
  1.000167652171756: 0.0001000000000000
  1.000162395837101: 0.0001000000000000
  1.000073879462783: 0.0001000000000000
  1.000004804855197: 0.0001000000000000
  mean = 4.87275
  HPD intervals:
  HPD interval (0.93): 1.00000480485520..8.93370606974302

  var : p
  Probabilities:
  false: 0.8895000000000000
  true: 0.1105000000000000
  mean = [false = 0.8895,true = 0.1105]



*/
go ?=>
  reset_store,
  run_model(10_000,$model,[show_probs_trunc,mean,
                           % show_percentiles,show_histogram,
                           show_hpd_intervals,hpd_intervals=[0.93]
                           % min_accepted_samples=1000,show_accepted_samples=true
                          ]),
  nl,
  % show_store_lengths,nl,
  % fail,
  nl.
go => true.
  
  
model() =>
  G = pareto1_dist(1,1.23),
  P = check(G > 6),
  
  add("g",G),
  add("p",P).