Skip to content

Commit

Permalink
no imath
Browse files Browse the repository at this point in the history
  • Loading branch information
thehowl committed Sep 18, 2023
1 parent e787c76 commit 86c3200
Show file tree
Hide file tree
Showing 12 changed files with 35 additions and 71 deletions.
10 changes: 4 additions & 6 deletions gnovm/stdlibs/math/all_test.gno
Original file line number Diff line number Diff line change
Expand Up @@ -7,8 +7,6 @@ package math
import (
"fmt"
"testing"

imath "internal/math"
)

var _vf = []float64{
Expand Down Expand Up @@ -1191,7 +1189,7 @@ var _vffdimSC = [][2]float64{
}

var (
nan = imath.Float64frombits(0xFFF8000000000000) // SSE2 DIVSD 0/0
nan = Float64frombits(0xFFF8000000000000) // SSE2 DIVSD 0/0
_vffdim2SC = [][2]float64{
{Inf(-1), Inf(-1)},
{Inf(-1), Inf(1)},
Expand Down Expand Up @@ -1982,7 +1980,7 @@ var _vfsqrtSC = []float64{
0,
Inf(1),
NaN(),
imath.Float64frombits(2), // subnormal; see https://golang.org/issue/13013
Float64frombits(2), // subnormal; see https://golang.org/issue/13013
}

var _sqrtSC = []float64{
Expand Down Expand Up @@ -3860,7 +3858,7 @@ func BenchmarkYn(b *testing.B) {
func BenchmarkFloat64bits(b *testing.B) {
y := uint64(0)
for i := 0; i < b.N; i++ {
y = imath.Float64bits(roundNeg)
y = Float64bits(roundNeg)
}
GlobalI = int(y)
}
Expand All @@ -3870,7 +3868,7 @@ var _roundUint64 = uint64(5)
func BenchmarkFloat64frombits(b *testing.B) {
x := 0.0
for i := 0; i < b.N; i++ {
x = imath.Float64frombits(_roundUint64)
x = Float64frombits(_roundUint64)
}
GlobalF = x
}
Expand Down
10 changes: 3 additions & 7 deletions gnovm/stdlibs/math/cbrt.gno
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,6 @@

package math

import (
imath "internal/math"
)

// The go code is a modified version of the original C code from
// http://www.netlib.org/fdlibm/s_cbrt.c and came with this notice.
//
Expand Down Expand Up @@ -55,12 +51,12 @@ func cbrt(x float64) float64 {
}

// rough cbrt to 5 bits
t := imath.Float64frombits(imath.Float64bits(x)/3 + B1<<32)
t := Float64frombits(Float64bits(x)/3 + B1<<32)
if x < SmallestNormal {
// subnormal number
t = float64(1 << 54) // set t= 2**54
t *= x
t = imath.Float64frombits(imath.Float64bits(t)/3 + B2<<32)
t = Float64frombits(Float64bits(t)/3 + B2<<32)
}

// new cbrt to 23 bits
Expand All @@ -69,7 +65,7 @@ func cbrt(x float64) float64 {
t *= G + F/(s+E+D/s)

// chop to 22 bits, make larger than cbrt(x)
t = imath.Float64frombits(imath.Float64bits(t)&(0xFFFFFFFFC<<28) + 1<<30)
t = Float64frombits(Float64bits(t)&(0xFFFFFFFFC<<28) + 1<<30)

// one step newton iteration to 53 bits with error less than 0.667ulps
s = t * t // t*t is exact
Expand Down
6 changes: 2 additions & 4 deletions gnovm/stdlibs/math/erf.gno
Original file line number Diff line number Diff line change
Expand Up @@ -4,8 +4,6 @@

package math

import imath "internal/math"

/*
Floating-point error function and complementary error function.
*/
Expand Down Expand Up @@ -255,7 +253,7 @@ func erf(x float64) float64 {
R = rb0 + s*(rb1+s*(rb2+s*(rb3+s*(rb4+s*(rb5+s*rb6)))))
S = 1 + s*(sb1+s*(sb2+s*(sb3+s*(sb4+s*(sb5+s*(sb6+s*sb7))))))
}
z := imath.Float64frombits(imath.Float64bits(x) & 0xffffffff00000000) // pseudo-single (20-bit) precision x
z := Float64frombits(Float64bits(x) & 0xffffffff00000000) // pseudo-single (20-bit) precision x
r := Exp(-z*z-0.5625) * Exp((z-x)*(z+x)+R/S)
if sign {
return r/x - 1
Expand Down Expand Up @@ -333,7 +331,7 @@ func erfc(x float64) float64 {
R = rb0 + s*(rb1+s*(rb2+s*(rb3+s*(rb4+s*(rb5+s*rb6)))))
S = 1 + s*(sb1+s*(sb2+s*(sb3+s*(sb4+s*(sb5+s*(sb6+s*sb7))))))
}
z := imath.Float64frombits(imath.Float64bits(x) & 0xffffffff00000000) // pseudo-single (20-bit) precision x
z := Float64frombits(Float64bits(x) & 0xffffffff00000000) // pseudo-single (20-bit) precision x
r := Exp(-z*z-0.5625) * Exp((z-x)*(z+x)+R/S)
if sign {
return 2 - r/x
Expand Down
12 changes: 5 additions & 7 deletions gnovm/stdlibs/math/expm1.gno
Original file line number Diff line number Diff line change
Expand Up @@ -4,8 +4,6 @@

package math

import imath "internal/math"

// The original C code, the long comment, and the constants
// below are from FreeBSD's /usr/src/lib/msun/src/s_expm1.c
// and came with this notice. The go code is a simplified
Expand Down Expand Up @@ -226,18 +224,18 @@ func expm1(x float64) float64 {
return 1 + 2*(x-e)
case k <= -2 || k > 56: // suffice to return exp(x)-1
y := 1 - (e - x)
y = imath.Float64frombits(imath.Float64bits(y) + uint64(k)<<52) // add k to y's exponent
y = Float64frombits(Float64bits(y) + uint64(k)<<52) // add k to y's exponent
return y - 1
}
if k < 20 {
t := imath.Float64frombits(0x3ff0000000000000 - (0x20000000000000 >> uint(k))) // t=1-2**-k
t := Float64frombits(0x3ff0000000000000 - (0x20000000000000 >> uint(k))) // t=1-2**-k
y := t - (e - x)
y = imath.Float64frombits(imath.Float64bits(y) + uint64(k)<<52) // add k to y's exponent
y = Float64frombits(Float64bits(y) + uint64(k)<<52) // add k to y's exponent
return y
}
t = imath.Float64frombits(uint64(0x3ff-k) << 52) // 2**-k
t = Float64frombits(uint64(0x3ff-k) << 52) // 2**-k
y := x - (e + t)
y += 1
y = imath.Float64frombits(imath.Float64bits(y) + uint64(k)<<52) // add k to y's exponent
y = Float64frombits(Float64bits(y) + uint64(k)<<52) // add k to y's exponent
return y
}
8 changes: 3 additions & 5 deletions gnovm/stdlibs/math/fma.gno
Original file line number Diff line number Diff line change
Expand Up @@ -6,8 +6,6 @@ package math

import (
"math/bits"

imath "internal/math"
)

func zero(x uint64) uint64 {
Expand Down Expand Up @@ -97,7 +95,7 @@ func split(b uint64) (sign uint32, exp int32, mantissa uint64) {
// FMA returns x * y + z, computed with only one rounding.
// (That is, FMA returns the fused multiply-add of x, y, and z.)
func FMA(x, y, z float64) float64 {
bx, by, bz := imath.Float64bits(x), imath.Float64bits(y), imath.Float64bits(z)
bx, by, bz := Float64bits(x), Float64bits(y), Float64bits(z)

// Inf or NaN or zero involved. At most one rounding will occur.
if x == 0.0 || y == 0.0 || z == 0.0 || bx&uvinf == uvinf || by&uvinf == uvinf {
Expand Down Expand Up @@ -161,7 +159,7 @@ func FMA(x, y, z float64) float64 {
// Round and break ties to even
if pe > 1022+bias || pe == 1022+bias && (m+1<<9)>>63 == 1 {
// rounded value overflows exponent range
return imath.Float64frombits(uint64(ps)<<63 | uvinf)
return Float64frombits(uint64(ps)<<63 | uvinf)
}
if pe < 0 {
n := uint(-pe)
Expand All @@ -170,5 +168,5 @@ func FMA(x, y, z float64) float64 {
}
m = ((m + 1<<9) >> 10) & ^zero((m&(1<<10-1))^1<<9)
pe &= -int32(nonzero(m))
return imath.Float64frombits(uint64(ps)<<63 + uint64(pe)<<52 + m)
return Float64frombits(uint64(ps)<<63 + uint64(pe)<<52 + m)
}
8 changes: 2 additions & 6 deletions gnovm/stdlibs/math/frexp.gno
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,6 @@

package math

import (
imath "internal/math"
)

// Frexp breaks f into a normalized fraction
// and an integral power of two.
// It returns frac and exp satisfying f == frac × 2**exp,
Expand All @@ -31,10 +27,10 @@ func frexp(f float64) (frac float64, exp int) {
return f, 0
}
f, exp = normalize(f)
x := imath.Float64bits(f)
x := Float64bits(f)
exp += int((x>>shift)&mask) - bias + 1
x &^= mask << shift
x |= (-1 + bias) << shift
frac = imath.Float64frombits(x)
frac = Float64frombits(x)
return
}
4 changes: 1 addition & 3 deletions gnovm/stdlibs/math/lgamma.gno
Original file line number Diff line number Diff line change
Expand Up @@ -4,8 +4,6 @@

package math

import imath "internal/math"

/*
Floating-point logarithm of the Gamma function.
*/
Expand Down Expand Up @@ -353,7 +351,7 @@ func sinPi(x float64) float64 {
if x < Two52 {
z = x + Two52 // exact
}
n = int(1 & imath.Float64bits(z))
n = int(1 & Float64bits(z))
x = float64(n)
n <<= 2
}
Expand Down
12 changes: 4 additions & 8 deletions gnovm/stdlibs/math/log1p.gno
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,6 @@

package math

import (
imath "internal/math"
)

// The original C code, the long comment, and the constants
// below are from FreeBSD's /usr/src/lib/msun/src/s_log1p.c
// and came with this notice. The go code is a simplified
Expand Down Expand Up @@ -153,7 +149,7 @@ func log1p(x float64) float64 {
var u float64
if absx < Two53 { // 1<<53
u = 1.0 + x
iu = imath.Float64bits(u)
iu = Float64bits(u)
k = int((iu >> 52) - 1023)
// correction term
if k > 0 {
Expand All @@ -164,16 +160,16 @@ func log1p(x float64) float64 {
c /= u
} else {
u = x
iu = imath.Float64bits(u)
iu = Float64bits(u)
k = int((iu >> 52) - 1023)
c = 0
}
iu &= 0x000fffffffffffff
if iu < 0x0006a09e667f3bcd { // mantissa of Sqrt(2)
u = imath.Float64frombits(iu | 0x3ff0000000000000) // normalize u
u = Float64frombits(iu | 0x3ff0000000000000) // normalize u
} else {
k++
u = imath.Float64frombits(iu | 0x3fe0000000000000) // normalize u/2
u = Float64frombits(iu | 0x3fe0000000000000) // normalize u/2
iu = (0x0010000000000000 - iu) >> 2
}
f = u - 1.0 // Sqrt(2)/2 < u < Sqrt(2)
Expand Down
6 changes: 1 addition & 5 deletions gnovm/stdlibs/math/logb.gno
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,6 @@

package math

import (
imath "internal/math"
)

// Logb returns the binary exponent of x.
//
// Special cases are:
Expand Down Expand Up @@ -52,5 +48,5 @@ func Ilogb(x float64) int {
// non-zero.
func ilogb(x float64) int {
x, exp := normalize(x)
return int((imath.Float64bits(x)>>shift)&mask) - bias + exp
return int((Float64bits(x)>>shift)&mask) - bias + exp
}
16 changes: 6 additions & 10 deletions gnovm/stdlibs/math/nextafter.gno
Original file line number Diff line number Diff line change
Expand Up @@ -4,10 +4,6 @@

package math

import (
imath "internal/math"
)

// Nextafter32 returns the next representable float32 value after x towards y.
//
// Special cases are:
Expand All @@ -22,11 +18,11 @@ func Nextafter32(x, y float32) (r float32) {
case x == y:
r = x
case x == 0:
r = float32(Copysign(float64(imath.Float32frombits(1)), float64(y)))
r = float32(Copysign(float64(Float32frombits(1)), float64(y)))
case (y > x) == (x > 0):
r = imath.Float32frombits(imath.Float32bits(x) + 1)
r = Float32frombits(Float32bits(x) + 1)
default:
r = imath.Float32frombits(imath.Float32bits(x) - 1)
r = Float32frombits(Float32bits(x) - 1)
}
return
}
Expand All @@ -45,11 +41,11 @@ func Nextafter(x, y float64) (r float64) {
case x == y:
r = x
case x == 0:
r = Copysign(imath.Float64frombits(1), y)
r = Copysign(Float64frombits(1), y)
case (y > x) == (x > 0):
r = imath.Float64frombits(imath.Float64bits(x) + 1)
r = Float64frombits(Float64bits(x) + 1)
default:
r = imath.Float64frombits(imath.Float64bits(x) - 1)
r = Float64frombits(Float64bits(x) - 1)
}
return
}
8 changes: 2 additions & 6 deletions gnovm/stdlibs/math/sqrt.gno
Original file line number Diff line number Diff line change
Expand Up @@ -5,10 +5,6 @@

package math

import (
imath "internal/math"
)

// The original C code and the long comment below are
// from FreeBSD's /usr/src/lib/msun/src/e_sqrt.c and
// came with this notice. The go code is a simplified
Expand Down Expand Up @@ -110,7 +106,7 @@ func sqrt(x float64) float64 {
case x < 0:
return NaN()
}
ix := imath.Float64bits(x)
ix := Float64bits(x)
// normalize x
exp := int((ix >> shift) & mask)
if exp == 0 { // subnormal x
Expand Down Expand Up @@ -146,5 +142,5 @@ func sqrt(x float64) float64 {
q += q & 1 // round according to extra bit
}
ix = q>>1 + uint64(exp-1+bias)<<shift // significand + biased exponent
return imath.Float64frombits(ix)
return Float64frombits(ix)
}
6 changes: 2 additions & 4 deletions gnovm/stdlibs/math/trig_reduce.gno
Original file line number Diff line number Diff line change
Expand Up @@ -6,8 +6,6 @@ package math

import (
"math/bits"

imath "internal/math"
)

// reduceThreshold is the maximum value of x where the reduction using Pi/4
Expand Down Expand Up @@ -37,7 +35,7 @@ func trigReduce(x float64) (j uint64, z float64) {
}
// Extract out the integer and exponent such that,
// x = ix * 2 ** exp.
ix := imath.Float64bits(x)
ix := Float64bits(x)
exp := int(ix>>shift&mask) - bias - shift
ix &^= mask << shift
ix |= 1 << shift
Expand Down Expand Up @@ -65,7 +63,7 @@ func trigReduce(x float64) (j uint64, z float64) {
hi >>= 64 - shift
// Include the exponent and convert to a float.
hi |= e << shift
z = imath.Float64frombits(hi)
z = Float64frombits(hi)
// Map zeros to origin.
if j&1 == 1 {
j++
Expand Down

0 comments on commit 86c3200

Please sign in to comment.