Robert Griesemer | 8d5fe68 | 2014-12-08 11:20:41 -0800 | [diff] [blame] | 1 | // skip |
| 2 | |
| 3 | // Copyright 2014 The Go Authors. All rights reserved. |
| 4 | // Use of this source code is governed by a BSD-style |
| 5 | // license that can be found in the LICENSE file. |
| 6 | |
| 7 | // Usage: |
Robert Griesemer | 7796821 | 2014-12-08 14:19:38 -0800 | [diff] [blame] | 8 | // fibo <n> compute fibonacci(n), n must be >= 0 |
| 9 | // fibo -bench benchmark fibonacci computation (takes about 1 min) |
Robert Griesemer | 8d5fe68 | 2014-12-08 11:20:41 -0800 | [diff] [blame] | 10 | // |
| 11 | // Additional flags: |
| 12 | // -half add values using two half-digit additions |
Robert Griesemer | 7796821 | 2014-12-08 14:19:38 -0800 | [diff] [blame] | 13 | // -opt optimize memory allocation through reuse |
| 14 | // -short only print the first 10 digits of very large fibonacci numbers |
Robert Griesemer | 8d5fe68 | 2014-12-08 11:20:41 -0800 | [diff] [blame] | 15 | |
| 16 | // Command fibo is a stand-alone test and benchmark to |
| 17 | // evaluate the performance of bignum arithmetic written |
| 18 | // entirely in Go. |
| 19 | package main |
| 20 | |
| 21 | import ( |
| 22 | "flag" |
| 23 | "fmt" |
| 24 | "math/big" // only used for printing |
| 25 | "os" |
| 26 | "strconv" |
| 27 | "testing" |
| 28 | "text/tabwriter" |
| 29 | "time" |
| 30 | ) |
| 31 | |
| 32 | var ( |
| 33 | bench = flag.Bool("bench", false, "run benchmarks") |
| 34 | half = flag.Bool("half", false, "use half-digit addition") |
| 35 | opt = flag.Bool("opt", false, "optimize memory usage") |
| 36 | short = flag.Bool("short", false, "only print first 10 digits of result") |
| 37 | ) |
| 38 | |
| 39 | // A large natural number is represented by a nat, each "digit" is |
| 40 | // a big.Word; the value zero corresponds to the empty nat slice. |
| 41 | type nat []big.Word |
| 42 | |
| 43 | const W = 1 << (5 + ^big.Word(0)>>63) // big.Word size in bits |
| 44 | |
| 45 | // The following methods are extracted from math/big to make this a |
| 46 | // stand-alone program that can easily be run without dependencies |
| 47 | // and compiled with different compilers. |
| 48 | |
| 49 | func (z nat) make(n int) nat { |
| 50 | if n <= cap(z) { |
| 51 | return z[:n] // reuse z |
| 52 | } |
| 53 | // Choosing a good value for e has significant performance impact |
| 54 | // because it increases the chance that a value can be reused. |
| 55 | const e = 4 // extra capacity |
| 56 | return make(nat, n, n+e) |
| 57 | } |
| 58 | |
| 59 | // z = x |
| 60 | func (z nat) set(x nat) nat { |
| 61 | z = z.make(len(x)) |
| 62 | copy(z, x) |
| 63 | return z |
| 64 | } |
| 65 | |
| 66 | // z = x + y |
| 67 | // (like add, but operating on half-digits at a time) |
| 68 | func (z nat) halfAdd(x, y nat) nat { |
| 69 | m := len(x) |
| 70 | n := len(y) |
| 71 | |
| 72 | switch { |
| 73 | case m < n: |
| 74 | return z.add(y, x) |
| 75 | case m == 0: |
| 76 | // n == 0 because m >= n; result is 0 |
| 77 | return z.make(0) |
| 78 | case n == 0: |
| 79 | // result is x |
| 80 | return z.set(x) |
| 81 | } |
| 82 | // m >= n > 0 |
| 83 | |
| 84 | const W2 = W / 2 // half-digit size in bits |
| 85 | const M2 = (1 << W2) - 1 // lower half-digit mask |
| 86 | |
| 87 | z = z.make(m + 1) |
| 88 | var c big.Word |
| 89 | for i := 0; i < n; i++ { |
| 90 | // lower half-digit |
| 91 | c += x[i]&M2 + y[i]&M2 |
| 92 | d := c & M2 |
| 93 | c >>= W2 |
| 94 | // upper half-digit |
| 95 | c += x[i]>>W2 + y[i]>>W2 |
| 96 | z[i] = c<<W2 | d |
| 97 | c >>= W2 |
| 98 | } |
| 99 | for i := n; i < m; i++ { |
| 100 | // lower half-digit |
| 101 | c += x[i] & M2 |
| 102 | d := c & M2 |
| 103 | c >>= W2 |
| 104 | // upper half-digit |
| 105 | c += x[i] >> W2 |
| 106 | z[i] = c<<W2 | d |
| 107 | c >>= W2 |
| 108 | } |
| 109 | if c != 0 { |
| 110 | z[m] = c |
| 111 | m++ |
| 112 | } |
| 113 | return z[:m] |
| 114 | } |
| 115 | |
| 116 | // z = x + y |
| 117 | func (z nat) add(x, y nat) nat { |
| 118 | m := len(x) |
| 119 | n := len(y) |
| 120 | |
| 121 | switch { |
| 122 | case m < n: |
| 123 | return z.add(y, x) |
| 124 | case m == 0: |
| 125 | // n == 0 because m >= n; result is 0 |
| 126 | return z.make(0) |
| 127 | case n == 0: |
| 128 | // result is x |
| 129 | return z.set(x) |
| 130 | } |
| 131 | // m >= n > 0 |
| 132 | |
| 133 | z = z.make(m + 1) |
| 134 | var c big.Word |
| 135 | |
| 136 | for i, xi := range x[:n] { |
| 137 | yi := y[i] |
| 138 | zi := xi + yi + c |
| 139 | z[i] = zi |
| 140 | // see "Hacker's Delight", section 2-12 (overflow detection) |
| 141 | c = ((xi & yi) | ((xi | yi) &^ zi)) >> (W - 1) |
| 142 | } |
| 143 | for i, xi := range x[n:] { |
| 144 | zi := xi + c |
| 145 | z[n+i] = zi |
| 146 | c = (xi &^ zi) >> (W - 1) |
| 147 | if c == 0 { |
| 148 | copy(z[n+i+1:], x[i+1:]) |
| 149 | break |
| 150 | } |
| 151 | } |
| 152 | if c != 0 { |
| 153 | z[m] = c |
| 154 | m++ |
| 155 | } |
| 156 | return z[:m] |
| 157 | } |
| 158 | |
| 159 | func bitlen(x big.Word) int { |
| 160 | n := 0 |
| 161 | for x > 0 { |
| 162 | x >>= 1 |
| 163 | n++ |
| 164 | } |
| 165 | return n |
| 166 | } |
| 167 | |
| 168 | func (x nat) bitlen() int { |
| 169 | if i := len(x); i > 0 { |
| 170 | return (i-1)*W + bitlen(x[i-1]) |
| 171 | } |
| 172 | return 0 |
| 173 | } |
| 174 | |
| 175 | func (x nat) String() string { |
| 176 | const shortLen = 10 |
| 177 | s := new(big.Int).SetBits(x).String() |
| 178 | if *short && len(s) > shortLen { |
| 179 | s = s[:shortLen] + "..." |
| 180 | } |
| 181 | return s |
| 182 | } |
| 183 | |
| 184 | func fibo(n int, half, opt bool) nat { |
| 185 | switch n { |
| 186 | case 0: |
| 187 | return nil |
| 188 | case 1: |
| 189 | return nat{1} |
| 190 | } |
| 191 | f0 := nat(nil) |
| 192 | f1 := nat{1} |
| 193 | if half { |
| 194 | if opt { |
| 195 | var f2 nat // reuse f2 |
| 196 | for i := 1; i < n; i++ { |
| 197 | f2 = f2.halfAdd(f1, f0) |
| 198 | f0, f1, f2 = f1, f2, f0 |
| 199 | } |
| 200 | } else { |
| 201 | for i := 1; i < n; i++ { |
| 202 | f2 := nat(nil).halfAdd(f1, f0) // allocate a new f2 each time |
| 203 | f0, f1 = f1, f2 |
| 204 | } |
| 205 | } |
| 206 | } else { |
| 207 | if opt { |
| 208 | var f2 nat // reuse f2 |
| 209 | for i := 1; i < n; i++ { |
| 210 | f2 = f2.add(f1, f0) |
| 211 | f0, f1, f2 = f1, f2, f0 |
| 212 | } |
| 213 | } else { |
| 214 | for i := 1; i < n; i++ { |
| 215 | f2 := nat(nil).add(f1, f0) // allocate a new f2 each time |
| 216 | f0, f1 = f1, f2 |
| 217 | } |
| 218 | } |
| 219 | } |
| 220 | return f1 // was f2 before shuffle |
| 221 | } |
| 222 | |
| 223 | var tests = []struct { |
| 224 | n int |
| 225 | want string |
| 226 | }{ |
| 227 | {0, "0"}, |
| 228 | {1, "1"}, |
| 229 | {2, "1"}, |
| 230 | {3, "2"}, |
| 231 | {4, "3"}, |
| 232 | {5, "5"}, |
| 233 | {6, "8"}, |
| 234 | {7, "13"}, |
| 235 | {8, "21"}, |
| 236 | {9, "34"}, |
| 237 | {10, "55"}, |
| 238 | {100, "354224848179261915075"}, |
| 239 | {1000, "43466557686937456435688527675040625802564660517371780402481729089536555417949051890403879840079255169295922593080322634775209689623239873322471161642996440906533187938298969649928516003704476137795166849228875"}, |
| 240 | } |
| 241 | |
| 242 | func test(half, opt bool) { |
| 243 | for _, test := range tests { |
| 244 | got := fibo(test.n, half, opt).String() |
| 245 | if got != test.want { |
| 246 | fmt.Printf("error: got std fibo(%d) = %s; want %s\n", test.n, got, test.want) |
| 247 | os.Exit(1) |
| 248 | } |
| 249 | } |
| 250 | } |
| 251 | |
| 252 | func selfTest() { |
| 253 | if W != 32 && W != 64 { |
| 254 | fmt.Printf("error: unexpected wordsize %d", W) |
| 255 | os.Exit(1) |
| 256 | } |
| 257 | for i := 0; i < 4; i++ { |
| 258 | test(i&2 == 0, i&1 != 0) |
| 259 | } |
| 260 | } |
| 261 | |
| 262 | func doFibo(n int) { |
| 263 | start := time.Now() |
| 264 | f := fibo(n, *half, *opt) |
| 265 | t := time.Since(start) |
| 266 | fmt.Printf("fibo(%d) = %s (%d bits, %s)\n", n, f, f.bitlen(), t) |
| 267 | } |
| 268 | |
| 269 | func benchFibo(b *testing.B, n int, half, opt bool) { |
| 270 | for i := 0; i < b.N; i++ { |
| 271 | fibo(n, half, opt) |
| 272 | } |
| 273 | } |
| 274 | |
| 275 | func doBench(half, opt bool) { |
| 276 | w := tabwriter.NewWriter(os.Stdout, 0, 8, 2, ' ', tabwriter.AlignRight) |
| 277 | fmt.Fprintf(w, "wordsize = %d, half = %v, opt = %v\n", W, half, opt) |
| 278 | fmt.Fprintf(w, "n\talloc count\talloc bytes\tns/op\ttime/op\t\n") |
| 279 | for n := 1; n <= 1e6; n *= 10 { |
| 280 | res := testing.Benchmark(func(b *testing.B) { benchFibo(b, n, half, opt) }) |
| 281 | fmt.Fprintf(w, "%d\t%d\t%d\t%d\t%s\t\n", n, res.AllocsPerOp(), res.AllocedBytesPerOp(), res.NsPerOp(), time.Duration(res.NsPerOp())) |
| 282 | } |
| 283 | fmt.Fprintln(w) |
| 284 | w.Flush() |
| 285 | } |
| 286 | |
| 287 | func main() { |
| 288 | selfTest() |
| 289 | flag.Parse() |
| 290 | |
| 291 | if args := flag.Args(); len(args) > 0 { |
| 292 | // command-line use |
| 293 | fmt.Printf("half = %v, opt = %v, wordsize = %d bits\n", *half, *opt, W) |
| 294 | for _, arg := range args { |
| 295 | n, err := strconv.Atoi(arg) |
| 296 | if err != nil || n < 0 { |
| 297 | fmt.Println("invalid argument", arg) |
| 298 | continue |
| 299 | } |
| 300 | doFibo(n) |
| 301 | } |
| 302 | return |
| 303 | } |
| 304 | |
| 305 | if *bench { |
| 306 | for i := 0; i < 4; i++ { |
| 307 | doBench(i&2 == 0, i&1 != 0) |
| 308 | } |
| 309 | } |
| 310 | } |