Minor cleanup to math constants

This commit is contained in:
gingerBill
2021-11-21 14:06:32 +00:00
parent 8a2c829e07
commit f40f12d480
4 changed files with 36 additions and 39 deletions
+14 -17
View File
@@ -607,9 +607,9 @@ floor_mod :: proc "contextless" (x, y: $T) -> T
} }
modf_f16 :: proc "contextless" (x: f16) -> (int: f16, frac: f16) { modf_f16 :: proc "contextless" (x: f16) -> (int: f16, frac: f16) {
shift :: 16 - 5 - 1 shift :: F16_SHIFT
mask :: 0x1f mask :: F16_MASK
bias :: 15 bias :: F16_BIAS
if x < 1 { if x < 1 {
switch { switch {
@@ -641,9 +641,9 @@ modf_f16be :: proc "contextless" (x: f16be) -> (int: f16be, frac: f16be) {
return f16be(i), f16be(f) return f16be(i), f16be(f)
} }
modf_f32 :: proc "contextless" (x: f32) -> (int: f32, frac: f32) { modf_f32 :: proc "contextless" (x: f32) -> (int: f32, frac: f32) {
shift :: 32 - 8 - 1 shift :: F32_SHIFT
mask :: 0xff mask :: F32_MASK
bias :: 127 bias :: F32_BIAS
if x < 1 { if x < 1 {
switch { switch {
@@ -674,10 +674,10 @@ modf_f32be :: proc "contextless" (x: f32be) -> (int: f32be, frac: f32be) {
i, f := #force_inline modf_f32(f32(x)) i, f := #force_inline modf_f32(f32(x))
return f32be(i), f32be(f) return f32be(i), f32be(f)
} }
modf_f64 :: proc "contextless" (x: f64) -> (int: f64, frac: f64) { modf_f64 :: proc "contextless" (x: f64) -> (int: f64, frac: f64) {
shift :: 64 - 11 - 1 shift :: F64_SHIFT
mask :: 0x7ff mask :: F64_MASK
bias :: 1023 bias :: F64_BIAS
if x < 1 { if x < 1 {
switch { switch {
@@ -708,7 +708,7 @@ modf_f64be :: proc "contextless" (x: f64be) -> (int: f64be, frac: f64be) {
i, f := #force_inline modf_f64(f64(x)) i, f := #force_inline modf_f64(f64(x))
return f64be(i), f64be(f) return f64be(i), f64be(f)
} }
modf :: proc{ modf :: proc{
modf_f16, modf_f16le, modf_f16be, modf_f16, modf_f16le, modf_f16be,
modf_f32, modf_f32le, modf_f32be, modf_f32, modf_f32le, modf_f32be,
modf_f64, modf_f64le, modf_f64be, modf_f64, modf_f64le, modf_f64be,
@@ -1127,13 +1127,11 @@ inf_f32be :: proc "contextless" (sign: int) -> f32be {
return f32be(inf_f64(sign)) return f32be(inf_f64(sign))
} }
inf_f64 :: proc "contextless" (sign: int) -> f64 { inf_f64 :: proc "contextless" (sign: int) -> f64 {
v: u64
if sign >= 0 { if sign >= 0 {
v = 0x7ff00000_00000000 return 0h7ff00000_00000000
} else { } else {
v = 0xfff00000_00000000 return 0hfff00000_00000000
} }
return transmute(f64)v
} }
inf_f64le :: proc "contextless" (sign: int) -> f64le { inf_f64le :: proc "contextless" (sign: int) -> f64le {
return f64le(inf_f64(sign)) return f64le(inf_f64(sign))
@@ -1161,8 +1159,7 @@ nan_f32be :: proc "contextless" () -> f32be {
return f32be(nan_f64()) return f32be(nan_f64())
} }
nan_f64 :: proc "contextless" () -> f64 { nan_f64 :: proc "contextless" () -> f64 {
v: u64 = 0x7ff80000_00000001 return 0h7ff80000_00000001
return transmute(f64)v
} }
nan_f64le :: proc "contextless" () -> f64le { nan_f64le :: proc "contextless" () -> f64le {
return f64le(nan_f64()) return f64le(nan_f64())
+9 -9
View File
@@ -68,17 +68,17 @@ package math
@(private="file") @(private="file")
stirling :: proc "contextless" (x: f64) -> (f64, f64) { stirling :: proc "contextless" (x: f64) -> (f64, f64) {
@(static) gamS := [?]f64{ @(static) gamS := [?]f64{
7.87311395793093628397e-04, +7.87311395793093628397e-04,
-2.29549961613378126380e-04, -2.29549961613378126380e-04,
-2.68132617805781232825e-03, -2.68132617805781232825e-03,
3.47222221605458667310e-03, +3.47222221605458667310e-03,
8.33333333333482257126e-02, +8.33333333333482257126e-02,
} }
if x > 200 { if x > 200 {
return inf_f64(1), 1 return inf_f64(1), 1
} }
SQRT_TWO_PI :: 2.506628274631000502417 SQRT_TWO_PI :: 0h40040d931ff62706 // 2.506628274631000502417
MAX_STIRLING :: 143.01608 MAX_STIRLING :: 143.01608
w := 1 / x w := 1 / x
w = 1 + w*((((gamS[0]*w+gamS[1])*w+gamS[2])*w+gamS[3])*w+gamS[4]) w = 1 + w*((((gamS[0]*w+gamS[1])*w+gamS[2])*w+gamS[3])*w+gamS[4])
@@ -113,13 +113,13 @@ gamma_f64 :: proc "contextless" (x: f64) -> f64 {
} }
@(static) gamQ := [?]f64{ @(static) gamQ := [?]f64{
-2.31581873324120129819e-05, -2.31581873324120129819e-05,
5.39605580493303397842e-04, +5.39605580493303397842e-04,
-4.45641913851797240494e-03, -4.45641913851797240494e-03,
1.18139785222060435552e-02, +1.18139785222060435552e-02,
3.58236398605498653373e-02, +3.58236398605498653373e-02,
-2.34591795718243348568e-01, -2.34591795718243348568e-01,
7.14304917030273074085e-02, +7.14304917030273074085e-02,
1.00000000000000000320e+00, +1.00000000000000000320e+00,
} }
+4 -4
View File
@@ -197,9 +197,9 @@ lgamma_f64 :: proc "contextless" (x: f64) -> (lgamma: f64, sign: int) {
} }
Y_MIN :: 1.461632144968362245 Y_MIN :: 0h3ff762d86356be3f // 1.461632144968362245
TWO_52 :: 0h4330000000000000 // ~4.5036e+15 TWO_52 :: 0h4330000000000000 // ~4.5036e+15
TWO_53 :: 0h4340000000000000 // ~9.0072e+15 TWO_53 :: 0h4340000000000000 // ~9.0072e+15
TWO_58 :: 0h4390000000000000 // ~2.8823e+17 TWO_58 :: 0h4390000000000000 // ~2.8823e+17
TINY :: 0h3b90000000000000 // ~8.47033e-22 TINY :: 0h3b90000000000000 // ~8.47033e-22
Tc :: 0h3FF762D86356BE3F Tc :: 0h3FF762D86356BE3F
@@ -345,8 +345,8 @@ lgamma_f64 :: proc "contextless" (x: f64) -> (lgamma: f64, sign: int) {
} }
lgamma_f16 :: proc "contextless" (x: f16) -> (lgamma: f16, sign: int) { r, s := lgamma_f64(f64(x)); return f16(r), s } lgamma_f16 :: proc "contextless" (x: f16) -> (lgamma: f16, sign: int) { r, s := lgamma_f64(f64(x)); return f16(r), s }
lgamma_f32 :: proc "contextless" (x: f32) -> (lgamma: f32, sign: int) { r, s := lgamma_f64(f64(x)); return f32(r), s } lgamma_f32 :: proc "contextless" (x: f32) -> (lgamma: f32, sign: int) { r, s := lgamma_f64(f64(x)); return f32(r), s }
lgamma_f16le :: proc "contextless" (x: f16le) -> (lgamma: f16le, sign: int) { r, s := lgamma_f64(f64(x)); return f16le(r), s } lgamma_f16le :: proc "contextless" (x: f16le) -> (lgamma: f16le, sign: int) { r, s := lgamma_f64(f64(x)); return f16le(r), s }
lgamma_f16be :: proc "contextless" (x: f16be) -> (lgamma: f16be, sign: int) { r, s := lgamma_f64(f64(x)); return f16be(r), s } lgamma_f16be :: proc "contextless" (x: f16be) -> (lgamma: f16be, sign: int) { r, s := lgamma_f64(f64(x)); return f16be(r), s }
lgamma_f32le :: proc "contextless" (x: f32le) -> (lgamma: f32le, sign: int) { r, s := lgamma_f64(f64(x)); return f32le(r), s } lgamma_f32le :: proc "contextless" (x: f32le) -> (lgamma: f32le, sign: int) { r, s := lgamma_f64(f64(x)); return f32le(r), s }
+9 -9
View File
@@ -100,11 +100,11 @@ log1p_f64le :: proc "contextless" (x: f64le) -> f64le { return f64le(log1p_f64(f
log1p_f64be :: proc "contextless" (x: f64be) -> f64be { return f64be(log1p_f64(f64(x))) } log1p_f64be :: proc "contextless" (x: f64be) -> f64be { return f64be(log1p_f64(f64(x))) }
log1p_f64 :: proc "contextless" (x: f64) -> f64 { log1p_f64 :: proc "contextless" (x: f64) -> f64 {
SQRT2_M1 :: 0h3fda827999fcef34 // Sqrt(2)-1 SQRT2_M1 :: 0h3fda827999fcef34 // sqrt(2)-1
SQRT2_HALF_M1 :: 0hbfd2bec333018866 // Sqrt(2)/2-1 SQRT2_HALF_M1 :: 0hbfd2bec333018866 // sqrt(2)/2-1
SMALL :: 0h3e20000000000000 // 2**-29 SMALL :: 0h3e20000000000000 // 2**-29
TINY :: 1.0 / (1 << 54) // 2**-54 TINY :: 0h3c90000000000000 // 2**-54
TWO53 :: 1 << 53 // 2**53 TWO53 :: 0h4340000000000000 // 2**53
LN2HI :: 0h3fe62e42fee00000 LN2HI :: 0h3fe62e42fee00000
LN2LO :: 0h3dea39ef35793c76 LN2LO :: 0h3dea39ef35793c76
LP1 :: 0h3FE5555555555593 LP1 :: 0h3FE5555555555593
@@ -128,15 +128,15 @@ log1p_f64 :: proc "contextless" (x: f64) -> f64 {
f: f64 f: f64
iu: u64 iu: u64
k := 1 k := 1
if absx < SQRT2_M1 { // |x| < Sqrt(2)-1 if absx < SQRT2_M1 { // |x| < sqrt(2)-1
if absx < SMALL { // |x| < 2**-29 if absx < SMALL { // |x| < 2**-29
if absx < TINY { // |x| < 2**-54 if absx < TINY { // |x| < 2**-54
return x return x
} }
return x - x*x*0.5 return x - x*x*0.5
} }
if x > SQRT2_HALF_M1 { // Sqrt(2)/2-1 < x if x > SQRT2_HALF_M1 { // sqrt(2)/2-1 < x
// (Sqrt(2)/2-1) < x < (Sqrt(2)-1) // (sqrt(2)/2-1) < x < (sqrt(2)-1)
k = 0 k = 0
f = x f = x
iu = 1 iu = 1
@@ -163,14 +163,14 @@ log1p_f64 :: proc "contextless" (x: f64) -> f64 {
c = 0 c = 0
} }
iu &= 0x000fffffffffffff iu &= 0x000fffffffffffff
if iu < 0x0006a09e667f3bcd { // mantissa of Sqrt(2) if iu < 0x0006a09e667f3bcd { // mantissa of sqrt(2)
u = transmute(f64)(iu | 0x3ff0000000000000) // normalize u u = transmute(f64)(iu | 0x3ff0000000000000) // normalize u
} else { } else {
k += 1 k += 1
u = transmute(f64)(iu | 0x3fe0000000000000) // normalize u/2 u = transmute(f64)(iu | 0x3fe0000000000000) // normalize u/2
iu = (0x0010000000000000 - iu) >> 2 iu = (0x0010000000000000 - iu) >> 2
} }
f = u - 1.0 // Sqrt(2)/2 < u < Sqrt(2) f = u - 1.0 // sqrt(2)/2 < u < sqrt(2)
} }
hfsq := 0.5 * f * f hfsq := 0.5 * f * f
s, R, z: f64 s, R, z: f64