libgo/go/math/big/arith.go

   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 func mulWW(x, y Word) (z1, z0 Word) { return mulWW_g(x, y) }
  56
  57 // Adapted from Warren, Hacker's Delight, p. 132.
  58 func mulWW_g(x, y Word) (z1, z0 Word) {
  59         x0 := x & _M2
  60         x1 := x >> _W2
  61         y0 := y & _M2
  62         y1 := y >> _W2
  63         w0 := x0 * y0
  64         t := x1*y0 + w0>>_W2
  65         w1 := t & _M2
  66         w2 := t >> _W2
  67         w1 += x0 * y1
  68         z1 = x1*y1 + w2 + w1>>_W2
  69         z0 = x * y
  70         return
  71 }
  72
  73 // z1<<_W + z0 = x*y + c
  74 func mulAddWWW_g(x, y, c Word) (z1, z0 Word) {
  75         z1, zz0 := mulWW(x, y)
  76         if z0 = zz0 + c; z0 < zz0 {
  77                 z1++
  78         }
  79         return
  80 }
  81
  82 // Length of x in bits.
  83 func bitLen(x Word) (n int) { return bitLen_g(x) }
  84 func bitLen_g(x Word) (n int) {
  85         for ; x >= 0x8000; x >>= 16 {
  86                 n += 16
  87         }
  88         if x >= 0x80 {
  89                 x >>= 8
  90                 n += 8
  91         }
  92         if x >= 0x8 {
  93                 x >>= 4
  94                 n += 4
  95         }
  96         if x >= 0x2 {
  97                 x >>= 2
  98                 n += 2
  99         }
 100         if x >= 0x1 {
 101                 n++
 102         }
 103         return
 104 }
 105
 106 // log2 computes the integer binary logarithm of x.
 107 // The result is the integer n for which 2^n <= x < 2^(n+1).
 108 // If x == 0, the result is -1.
 109 func log2(x Word) int {
 110         return bitLen(x) - 1
 111 }
 112
 113 // Number of leading zeros in x.
 114 func leadingZeros(x Word) uint {
 115         return uint(_W - bitLen(x))
 116 }
 117
 118 // q = (u1<<_W + u0 - r)/y
 119 func divWW(x1, x0, y Word) (q, r Word) { return divWW_g(x1, x0, y) }
 120
 121 // Adapted from Warren, Hacker's Delight, p. 152.
 122 func divWW_g(u1, u0, v Word) (q, r Word) {
 123         if u1 >= v {
 124                 return 1<<_W - 1, 1<<_W - 1
 125         }
 126
 127         s := leadingZeros(v)
 128         v <<= s
 129
 130         vn1 := v >> _W2
 131         vn0 := v & _M2
 132         un32 := u1<<s | u0>>(_W-s)
 133         un10 := u0 << s
 134         un1 := un10 >> _W2
 135         un0 := un10 & _M2
 136         q1 := un32 / vn1
 137         rhat := un32 - q1*vn1
 138
 139 again1:
 140         if q1 >= _B2 || q1*vn0 > _B2*rhat+un1 {
 141                 q1--
 142                 rhat += vn1
 143                 if rhat < _B2 {
 144                         goto again1
 145                 }
 146         }
 147
 148         un21 := un32*_B2 + un1 - q1*v
 149         q0 := un21 / vn1
 150         rhat = un21 - q0*vn1
 151
 152 again2:
 153         if q0 >= _B2 || q0*vn0 > _B2*rhat+un0 {
 154                 q0--
 155                 rhat += vn1
 156                 if rhat < _B2 {
 157                         goto again2
 158                 }
 159         }
 160
 161         return q1*_B2 + q0, (un21*_B2 + un0 - q0*v) >> s
 162 }
 163
 164 func addVV(z, x, y []Word) (c Word) { return addVV_g(z, x, y) }
 165 func addVV_g(z, x, y []Word) (c Word) {
 166         for i := range z {
 167                 c, z[i] = addWW_g(x[i], y[i], c)
 168         }
 169         return
 170 }
 171
 172 func subVV(z, x, y []Word) (c Word) { return subVV_g(z, x, y) }
 173 func subVV_g(z, x, y []Word) (c Word) {
 174         for i := range z {
 175                 c, z[i] = subWW_g(x[i], y[i], c)
 176         }
 177         return
 178 }
 179
 180 func addVW(z, x []Word, y Word) (c Word) { return addVW_g(z, x, y) }
 181 func addVW_g(z, x []Word, y Word) (c Word) {
 182         c = y
 183         for i := range z {
 184                 c, z[i] = addWW_g(x[i], c, 0)
 185         }
 186         return
 187 }
 188
 189 func subVW(z, x []Word, y Word) (c Word) { return subVW_g(z, x, y) }
 190 func subVW_g(z, x []Word, y Word) (c Word) {
 191         c = y
 192         for i := range z {
 193                 c, z[i] = subWW_g(x[i], c, 0)
 194         }
 195         return
 196 }
 197
 198 func shlVU(z, x []Word, s uint) (c Word) { return shlVU_g(z, x, s) }
 199 func shlVU_g(z, x []Word, s uint) (c Word) {
 200         if n := len(z); n > 0 {
 201                 ŝ := _W - s
 202                 w1 := x[n-1]
 203                 c = w1 >> ŝ
 204                 for i := n - 1; i > 0; i-- {
 205                         w := w1
 206                         w1 = x[i-1]
 207                         z[i] = w<<s | w1>>ŝ
 208                 }
 209                 z[0] = w1 << s
 210         }
 211         return
 212 }
 213
 214 func shrVU(z, x []Word, s uint) (c Word) { return shrVU_g(z, x, s) }
 215 func shrVU_g(z, x []Word, s uint) (c Word) {
 216         if n := len(z); n > 0 {
 217                 ŝ := _W - s
 218                 w1 := x[0]
 219                 c = w1 << ŝ
 220                 for i := 0; i < n-1; i++ {
 221                         w := w1
 222                         w1 = x[i+1]
 223                         z[i] = w>>s | w1<<ŝ
 224                 }
 225                 z[n-1] = w1 >> s
 226         }
 227         return
 228 }
 229
 230 func mulAddVWW(z, x []Word, y, r Word) (c Word) { return mulAddVWW_g(z, x, y, r) }
 231 func mulAddVWW_g(z, x []Word, y, r Word) (c Word) {
 232         c = r
 233         for i := range z {
 234                 c, z[i] = mulAddWWW_g(x[i], y, c)
 235         }
 236         return
 237 }
 238
 239 func addMulVVW(z, x []Word, y Word) (c Word) { return addMulVVW_g(z, x, y) }
 240 func addMulVVW_g(z, x []Word, y Word) (c Word) {
 241         for i := range z {
 242                 z1, z0 := mulAddWWW_g(x[i], y, z[i])
 243                 c, z[i] = addWW_g(z0, c, 0)
 244                 c += z1
 245         }
 246         return
 247 }
 248
 249 func divWVW(z []Word, xn Word, x []Word, y Word) (r Word) { return divWVW_g(z, xn, x, y) }
 250 func divWVW_g(z []Word, xn Word, x []Word, y Word) (r Word) {
 251         r = xn
 252         for i := len(z) - 1; i >= 0; i-- {
 253                 z[i], r = divWW_g(r, x[i], y)
 254         }
 255         return
 256 }