github.com/bgentry/go@v0.0.0-20150121062915-6cf5a733d54d/src/math/big/arith.go (about) 1 // Copyright 2009 The Go Authors. All rights reserved. 2 // Use of this source code is governed by a BSD-style 3 // license that can be found in the LICENSE file. 4 5 // This file provides Go implementations of elementary multi-precision 6 // arithmetic operations on word vectors. Needed for platforms without 7 // assembly implementations of these routines. 8 9 package big 10 11 // A Word represents a single digit of a multi-precision unsigned integer. 12 type Word uintptr 13 14 const ( 15 // Compute the size _S of a Word in bytes. 16 _m = ^Word(0) 17 _logS = _m>>8&1 + _m>>16&1 + _m>>32&1 18 _S = 1 << _logS 19 20 _W = _S << 3 // word size in bits 21 _B = 1 << _W // digit base 22 _M = _B - 1 // digit mask 23 24 _W2 = _W / 2 // half word size in bits 25 _B2 = 1 << _W2 // half digit base 26 _M2 = _B2 - 1 // half digit mask 27 ) 28 29 // ---------------------------------------------------------------------------- 30 // Elementary operations on words 31 // 32 // These operations are used by the vector operations below. 33 34 // z1<<_W + z0 = x+y+c, with c == 0 or 1 35 func addWW_g(x, y, c Word) (z1, z0 Word) { 36 yc := y + c 37 z0 = x + yc 38 if z0 < x || yc < y { 39 z1 = 1 40 } 41 return 42 } 43 44 // z1<<_W + z0 = x-y-c, with c == 0 or 1 45 func subWW_g(x, y, c Word) (z1, z0 Word) { 46 yc := y + c 47 z0 = x - yc 48 if z0 > x || yc < y { 49 z1 = 1 50 } 51 return 52 } 53 54 // z1<<_W + z0 = x*y 55 // Adapted from Warren, Hacker's Delight, p. 132. 56 func mulWW_g(x, y Word) (z1, z0 Word) { 57 x0 := x & _M2 58 x1 := x >> _W2 59 y0 := y & _M2 60 y1 := y >> _W2 61 w0 := x0 * y0 62 t := x1*y0 + w0>>_W2 63 w1 := t & _M2 64 w2 := t >> _W2 65 w1 += x0 * y1 66 z1 = x1*y1 + w2 + w1>>_W2 67 z0 = x * y 68 return 69 } 70 71 // z1<<_W + z0 = x*y + c 72 func mulAddWWW_g(x, y, c Word) (z1, z0 Word) { 73 z1, zz0 := mulWW_g(x, y) 74 if z0 = zz0 + c; z0 < zz0 { 75 z1++ 76 } 77 return 78 } 79 80 // Length of x in bits. 81 func bitLen_g(x Word) (n int) { 82 for ; x >= 0x8000; x >>= 16 { 83 n += 16 84 } 85 if x >= 0x80 { 86 x >>= 8 87 n += 8 88 } 89 if x >= 0x8 { 90 x >>= 4 91 n += 4 92 } 93 if x >= 0x2 { 94 x >>= 2 95 n += 2 96 } 97 if x >= 0x1 { 98 n++ 99 } 100 return 101 } 102 103 // log2 computes the integer binary logarithm of x. 104 // The result is the integer n for which 2^n <= x < 2^(n+1). 105 // If x == 0, the result is -1. 106 func log2(x Word) int { 107 return bitLen(x) - 1 108 } 109 110 // Number of leading zeros in x. 111 func leadingZeros(x Word) uint { 112 return uint(_W - bitLen(x)) 113 } 114 115 // q = (u1<<_W + u0 - r)/y 116 // Adapted from Warren, Hacker's Delight, p. 152. 117 func divWW_g(u1, u0, v Word) (q, r Word) { 118 if u1 >= v { 119 return 1<<_W - 1, 1<<_W - 1 120 } 121 122 s := leadingZeros(v) 123 v <<= s 124 125 vn1 := v >> _W2 126 vn0 := v & _M2 127 un32 := u1<<s | u0>>(_W-s) 128 un10 := u0 << s 129 un1 := un10 >> _W2 130 un0 := un10 & _M2 131 q1 := un32 / vn1 132 rhat := un32 - q1*vn1 133 134 for q1 >= _B2 || q1*vn0 > _B2*rhat+un1 { 135 q1-- 136 rhat += vn1 137 if rhat >= _B2 { 138 break 139 } 140 } 141 142 un21 := un32*_B2 + un1 - q1*v 143 q0 := un21 / vn1 144 rhat = un21 - q0*vn1 145 146 for q0 >= _B2 || q0*vn0 > _B2*rhat+un0 { 147 q0-- 148 rhat += vn1 149 if rhat >= _B2 { 150 break 151 } 152 } 153 154 return q1*_B2 + q0, (un21*_B2 + un0 - q0*v) >> s 155 } 156 157 // Keep for performance debugging. 158 // Using addWW_g is likely slower. 159 const use_addWW_g = false 160 161 // The resulting carry c is either 0 or 1. 162 func addVV_g(z, x, y []Word) (c Word) { 163 if use_addWW_g { 164 for i := range z { 165 c, z[i] = addWW_g(x[i], y[i], c) 166 } 167 return 168 } 169 170 for i, xi := range x[:len(z)] { 171 yi := y[i] 172 zi := xi + yi + c 173 z[i] = zi 174 // see "Hacker's Delight", section 2-12 (overflow detection) 175 c = (xi&yi | (xi|yi)&^zi) >> (_W - 1) 176 } 177 return 178 } 179 180 // The resulting carry c is either 0 or 1. 181 func subVV_g(z, x, y []Word) (c Word) { 182 if use_addWW_g { 183 for i := range z { 184 c, z[i] = subWW_g(x[i], y[i], c) 185 } 186 return 187 } 188 189 for i, xi := range x[:len(z)] { 190 yi := y[i] 191 zi := xi - yi - c 192 z[i] = zi 193 // see "Hacker's Delight", section 2-12 (overflow detection) 194 c = (yi&^xi | (yi|^xi)&zi) >> (_W - 1) 195 } 196 return 197 } 198 199 // Argument y must be either 0 or 1. 200 // The resulting carry c is either 0 or 1. 201 func addVW_g(z, x []Word, y Word) (c Word) { 202 if use_addWW_g { 203 c = y 204 for i := range z { 205 c, z[i] = addWW_g(x[i], c, 0) 206 } 207 return 208 } 209 210 c = y 211 for i, xi := range x[:len(z)] { 212 zi := xi + c 213 z[i] = zi 214 c = xi &^ zi >> (_W - 1) 215 } 216 return 217 } 218 219 func subVW_g(z, x []Word, y Word) (c Word) { 220 if use_addWW_g { 221 c = y 222 for i := range z { 223 c, z[i] = subWW_g(x[i], c, 0) 224 } 225 return 226 } 227 228 c = y 229 for i, xi := range x[:len(z)] { 230 zi := xi - c 231 z[i] = zi 232 c = (zi &^ xi) >> (_W - 1) 233 } 234 return 235 } 236 237 func shlVU_g(z, x []Word, s uint) (c Word) { 238 if n := len(z); n > 0 { 239 ŝ := _W - s 240 w1 := x[n-1] 241 c = w1 >> ŝ 242 for i := n - 1; i > 0; i-- { 243 w := w1 244 w1 = x[i-1] 245 z[i] = w<<s | w1>>ŝ 246 } 247 z[0] = w1 << s 248 } 249 return 250 } 251 252 func shrVU_g(z, x []Word, s uint) (c Word) { 253 if n := len(z); n > 0 { 254 ŝ := _W - s 255 w1 := x[0] 256 c = w1 << ŝ 257 for i := 0; i < n-1; i++ { 258 w := w1 259 w1 = x[i+1] 260 z[i] = w>>s | w1<<ŝ 261 } 262 z[n-1] = w1 >> s 263 } 264 return 265 } 266 267 func mulAddVWW_g(z, x []Word, y, r Word) (c Word) { 268 c = r 269 for i := range z { 270 c, z[i] = mulAddWWW_g(x[i], y, c) 271 } 272 return 273 } 274 275 // TODO(gri) Remove use of addWW_g here and then we can remove addWW_g and subWW_g. 276 func addMulVVW_g(z, x []Word, y Word) (c Word) { 277 for i := range z { 278 z1, z0 := mulAddWWW_g(x[i], y, z[i]) 279 c, z[i] = addWW_g(z0, c, 0) 280 c += z1 281 } 282 return 283 } 284 285 func divWVW_g(z []Word, xn Word, x []Word, y Word) (r Word) { 286 r = xn 287 for i := len(z) - 1; i >= 0; i-- { 288 z[i], r = divWW_g(r, x[i], y) 289 } 290 return 291 }