Math 344: Lecture 13: 
Mathematica tutorial: Questions? Email sjm1@williams.edu

Print[Hyperlink["https://www.youtube.com/watch?v=g1oj7CIqGM8"]]
Print[Hyperlink["http://web.williams.edu/Mathematics/sjmiller/public_html/math/handouts/latex.htm"]]

Lecture on problem solving and coding

Print[Hyperlink["https://web.williams.edu/Mathematics/sjmiller/public_html/341Fa21/mathematicaprograms/Math341_DiplomaProblems.nb"]]
Print[Hyperlink["https://web.williams.edu/Mathematics/sjmiller/public_html/341Fa21/mathematicaprograms/Math341_DiplomaProblems.pdf"]]

(* Want the probability of at least 7 in same suit in a hand of 13 cards*)
(* create a simple deck good enough for this problem *)
deck = {};  (* initializes deck to be empty *)
For[i = 1, i <= 4, i++,
For[j = 1, j <= 13, j++, deck = AppendTo[deck, i]];
] (* end of the i loop *)
Print[deck]

weightedatleastseven[numdo_] := Module[{},
success = 0; (* counts how often happens *)
For[n = 1, n <= numdo, n++,
{
hand = RandomSample[deck,13]; (* chooses a hand of 13 *)
For[s = 0, s <= 5, s++, suit[s] = 0]; (* initialize each suit count to be 0 *)
For[c = 1, c <= 13, c++, 
{
x = hand[[c]]; (*this calls in the c-th element*)
suit[x] = suit[x] + 1;
For[s = 1, s <= 4, s++, If[suit[s] >= 7, success = success+1]];
(* this gives a WEIGHTED count - get one point for 7 in a row, 2 for 8 in a row, .... Even worse as can get multiple credit for 7 in a row! Exit the loop once you find 7 in a row or better put this outside the c loop!*)
}];
(* this is where should put the s For loop *)
}]; (* end of the n loop *)
Print["Observed success rate is ", 100.0 success / numdo, "%."];
theory = 4.0 Sum[ Binomial[13,k] Binomial[39,13-k] ,{k,7,13}]/Binomial[52,13];
Print["Theoretical success rate is about ", 100.0 theory, "%."];
] (* end of the module *)

Timing[weightedatleastseven[10000]]

atleastseven[numdo_] := Module[{},
success = 0; (* counts how often happens *)
For[n = 1, n <= numdo, n++,
{
hand = RandomSample[deck,13]; (* chooses a hand of 13 *)
For[s = 0, s <= 5, s++, suit[s] = 0]; (* initialize each suit count to be 0 *)
For[c = 1, c <= 13, c++, 
{
x = hand[[c]]; (*this calls in the c-th element*)
suit[x] = suit[x] + 1;
}];
For[s = 1, s <= 4, s++, If[suit[s] >= 7, success = success+1]];
}
]; (* end of the n loop *)
Print["Observed success rate is ", 100.0 success / numdo, "%."];
theory = 4.0 Sum[ Binomial[13,k] Binomial[39,13-k] ,{k,7,13}]/Binomial[52,13];
Print["Theoretical success rate is about ", 100.0 theory, "%."];
] (* end of the module *)

atleastseven[100000]

(* darth vader problem - fail with prob p, die after first fail *)
darthvader[numdo_, p_, numfailsallowed_] := Module[{},
(* p is the prob of failing, numdo is number of iterations *)
servicetime = {}; (* how long live *)
For[n = 1, n <= numdo, n++,
{
havenotfailed =  numfailsallowed;
tasksdone = 0;
While[havenotfailed  > 0,
{
x = Random[]; (* if x <= p, failed, else live *)
tasksdone = tasksdone + 1; 
If[x <= p, havenotfailed = havenotfailed - 1];
}]; (* end of while loop *)
servicetime = AppendTo[servicetime, tasksdone];
}]; (* end of the n loop *)
Print[Histogram[servicetime]];
Print[1.0 Mean[servicetime]];
Print[1.0 StandardDeviation[servicetime]];
] (* end of module *)

darthvader[10000,.25, 2]