// compute probability of landing on tile after certain number of steps of brownian package main import ( "flag" "fmt" "math/big" "math/rand" "os" "time" ) // Tile is the fundamental data object for this puzzle // in the puzzle, it's a hex piece, but this is more general because that's just the way I thought to do it type Tile struct { index int nbrs map[int]bool } var ( numIters = flag.Int("num_iters", 19, "number of iterations of brownian motion to do") numTiles = flag.Int("num_tiles", 74, "number of tiles to put in the grid") pcxn = flag.Float64("pcxn", 0, "if non-zero, generates a random grid with pcxn probability of being connected to another neighbor") verbose = flag.Bool("verbose", false, "show verbose log messages") printAsRat = flag.Bool("print_as_rat", false, "prints the final probability as a floating point number") printTop = flag.Int("print_top", 1, "number of tiles to print ordered from most probability, descending.") ) var tile = make([]*Tile, 0) var p = make([]*big.Rat, 0) // Reset the state of the Tile Generator for hermetic Tests func Reset() { tile = make([]*Tile, 0) } // AddNeighbors to the associated Tile func (t *Tile) AddNeighbors(nbrs ...*Tile) { for _, n := range nbrs { t.nbrs[n.index] = true n.nbrs[t.index] = true } } // AddNeighborsByIndex instead of their object func (t *Tile) AddNeighborsByIndex(idx ...int) { for _, i := range idx { t.nbrs[i] = true tile[i].nbrs[t.index] = true } } // NumNeighbors returns the number of neighbors for the associated Tile func (t *Tile) NumNeighbors() int { return len(t.nbrs) } // GetNeighbors of the associated Tile func (t *Tile) GetNeighbors() []*Tile { ret := []*Tile{} for idx := range t.nbrs { ret = append(ret, tile[idx]) } return ret } // HasNeighbor returns true iff the the queried tile is a neighbor of the associated Tile func (t *Tile) HasNeighbor(n *Tile) bool { _, ok := t.nbrs[n.index] return ok } func (t *Tile) String() string { return fmt.Sprintf("{index: %d, neighbors: %v}", t.index, t.nbrs) } // NewTile creates a new Tile and returns it. index will increment every time it's called. func NewTile() *Tile { ret := &Tile{ index: len(tile), nbrs: make(map[int]bool), } tile = append(tile, ret) return ret } // NewTiles creates a number of new tiles based on the parameter, n. func NewTiles(n int) []*Tile { ret := []*Tile{} for i := 0; i < n; i++ { ret = append(ret, NewTile()) } return ret } func main() { flag.Parse() NewTiles(*numTiles) // set up neighbors if *pcxn == 0 { tile[1].AddNeighborsByIndex(3, 4) tile[2].AddNeighborsByIndex(3, 5) tile[3].AddNeighborsByIndex(1, 2, 4, 5, 7) tile[4].AddNeighborsByIndex(1, 3, 6, 7, 9) tile[5].AddNeighborsByIndex(2, 3, 7, 10) tile[6].AddNeighborsByIndex(4, 9, 12) tile[7].AddNeighborsByIndex(3, 4, 5, 9, 10, 13) tile[8].AddNeighborsByIndex(11, 14) tile[9].AddNeighborsByIndex(4, 6, 7, 12, 13, 17) tile[10].AddNeighborsByIndex(5, 7, 13, 18) tile[11].AddNeighborsByIndex(8, 14, 15, 19) tile[12].AddNeighborsByIndex(6, 9, 16, 17, 22) tile[13].AddNeighborsByIndex(7, 9, 10, 17, 18, 23) tile[14].AddNeighborsByIndex(8, 11, 19) tile[15].AddNeighborsByIndex(11, 19, 20, 25) tile[16].AddNeighborsByIndex(12, 21, 22, 28) tile[17].AddNeighborsByIndex(9, 12, 13, 22, 23, 29) tile[18].AddNeighborsByIndex(10, 13, 23, 24, 30) tile[19].AddNeighborsByIndex(11, 14, 15, 25, 31) tile[20].AddNeighborsByIndex(15, 25, 26, 32) tile[21].AddNeighborsByIndex(16, 27, 28, 34) tile[22].AddNeighborsByIndex(12, 16, 17, 28, 29, 35) tile[23].AddNeighborsByIndex(13, 17, 18, 29, 30, 36) tile[24].AddNeighborsByIndex(18, 30) tile[25].AddNeighborsByIndex(15, 19, 20, 31, 32, 38) tile[26].AddNeighborsByIndex(20, 32, 33, 39) tile[27].AddNeighborsByIndex(21, 33, 34, 40) tile[28].AddNeighborsByIndex(16, 21, 22, 34, 35, 41) tile[29].AddNeighborsByIndex(17, 22, 23, 35, 36, 42) tile[30].AddNeighborsByIndex(18, 23, 24, 36) tile[31].AddNeighborsByIndex(19, 25, 37, 38, 43) tile[32].AddNeighborsByIndex(20, 25, 26, 38, 39, 44) tile[33].AddNeighborsByIndex(26, 27, 39, 40, 45) tile[34].AddNeighborsByIndex(21, 27, 28, 40, 41, 46) tile[35].AddNeighborsByIndex(22, 28, 29, 41, 42, 47) tile[36].AddNeighborsByIndex(23, 29, 30, 42) tile[37].AddNeighborsByIndex(31, 43, 48) tile[38].AddNeighborsByIndex(25, 31, 32, 43, 44, 49) tile[39].AddNeighborsByIndex(26, 32, 33, 44, 45, 50) tile[40].AddNeighborsByIndex(27, 33, 34, 45, 46, 51) tile[41].AddNeighborsByIndex(28, 34, 35, 46, 47, 52) tile[42].AddNeighborsByIndex(29, 35, 36, 47, 53) tile[43].AddNeighborsByIndex(31, 37, 38, 48, 49) tile[44].AddNeighborsByIndex(32, 38, 39, 49, 50, 55) tile[45].AddNeighborsByIndex(33, 39, 40, 50, 51, 56) tile[46].AddNeighborsByIndex(34, 40, 41, 51, 52, 57) tile[47].AddNeighborsByIndex(35, 41, 42, 52, 53, 58) tile[48].AddNeighborsByIndex(37, 43, 54) tile[49].AddNeighborsByIndex(38, 43, 44, 55) tile[50].AddNeighborsByIndex(39, 44, 45, 55, 56) tile[51].AddNeighborsByIndex(40, 45, 46, 56, 57, 59) tile[52].AddNeighborsByIndex(41, 46, 47, 57, 58, 60) tile[53].AddNeighborsByIndex(42, 47, 58, 61) tile[54].AddNeighborsByIndex(48) tile[55].AddNeighborsByIndex(44, 49, 50) tile[56].AddNeighborsByIndex(45, 50, 51, 59) tile[57].AddNeighborsByIndex(46, 51, 52, 59, 60, 62) tile[58].AddNeighborsByIndex(47, 52, 53, 60, 61, 63) tile[59].AddNeighborsByIndex(51, 56, 57, 62) tile[60].AddNeighborsByIndex(52, 57, 58, 62, 63, 64) tile[61].AddNeighborsByIndex(53, 58, 63) tile[62].AddNeighborsByIndex(57, 59, 60, 64, 65) tile[63].AddNeighborsByIndex(58, 60, 61, 64, 66) tile[64].AddNeighborsByIndex(60, 62, 63, 65, 66, 67) tile[65].AddNeighborsByIndex(62, 64, 67, 68) tile[66].AddNeighborsByIndex(63, 64, 67, 69) tile[67].AddNeighborsByIndex(64, 65, 66, 68, 69, 70) tile[68].AddNeighborsByIndex(65, 67, 70, 71) tile[69].AddNeighborsByIndex(66, 67, 70, 72) tile[70].AddNeighborsByIndex(67, 68, 69, 71, 72, 73) tile[71].AddNeighborsByIndex(68, 70, 73) tile[72].AddNeighborsByIndex(69, 70, 73) tile[73].AddNeighborsByIndex(70, 71, 72) } else { r := rand.New(rand.NewSource(time.Now().UnixNano())) for i := 0; i < len(tile); i++ { for j := i + 1; j < len(tile); j++ { if r.Float64() < *pcxn { tile[i].AddNeighborsByIndex(j) } } } } // set up probs pLast := make([]*big.Rat, len(tile)) for i := 0; i < len(pLast); i++ { pLast[i] = big.NewRat(0, 1) } if *pcxn > 0 { pLast[0] = big.NewRat(1, 1) } else { pLast[1] = big.NewRat(1, 1) } pNext := make([]*big.Rat, len(tile)) for i := 0; i < len(pNext); i++ { pNext[i] = big.NewRat(0, 1) } // loop through iterations and perform brownian motion, storing // probabilities at each iteration if *verbose { fmt.Println("iteration 0") fmt.Println(tile) fmt.Println(pLast) fmt.Println() } for t := 0; t < *numIters; t++ { for tIdx, t := range tile { pNext[tIdx] = big.NewRat(0, 1) for _, n := range t.GetNeighbors() { pOfN := big.NewRat(0, 1).Set(pLast[n.index]) pNext[tIdx] = pNext[tIdx].Add( pNext[tIdx], pOfN.Mul(pOfN, big.NewRat(1, int64(n.NumNeighbors())))) } } // swap pLast and pNext for i := range pLast { t := pLast[i] pLast[i] = pNext[i] pNext[i] = t } if *verbose { fmt.Println("iteration ", t+1) fmt.Println(tile) fmt.Println(pLast) fmt.Println() } } // sanity check that all probs add up to 1 s := big.NewRat(0, 1) for idx := range pLast { s.Add(s, pLast[idx]) } if s.Cmp(big.NewRat(1, 1)) != 0 { fmt.Fprintf(os.Stderr, "All probabilities do not add up to 1. Their sum is %v.", s) return } printed := make(map[int]bool) for i := 0; i < *printTop; i++ { // find most likely ending location maxP := big.NewRat(0, 1) maxI := -1 for idx := range pLast { if maxP.Cmp(pLast[idx]) < 0 { if _, ok := printed[idx]; !ok { maxP = pLast[idx] maxI = idx } } } printed[maxI] = true if *printAsRat { fmt.Printf("Most Likely #%d: index %d with probability %v\n", i+1, maxI, maxP) } else { fmt.Printf("Most Likely #%d: index %d with probability %v\n", i+1, maxI, big.NewFloat(0).SetRat(maxP)) } } }