save final version, which goes on for about 1h

makefile

@@ -1,2 +1,6 @@
 dev:
 	go run probppl.go
+
+prod:
+	go build -o probppl
+	./probppl

probppl.go

@@ -8,24 +8,24 @@ import (
 )
 
 type src = *rand.Rand
-type IntProbability struct {
+type IntProb struct {
 	N int64
 	p float64
 }
-type IntProbabilities = []IntProbability
-type IntProbabilitiesWeights struct {
-	IntProb IntProbabilities
-	w       int64
+type IntProbs = []IntProb
+type IntProbsWeights struct {
+	IntProbs IntProbs
+	w        int64
 }
 
-func generatePeopleKnownDistribution(r src) IntProbabilities {
+func generatePeopleKnownDistribution(r src) IntProbs {
 
-	var probabilities IntProbabilities
+	var probabilities IntProbs
 	sum := 0.0
 
 	for i := 16; i <= 2048; i *= 2 {
 		p := r.Float64()
-		probabilities = append(probabilities, IntProbability{N: int64(i), p: p})
+		probabilities = append(probabilities, IntProb{N: int64(i), p: p})
 		sum += p
 	}
 
@@ -74,7 +74,7 @@ func getMatchesDrawGivenNPeopleKnown(n int64, r src) int64 {
 	≥3: 9.5% | 14
 */
 
-func drawFromDistributionWithReplacement(d IntProbabilities, r src) int64 {
+func drawFromDistributionWithReplacement(d IntProbs, r src) int64 {
 	pp := r.Float64()
 	sum := 0.0
 	for i := range d {
@@ -94,7 +94,7 @@ func aboutEq(a int64, b int64) bool {
 	return ((-h) <= (a - b)) && ((a - b) <= h)
 }
 
-func draw148PplFromDistributionAndCheck(d IntProbabilities, r src, show bool) int64 {
+func draw148PplFromDistributionAndCheck(d IntProbs, r src, show bool) int64 {
 
 	count := make(map[int64]int64)
 	count[0] = 0
@@ -117,9 +117,9 @@ func draw148PplFromDistributionAndCheck(d IntProbabilities, r src, show bool) in
 	}
 }
 
-func getUnnormalizedBayesianUpdateForDistribution(d IntProbabilities, r src) int64 {
+func getUnnormalizedBayesianUpdateForDistribution(d IntProbs, r src) int64 {
 	var sum int64 = 0
-	n := 1000
+	n := 10_000
 	for i := 0; i < n; i++ {
 		/* if i%1000 == 0 {
 			fmt.Println(i)
@@ -135,24 +135,37 @@ func main() {
 
 	var r = rand.New(rand.NewPCG(uint64(1), uint64(2)))
 
-	var distribs []IntProbabilitiesWeights
+	var dists []IntProbsWeights
 
 	sum_weights := int64(0)
-	for i := 0; i < 100; i++ {
+	for i := 0; i < 10_000; i++ {
 
 		people_known_distribution := generatePeopleKnownDistribution(r)
 		// fmt.Println(people_known_distribution)
 		result := getUnnormalizedBayesianUpdateForDistribution(people_known_distribution, r)
 		// fmt.Println(i)
 		if result > 0 {
-			fmt.Println(people_known_distribution)
-			fmt.Println(result)
-			distribs = append(distribs, IntProbabilitiesWeights{IntProb: people_known_distribution, w: result})
+			// fmt.Println(people_known_distribution)
+			// fmt.Println(result)
+			dists = append(dists, IntProbsWeights{IntProbs: people_known_distribution, w: result})
 		}
 		sum_weights += result
 		// fmt.Println(result)
 	}
-	// fmt.Println(distribs)
+	// fmt.Println(dists)
 
 	// Now calculate the posterior
+
+	for i := int64(16); i <= 2048; i *= 2 {
+		p := 0.0
+		for _, dist := range dists {
+			for _, int_prob := range dist.IntProbs {
+				if int_prob.N == i {
+					p += float64(dist.w) * int_prob.p
+				}
+			}
+		}
+		p = p / float64(sum_weights)
+		fmt.Printf("%d: %f\n", i, p)
+	}
 }