magic.go

     1  // Copyright 2016 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  package ssacompile
     6  
     7  // So you want to compute x / c for some constant c?
     8  // Machine division instructions are slow, so we try to
     9  // compute this division with a multiplication + a few
    10  // other cheap instructions instead.
    11  // (We assume here that c != 0, +/- 1, or +/- 2^i.  Those
    12  // cases are easy to handle in different ways).
    13  
    14  // Technique from https://gmplib.org/~tege/divcnst-pldi94.pdf
    15  
    16  // First consider unsigned division.
    17  // Our strategy is to precompute 1/c then do
    18  //   ⎣x / c⎦ = ⎣x * (1/c)⎦.
    19  // 1/c is less than 1, so we can't compute it directly in
    20  // integer arithmetic.  Let's instead compute 2^e/c
    21  // for a value of e TBD (^ = exponentiation).  Then
    22  //   ⎣x / c⎦ = ⎣x * (2^e/c) / 2^e⎦.
    23  // Dividing by 2^e is easy.  2^e/c isn't an integer, unfortunately.
    24  // So we must approximate it.  Let's call its approximation m.
    25  // We'll then compute
    26  //   ⎣x * m / 2^e⎦
    27  // Which we want to be equal to ⎣x / c⎦ for 0 <= x < 2^n-1
    28  // where n is the word size.
    29  // Setting x = c gives us c * m >= 2^e.
    30  // We'll chose m = ⎡2^e/c⎤ to satisfy that equation.
    31  // What remains is to choose e.
    32  // Let m = 2^e/c + delta, 0 <= delta < 1
    33  //   ⎣x * (2^e/c + delta) / 2^e⎦
    34  //   ⎣x / c + x * delta / 2^e⎦
    35  // We must have x * delta / 2^e < 1/c so that this
    36  // additional term never rounds differently than ⎣x / c⎦ does.
    37  // Rearranging,
    38  //   2^e > x * delta * c
    39  // x can be at most 2^n-1 and delta can be at most 1.
    40  // So it is sufficient to have 2^e >= 2^n*c.
    41  // So we'll choose e = n + s, with s = ⎡log2(c)⎤.
    42  //
    43  // An additional complication arises because m has n+1 bits in it.
    44  // Hardware restricts us to n bit by n bit multiplies.
    45  // We divide into 3 cases:
    46  //
    47  // Case 1: m is even.
    48  //   ⎣x / c⎦ = ⎣x * m / 2^(n+s)⎦
    49  //   ⎣x / c⎦ = ⎣x * (m/2) / 2^(n+s-1)⎦
    50  //   ⎣x / c⎦ = ⎣x * (m/2) / 2^n / 2^(s-1)⎦
    51  //   ⎣x / c⎦ = ⎣⎣x * (m/2) / 2^n⎦ / 2^(s-1)⎦
    52  //   multiply + shift
    53  //
    54  // Case 2: c is even.
    55  //   ⎣x / c⎦ = ⎣(x/2) / (c/2)⎦
    56  //   ⎣x / c⎦ = ⎣⎣x/2⎦ / (c/2)⎦
    57  //     This is just the original problem, with x' = ⎣x/2⎦, c' = c/2, n' = n-1.
    58  //       s' = s-1
    59  //       m' = ⎡2^(n'+s')/c'⎤
    60  //          = ⎡2^(n+s-1)/c⎤
    61  //          = ⎡m/2⎤
    62  //   ⎣x / c⎦ = ⎣x' * m' / 2^(n'+s')⎦
    63  //   ⎣x / c⎦ = ⎣⎣x/2⎦ * ⎡m/2⎤ / 2^(n+s-2)⎦
    64  //   ⎣x / c⎦ = ⎣⎣⎣x/2⎦ * ⎡m/2⎤ / 2^n⎦ / 2^(s-2)⎦
    65  //   shift + multiply + shift
    66  //
    67  // Case 3: everything else
    68  //   let k = m - 2^n. k fits in n bits.
    69  //   ⎣x / c⎦ = ⎣x * m / 2^(n+s)⎦
    70  //   ⎣x / c⎦ = ⎣x * (2^n + k) / 2^(n+s)⎦
    71  //   ⎣x / c⎦ = ⎣(x + x * k / 2^n) / 2^s⎦
    72  //   ⎣x / c⎦ = ⎣(x + ⎣x * k / 2^n⎦) / 2^s⎦
    73  //   ⎣x / c⎦ = ⎣(x + ⎣x * k / 2^n⎦) / 2^s⎦
    74  //   ⎣x / c⎦ = ⎣⎣(x + ⎣x * k / 2^n⎦) / 2⎦ / 2^(s-1)⎦
    75  //   multiply + avg + shift
    76  //
    77  // These can be implemented in hardware using:
    78  //  ⎣a * b / 2^n⎦ - aka high n bits of an n-bit by n-bit multiply.
    79  //  ⎣(a+b) / 2⎦   - aka "average" of two n-bit numbers.
    80  //                  (Not just a regular add & shift because the intermediate result
    81  //                   a+b has n+1 bits in it.  Nevertheless, can be done
    82  //                   in 2 instructions on x86.)
    83  
    84  // umagicOK reports whether we should strength reduce a n-bit divide by c.
    85  func umagicOK(n uint, c int64) bool {
    86  	// Convert from ConstX auxint values to the real uint64 constant they represent.
    87  	d := uint64(c) << (64 - n) >> (64 - n)
    88  
    89  	// Doesn't work for 0.
    90  	// Don't use for powers of 2.
    91  	return d&(d-1) != 0
    92  }
    93
View as plain text