[WASM] Add kotlin.math implementation

This commit is contained in:
Igor Yakovlev
2021-12-15 16:21:13 +01:00
committed by igoriakovlev
parent 5e42be9849
commit adee0f1cb0
37 changed files with 3743 additions and 102 deletions
@@ -246,6 +246,38 @@ public external fun wasm_f64_min(a: Double, b: Double): Double
@WasmOp(WasmOp.F64_MAX)
public external fun wasm_f64_max(a: Double, b: Double): Double
@WasmOp(WasmOp.F64_SQRT)
public external fun wasm_f64_sqrt(a: Double): Double
@WasmOp(WasmOp.F64_CEIL)
public external fun wasm_f64_ceil(a: Double): Double
@WasmOp(WasmOp.F64_FLOOR)
public external fun wasm_f64_floor(a: Double): Double
@WasmOp(WasmOp.F64_TRUNC)
public external fun wasm_f64_truncate(a: Double): Double
@WasmOp(WasmOp.F64_COPYSIGN)
public external fun wasm_f64_copysign(a: Double, b: Double): Double
@WasmOp(WasmOp.F64_ABS)
public external fun wasm_f64_abs(a: Double): Double
@WasmOp(WasmOp.F32_SQRT)
public external fun wasm_f32_sqrt(a: Float): Float
@WasmOp(WasmOp.F32_CEIL)
public external fun wasm_f32_ceil(a: Float): Float
@WasmOp(WasmOp.F32_FLOOR)
public external fun wasm_f32_floor(a: Float): Float
@WasmOp(WasmOp.F32_TRUNC)
public external fun wasm_f32_truncate(a: Float): Float
@WasmOp(WasmOp.F32_ABS)
public external fun wasm_f32_abs(a: Float): Float
@WasmOp(WasmOp.REF_IS_NULL)
public external fun wasm_ref_is_null(a: Any?): Boolean
@@ -0,0 +1,103 @@
/* @(#)e_acos.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_acos(x)
* Method :
* acos(x) = pi/2 - asin(x)
* acos(-x) = pi/2 + asin(x)
* For |x|<=0.5
* acos(x) = pi/2 - (x + x*x^2*R(x^2)) (see asin.c)
* For x>0.5
* acos(x) = pi/2 - (pi/2 - 2asin(sqrt((1-x)/2)))
* = 2asin(sqrt((1-x)/2))
* = 2s + 2s*z*R(z) ...z=(1-x)/2, s=sqrt(z)
* = 2f + (2c + 2s*z*R(z))
* where f=hi part of s, and c = (z-f*f)/(s+f) is the correction term
* for f so that f+c ~ sqrt(z).
* For x<-0.5
* acos(x) = pi - 2asin(sqrt((1-|x|)/2))
* = pi - 0.5*(s+s*z*R(z)), where z=(1-|x|)/2,s=sqrt(z)
*
* Special cases:
* if x is NaN, return x itself;
* if |x|>1, return NaN with invalid signal.
*
* Function needed: sqrt
*/
package kotlin.math.fdlibm
import kotlin.wasm.internal.wasm_f64_sqrt as sqrt
private const val one = 1.00000000000000000000e+00 /* 0x3FF00000, 0x00000000 */
private const val pi = 3.14159265358979311600e+00 /* 0x400921FB, 0x54442D18 */
private const val pio2_hi = 1.57079632679489655800e+00 /* 0x3FF921FB, 0x54442D18 */
private const val pio2_lo = 6.12323399573676603587e-17 /* 0x3C91A626, 0x33145C07 */
private const val pS0 = 1.66666666666666657415e-01 /* 0x3FC55555, 0x55555555 */
private const val pS1 = -3.25565818622400915405e-01 /* 0xBFD4D612, 0x03EB6F7D */
private const val pS2 = 2.01212532134862925881e-01 /* 0x3FC9C155, 0x0E884455 */
private const val pS3 = -4.00555345006794114027e-02 /* 0xBFA48228, 0xB5688F3B */
private const val pS4 = 7.91534994289814532176e-04 /* 0x3F49EFE0, 0x7501B288 */
private const val pS5 = 3.47933107596021167570e-05 /* 0x3F023DE1, 0x0DFDF709 */
private const val qS1 = -2.40339491173441421878e+00 /* 0xC0033A27, 0x1C8A2D4B */
private const val qS2 = 2.02094576023350569471e+00 /* 0x40002AE5, 0x9C598AC8 */
private const val qS3 = -6.88283971605453293030e-01 /* 0xBFE6066C, 0x1B8D0159 */
private const val qS4 = 7.70381505559019352791e-02 /* 0x3FB3B8C5, 0xB12E9282 */
internal fun __ieee754_acos(x: Double): Double {
var z: Double = 0.0
var p: Double = 0.0
var q: Double = 0.0
var r: Double = 0.0
var w: Double = 0.0
var s: Double = 0.0
var c: Double = 0.0
var df: Double = 0.0
var hx: Int = 0
var ix: Int = 0
hx = __HI(x)
ix = hx and 0x7fffffff
if (ix >= 0x3ff00000) { /* |x| >= 1 */
if (((ix - 0x3ff00000) or __LO(x)) == 0) { /* |x|==1 */
if (hx > 0) return 0.0 /* acos(1) = 0 */
else return pi + 2.0 * pio2_lo /* acos(-1)= pi */
}
return (x - x) / (x - x) /* acos(|x|>1) is NaN */
}
if (ix < 0x3fe00000) { /* |x| < 0.5 */
if (ix <= 0x3c600000) return pio2_hi + pio2_lo/*if|x|<2**-57*/
z = x * x
p = z * (pS0 + z * (pS1 + z * (pS2 + z * (pS3 + z * (pS4 + z * pS5)))))
q = one + z * (qS1 + z * (qS2 + z * (qS3 + z * qS4)))
r = p / q
return pio2_hi - (x - (pio2_lo - x * r))
} else if (hx < 0) { /* x < -0.5 */
z = (one + x) * 0.5
p = z * (pS0 + z * (pS1 + z * (pS2 + z * (pS3 + z * (pS4 + z * pS5)))))
q = one + z * (qS1 + z * (qS2 + z * (qS3 + z * qS4)))
s = sqrt(z)
r = p / q
w = r * s - pio2_lo
return pi - 2.0 * (s + w)
} else { /* x > 0.5 */
z = (one - x) * 0.5
s = sqrt(z)
df = s
df = doubleSetWord(d = df, lo = 0)
c = (z - df * df) / (s + df)
p = z * (pS0 + z * (pS1 + z * (pS2 + z * (pS3 + z * (pS4 + z * pS5)))))
q = one + z * (qS1 + z * (qS2 + z * (qS3 + z * qS4)))
r = p / q
w = r * s + c
return 2.0 * (df + w)
}
}
@@ -0,0 +1,55 @@
/* @(#)e_acosh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
/* __ieee754_acosh(x)
* Method :
* Based on
* acosh(x) = log [ x + sqrt(x*x-1) ]
* we have
* acosh(x) := log(x)+ln2, if x is large; else
* acosh(x) := log(2x-1/(sqrt(x*x-1)+x)) if x>2; else
* acosh(x) := log1p(t+sqrt(2.0*t+t*t)); where t=x-1.
*
* Special cases:
* acosh(x) is NaN with signal if x<1.
* acosh(NaN) is NaN without signal.
*/
package kotlin.math.fdlibm
import kotlin.wasm.internal.wasm_f64_sqrt as sqrt
private const val one = 1.0
private const val ln2 = 6.93147180559945286227e-01 /* 0x3FE62E42, 0xFEFA39EF */
internal fun __ieee754_acosh(x: Double): Double {
var t: Double = 0.0
var hx: Int = 0
hx = __HI(x)
if (hx < 0x3ff00000) { /* x < 1 */
return (x - x) / (x - x)
} else if (hx >= 0x41b00000) { /* x > 2**28 */
if (hx >= 0x7ff00000) { /* x is inf of NaN */
return x + x
} else
return __ieee754_log(x) + ln2 /* acosh(huge)=log(2x) */
} else if (((hx - 0x3ff00000) or __LO(x)) == 0) {
return 0.0 /* acosh(1) = 0 */
} else if (hx > 0x40000000) { /* 2**28 > x > 2 */
t = x * x
return __ieee754_log(2.0 * x - one / (x + sqrt(t - one)))
} else { /* 1<x<2 */
t = x - one
return log1p(t + sqrt(2.0 * t + t * t))
}
}
@@ -0,0 +1,113 @@
/* @(#)e_asin.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_asin(x)
* Method :
* Since asin(x) = x + x^3/6 + x^5*3/40 + x^7*15/336 + ...
* we approximate asin(x) on [0,0.5] by
* asin(x) = x + x*x^2*R(x^2)
* where
* R(x^2) is a rational approximation of (asin(x)-x)/x^3
* and its remez error is bounded by
* |(asin(x)-x)/x^3 - R(x^2)| < 2^(-58.75)
*
* For x in [0.5,1]
* asin(x) = pi/2-2*asin(sqrt((1-x)/2))
* Let y = (1-x), z = y/2, s := sqrt(z), and pio2_hi+pio2_lo=pi/2;
* then for x>0.98
* asin(x) = pi/2 - 2*(s+s*z*R(z))
* = pio2_hi - (2*(s+s*z*R(z)) - pio2_lo)
* For x<=0.98, let pio4_hi = pio2_hi/2, then
* f = hi part of s;
* c = sqrt(z) - f = (z-f*f)/(s+f) ...f+c=sqrt(z)
* and
* asin(x) = pi/2 - 2*(s+s*z*R(z))
* = pio4_hi+(pio4-2s)-(2s*z*R(z)-pio2_lo)
* = pio4_hi+(pio4-2f)-(2s*z*R(z)-(pio2_lo+2c))
*
* Special cases:
* if x is NaN, return x itself;
* if |x|>1, return NaN with invalid signal.
*
*/
package kotlin.math.fdlibm
import kotlin.wasm.internal.wasm_f64_sqrt as sqrt
private const val one = 1.00000000000000000000e+00 /* 0x3FF00000, 0x00000000 */
private const val huge = 1.000e+300
private const val pio2_hi = 1.57079632679489655800e+00 /* 0x3FF921FB, 0x54442D18 */
private const val pio2_lo = 6.12323399573676603587e-17 /* 0x3C91A626, 0x33145C07 */
private const val pio4_hi = 7.85398163397448278999e-01 /* 0x3FE921FB, 0x54442D18 */
/* coefficient for R(x^2) */
private const val pS0 = 1.66666666666666657415e-01 /* 0x3FC55555, 0x55555555 */
private const val pS1 = -3.25565818622400915405e-01 /* 0xBFD4D612, 0x03EB6F7D */
private const val pS2 = 2.01212532134862925881e-01 /* 0x3FC9C155, 0x0E884455 */
private const val pS3 = -4.00555345006794114027e-02 /* 0xBFA48228, 0xB5688F3B */
private const val pS4 = 7.91534994289814532176e-04 /* 0x3F49EFE0, 0x7501B288 */
private const val pS5 = 3.47933107596021167570e-05 /* 0x3F023DE1, 0x0DFDF709 */
private const val qS1 = -2.40339491173441421878e+00 /* 0xC0033A27, 0x1C8A2D4B */
private const val qS2 = 2.02094576023350569471e+00 /* 0x40002AE5, 0x9C598AC8 */
private const val qS3 = -6.88283971605453293030e-01 /* 0xBFE6066C, 0x1B8D0159 */
private const val qS4 = 7.70381505559019352791e-02 /* 0x3FB3B8C5, 0xB12E9282 */
internal fun __ieee754_asin(x: Double): Double {
var t: Double = 0.0
var w: Double = 0.0
var p: Double = 0.0
var q: Double = 0.0
var c: Double = 0.0
var r: Double = 0.0
var s: Double = 0.0
var hx: Int = 0
var ix: Int = 0
hx = __HI(x)
ix = hx and 0x7fffffff
if (ix >= 0x3ff00000) { /* |x|>= 1 */
if (((ix - 0x3ff00000) or __LO(x)) == 0)
/* asin(1)=+-pi/2 with inexact */
return x * pio2_hi + x * pio2_lo
return (x - x) / (x - x) /* asin(|x|>1) is NaN */
} else if (ix < 0x3fe00000) { /* |x|<0.5 */
if (ix < 0x3e400000) { /* if |x| < 2**-27 */
if (huge + x > one) return x/* return x with inexact if x!=0*/
} else
t = x * x
p = t * (pS0 + t * (pS1 + t * (pS2 + t * (pS3 + t * (pS4 + t * pS5)))))
q = one + t * (qS1 + t * (qS2 + t * (qS3 + t * qS4)))
w = p / q
return x + x * w
}
/* 1> |x|>= 0.5 */
w = one - fabs(x)
t = w * 0.5
p = t * (pS0 + t * (pS1 + t * (pS2 + t * (pS3 + t * (pS4 + t * pS5)))))
q = one + t * (qS1 + t * (qS2 + t * (qS3 + t * qS4)))
s = sqrt(t)
if (ix >= 0x3FEF3333) { /* if |x| > 0.975 */
w = p / q
t = pio2_hi - (2.0 * (s + s * w) - pio2_lo)
} else {
w = s
w = doubleSetWord(d = w, lo = 0)
c = (t - w * w) / (s + w)
r = p / q
p = 2.0 * s * r - (pio2_lo - 2.0 * c)
q = pio4_hi - 2.0 * w
t = pio4_hi - (p - q)
}
if (hx > 0) return t; else return -t
}
@@ -0,0 +1,119 @@
/* @(#)e_atan2.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
/* __ieee754_atan2(y,x)
* Method :
* 1. Reduce y to positive by atan2(y,x)=-atan2(-y,x).
* 2. Reduce x to positive by (if x and y are unexceptional):
* ARG (x+iy) = arctan(y/x) ... if x > 0,
* ARG (x+iy) = pi - arctan[y/(-x)] ... if x < 0,
*
* Special cases:
*
* ATAN2((anything), NaN ) is NaN;
* ATAN2(NAN , (anything) ) is NaN;
* ATAN2(+-0, +(anything but NaN)) is +-0 ;
* ATAN2(+-0, -(anything but NaN)) is +-pi ;
* ATAN2(+-(anything but 0 and NaN), 0) is +-pi/2;
* ATAN2(+-(anything but INF and NaN), +INF) is +-0 ;
* ATAN2(+-(anything but INF and NaN), -INF) is +-pi;
* ATAN2(+-INF,+INF ) is +-pi/4 ;
* ATAN2(+-INF,-INF ) is +-3pi/4;
* ATAN2(+-INF, (anything but,0,NaN, and INF)) is +-pi/2;
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
package kotlin.math.fdlibm
private const val tiny = 1.0e-300
private const val zero = 0.0
private const val pi_o_4 = 7.8539816339744827900E-01 /* 0x3FE921FB, 0x54442D18 */
private const val pi_o_2 = 1.5707963267948965580E+00 /* 0x3FF921FB, 0x54442D18 */
private const val pi = 3.1415926535897931160E+00 /* 0x400921FB, 0x54442D18 */
private const val pi_lo = 1.2246467991473531772E-16 /* 0x3CA1A626, 0x33145C07 */
internal fun __ieee754_atan2(y: Double, x: Double): Double {
var z: Double = 0.0
var k: Int = 0
var m: Int = 0
var hx: Int = 0
var hy: Int = 0
var ix: Int = 0
var iy: Int = 0
var lx: UInt = 0U
var ly: UInt = 0U
hx = __HI(x); ix = hx and 0x7fffffff
lx = __LOu(x)
hy = __HI(y); iy = hy and 0x7fffffff
ly = __LOu(y)
if (((ix or ((lx or lx.negate()) shr 31).toInt()) > 0x7ff00000) ||
((iy or ((ly or ly.negate()) shr 31).toInt()) > 0x7ff00000)
) /* x or y is NaN */
return x + y
if (((hx - 0x3ff00000) or lx.toInt()) == 0) return atan(y) /* x=1.0 */
m = ((hy shr 31) and 1) or ((hx shr 30) and 2) /* 2*sign(x)+sign(y) */
/* when y = 0 */
if ((iy or ly.toInt()) == 0) {
when (m) {
0, 1 -> return y /* atan(+-0,+anything)=+-0 */
2 -> return pi + tiny/* atan(+0,-anything) = pi */
3 -> return -pi - tiny/* atan(-0,-anything) =-pi */
}
}
/* when x = 0 */
if ((ix or lx.toInt()) == 0) return if (hy < 0) -pi_o_2 - tiny else pi_o_2 + tiny
/* when x is INF */
if (ix == 0x7ff00000) {
if (iy == 0x7ff00000) {
when (m) {
0 -> return pi_o_4 + tiny/* atan(+INF,+INF) */
1 -> return -pi_o_4 - tiny/* atan(-INF,+INF) */
2 -> return 3.0 * pi_o_4 + tiny/*atan(+INF,-INF)*/
3 -> return -3.0 * pi_o_4 - tiny/*atan(-INF,-INF)*/
}
} else {
when (m) {
0 -> return zero /* atan(+...,+INF) */
1 -> return -zero /* atan(-...,+INF) */
2 -> return pi + tiny /* atan(+...,-INF) */
3 -> return -pi - tiny /* atan(-...,-INF) */
}
}
}
/* when y is INF */
if (iy == 0x7ff00000) return if (hy < 0) -pi_o_2 - tiny else pi_o_2 + tiny
/* compute y/x */
k = (iy - ix) shr 20
if (k > 60) z = pi_o_2 + 0.5 * pi_lo /* |y/x| > 2**60 */
else if (hx < 0 && k < -60) z = 0.0 /* |y|/x < -2**60 */
else z = atan(fabs(y / x)) /* safe to do y/x */
when (m) {
0 -> return z /* atan(+,+) */
1 -> {
z = doubleSetWord(d = z, hi = __HI(z) xor Int.MIN_VALUE)
return z /* atan(-,+) */
}
2 -> return pi - (z - pi_lo)/* atan(+,-) */
else -> /* case 3 */ return (z - pi_lo) - pi/* atan(-,-) */
}
}
@@ -0,0 +1,59 @@
/* @(#)e_atanh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
/* __ieee754_atanh(x)
* Method :
* 1.Reduced x to positive by atanh(-x) = -atanh(x)
* 2.For x>=0.5
* 1 2x x
* atanh(x) = --- * log(1 + -------) = 0.5 * log1p(2 * --------)
* 2 1 - x 1 - x
*
* For x<0.5
* atanh(x) = 0.5*log1p(2x+2x*x/(1-x))
*
* Special cases:
* atanh(x) is NaN if |x| > 1 with signal;
* atanh(NaN) is that NaN with no signal;
* atanh(+-1) is +-INF with signal.
*
*/
package kotlin.math.fdlibm
private const val one = 1.0
private const val huge = 1e300
private const val zero = 0.0
internal fun __ieee754_atanh(x: Double): Double {
var x: Double = x
var t: Double = 0.0
var hx: Int = 0
var ix: Int = 0
var lx: UInt = 0U
hx = __HI(x) /* high word */
lx = __LOu(x) /* low word */
ix = hx and 0x7fffffff
if ((ix or ((lx or lx.negate()) shr 31).toInt()) > 0x3ff00000) /* |x|>1 */
return (x - x) / (x - x)
if (ix == 0x3ff00000)
return x / zero
if (ix < 0x3e300000 && (huge + x) > zero) return x /* x<2**-28 */
x = doubleSetWord(d = x, hi = ix) /* x <- |x| */
if (ix < 0x3fe00000) { /* x < 0.5 */
t = x + x
t = 0.5 * log1p(t + t * x / (one - x))
} else
t = 0.5 * log1p((x + x) / (one - x))
if (hx >= 0) return t; else return -t
}
@@ -0,0 +1,83 @@
/* @(#)e_cosh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_cosh(x)
* Method :
* mathematically cosh(x) if defined to be (exp(x)+exp(-x))/2
* 1. Replace x by |x| (cosh(x) = cosh(-x)).
* 2.
* [ exp(x) - 1 ]^2
* 0 <= x <= ln2/2 : cosh(x) := 1 + -------------------
* 2*exp(x)
*
* exp(x) + 1/exp(x)
* ln2/2 <= x <= 22 : cosh(x) := -------------------
* 2
* 22 <= x <= lnovft : cosh(x) := exp(x)/2
* lnovft <= x <= ln2ovft: cosh(x) := exp(x/2)/2 * exp(x/2)
* ln2ovft < x : cosh(x) := huge*huge (overflow)
*
* Special cases:
* cosh(x) is |x| if x is +INF, -INF, or NaN.
* only cosh(0)=1 is exact for finite x.
*/
package kotlin.math.fdlibm
private const val one = 1.0
private const val half = 0.5
private const val huge = 1.0e300
internal fun __ieee754_cosh(x: Double): Double {
var t: Double = 0.0
var w: Double = 0.0
var ix: Int = 0
var lx: UInt = 0U
/* High word of |x|. */
ix = __HI(x)
ix = ix and 0x7fffffff
/* x is INF or NaN */
if (ix >= 0x7ff00000) return x * x
/* |x| in [0,0.5*ln2], return 1+expm1(|x|)^2/(2*exp(|x|)) */
if (ix < 0x3fd62e43) {
t = expm1(fabs(x))
w = one + t
if (ix < 0x3c800000) return w /* cosh(tiny) = 1 */
return one + (t * t) / (w + w)
}
/* |x| in [0.5*ln2,22], return (exp(|x|)+1/exp(|x|)/2; */
if (ix < 0x40360000) {
t = __ieee754_exp(fabs(x))
return half * t + half / t
}
/* |x| in [22, log(maxdouble)] return half*exp(|x|) */
if (ix < 0x40862E42) return half * __ieee754_exp(fabs(x))
/* |x| in [log(maxdouble), overflowthresold] */
//lx = *( (((*(unsigned*)&one) shr 29)) + (unsigned*)&x);
lx = __LOu(x)
if (ix < 0x408633CE ||
(ix == 0x408633ce) && (lx <= 0x8fb9f87d.toUInt())
) {
w = __ieee754_exp(half * fabs(x))
t = half * w
return t * w
}
/* |x| > overflowthresold, cosh(x) overflow */
return huge * huge
}
@@ -0,0 +1,153 @@
/* @(#)e_exp.c 1.6 04/04/22 */
/*
* ====================================================
* Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
*
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_exp(x)
* Returns the exponential of x.
*
* Method
* 1. Argument reduction:
* Reduce x to an r so that |r| <= 0.5*ln2 ~ 0.34658.
* Given x, find r and integer k such that
*
* x = k*ln2 + r, |r| <= 0.5*ln2.
*
* Here r will be represented as r = hi-lo for better
* accuracy.
*
* 2. Approximation of exp(r) by a special rational function on
* the interval [0,0.34658]:
* Write
* R(r**2) = r*(exp(r)+1)/(exp(r)-1) = 2 + r*r/6 - r**4/360 + ...
* We use a special Remes algorithm on [0,0.34658] to generate
* a polynomial of degree 5 to approximate R. The maximum error
* of this polynomial approximation is bounded by 2**-59. In
* other words,
* R(z) ~ 2.0 + P1*z + P2*z**2 + P3*z**3 + P4*z**4 + P5*z**5
* (where z=r*r, and the values of P1 to P5 are listed below)
* and
* | 5 | -59
* | 2.0+P1*z+...+P5*z - R(z) | <= 2
* | |
* The computation of exp(r) thus becomes
* 2*r
* exp(r) = 1 + -------
* R - r
* r*R1(r)
* = 1 + r + ----------- (for better accuracy)
* 2 - R1(r)
* where
* 2 4 10
* R1(r) = r - (P1*r + P2*r + ... + P5*r ).
*
* 3. Scale back to obtain exp(x):
* From step 1, we have
* exp(x) = 2^k * exp(r)
*
* Special cases:
* exp(INF) is INF, exp(NaN) is NaN;
* exp(-INF) is 0, and
* for finite argument, only exp(0)=1 is exact.
*
* Accuracy:
* according to an error analysis, the error is always less than
* 1 ulp (unit in the last place).
*
* Misc. info.
* For IEEE double
* if x > 7.09782712893383973096e+02 then exp(x) overflow
* if x < -7.45133219101941108420e+02 then exp(x) underflow
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
package kotlin.math.fdlibm
private const val one = 1.0
private val halF = doubleArrayOf(0.5, -0.5)
private const val huge = 1.0e+300
private const val twom1000 = 9.33263618503218878990e-302 /* 2**-1000=0x01700000,0*/
private const val o_threshold = 7.09782712893383973096e+02 /* 0x40862E42, 0xFEFA39EF */
private const val u_threshold = -7.45133219101941108420e+02 /* 0xc0874910, 0xD52D3051 */
private val ln2HI = doubleArrayOf(
6.93147180369123816490e-01, /* 0x3fe62e42, 0xfee00000 */
-6.93147180369123816490e-01,
)/* 0xbfe62e42, 0xfee00000 */
private val ln2LO = doubleArrayOf(
1.90821492927058770002e-10, /* 0x3dea39ef, 0x35793c76 */
-1.90821492927058770002e-10,
)/* 0xbdea39ef, 0x35793c76 */
private const val invln2 = 1.44269504088896338700e+00 /* 0x3ff71547, 0x652b82fe */
private const val P1 = 1.66666666666666019037e-01 /* 0x3FC55555, 0x5555553E */
private const val P2 = -2.77777777770155933842e-03 /* 0xBF66C16C, 0x16BEBD93 */
private const val P3 = 6.61375632143793436117e-05 /* 0x3F11566A, 0xAF25DE2C */
private const val P4 = -1.65339022054652515390e-06 /* 0xBEBBBD41, 0xC5D26BF1 */
private const val P5 = 4.13813679705723846039e-08 /* 0x3E663769, 0x72BEA4D0 */
internal fun __ieee754_exp(x: Double): Double /* default IEEE double exp */ {
var x: Double = x
var y: Double = 0.0
var hi: Double = 0.0
var lo: Double = 0.0
var c: Double = 0.0
var t: Double = 0.0
var k: Int = 0
var xsb: Int = 0
var hx: UInt = 0U
hx = __HIu(x) /* high word of x */
xsb = ((hx shr 31) and 1U).toInt() /* sign bit of x */
hx = hx and 0x7fffffffU /* high word of |x| */
/* filter out non-finite argument */
if (hx >= 0x40862E42U) { /* if |x|>=709.78... */
if (hx >= 0x7ff00000U) {
if (((hx and 0xfffffU) or __LOu(x)) != 0U)
return x + x /* NaN */
else return if (xsb == 0) x else 0.0 /* exp(+-inf)={inf,0} */
}
if (x > o_threshold) return huge * huge /* overflow */
if (x < u_threshold) return twom1000 * twom1000 /* underflow */
}
/* argument reduction */
if (hx > 0x3fd62e42U) { /* if |x| > 0.5 ln2 */
if (hx < 0x3FF0A2B2U) { /* and |x| < 1.5 ln2 */
hi = x - ln2HI[xsb]; lo = ln2LO[xsb]; k = 1 - xsb - xsb
} else {
k = (invln2 * x + halF[xsb]).toInt()
t = k.toDouble()
hi = x - t * ln2HI[0] /* t*ln2HI is exact here */
lo = t * ln2LO[0]
}
x = hi - lo
} else if (hx < 0x3e300000U) { /* when |x|<2**-28 */
if (huge + x > one) return one + x/* trigger inexact */
} else k = 0
/* x is now in primary range */
t = x * x
c = x - t * (P1 + t * (P2 + t * (P3 + t * (P4 + t * P5))))
if (k == 0) return one - ((x * c) / (c - 2.0) - x)
else y = one - ((lo - (x * c) / (2.0 - c)) - hi)
if (k >= -1021) {
y = doubleSetWord(d = y, hi = __HI(y) + (k shl 20)) /* add k to y's exponent */
return y
} else {
y = doubleSetWord(d = y, hi = __HI(y) + ((k + 1000) shl 20))/* add k to y's exponent */
return y * twom1000
}
}
@@ -0,0 +1,125 @@
/* @(#)e_hypot.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_hypot(x,y)
*
* Method :
* If (assume round-to-nearest) z=x*x+y*y
* has error less than sqrt(2)/2 ulp, than
* sqrt(z) has error less than 1 ulp (exercise).
*
* So, compute sqrt(x*x+y*y) with some care as
* follows to get the error below 1 ulp:
*
* Assume x>y>0;
* (if possible, set rounding to round-to-nearest)
* 1. if x > 2y use
* x1*x1+(y*y+(x2*(x+x1))) for x*x+y*y
* where x1 = x with lower 32 bits cleared, x2 = x-x1; else
* 2. if x <= 2y use
* t1*y1+((x-y)*(x-y)+(t1*y2+t2*y))
* where t1 = 2x with lower 32 bits cleared, t2 = 2x-t1,
* y1= y with lower 32 bits chopped, y2 = y-y1.
*
* NOTE: scaling may be necessary if some argument is too
* large or too tiny
*
* Special cases:
* hypot(x,y) is INF if x or y is +INF or -INF; else
* hypot(x,y) is NAN if x or y is NAN.
*
* Accuracy:
* hypot(x,y) returns sqrt(x^2+y^2) with error less
* than 1 ulps (units in the last place)
*/
package kotlin.math.fdlibm
import kotlin.wasm.internal.wasm_f64_sqrt as sqrt
internal fun __ieee754_hypot(x: Double, y: Double): Double {
var a = x
var b = y
var t1: Double = 0.0
var t2: Double = 0.0
var y1: Double = 0.0
var y2: Double = 0.0
var w: Double = 0.0
var j: Int = 0
var k: Int = 0
var ha: Int = 0
var hb: Int = 0
ha = __HI(x) and 0x7fffffff /* high word of x */
hb = __HI(y) and 0x7fffffff /* high word of y */
if (hb > ha) {
a = y;b = x;j = ha; ha = hb;hb = j
} else {
a = x;b = y
}
a = doubleSetWord(d = a, hi = ha) /* a <- |a| */
b = doubleSetWord(d = b, hi = hb) /* b <- |b| */
if ((ha - hb) > 0x3c00000) {
return a + b
} /* x/y > 2**60 */
k = 0
if (ha > 0x5f300000) { /* a>2**500 */
if (ha >= 0x7ff00000) { /* Inf or NaN */
w = a + b /* for sNaN */
if (((ha and 0xfffff) or __LO(a)) == 0) w = a
if (((hb xor 0x7ff00000) or __LO(b)) == 0) w = b
return w
}
/* scale a and b by 2**-600 */
ha -= 0x25800000; hb -= 0x25800000; k += 600
a = doubleSetWord(d = a, hi = ha)
b = doubleSetWord(d = b, hi = hb)
}
if (hb < 0x20b00000) { /* b < 2**-500 */
if (hb <= 0x000fffff) { /* subnormal b or 0 */
if ((hb or (__LO(b))) == 0) return a
t1 = 0.0
t1 = doubleSetWord(d = t1, hi = 0x7fd00000) /* t1=2^1022 */
b *= t1
a *= t1
k -= 1022
} else { /* scale a and b by 2^600 */
ha += 0x25800000 /* a *= 2^600 */
hb += 0x25800000 /* b *= 2^600 */
k -= 600
a = doubleSetWord(d = a, hi = ha)
b = doubleSetWord(d = b, hi = hb)
}
}
/* medium size a and b */
w = a - b
if (w > b) {
t1 = 0.0
t1 = doubleSetWord(d = t1, hi = ha)
t2 = a - t1
w = sqrt(t1 * t1 - (b * (-b) - t2 * (a + t1)))
} else {
a = a + a
y1 = 0.0
y1 = doubleSetWord(d = y1, hi = hb)
y2 = b - y1
t1 = 0.0
t1 = doubleSetWord(d = t1, hi = ha + 0x00100000)
t2 = a - t1
w = sqrt(t1 * y1 - (w * (-w) - (t1 * y2 + t2 * b)))
}
if (k != 0) {
t1 = 1.0
t1 = doubleSetWord(d = t1, hi = __HI(t1) + (k shl 20))
return t1 * w
} else return w
}
@@ -0,0 +1,144 @@
/* @(#)e_log.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_log(x)
* Return the logrithm of x
*
* Method :
* 1. Argument Reduction: find k and f such that
* x = 2^k * (1+f),
* where sqrt(2)/2 < 1+f < sqrt(2) .
*
* 2. Approximation of log(1+f).
* Let s = f/(2+f) ; based on log(1+f) = log(1+s) - log(1-s)
* = 2s + 2/3 s**3 + 2/5 s**5 + .....,
* = 2s + s*R
* We use a special Reme algorithm on [0,0.1716] to generate
* a polynomial of degree 14 to approximate R The maximum error
* of this polynomial approximation is bounded by 2**-58.45. In
* other words,
* 2 4 6 8 10 12 14
* R(z) ~ Lg1*s +Lg2*s +Lg3*s +Lg4*s +Lg5*s +Lg6*s +Lg7*s
* (the values of Lg1 to Lg7 are listed in the program)
* and
* | 2 14 | -58.45
* | Lg1*s +...+Lg7*s - R(z) | <= 2
* | |
* Note that 2s = f - s*f = f - hfsq + s*hfsq, where hfsq = f*f/2.
* In order to guarantee error in log below 1ulp, we compute log
* by
* log(1+f) = f - s*(f - R) (if f is not too large)
* log(1+f) = f - (hfsq - s*(hfsq+R)). (better accuracy)
*
* 3. Finally, log(x) = k*ln2 + log(1+f).
* = k*ln2_hi+(f-(hfsq-(s*(hfsq+R)+k*ln2_lo)))
* Here ln2 is split into two floating point number:
* ln2_hi + ln2_lo,
* where n*ln2_hi is always exact for |n| < 2000.
*
* Special cases:
* log(x) is NaN with signal if x < 0 (including -INF) ;
* log(+INF) is +INF; log(0) is -INF with signal;
* log(NaN) is that NaN with no signal.
*
* Accuracy:
* according to an error analysis, the error is always less than
* 1 ulp (unit in the last place).
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
package kotlin.math.fdlibm
private const val ln2_hi = 6.93147180369123816490e-01 /* 3fe62e42 fee00000 */
private const val ln2_lo = 1.90821492927058770002e-10 /* 3dea39ef 35793c76 */
private const val two54 = 1.80143985094819840000e+16 /* 43500000 00000000 */
private const val Lg1 = 6.666666666666735130e-01 /* 3FE55555 55555593 */
private const val Lg2 = 3.999999999940941908e-01 /* 3FD99999 9997FA04 */
private const val Lg3 = 2.857142874366239149e-01 /* 3FD24924 94229359 */
private const val Lg4 = 2.222219843214978396e-01 /* 3FCC71C5 1D8E78AF */
private const val Lg5 = 1.818357216161805012e-01 /* 3FC74664 96CB03DE */
private const val Lg6 = 1.531383769920937332e-01 /* 3FC39A09 D078C69F */
private const val Lg7 = 1.479819860511658591e-01 /* 3FC2F112 DF3E5244 */
private const val zero = 0.0
internal fun __ieee754_log(x: Double): Double {
var x: Double = x
var hfsq: Double = 0.0
var f: Double = 0.0
var s: Double = 0.0
var z: Double = 0.0
var R: Double = 0.0
var w: Double = 0.0
var t1: Double = 0.0
var t2: Double = 0.0
var dk: Double = 0.0
var k: Int = 0
var hx: Int = 0
var i: Int = 0
var j: Int = 0
var lx: UInt = 0U
hx = __HI(x) /* high word of x */
lx = __LOu(x) /* low word of x */
k = 0
if (hx < 0x00100000) { /* x < 2**-1022 */
if (((hx and 0x7fffffff) or lx.toInt()) == 0)
return -two54 / zero /* log(+-0)=-inf */
if (hx < 0) return (x - x) / zero /* log(-#) = NaN */
k -= 54; x *= two54 /* subnormal number, scale up x */
hx = __HI(x) /* high word of x */
}
if (hx >= 0x7ff00000) return x + x
k += (hx shr 20) - 1023
hx = hx and 0x000fffff
i = (hx + 0x95f64) and 0x100000
x = doubleSetWord(d = x, hi = hx or (i xor 0x3ff00000)) /* normalize x or x/2 */
k += (i shr 20)
f = x - 1.0
if ((0x000fffff and (2 + hx)) < 3) { /* |f| < 2**-20 */
if (f == zero) if (k == 0) return zero; else {
dk = k.toDouble()
return dk * ln2_hi + dk * ln2_lo
}
R = f * f * (0.5 - 0.33333333333333333 * f)
if (k == 0) return f - R; else {
dk = k.toDouble()
return dk * ln2_hi - ((R - dk * ln2_lo) - f)
}
}
s = f / (2.0 + f)
dk = k.toDouble()
z = s * s
i = hx - 0x6147a
w = z * z
j = 0x6b851 - hx
t1 = w * (Lg2 + w * (Lg4 + w * Lg6))
t2 = z * (Lg1 + w * (Lg3 + w * (Lg5 + w * Lg7)))
i = i or j
R = t2 + t1
if (i > 0) {
hfsq = 0.5 * f * f
if (k == 0) return f - (hfsq - s * (hfsq + R)); else
return dk * ln2_hi - ((hfsq - (s * (hfsq + R) + dk * ln2_lo)) - f)
} else {
if (k == 0) return f - s * (f - R); else
return dk * ln2_hi - ((s * (f - R) - dk * ln2_lo) - f)
}
}
@@ -0,0 +1,83 @@
/* @(#)e_log10.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_log10(x)
* Return the base 10 logarithm of x
*
* Method :
* Let log10_2hi = leading 40 bits of log10(2) and
* log10_2lo = log10(2) - log10_2hi,
* ivln10 = 1/log(10) rounded.
* Then
* n = ilogb(x),
* if(n<0) n = n+1;
* x = scalbn(x,-n);
* log10(x) := n*log10_2hi + (n*log10_2lo + ivln10*log(x))
*
* Note 1:
* To guarantee log10(10**n)=n, where 10**n is normal, the rounding
* mode must set to Round-to-Nearest.
* Note 2:
* [1/log(10)] rounded to 53 bits has error .198 ulps;
* log10 is monotonic at all binary break points.
*
* Special cases:
* log10(x) is NaN with signal if x < 0;
* log10(+INF) is +INF with no signal; log10(0) is -INF with signal;
* log10(NaN) is that NaN with no signal;
* log10(10**N) = N for N=0,1,...,22.
*
* Constants:
* The hexadecimal values are the intended ones for the following constants.
* The decimal values may be used, provided that the compiler will convert
* from decimal to binary accurately enough to produce the hexadecimal values
* shown.
*/
package kotlin.math.fdlibm
private const val two54 = 1.80143985094819840000e+16 /* 0x43500000, 0x00000000 */
private const val ivln10 = 4.34294481903251816668e-01 /* 0x3FDBCB7B, 0x1526E50E */
private const val log10_2hi = 3.01029995663611771306e-01 /* 0x3FD34413, 0x509F6000 */
private const val log10_2lo = 3.69423907715893078616e-13 /* 0x3D59FEF3, 0x11F12B36 */
private const val zero = 0.0
internal fun __ieee754_log10(x: Double): Double {
var x: Double = x
var y: Double = 0.0
var z: Double = 0.0
var i: Int = 0
var k: Int = 0
var hx: Int = 0
var lx: UInt = 0U
hx = __HI(x) /* high word of x */
lx = __LOu(x) /* low word of x */
k = 0
if (hx < 0x00100000) { /* x < 2**-1022 */
if (((hx and 0x7fffffff) or lx.toInt()) == 0)
return -two54 / zero /* log(+-0)=-inf */
if (hx < 0) return (x - x) / zero /* log(-#) = NaN */
k -= 54; x *= two54 /* subnormal number, scale up x */
hx = __HI(x) /* high word of x */
}
if (hx >= 0x7ff00000) return x + x
k += (hx shr 20) - 1023
i = ((k.toUInt() and Int.MIN_VALUE.toUInt()) shr 31).toInt()
hx = (hx and 0x000fffff) or ((0x3ff - i) shl 20)
y = (k + i).toDouble()
x = doubleSetWord(d = x, hi = hx)
z = y * log10_2lo + ivln10 * __ieee754_log(x)
return z + y * log10_2hi
}
@@ -0,0 +1,72 @@
/* @(#)e_log2.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_log2(x)
* Return the base 2 logarithm of x
*
* Method :
* Let ivln2 = 1/log(2) rounded.
* Then
* n = ilogb(x),
* if(n<0) n = n+1;
* x = scalbn(x,-n);
* log2(x) := n + ivln2*log(x)
*
* Special cases:
* log2(x) is NaN with signal if x < 0;
* log2(+INF) is +INF with no signal; log2(0) is -INF with signal;
* log2(NaN) is that NaN with no signal;
* log2(2**N) = N for N=1022 to +1023.
*
* Constants:
* The hexadecimal values are the intended ones for the following constants.
* The decimal values may be used, provided that the compiler will convert
* from decimal to binary accurately enough to produce the hexadecimal values
* shown.
*/
package kotlin.math.fdlibm
private const val two54 = 1.80143985094819840000e+16 /* 0x43500000, 0x00000000 */
private const val ivln2 = 0.14426950408889634073e+01 /* 0x3ff71547, 0x652b82fe */
private const val zero = 0.0
internal fun __ieee754_log2(x: Double): Double {
var x: Double = x
var y: Double = 0.0
var z: Double = 0.0
var i: Int = 0
var k: Int = 0
var hx: Int = 0
var lx: UInt = 0U
hx = __HI(x) /* high word of x */
lx = __LOu(x) /* low word of x */
k = 0
if (hx < 0x00100000) { /* x < 2**-1022 */
if (((hx and 0x7fffffff) or lx.toInt()) == 0)
return -two54 / zero /* log(+-0)=-inf */
if (hx < 0) return (x - x) / zero /* log(-#) = NaN */
k -= 54; x *= two54 /* subnormal number, scale up x */
hx = __HI(x) /* high word of x */
}
if (hx >= 0x7ff00000) return x + x
k += (hx shr 20) - 1023
i = ((k.toUInt() and Int.MIN_VALUE.toUInt()) shr 31).toInt()
hx = (hx and 0x000fffff) or ((0x3ff - i) shl 20)
y = (k + i).toDouble()
x = doubleSetWord(d = x, hi = hx)
z = y + ivln2 * __ieee754_log(x)
return z
}
@@ -0,0 +1,333 @@
//#ifndef lint
//static char sccsid[] = "@(#)e_pow.c 1.5 04/04/22 SMI";
//#endif
/*
* ====================================================
* Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
*
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_pow(x,y) return x**y
*
* n
* Method: Let x = 2 * (1+f)
* 1. Compute and return log2(x) in two pieces:
* log2(x) = w1 + w2,
* where w1 has 53-24 = 29 bit trailing zeros.
* 2. Perform y*log2(x) = n+y' by simulating muti-precision
* arithmetic, where |y'|<=0.5.
* 3. Return x**y = 2**n*exp(y'*log2)
*
* Special cases:
* 1. (anything) ** 0 is 1
* 2. (anything) ** 1 is itself
* 3. (anything) ** NAN is NAN
* 4. NAN ** (anything except 0) is NAN
* 5. +-(|x| > 1) ** +INF is +INF
* 6. +-(|x| > 1) ** -INF is +0
* 7. +-(|x| < 1) ** +INF is +0
* 8. +-(|x| < 1) ** -INF is +INF
* 9. +-1 ** +-INF is NAN
* 10. +0 ** (+anything except 0, NAN) is +0
* 11. -0 ** (+anything except 0, NAN, odd integer) is +0
* 12. +0 ** (-anything except 0, NAN) is +INF
* 13. -0 ** (-anything except 0, NAN, odd integer) is +INF
* 14. -0 ** (odd integer) = -( +0 ** (odd integer) )
* 15. +INF ** (+anything except 0,NAN) is +INF
* 16. +INF ** (-anything except 0,NAN) is +0
* 17. -INF ** (anything) = -0 ** (-anything)
* 18. (-anything) ** (integer) is (-1)**(integer)*(+anything**integer)
* 19. (-anything except 0 and inf) ** (non-integer) is NAN
*
* Accuracy:
* pow(x,y) returns x**y nearly rounded. In particular
* pow(integer,integer)
* always returns the correct integer provided it is
* representable.
*
* Constants :
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
package kotlin.math.fdlibm
import kotlin.wasm.internal.wasm_f64_sqrt as sqrt
private val bp = doubleArrayOf(1.0, 1.5)
private val dp_h = doubleArrayOf(0.0, 5.84962487220764160156e-01) /* 0x3FE2B803, 0x40000000 */
private val dp_l = doubleArrayOf(0.0, 1.35003920212974897128e-08) /* 0x3E4CFDEB, 0x43CFD006 */
private const val zero = 0.0
private const val one = 1.0
private const val two = 2.0
private const val two53 = 9007199254740992.0 /* 0x43400000, 0x00000000 */
private const val huge = 1.0e300
private const val tiny = 1.0e-300
/* poly coefs for (3/2)*(log(x)-2s-2/3*s**3 */
private const val L1 = 5.99999999999994648725e-01 /* 0x3FE33333, 0x33333303 */
private const val L2 = 4.28571428578550184252e-01 /* 0x3FDB6DB6, 0xDB6FABFF */
private const val L3 = 3.33333329818377432918e-01 /* 0x3FD55555, 0x518F264D */
private const val L4 = 2.72728123808534006489e-01 /* 0x3FD17460, 0xA91D4101 */
private const val L5 = 2.30660745775561754067e-01 /* 0x3FCD864A, 0x93C9DB65 */
private const val L6 = 2.06975017800338417784e-01 /* 0x3FCA7E28, 0x4A454EEF */
private const val P1 = 1.66666666666666019037e-01 /* 0x3FC55555, 0x5555553E */
private const val P2 = -2.77777777770155933842e-03 /* 0xBF66C16C, 0x16BEBD93 */
private const val P3 = 6.61375632143793436117e-05 /* 0x3F11566A, 0xAF25DE2C */
private const val P4 = -1.65339022054652515390e-06 /* 0xBEBBBD41, 0xC5D26BF1 */
private const val P5 = 4.13813679705723846039e-08 /* 0x3E663769, 0x72BEA4D0 */
private const val lg2 = 6.93147180559945286227e-01 /* 0x3FE62E42, 0xFEFA39EF */
private const val lg2_h = 6.93147182464599609375e-01 /* 0x3FE62E43, 0x00000000 */
private const val lg2_l = -1.90465429995776804525e-09 /* 0xBE205C61, 0x0CA86C39 */
private const val ovt = 8.0085662595372944372e-0017 /* -(1024-log2(ovfl+.5ulp)) */
private const val cp = 9.61796693925975554329e-01 /* 0x3FEEC709, 0xDC3A03FD =2/(3ln2) */
private const val cp_h = 9.61796700954437255859e-01 /* 0x3FEEC709, 0xE0000000 =(float)cp */
private const val cp_l = -7.02846165095275826516e-09 /* 0xBE3E2FE0, 0x145B01F5 =tail of cp_h*/
private const val ivln2 = 1.44269504088896338700e+00 /* 0x3FF71547, 0x652B82FE =1/ln2 */
private const val ivln2_h = 1.44269502162933349609e+00 /* 0x3FF71547, 0x60000000 =24b 1/ln2*/
private const val ivln2_l = 1.92596299112661746887e-08 /* 0x3E54AE0B, 0xF85DDF44 =1/ln2 tail*/
internal fun __ieee754_pow(x: Double, y: Double): Double {
var z: Double = 0.0
var ax: Double = 0.0
var z_h: Double = 0.0
var z_l: Double = 0.0
var p_h: Double = 0.0
var p_l: Double = 0.0
var y1: Double = 0.0
var t1: Double = 0.0
var t2: Double = 0.0
var r: Double = 0.0
var s: Double = 0.0
var t: Double = 0.0
var u: Double = 0.0
var v: Double = 0.0
var w: Double = 0.0
var i0: Int = 0
var i1: Int = 0
var i: Int = 0
var j: Int = 0
var k: Int = 0
var yisint: Int = 0
var n = 0
var hx: Int = 0
var hy: Int = 0
var ix: Int = 0
var iy: Int = 0
var lx: UInt = 0U
var ly: UInt = 0U
//i0 = ((*(int*)&one)>>29)^1
i0 = 1
i1 = 1 - i0
hx = __HI(x); lx = __LO(x).toUInt()
hy = __HI(y); ly = __LO(y).toUInt()
ix = hx and 0x7fffffff; iy = hy and 0x7fffffff
/* y==zero: x**0 = 1 */
if ((iy or ly.toInt()) == 0) return one
/* +-NaN return x+y */
if (ix > 0x7ff00000 || ((ix == 0x7ff00000) && (lx != 0U)) ||
iy > 0x7ff00000 || ((iy == 0x7ff00000) && (ly != 0U))
)
return x + y
/* determine if y is an odd int when x < 0
* yisint = 0 ... y is not an integer
* yisint = 1 ... y is an odd int
* yisint = 2 ... y is an even int
*/
yisint = 0
if (hx < 0) {
if (iy >= 0x43400000) yisint = 2 /* even integer y */
else if (iy >= 0x3ff00000) {
k = (iy shr 20) - 0x3ff /* exponent */
if (k > 20) {
j = (ly shr (52 - k)).toInt()
if ((j shl (52 - k)) == ly.toInt()) yisint = 2 - (j and 1)
} else if (ly == 0U) {
j = iy shr (20 - k)
if ((j shl (20 - k)) == iy) yisint = 2 - (j and 1)
}
}
}
/* special value of y */
if (ly == 0U) {
if (iy == 0x7ff00000) { /* y is +-inf */
if (((ix - 0x3ff00000) or lx.toInt()) == 0)
return y - y /* inf**+-1 is NaN */
else if (ix >= 0x3ff00000)/* (|x|>1)**+-inf = inf,0 */
return if (hy >= 0) y else zero
else /* (|x|<1)**-,+inf = inf,0 */
return if (hy < 0) -y else zero
}
if (iy == 0x3ff00000) { /* y is +-1 */
if (hy < 0) return one / x; else return x
}
if (hy == 0x40000000) return x * x /* y is 2 */
if (hy == 0x3fe00000) { /* y is 0.5 */
if (hx >= 0) /* x >= +0 */
return sqrt(x)
}
}
ax = fabs(x)
/* special value of x */
if (lx == 0U) {
if (ix == 0x7ff00000 || ix == 0 || ix == 0x3ff00000) {
z = ax /*x is +-0,+-inf,+-1*/
if (hy < 0) z = one / z /* z = (1/|x|) */
if (hx < 0) {
if (((ix - 0x3ff00000) or yisint) == 0) {
z = (z - z) / (z - z) /* (-1)**non-int is NaN */
} else if (yisint == 1)
z = -z /* (x<0)**odd = -(|x|**odd) */
}
return z
}
}
n = (hx shr 31) + 1
/* (x<0)**(non-int) is NaN */
if ((n or yisint) == 0) return (x - x) / (x - x)
s = one /* s (sign of result -ve**odd) = -1 else = 1 */
if ((n or (yisint - 1)) == 0) s = -one/* (-ve)**(odd int) */
/* |y| is huge */
if (iy > 0x41e00000) { /* if |y| > 2**31 */
if (iy > 0x43f00000) { /* if |y| > 2**64, must o/uflow */
if (ix <= 0x3fefffff) return if (hy < 0) huge * huge else tiny * tiny
if (ix >= 0x3ff00000) return if (hy > 0) huge * huge else tiny * tiny
}
/* over/underflow if x is not close to one */
if (ix < 0x3fefffff) return if (hy < 0) s * huge * huge else s * tiny * tiny
if (ix > 0x3ff00000) return if (hy > 0) s * huge * huge else s * tiny * tiny
/* now |1-x| is tiny <= 2**-20, suffice to compute
log(x) by x-x^2/2+x^3/3-x^4/4 */
t = ax - one /* t has 20 trailing zeros */
w = (t * t) * (0.5 - t * (0.3333333333333333333333 - t * 0.25))
u = ivln2_h * t /* ivln2_h has 21 sig. bits */
v = t * ivln2_l - w * ivln2
t1 = u + v
t1 = doubleSetWord(d = t1, lo = 0)
t2 = v - (t1 - u)
} else {
var ss: Double = 0.0
var s2: Double = 0.0
var s_h: Double = 0.0
var s_l: Double = 0.0
var t_h: Double = 0.0
var t_l: Double = 0.0
n = 0
/* take care subnormal number */
if (ix < 0x00100000) {
ax *= two53; n -= 53; ix = __HI(ax); }
n += ((ix) shr 20) - 0x3ff
j = ix and 0x000fffff
/* determine interval */
ix = j or 0x3ff00000 /* normalize ix */
if (j <= 0x3988E) k = 0 /* |x|<sqrt(3/2) */
else if (j < 0xBB67A) k = 1 /* |x|<sqrt(3) */
else {
k = 0;n += 1;ix -= 0x00100000;
}
ax = doubleSetWord(d = ax, hi = ix)
/* compute ss = s_h+s_l = (x-1)/(x+1) or (x-1.5)/(x+1.5) */
u = ax - bp[k] /* bp[0]=1.0, bp[1]=1.5 */
v = one / (ax + bp[k])
ss = u * v
s_h = ss
s_h = doubleSetWord(d = s_h, lo = 0)
/* t_h=ax+bp[k] High */
t_h = zero
t_h = doubleSetWord(d = t_h, hi = ((ix shr 1) or 0x20000000) + 0x00080000 + (k shl 18))
t_l = ax - (t_h - bp[k])
s_l = v * ((u - s_h * t_h) - s_h * t_l)
/* compute log(ax) */
s2 = ss * ss
r = s2 * s2 * (L1 + s2 * (L2 + s2 * (L3 + s2 * (L4 + s2 * (L5 + s2 * L6)))))
r += s_l * (s_h + ss)
s2 = s_h * s_h
t_h = 3.0 + s2 + r
t_h = doubleSetWord(d = t_h, lo = 0)
t_l = r - ((t_h - 3.0) - s2)
/* u+v = ss*(1+...) */
u = s_h * t_h
v = s_l * t_h + t_l * ss
/* 2/(3log2)*(ss+...) */
p_h = u + v
p_h = doubleSetWord(d = p_h, lo = 0)
p_l = v - (p_h - u)
z_h = cp_h * p_h /* cp_h+cp_l = 2/(3*log2) */
z_l = cp_l * p_h + p_l * cp + dp_l[k]
/* log2(ax) = (ss+..)*2/(3*log2) = n + dp_h + z_h + z_l */
t = n.toDouble()
t1 = (((z_h + z_l) + dp_h[k]) + t)
t1 = doubleSetWord(d = t1, lo = 0)
t2 = z_l - (((t1 - t) - dp_h[k]) - z_h)
}
/* split up y into y1+y2 and compute (y1+y2)*(t1+t2) */
y1 = y
y1 = doubleSetWord(d = y1, lo = 0)
p_l = (y - y1) * t1 + y * t2
p_h = y1 * t1
z = p_l + p_h
j = __HI(z)
i = __LO(z)
if (j >= 0x40900000) { /* z >= 1024 */
if (((j - 0x40900000) or i) != 0) /* if z > 1024 */
return s * huge * huge /* overflow */
else {
if (p_l + ovt > z - p_h) return s * huge * huge /* overflow */
}
} else if ((j and 0x7fffffff) >= 0x4090cc00) { /* z <= -1075 */
if (((j - 0xc090cc00) or i.toLong()) != 0L) /* z < -1075 */
return s * tiny * tiny /* underflow */
else {
if (p_l <= z - p_h) return s * tiny * tiny /* underflow */
}
}
/*
* compute 2**(p_h+p_l)
*/
i = j and 0x7fffffff
k = (i shr 20) - 0x3ff
n = 0
if (i > 0x3fe00000) { /* if |z| > 0.5, set n = [z+0.5] */
n = j + (0x00100000 shr (k + 1))
k = ((n and 0x7fffffff) shr 20) - 0x3ff /* new k for n */
t = zero
t = doubleSetWord(d = t, hi = (n and (0x000fffff shr k).inv()))
n = ((n and 0x000fffff) or 0x00100000) shr (20 - k)
if (j < 0) n = -n
p_h -= t
}
t = p_l + p_h
t = doubleSetWord(d = t, lo = 0)
u = t * lg2_h
v = (p_l - (t - p_h)) * lg2 + t * lg2_l
z = u + v
w = v - (z - u)
t = z * z
t1 = z - t * (P1 + t * (P2 + t * (P3 + t * (P4 + t * P5))))
r = (z * t1) / (t1 - two) - (w + z * w)
z = one - (r - z)
j = __HI(z)
j += (n shl 20)
if ((j shr 20) <= 0) z = scalbn(z, n) /* subnormal output */
else z = doubleSetWord(d = z, hi = __HI(z) + (n shl 20))
return s * z
}
@@ -0,0 +1,175 @@
/* @(#)e_rem_pio2.c 1.4 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
/* __ieee754_rem_pio2(x,y)
*
* return the remainder of x rem pi/2 in y[0]+y[1]
* use __kernel_rem_pio2()
*/
package kotlin.math.fdlibm
/*
* Table of constants for 2/pi, 396 Hex digits (476 decimal) of 2/pi
*/
private val two_over_pi = intArrayOf(
0xA2F983, 0x6E4E44, 0x1529FC, 0x2757D1, 0xF534DD, 0xC0DB62,
0x95993C, 0x439041, 0xFE5163, 0xABDEBB, 0xC561B7, 0x246E3A,
0x424DD2, 0xE00649, 0x2EEA09, 0xD1921C, 0xFE1DEB, 0x1CB129,
0xA73EE8, 0x8235F5, 0x2EBB44, 0x84E99C, 0x7026B4, 0x5F7E41,
0x3991D6, 0x398353, 0x39F49C, 0x845F8B, 0xBDF928, 0x3B1FF8,
0x97FFDE, 0x05980F, 0xEF2F11, 0x8B5A0A, 0x6D1F6D, 0x367ECF,
0x27CB09, 0xB74F46, 0x3F669E, 0x5FEA2D, 0x7527BA, 0xC7EBE5,
0xF17B3D, 0x0739F7, 0x8A5292, 0xEA6BFB, 0x5FB11F, 0x8D5D08,
0x560330, 0x46FC7B, 0x6BABF0, 0xCFBC20, 0x9AF436, 0x1DA9E3,
0x91615E, 0xE61B08, 0x659985, 0x5F14A0, 0x68408D, 0xFFD880,
0x4D7327, 0x310606, 0x1556CA, 0x73A8C9, 0x60E27B, 0xC08C6B,
)
private val npio2_hw = intArrayOf(
0x3FF921FB, 0x400921FB, 0x4012D97C, 0x401921FB, 0x401F6A7A, 0x4022D97C,
0x4025FDBB, 0x402921FB, 0x402C463A, 0x402F6A7A, 0x4031475C, 0x4032D97C,
0x40346B9C, 0x4035FDBB, 0x40378FDB, 0x403921FB, 0x403AB41B, 0x403C463A,
0x403DD85A, 0x403F6A7A, 0x40407E4C, 0x4041475C, 0x4042106C, 0x4042D97C,
0x4043A28C, 0x40446B9C, 0x404534AC, 0x4045FDBB, 0x4046C6CB, 0x40478FDB,
0x404858EB, 0x404921FB,
)
/*
* invpio2: 53 bits of 2/pi
* pio2_1: first 33 bit of pi/2
* pio2_1t: pi/2 - pio2_1
* pio2_2: second 33 bit of pi/2
* pio2_2t: pi/2 - (pio2_1+pio2_2)
* pio2_3: third 33 bit of pi/2
* pio2_3t: pi/2 - (pio2_1+pio2_2+pio2_3)
*/
private const val zero = 0.00000000000000000000e+00 /* 0x00000000, 0x00000000 */
private const val half = 5.00000000000000000000e-01 /* 0x3FE00000, 0x00000000 */
private const val two24 = 1.67772160000000000000e+07 /* 0x41700000, 0x00000000 */
private const val invpio2 = 6.36619772367581382433e-01 /* 0x3FE45F30, 0x6DC9C883 */
private const val pio2_1 = 1.57079632673412561417e+00 /* 0x3FF921FB, 0x54400000 */
private const val pio2_1t = 6.07710050650619224932e-11 /* 0x3DD0B461, 0x1A626331 */
private const val pio2_2 = 6.07710050630396597660e-11 /* 0x3DD0B461, 0x1A600000 */
private const val pio2_2t = 2.02226624879595063154e-21 /* 0x3BA3198A, 0x2E037073 */
private const val pio2_3 = 2.02226624871116645580e-21 /* 0x3BA3198A, 0x2E000000 */
private const val pio2_3t = 8.47842766036889956997e-32 /* 0x397B839A, 0x252049C1 */
internal fun __ieee754_rem_pio2(x: Double, y: DoubleArray): Int {
var z: Double = 0.0
var w: Double = 0.0
var t: Double = 0.0
var r: Double = 0.0
var fn: Double = 0.0
val tx: DoubleArray = DoubleArray(3)
var e0: Int = 0
var i: Int = 0
var j: Int = 0
var nx: Int = 0
var n: Int = 0
var ix: Int = 0
var hx: Int = 0
hx = __HI(x) /* high word of x */
ix = hx and 0x7fffffff
if (ix <= 0x3fe921fb) /* |x| ~<= pi/4 , no need for reduction */ {
y[0] = x; y[1] = 0.0; return 0
}
if (ix < 0x4002d97c) { /* |x| < 3pi/4, special case with n=+-1 */
if (hx > 0) {
z = x - pio2_1
if (ix != 0x3ff921fb) { /* 33+53 bit pi is good enough */
y[0] = z - pio2_1t
y[1] = (z - y[0]) - pio2_1t
} else { /* near pi/2, use 33+33+53 bit pi */
z -= pio2_2
y[0] = z - pio2_2t
y[1] = (z - y[0]) - pio2_2t
}
return 1
} else { /* negative x */
z = x + pio2_1
if (ix != 0x3ff921fb) { /* 33+53 bit pi is good enough */
y[0] = z + pio2_1t
y[1] = (z - y[0]) + pio2_1t
} else { /* near pi/2, use 33+33+53 bit pi */
z += pio2_2
y[0] = z + pio2_2t
y[1] = (z - y[0]) + pio2_2t
}
return -1
}
}
if (ix <= 0x413921fb) { /* |x| ~<= 2^19*(pi/2), medium size */
t = fabs(x)
n = (t * invpio2 + half).toInt()
fn = n.toDouble()
r = t - fn * pio2_1
w = fn * pio2_1t /* 1st round good to 85 bit */
if (n < 32 && ix != npio2_hw[n - 1]) {
y[0] = r - w /* quick check no cancellation */
} else {
j = ix shr 20
y[0] = r - w
i = j - (((__HI(y[0])) shr 20) and 0x7ff)
if (i > 16) { /* 2nd iteration needed, good to 118 */
t = r
w = fn * pio2_2
r = t - w
w = fn * pio2_2t - ((t - r) - w)
y[0] = r - w
i = j - (((__HI(y[0])) shr 20) and 0x7ff)
if (i > 49) { /* 3rd iteration need, 151 bits acc */
t = r /* will cover all possible cases */
w = fn * pio2_3
r = t - w
w = fn * pio2_3t - ((t - r) - w)
y[0] = r - w
}
}
}
y[1] = (r - y[0]) - w
if (hx < 0) {
y[0] = -y[0]; y[1] = -y[1]; return -n
} else return n
}
/*
* all other (large) arguments
*/
if (ix >= 0x7ff00000) { /* x is inf or NaN */
y[1] = x - x
y[0] = y[1]; return 0
}
/* set z = scalbn(|x|,ilogb(x)-23) */
z = doubleSetWord(d = z, lo = __LO(x))
e0 = (ix shr 20) - 1046 /* e0 = ilogb(z)-23; */
z = doubleSetWord(d = z, hi = ix - (e0 shl 20))
//for(i=0;i<2;i++) {
i = 0
while (i < 2) {
tx[i] = (z.toInt()).toDouble()
z = (z - tx[i]) * two24
//--
i++
}
tx[2] = z
nx = 3
while (tx[nx - 1] == zero) nx-- /* skip zero term */
n = __kernel_rem_pio2(tx, y, e0, nx, 2, two_over_pi)
if (hx < 0) {
y[0] = -y[0]; y[1] = -y[1]; return -n
}
return n
}
@@ -0,0 +1,75 @@
/* @(#)e_sinh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __ieee754_sinh(x)
* Method :
* mathematically sinh(x) if defined to be (exp(x)-exp(-x))/2
* 1. Replace x by |x| (sinh(-x) = -sinh(x)).
* 2.
* E + E/(E+1)
* 0 <= x <= 22 : sinh(x) := --------------, E=expm1(x)
* 2
*
* 22 <= x <= lnovft : sinh(x) := exp(x)/2
* lnovft <= x <= ln2ovft: sinh(x) := exp(x/2)/2 * exp(x/2)
* ln2ovft < x : sinh(x) := x*shuge (overflow)
*
* Special cases:
* sinh(x) is |x| if x is +INF, -INF, or NaN.
* only sinh(0)=0 is exact for finite x.
*/
package kotlin.math.fdlibm
private const val one = 1.0
private const val shuge = 1.0e307
internal fun __ieee754_sinh(x: Double): Double {
var t: Double = 0.0
var w: Double = 0.0
var h: Double = 0.0
var ix: Int = 0
var jx: Int = 0
var lx: UInt = 0U
/* High word of |x|. */
jx = __HI(x)
ix = jx and 0x7fffffff
/* x is INF or NaN */
if (ix >= 0x7ff00000) return x + x
h = 0.5
if (jx < 0) h = -h
/* |x| in [0,22], return sign(x)*0.5*(E+E/(E+1))) */
if (ix < 0x40360000) { /* |x|<22 */
if (ix < 0x3e300000) /* |x|<2**-28 */
if (shuge + x > one) return x/* sinh(tiny) = tiny with inexact */
t = expm1(fabs(x))
if (ix < 0x3ff00000) return h * (2.0 * t - t * t / (t + one))
return h * (t + t / (t + one))
}
/* |x| in [22, log(maxdouble)] return 0.5*exp(|x|) */
if (ix < 0x40862E42) return h * __ieee754_exp(fabs(x))
/* |x| in [log(maxdouble), overflowthresold] */
//lx = *( (((*(unsigned*)&one) shr 29)) + (unsigned*)&x);
lx = (x.toBits() and 0xFFFFFFFF).toUInt()
if (ix < 0x408633CE || (ix == 0x408633ce) && (lx <= 0x8fb9f87d.toUInt())) {
w = __ieee754_exp(0.5 * fabs(x))
t = h * w
return t * w
}
/* |x| > overflowthresold, sinh(x) overflow */
return x * shuge
}
@@ -0,0 +1,85 @@
/* @(#)k_cos.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/*
* __kernel_cos( x, y )
* kernel cos function on [-pi/4, pi/4], pi/4 ~ 0.785398164
* Input x is assumed to be bounded by ~pi/4 in magnitude.
* Input y is the tail of x.
*
* Algorithm
* 1. Since cos(-x) = cos(x), we need only to consider positive x.
* 2. if x < 2^-27 (hx<0x3e400000 0), return 1 with inexact if x!=0.
* 3. cos(x) is approximated by a polynomial of degree 14 on
* [0,pi/4]
* 4 14
* cos(x) ~ 1 - x*x/2 + C1*x + ... + C6*x
* where the remez error is
*
* | 2 4 6 8 10 12 14 | -58
* |cos(x)-(1-.5*x +C1*x +C2*x +C3*x +C4*x +C5*x +C6*x )| <= 2
* | |
*
* 4 6 8 10 12 14
* 4. let r = C1*x +C2*x +C3*x +C4*x +C5*x +C6*x , then
* cos(x) = 1 - x*x/2 + r
* since cos(x+y) ~ cos(x) - sin(x)*y
* ~ cos(x) - x*y,
* a correction term is necessary in cos(x) and hence
* cos(x+y) = 1 - (x*x/2 - (r - x*y))
* For better accuracy when x > 0.3, let qx = |x|/4 with
* the last 32 bits mask off, and if x > 0.78125, let qx = 0.28125.
* Then
* cos(x+y) = (1-qx) - ((x*x/2-qx) - (r-x*y)).
* Note that 1-qx and (x*x/2-qx) is EXACT here, and the
* magnitude of the latter is at least a quarter of x*x/2,
* thus, reducing the rounding error in the subtraction.
*/
package kotlin.math.fdlibm
private const val one = 1.00000000000000000000e+00 /* 0x3FF00000, 0x00000000 */
private const val C1 = 4.16666666666666019037e-02 /* 0x3FA55555, 0x5555554C */
private const val C2 = -1.38888888888741095749e-03 /* 0xBF56C16C, 0x16C15177 */
private const val C3 = 2.48015872894767294178e-05 /* 0x3EFA01A0, 0x19CB1590 */
private const val C4 = -2.75573143513906633035e-07 /* 0xBE927E4F, 0x809C52AD */
private const val C5 = 2.08757232129817482790e-09 /* 0x3E21EE9E, 0xBDB4B1C4 */
private const val C6 = -1.13596475577881948265e-11 /* 0xBDA8FAE9, 0xBE8838D4 */
internal fun __kernel_cos(x: Double, y: Double): Double {
var a: Double = 0.0
var hz: Double = 0.0
var z: Double = 0.0
var r: Double = 0.0
var qx: Double = 0.0
var ix: Int = 0
ix = __HI(x) and 0x7fffffff /* ix = |x|'s high word*/
if (ix < 0x3e400000) { /* if x < 2**27 */
if ((x.toInt()) == 0) return one /* generate inexact */
}
z = x * x
r = z * (C1 + z * (C2 + z * (C3 + z * (C4 + z * (C5 + z * C6)))))
if (ix < 0x3FD33333) /* if |x| < 0.3 */
return one - (0.5 * z - (z * r - x * y))
else {
if (ix > 0x3fe90000) { /* x > 0.78125 */
qx = 0.28125
} else {
qx = doubleSetWord(d = qx, hi = ix - 0x00200000) /* x/4 */
qx = doubleSetWord(d = qx, lo = 0)
}
hz = 0.5 * z - qx
a = one - qx
return a - (hz - (z * r - x * y))
}
}
@@ -0,0 +1,408 @@
/* @(#)k_rem_pio2.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/*
* __kernel_rem_pio2(x,y,e0,nx,prec,ipio2)
* double x[],y[]; int e0,nx,prec; int ipio2[];
*
* __kernel_rem_pio2 return the last three digits of N with
* y = x - N*pi/2
* so that |y| < pi/2.
*
* The method is to compute the integer (mod 8) and fraction parts of
* (2/pi)*x without doing the full multiplication. In general we
* skip the part of the product that are known to be a huge integer (
* more accurately, = 0 mod 8 ). Thus the number of operations are
* independent of the exponent of the input.
*
* (2/pi) is represented by an array of 24-bit integers in ipio2[].
*
* Input parameters:
* x[] The input value (must be positive) is broken into nx
* pieces of 24-bit integers in double precision format.
* x[i] will be the i-th 24 bit of x. The scaled exponent
* of x[0] is given in input parameter e0 (i.e., x[0]*2^e0
* match x's up to 24 bits.
*
* Example of breaking a double positive z into x[0]+x[1]+x[2]:
* e0 = ilogb(z)-23
* z = scalbn(z,-e0)
* for i = 0,1,2
* x[i] = floor(z)
* z = (z-x[i])*2**24
*
*
* y[] ouput result in an array of double precision numbers.
* The dimension of y[] is:
* 24-bit precision 1
* 53-bit precision 2
* 64-bit precision 2
* 113-bit precision 3
* The actual value is the sum of them. Thus for 113-bit
* precison, one may have to do something like:
*
* long double t,w,r_head, r_tail;
* t = (long double)y[2] + (long double)y[1];
* w = (long double)y[0];
* r_head = t+w;
* r_tail = w - (r_head - t);
*
* e0 The exponent of x[0]
*
* nx dimension of x[]
*
* prec an integer indicating the precision:
* 0 24 bits (single)
* 1 53 bits (double)
* 2 64 bits (extended)
* 3 113 bits (quad)
*
* ipio2[]
* integer array, contains the (24*i)-th to (24*i+23)-th
* bit of 2/pi after binary point. The corresponding
* floating value is
*
* ipio2[i] * 2^(-24(i+1)).
*
* External function:
* double scalbn(), floor();
*
*
* Here is the description of some local variables:
*
* jk jk+1 is the initial number of terms of ipio2[] needed
* in the computation. The recommended value is 2,3,4,
* 6 for single, double, extended,and quad.
*
* jz local integer variable indicating the number of
* terms of ipio2[] used.
*
* jx nx - 1
*
* jv index for pointing to the suitable ipio2[] for the
* computation. In general, we want
* ( 2^e0*x[0] * ipio2[jv-1]*2^(-24jv) )/8
* is an integer. Thus
* e0-3-24*jv >= 0 or (e0-3)/24 >= jv
* Hence jv = max(0,(e0-3)/24).
*
* jp jp+1 is the number of terms in PIo2[] needed, jp = jk.
*
* q[] double array with integral value, representing the
* 24-bits chunk of the product of x and 2/pi.
*
* q0 the corresponding exponent of q[0]. Note that the
* exponent for q[i] would be q0-24*i.
*
* PIo2[] double precision array, obtained by cutting pi/2
* into 24 bits chunks.
*
* f[] ipio2[] in floating point
*
* iq[] integer array by breaking up q[] in 24-bits chunk.
*
* fq[] final product of x*(2/pi) in fq[0],..,fq[jk]
*
* ih integer. If >0 it indicates q[] is >= 0.5, hence
* it also indicates the *sign* of the result.
*
*/
/*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
package kotlin.math.fdlibm
import kotlin.wasm.internal.wasm_f64_floor as floor
private val init_jk = intArrayOf(2, 3, 4, 6) /* initial value for jk */
private val PIo2 = doubleArrayOf(
1.57079625129699707031e+00, /* 0x3FF921FB, 0x40000000 */
7.54978941586159635335e-08, /* 0x3E74442D, 0x00000000 */
5.39030252995776476554e-15, /* 0x3CF84698, Int.MIN_VALUE */
3.28200341580791294123e-22, /* 0x3B78CC51, 0x60000000 */
1.27065575308067607349e-29, /* 0x39F01B83, Int.MIN_VALUE */
1.22933308981111328932e-36, /* 0x387A2520, 0x40000000 */
2.73370053816464559624e-44, /* 0x36E38222, Int.MIN_VALUE */
2.16741683877804819444e-51, /* 0x3569F31D, 0x00000000 */
)
private const val zero = 0.0
private const val one = 1.0
private const val two24 = 1.67772160000000000000e+07 /* 0x41700000, 0x00000000 */
private const val twon24 = 5.96046447753906250000e-08 /* 0x3E700000, 0x00000000 */
internal fun __kernel_rem_pio2(x: DoubleArray, y: DoubleArray, e0: Int, nx: Int, prec: Int, ipio2: IntArray): Int {
var jz: Int = 0
var jx: Int = 0
var jv: Int = 0
var jp: Int = 0
var jk: Int = 0
var carry: Int = 0
var n: Int = 0
var iq: IntArray = IntArray(20)
var i: Int = 0
var j: Int = 0
var k: Int = 0
var m: Int = 0
var q0: Int = 0
var ih: Int = 0
var z: Double = 0.0
var fw: Double = 0.0
var f: DoubleArray = DoubleArray(20)
var fq: DoubleArray = DoubleArray(20)
var q: DoubleArray = DoubleArray(20)
/* initialize jk*/
jk = init_jk[prec]
jp = jk
/* determine jx,jv,q0, note that 3>q0 */
jx = nx - 1
jv = (e0 - 3) / 24; if (jv < 0) jv = 0
q0 = e0 - 24 * (jv + 1)
/* set up f[0] to f[jx+jk] where f[jx+jk] = ipio2[jv+jk] */
j = jv - jx; m = jx + jk
//for(i=0;i<=m;i++,j++) {
i = 0
while (i <= m) {
f[i] = if (j < 0) zero else ipio2[j].toDouble()
//--
i++; j++
}
/* compute q[0],q[1],...q[jk] */
//for (i=0;i<=jk;i++) {
i = 0
while (i <= jk) {
j = 0; fw = 0.0
while (j <= jx) {
fw += x[j] * f[jx + i - j]; q[i] = fw
//--
j++
}
//--
i++
}
jz = jk
goto@ while (true) {
/* distill q[] into iq[] reversingly */
//for(i=0,j=jz,z=q[jz];j>0;i++,j--) {
i = 0; j = jz; z = q[jz]
while (j > 0) {
fw = ((twon24 * z).toInt()).toDouble()
iq[i] = (z - two24 * fw).toInt()
z = q[j - 1] + fw
//--
i++; j--
}
/* compute n */
z = scalbn(z, q0) /* actual value of z */
z -= 8.0 * floor(z * 0.125) /* trim off integer >= 8 */
n = z.toInt()
z -= n.toDouble()
ih = 0
if (q0 > 0) { /* need iq[jz-1] to determine n */
i = (iq[jz - 1] shr (24 - q0)); n += i
iq[jz - 1] -= i shl (24 - q0)
ih = iq[jz - 1] shr (23 - q0)
} else if (q0 == 0) ih = iq[jz - 1] shr 23
else if (z >= 0.5) ih = 2
if (ih > 0) { /* q > 0.5 */
n += 1; carry = 0
//for(i=0;i<jz ;i++) { /* compute 1-q */
i = 0
while (i < jz) { /* compute 1-q */
j = iq[i]
if (carry == 0) {
if (j != 0) {
carry = 1; iq[i] = 0x1000000 - j
}
} else iq[i] = 0xffffff - j
//--
i++
}
if (q0 > 0) { /* rare case: chance is 1 in 12 */
when (q0) {
1 ->
iq[jz - 1] = iq[jz - 1] and 0x7fffff
2 ->
iq[jz - 1] = iq[jz - 1] and 0x3fffff
}
}
if (ih == 2) {
z = one - z
if (carry != 0) z -= scalbn(one, q0)
}
}
/* check if recomputation is needed */
if (z == zero) {
j = 0
//for (i=jz-1;i>=jk;i--) {
i = jz - 1
while (i >= jk) {
j = j or iq[i]
//--
i--
}
if (j == 0) { /* need recomputation */
//for(k=1;iq[jk-k]==0;k++); /* k = no. of terms needed */
k = 1
while (iq[jk - k] == 0) { /* k = no. of terms needed */
//--
k++
}
//for(i=jz+1;i<=jz+k;i++) { /* add q[jz+1] to q[jz+k] */
i = jz + 1
while (i <= jz + k) { /* add q[jz+1] to q[jz+k] */
f[jx + i] = ipio2[jv + i].toDouble()
//for(j=0,fw=0.0;j<=jx;j++) {
j = 0; fw = 0.0
while (j <= jx) {
fw += x[j] * f[jx + i - j]
//--
j++
}
q[i] = fw
//--
i++
}
jz += k
continue@goto
}
}
break@goto
}
/* chop off zero terms */
if (z == 0.0) {
jz -= 1; q0 -= 24
while (iq[jz] == 0) {
jz--; q0 -= 24
}
} else { /* break z into 24-bit if necessary */
z = scalbn(z, -q0)
if (z >= two24) {
fw = ((twon24 * z).toInt()).toDouble()
iq[jz] = (z - two24 * fw).toInt()
jz += 1; q0 += 24
iq[jz] = fw.toInt()
} else iq[jz] = z.toInt()
}
/* convert integer "bit" chunk to floating-point value */
fw = scalbn(one, q0)
//for(i=jz;i>=0;i--) {
i = jz
while (i >= 0) {
q[i] = fw * iq[i].toDouble(); fw *= twon24
//--
i--
}
/* compute PIo2[0,...,jp]*q[jz,...,0] */
//for(i=jz;i>=0;i--) {
i = jz
while (i >= 0) {
//for(fw=0.0,k=0;k<=jp&&k<=jz-i;k++) {
fw = 0.0; k = 0
while (k <= jp && k <= jz - i) {
fw += PIo2[k] * q[i + k]
//--
k++
}
fq[jz - i] = fw
//--
i--
}
/* compress fq[] into y[] */
when (prec) {
0 -> {
fw = 0.0
//for (i= jz;i >= 0;i--) {
i = jz
while (i >= 0) {
fw += fq[i]
//--
i--
}
y[0] = if (ih == 0) fw else -fw
}
1, 2 -> {
fw = 0.0
//for (i= jz;i >= 0;i--) {
i = jz
while (i >= 0) {
fw += fq[i]
//--
i--
}
y[0] = if (ih == 0) fw else -fw
fw = fq[0] - fw
//for (i= 1;i <= jz;i++) {
i = 1
while (i <= jz) {
fw += fq[i]
//--
i++
}
y[1] = if (ih == 0) fw else -fw
}
3 -> { /* painful */
//for (i= jz;i > 0;i--) {
i = jz
while (i > 0) {
fw = fq[i - 1] + fq[i]
fq[i] += fq[i - 1] - fw
fq[i - 1] = fw
//--
i--
}
//for (i= jz;i > 1;i--) {
i = jz
while (i > 1) {
fw = fq[i - 1] + fq[i]
fq[i] += fq[i - 1] - fw
fq[i - 1] = fw
//--
i--
}
//for (fw= 0.0, i = jz;i >= 2;i--) {
fw = 0.0; i = jz
while (i >= 2) {
fw += fq[i]
//--
i--
}
if (ih == 0) {
y[0] = fq[0]; y[1] = fq[1]; y[2] = fw
} else {
y[0] = -fq[0]; y[1] = -fq[1]; y[2] = -fw
}
}
}
return n and 7
}
@@ -0,0 +1,66 @@
/* @(#)k_sin.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* __kernel_sin( x, y, iy)
* kernel sin function on [-pi/4, pi/4], pi/4 ~ 0.7854
* Input x is assumed to be bounded by ~pi/4 in magnitude.
* Input y is the tail of x.
* Input iy indicates whether y is 0. (if iy=0, y assume to be 0).
*
* Algorithm
* 1. Since sin(-x) = -sin(x), we need only to consider positive x.
* 2. if x < 2^-27 (hx<0x3e400000 0), return x with inexact if x!=0.
* 3. sin(x) is approximated by a polynomial of degree 13 on
* [0,pi/4]
* 3 13
* sin(x) ~ x + S1*x + ... + S6*x
* where
*
* |sin(x) 2 4 6 8 10 12 | -58
* |----- - (1+S1*x +S2*x +S3*x +S4*x +S5*x +S6*x )| <= 2
* | x |
*
* 4. sin(x+y) = sin(x) + sin'(x')*y
* ~ sin(x) + (1-x*x/2)*y
* For better accuracy, let
* 3 2 2 2 2
* r = x *(S2+x *(S3+x *(S4+x *(S5+x *S6))))
* then 3 2
* sin(x) = x + (S1*x + (x *(r-y/2)+y))
*/
package kotlin.math.fdlibm
private const val half = 5.00000000000000000000e-01 /* 0x3FE00000, 0x00000000 */
private const val S1 = -1.66666666666666324348e-01 /* 0xBFC55555, 0x55555549 */
private const val S2 = 8.33333333332248946124e-03 /* 0x3F811111, 0x1110F8A6 */
private const val S3 = -1.98412698298579493134e-04 /* 0xBF2A01A0, 0x19C161D5 */
private const val S4 = 2.75573137070700676789e-06 /* 0x3EC71DE3, 0x57B1FE7D */
private const val S5 = -2.50507602534068634195e-08 /* 0xBE5AE5E6, 0x8A2B9CEB */
private const val S6 = 1.58969099521155010221e-10 /* 0x3DE5D93A, 0x5ACFD57C */
internal fun __kernel_sin(x: Double, y: Double, iy: Int): Double {
var z: Double = 0.0
var r: Double = 0.0
var v: Double = 0.0
var ix: Int = 0
ix = __HI(x) and 0x7fffffff /* high word of x */
if (ix < 0x3e400000) /* |x| < 2**-27 */ {
if (x.toInt() == 0) return x
} /* generate inexact */
z = x * x
v = z * x
r = S2 + z * (S3 + z * (S4 + z * (S5 + z * S6)))
if (iy == 0) return x + v * (S1 + z * r)
else return x - ((z * (half * y - v * r) - y) - v * S1)
}
@@ -0,0 +1,158 @@
//#pragma ident "@(#)k_tan.c 1.5 04/04/22 SMI"
/*
* ====================================================
* Copyright 2004 Sun Microsystems, Inc. All Rights Reserved.
*
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* INDENT OFF */
/* __kernel_tan( x, y, k )
* kernel tan function on [-pi/4, pi/4], pi/4 ~ 0.7854
* Input x is assumed to be bounded by ~pi/4 in magnitude.
* Input y is the tail of x.
* Input k indicates whether tan (if k = 1) or -1/tan (if k = -1) is returned.
*
* Algorithm
* 1. Since tan(-x) = -tan(x), we need only to consider positive x.
* 2. if x < 2^-28 (hx<0x3e300000 0), return x with inexact if x!=0.
* 3. tan(x) is approximated by a odd polynomial of degree 27 on
* [0,0.67434]
* 3 27
* tan(x) ~ x + T1*x + ... + T13*x
* where
*
* |tan(x) 2 4 26 | -59.2
* |----- - (1+T1*x +T2*x +.... +T13*x )| <= 2
* | x |
*
* Note: tan(x+y) = tan(x) + tan'(x)*y
* ~ tan(x) + (1+x*x)*y
* Therefore, for better accuracy in computing tan(x+y), let
* 3 2 2 2 2
* r = x *(T2+x *(T3+x *(...+x *(T12+x *T13))))
* then
* 3 2
* tan(x+y) = x + (T1*x + (x *(r+y)+y))
*
* 4. For x in [0.67434,pi/4], let y = pi/4 - x, then
* tan(x) = tan(pi/4-y) = (1-tan(y))/(1+tan(y))
* = 1 - 2*(tan(y) - (tan(y)^2)/(1+tan(y)))
*/
package kotlin.math.fdlibm
private val xxx = doubleArrayOf(
3.33333333333334091986e-01, /* 3FD55555, 55555563 */
1.33333333333201242699e-01, /* 3FC11111, 1110FE7A */
5.39682539762260521377e-02, /* 3FABA1BA, 1BB341FE */
2.18694882948595424599e-02, /* 3F9664F4, 8406D637 */
8.86323982359930005737e-03, /* 3F8226E3, E96E8493 */
3.59207910759131235356e-03, /* 3F6D6D22, C9560328 */
1.45620945432529025516e-03, /* 3F57DBC8, FEE08315 */
5.88041240820264096874e-04, /* 3F4344D8, F2F26501 */
2.46463134818469906812e-04, /* 3F3026F7, 1A8D1068 */
7.81794442939557092300e-05, /* 3F147E88, A03792A6 */
7.14072491382608190305e-05, /* 3F12B80F, 32F0A7E9 */
-1.85586374855275456654e-05, /* BEF375CB, DB605373 */
2.59073051863633712884e-05, /* 3EFB2A70, 74BF7AD4 */
/* one */ 1.00000000000000000000e+00, /* 3FF00000, 00000000 */
/* pio4 */ 7.85398163397448278999e-01, /* 3FE921FB, 54442D18 */
/* pio4lo */ 3.06161699786838301793e-17 /* 3C81A626, 33145C07 */
)
private val one = xxx[13]
private val pio4 = xxx[14]
private val pio4lo = xxx[15]
private val T = xxx
/* INDENT ON */
internal fun __kernel_tan(x: Double, y: Double, iy: Int): Double {
var x: Double = x
var y: Double = y
var z: Double = 0.0
var r: Double = 0.0
var v: Double = 0.0
var w: Double = 0.0
var s: Double = 0.0
var ix: Int = 0
var hx: Int = 0
hx = __HI(x) /* high word of x */
ix = hx and 0x7fffffff /* high word of |x| */
if (ix < 0x3e300000) { /* x < 2**-28 */
if (x.toInt() == 0) { /* generate inexact */
if (((ix or __LO(x)) or (iy + 1)) == 0)
return one / fabs(x)
else {
if (iy == 1) {
return x
} else { /* compute -1 / (x+y) carefully */
var a: Double = 0.0
var t: Double = 0.0
w = x + y
z = w
z = doubleSetWord(d = z, lo = 0)
v = y - (z - x)
a = -one / w
t = a
t = doubleSetWord(d = t, lo = 0)
s = one + t * z
return t + a * (s + t * v)
}
}
}
}
if (ix >= 0x3FE59428) { /* |x| >= 0.6744 */
if (hx < 0) {
x = -x
y = -y
}
z = pio4 - x
w = pio4lo - y
x = z + w
y = 0.0
}
z = x * x
w = z * z
/*
* Break x^5*(T[1]+x^2*T[2]+...) into
* x^5(T[1]+x^4*T[3]+...+x^20*T[11]) +
* x^5(x^2*(T[2]+x^4*T[4]+...+x^22*[T12]))
*/
r = T[1] + w * (T[3] + w * (T[5] + w * (T[7] + w * (T[9] +
w * T[11]))))
v = z * (T[2] + w * (T[4] + w * (T[6] + w * (T[8] + w * (T[10] +
w * T[12])))))
s = z * x
r = y + z * (s * (r + v) + y)
r += T[0] * s
w = x + r
if (ix >= 0x3FE59428) {
v = iy.toDouble()
return (1 - ((hx shr 30) and 2)).toDouble() *
(v - 2.0 * (x - (w * w / (w + v) - r)))
}
if (iy == 1)
return w
else {
/*
* if allow error up to 2 ulp, simply return
* -1.0 / (x+r) here
*/
/* compute -1.0 / (x+r) accurately */
var a: Double = 0.0
var t: Double = 0.0
z = w
z = doubleSetWord(d = z, lo = 0)
v = r - (z - x) /* z+v = r+x */
a = -1.0 / w
t = a /* a = -1.0/w */
t = doubleSetWord(d = t, lo = 0)
s = 1.0 + t * z
return t + a * (s + t * v)
}
}
@@ -0,0 +1,53 @@
/* @(#)s_asinh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* asinh(x)
* Method :
* Based on
* asinh(x) = sign(x) * log [ |x| + sqrt(x*x+1) ]
* we have
* asinh(x) := x if 1+x*x=1,
* := sign(x)*(log(x)+ln2)) for large |x|, else
* := sign(x)*log(2|x|+1/(|x|+sqrt(x*x+1))) if|x|>2, else
* := sign(x)*log1p(|x| + x^2/(1 + sqrt(1+x^2)))
*/
package kotlin.math.fdlibm
import kotlin.wasm.internal.wasm_f64_sqrt as sqrt
private const val one = 1.00000000000000000000e+00 /* 0x3FF00000, 0x00000000 */
private const val ln2 = 6.93147180559945286227e-01 /* 0x3FE62E42, 0xFEFA39EF */
private const val huge = 1.00000000000000000000e+300
internal fun asinh(x: Double): Double {
var t: Double = 0.0
var w: Double = 0.0
var hx: Int = 0
var ix: Int = 0
hx = __HI(x)
ix = hx and 0x7fffffff
if (ix >= 0x7ff00000) return x + x /* x is inf or NaN */
if (ix < 0x3e300000) { /* |x|<2**-28 */
if (huge + x > one) return x /* return x inexact except 0 */
}
if (ix > 0x41b00000) { /* |x| > 2**28 */
w = __ieee754_log(fabs(x)) + ln2
} else if (ix > 0x40000000) { /* 2**28 > |x| > 2.0 */
t = fabs(x)
w = __ieee754_log(2.0 * t + one / (sqrt(x * x + one) + t))
} else { /* 2.0 > |x| > 2**-28 */
t = x * x
w = log1p(fabs(x) + t / (one + sqrt(one + t)))
}
if (hx > 0) return w; else return -w
}
@@ -0,0 +1,120 @@
/* @(#)s_atan.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
*/
/* atan(x)
* Method
* 1. Reduce x to positive by atan(x) = -atan(-x).
* 2. According to the integer k=4t+0.25 chopped, t=x, the argument
* is further reduced to one of the following intervals and the
* arctangent of t is evaluated by the corresponding formula:
*
* [0,7/16] atan(x) = t-t^3*(a1+t^2*(a2+...(a10+t^2*a11)...)
* [7/16,11/16] atan(x) = atan(1/2) + atan( (t-0.5)/(1+t/2) )
* [11/16.19/16] atan(x) = atan( 1 ) + atan( (t-1)/(1+t) )
* [19/16,39/16] atan(x) = atan(3/2) + atan( (t-1.5)/(1+1.5t) )
* [39/16,INF] atan(x) = atan(INF) + atan( -1/t )
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
package kotlin.math.fdlibm
private val atanhi = doubleArrayOf(
4.63647609000806093515e-01, /* atan(0.5)hi 0x3FDDAC67, 0x0561BB4F */
7.85398163397448278999e-01, /* atan(1.0)hi 0x3FE921FB, 0x54442D18 */
9.82793723247329054082e-01, /* atan(1.5)hi 0x3FEF730B, 0xD281F69B */
1.57079632679489655800e+00, /* atan(inf)hi 0x3FF921FB, 0x54442D18 */
)
private val atanlo = doubleArrayOf(
2.26987774529616870924e-17, /* atan(0.5)lo 0x3C7A2B7F, 0x222F65E2 */
3.06161699786838301793e-17, /* atan(1.0)lo 0x3C81A626, 0x33145C07 */
1.39033110312309984516e-17, /* atan(1.5)lo 0x3C700788, 0x7AF0CBBD */
6.12323399573676603587e-17, /* atan(inf)lo 0x3C91A626, 0x33145C07 */
)
private val aT = doubleArrayOf(
3.33333333333329318027e-01, /* 0x3FD55555, 0x5555550D */
-1.99999999998764832476e-01, /* 0xBFC99999, 0x9998EBC4 */
1.42857142725034663711e-01, /* 0x3FC24924, 0x920083FF */
-1.11111104054623557880e-01, /* 0xBFBC71C6, 0xFE231671 */
9.09088713343650656196e-02, /* 0x3FB745CD, 0xC54C206E */
-7.69187620504482999495e-02, /* 0xBFB3B0F2, 0xAF749A6D */
6.66107313738753120669e-02, /* 0x3FB10D66, 0xA0D03D51 */
-5.83357013379057348645e-02, /* 0xBFADDE2D, 0x52DEFD9A */
4.97687799461593236017e-02, /* 0x3FA97B4B, 0x24760DEB */
-3.65315727442169155270e-02, /* 0xBFA2B444, 0x2C6A6C2F */
1.62858201153657823623e-02, /* 0x3F90AD3A, 0xE322DA11 */
)
private const val one = 1.0
private const val huge = 1.0e300
internal fun atan(x: Double): Double {
var x: Double = x
var w: Double = 0.0
var s1: Double = 0.0
var s2: Double = 0.0
var z: Double = 0.0
var ix: Int = 0
var hx: Int = 0
var id: Int = 0
hx = __HI(x)
ix = hx and 0x7fffffff
if (ix >= 0x44100000) { /* if |x| >= 2^66 */
if (ix > 0x7ff00000 ||
(ix == 0x7ff00000 && (__LO(x) != 0))
)
return x + x /* NaN */
if (hx > 0) return atanhi[3] + atanlo[3]
else return -atanhi[3] - atanlo[3]
}
if (ix < 0x3fdc0000) { /* |x| < 0.4375 */
if (ix < 0x3e200000) { /* |x| < 2^-29 */
if (huge + x > one) return x /* raise inexact */
}
id = -1
} else {
x = fabs(x)
if (ix < 0x3ff30000) { /* |x| < 1.1875 */
if (ix < 0x3fe60000) { /* 7/16 <=|x|<11/16 */
id = 0; x = (2.0 * x - one) / (2.0 + x)
} else { /* 11/16<=|x|< 19/16 */
id = 1; x = (x - one) / (x + one)
}
} else {
if (ix < 0x40038000) { /* |x| < 2.4375 */
id = 2; x = (x - 1.5) / (one + 1.5 * x)
} else { /* 2.4375 <= |x| < 2^66 */
id = 3; x = -1.0 / x
}
}
}
/* end of argument reduction */
z = x * x
w = z * z
/* break sum from i=0 to 10 aT[i]z**(i+1) into odd and even poly */
s1 = z * (aT[0] + w * (aT[2] + w * (aT[4] + w * (aT[6] + w * (aT[8] + w * aT[10])))))
s2 = w * (aT[1] + w * (aT[3] + w * (aT[5] + w * (aT[7] + w * aT[9]))))
if (id < 0) return x - x * (s1 + s2)
else {
z = atanhi[id] - ((x * (s1 + s2) - atanlo[id]) - x)
return if (hx < 0) -z else z
}
}
@@ -0,0 +1,72 @@
/* @(#)s_cos.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* cos(x)
* Return cosine function of x.
*
* kernel function:
* __kernel_sin ... sine function on [-pi/4,pi/4]
* __kernel_cos ... cosine function on [-pi/4,pi/4]
* __ieee754_rem_pio2 ... argument reduction routine
*
* Method.
* Let S,C and T denote the sin, cos and tan respectively on
* [-PI/4, +PI/4]. Reduce the argument x to y1+y2 = x-k*pi/2
* in [-pi/4 , +pi/4], and let n = k mod 4.
* We have
*
* n sin(x) cos(x) tan(x)
* ----------------------------------------------------------
* 0 S C T
* 1 C -S -1/T
* 2 -S -C T
* 3 -C S -1/T
* ----------------------------------------------------------
*
* Special cases:
* Let trig be any of sin, cos, or tan.
* trig(+-INF) is NaN, with signals;
* trig(NaN) is that NaN;
*
* Accuracy:
* TRIG(x) returns trig(x) nearly rounded
*/
package kotlin.math.fdlibm
internal fun cos(x: Double): Double {
val y: DoubleArray = DoubleArray(2)
var z: Double = 0.0
var n: Int = 0
var ix: Int = 0
/* High word of x. */
ix = __HI(x)
/* |x| ~< pi/4 */
ix = ix and 0x7fffffff
if (ix <= 0x3fe921fb) return __kernel_cos(x, z)
/* cos(Inf or NaN) is NaN */
else if (ix >= 0x7ff00000) return x - x
/* argument reduction needed */
else {
n = __ieee754_rem_pio2(x, y)
when (n and 3) {
0 -> return __kernel_cos(y[0], y[1])
1 -> return -__kernel_sin(y[0], y[1], 1)
2 -> return -__kernel_cos(y[0], y[1])
else -> return __kernel_sin(y[0], y[1], 1)
}
}
}
@@ -0,0 +1,214 @@
/* @(#)s_expm1.c 1.5 04/04/22 */
/*
* ====================================================
* Copyright (C) 2004 by Sun Microsystems, Inc. All rights reserved.
*
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* expm1(x)
* Returns exp(x)-1, the exponential of x minus 1.
*
* Method
* 1. Argument reduction:
* Given x, find r and integer k such that
*
* x = k*ln2 + r, |r| <= 0.5*ln2 ~ 0.34658
*
* Here a correction term c will be computed to compensate
* the error in r when rounded to a floating-point number.
*
* 2. Approximating expm1(r) by a special rational function on
* the interval [0,0.34658]:
* Since
* r*(exp(r)+1)/(exp(r)-1) = 2+ r^2/6 - r^4/360 + ...
* we define R1(r*r) by
* r*(exp(r)+1)/(exp(r)-1) = 2+ r^2/6 * R1(r*r)
* That is,
* R1(r**2) = 6/r *((exp(r)+1)/(exp(r)-1) - 2/r)
* = 6/r * ( 1 + 2.0*(1/(exp(r)-1) - 1/r))
* = 1 - r^2/60 + r^4/2520 - r^6/100800 + ...
* We use a special Remes algorithm on [0,0.347] to generate
* a polynomial of degree 5 in r*r to approximate R1. The
* maximum error of this polynomial approximation is bounded
* by 2**-61. In other words,
* R1(z) ~ 1.0 + Q1*z + Q2*z**2 + Q3*z**3 + Q4*z**4 + Q5*z**5
* where Q1 = -1.6666666666666567384E-2,
* Q2 = 3.9682539681370365873E-4,
* Q3 = -9.9206344733435987357E-6,
* Q4 = 2.5051361420808517002E-7,
* Q5 = -6.2843505682382617102E-9;
* (where z=r*r, and the values of Q1 to Q5 are listed below)
* with error bounded by
* | 5 | -61
* | 1.0+Q1*z+...+Q5*z - R1(z) | <= 2
* | |
*
* expm1(r) = exp(r)-1 is then computed by the following
* specific way which minimize the accumulation rounding error:
* 2 3
* r r [ 3 - (R1 + R1*r/2) ]
* expm1(r) = r + --- + --- * [--------------------]
* 2 2 [ 6 - r*(3 - R1*r/2) ]
*
* To compensate the error in the argument reduction, we use
* expm1(r+c) = expm1(r) + c + expm1(r)*c
* ~ expm1(r) + c + r*c
* Thus c+r*c will be added in as the correction terms for
* expm1(r+c). Now rearrange the term to avoid optimization
* screw up:
* ( 2 2 )
* ({ ( r [ R1 - (3 - R1*r/2) ] ) } r )
* expm1(r+c)~r - ({r*(--- * [--------------------]-c)-c} - --- )
* ({ ( 2 [ 6 - r*(3 - R1*r/2) ] ) } 2 )
* ( )
*
* = r - E
* 3. Scale back to obtain expm1(x):
* From step 1, we have
* expm1(x) = either 2^k*[expm1(r)+1] - 1
* = or 2^k*[expm1(r) + (1-2^-k)]
* 4. Implementation notes:
* (A). To save one multiplication, we scale the coefficient Qi
* to Qi*2^i, and replace z by (x^2)/2.
* (B). To achieve maximum accuracy, we compute expm1(x) by
* (i) if x < -56*ln2, return -1.0, (raise inexact if x!=inf)
* (ii) if k=0, return r-E
* (iii) if k=-1, return 0.5*(r-E)-0.5
* (iv) if k=1 if r < -0.25, return 2*((r+0.5)- E)
* else return 1.0+2.0*(r-E);
* (v) if (k<-2||k>56) return 2^k(1-(E-r)) - 1 (or exp(x)-1)
* (vi) if k <= 20, return 2^k((1-2^-k)-(E-r)), else
* (vii) return 2^k(1-((E+2^-k)-r))
*
* Special cases:
* expm1(INF) is INF, expm1(NaN) is NaN;
* expm1(-INF) is -1, and
* for finite argument, only expm1(0)=0 is exact.
*
* Accuracy:
* according to an error analysis, the error is always less than
* 1 ulp (unit in the last place).
*
* Misc. info.
* For IEEE double
* if x > 7.09782712893383973096e+02 then expm1(x) overflow
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*/
package kotlin.math.fdlibm
private const val one = 1.0
private const val huge = 1.0e+300
private const val tiny = 1.0e-300
private const val o_threshold = 7.09782712893383973096e+02 /* 0x40862E42, 0xFEFA39EF */
private const val ln2_hi = 6.93147180369123816490e-01 /* 0x3fe62e42, 0xfee00000 */
private const val ln2_lo = 1.90821492927058770002e-10 /* 0x3dea39ef, 0x35793c76 */
private const val invln2 = 1.44269504088896338700e+00 /* 0x3ff71547, 0x652b82fe */
/* scaled coefficients related to expm1 */
private const val Q1 = -3.33333333333331316428e-02 /* BFA11111 111110F4 */
private const val Q2 = 1.58730158725481460165e-03 /* 3F5A01A0 19FE5585 */
private const val Q3 = -7.93650757867487942473e-05 /* BF14CE19 9EAADBB7 */
private const val Q4 = 4.00821782732936239552e-06 /* 3ED0CFCA 86E65239 */
private const val Q5 = -2.01099218183624371326e-07 /* BE8AFDB7 6E09C32D */
internal fun expm1(x: Double): Double {
var x: Double = x
var y: Double = 0.0
var hi: Double = 0.0
var lo: Double = 0.0
var c: Double = 0.0
var t: Double = 0.0
var e: Double = 0.0
var hxs: Double = 0.0
var hfx: Double = 0.0
var r1: Double = 0.0
var k: Int = 0
var xsb: Int = 0
var hx: UInt = 0U
hx = __HIu(x) /* high word of x */
xsb = (hx and Int.MIN_VALUE.toUInt()).toInt() /* sign bit of x */
if (xsb == 0) y = x; else y = -x /* y = |x| */
hx = (hx and 0x7fffffffU) /* high word of |x| */
/* filter out huge and non-finite argument */
if (hx >= 0x4043687AU) { /* if |x|>=56*ln2 */
if (hx >= 0x40862E42U) { /* if |x|>=709.78... */
if (hx >= 0x7ff00000U) {
if (((hx and 0xfffffU) or __LOu(x)) != 0U)
return x + x /* NaN */
else return if (xsb == 0) x else -1.0/* exp(+-inf)={inf,-1} */
}
if (x > o_threshold) return huge * huge /* overflow */
}
if (xsb != 0) { /* x < -56*ln2, return -1.0 with inexact */
if (x + tiny < 0.0) /* raise inexact */
return tiny - one /* return -1 */
}
}
/* argument reduction */
if (hx > 0x3fd62e42U) { /* if |x| > 0.5 ln2 */
if (hx < 0x3FF0A2B2U) { /* and |x| < 1.5 ln2 */
if (xsb == 0) {
hi = x - ln2_hi; lo = ln2_lo; k = 1
} else {
hi = x + ln2_hi; lo = -ln2_lo; k = -1
}
} else {
k = (invln2 * x + (if (xsb == 0) 0.5 else -0.5)).toInt()
t = k.toDouble()
hi = x - t * ln2_hi /* t*ln2_hi is exact here */
lo = t * ln2_lo
}
x = hi - lo
c = (hi - x) - lo
} else if (hx < 0x3c900000U) { /* when |x|<2**-54, return x */
t = huge + x /* return x with inexact flags when x!=0 */
return x - (t - (huge + x))
} else k = 0
/* x is now in primary range */
hfx = 0.5 * x
hxs = x * hfx
r1 = one + hxs * (Q1 + hxs * (Q2 + hxs * (Q3 + hxs * (Q4 + hxs * Q5))))
t = 3.0 - r1 * hfx
e = hxs * ((r1 - t) / (6.0 - x * t))
if (k == 0) return x - (x * e - hxs) /* c is 0 */
else {
e = (x * (e - c) - c)
e -= hxs
if (k == -1) return 0.5 * (x - e) - 0.5
if (k == 1)
if (x < -0.25) return -2.0 * (e - (x + 0.5))
else return one + 2.0 * (x - e)
if (k <= -2 || k > 56) { /* suffice to return exp(x)-1 */
y = one - (e - x)
y = doubleSetWord(d = y, hi = __HI(y) + (k shl 20)) /* add k to y's exponent */
return y - one
}
t = one
if (k < 20) {
t = doubleSetWord(d = t, hi = 0x3ff00000 - (0x200000 shr k)) /* t=1-2^-k */
y = t - (e - x)
y = doubleSetWord(d = y, hi = __HI(y) + (k shl 20)) /* add k to y's exponent */
} else {
t = doubleSetWord(d = t, hi = ((0x3ff - k) shl 20)) /* 2^-k */
y = x - (e + t)
y += one
y = doubleSetWord(d = y, hi = __HI(y) + (k shl 20)) /* add k to y's exponent */
}
}
return y
}
@@ -0,0 +1,21 @@
/* @(#)s_fabs.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/*
* fabs(x) returns the absolute value of x.
*/
package kotlin.math.fdlibm
internal fun fabs(x: Double): Double {
return doubleSetWord(d = x, hi = __HI(x) and 0x7fffffff) /* add k to y's exponent */
}
@@ -0,0 +1,52 @@
/* @(#)s_ilogb.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* ilogb(x: Double)
* return the binary exponent of non-zero x
* ilogb(0) = 0x80000001
* ilogb(inf/NaN) = 0x7fffffff (no signal is raised)
*/
package kotlin.math.fdlibm
internal fun ilogb(x: Double): Int {
var hx: Int = 0
var lx: Int = 0
var ix: Int = 0;
hx = (__HI(x)) and 0x7fffffff /* high word of x */
if (hx < 0x00100000) {
lx = __LO(x)
if ((hx or lx) == 0)
return 0x80000001.toInt() /* ilogb(0) = 0x80000001 */
else /* subnormal x */
if (hx == 0) {
//for (ix = -1043; lx>0; lx<<=1) {
ix = -1043
while (lx > 0) {
ix -= 1
//--
lx = lx shl 1
}
} else {
//for (ix = -1022,hx<<=11; hx>0; hx<<=1) {
ix = -1022; hx = hx shl 11
while (hx > 0) {
ix -= 1
//--
hx = hx shl 1
}
}
return ix
} else if (hx < 0x7ff00000) return (hx shr 20) - 1023
else return 0x7fffffff
}
@@ -0,0 +1,167 @@
/* @(#)s_log1p.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* double log1p(x: Double)
*
* Method :
* 1. Argument Reduction: find k and f such that
* 1+x = 2^k * (1+f),
* where sqrt(2)/2 < 1+f < sqrt(2) .
*
* Note. If k=0, then f=x is exact. However, if k!=0, then f
* may not be representable exactly. In that case, a correction
* term is need. Let u=1+x rounded. Let c = (1+x)-u, then
* log(1+x) - log(u) ~ c/u. Thus, we proceed to compute log(u),
* and add back the correction term c/u.
* (Note: when x > 2**53, one can simply return log(x))
*
* 2. Approximation of log1p(f).
* Let s = f/(2+f) ; based on log(1+f) = log(1+s) - log(1-s)
* = 2s + 2/3 s**3 + 2/5 s**5 + .....,
* = 2s + s*R
* We use a special Reme algorithm on [0,0.1716] to generate
* a polynomial of degree 14 to approximate R The maximum error
* of this polynomial approximation is bounded by 2**-58.45. In
* other words,
* 2 4 6 8 10 12 14
* R(z) ~ Lp1*s +Lp2*s +Lp3*s +Lp4*s +Lp5*s +Lp6*s +Lp7*s
* (the values of Lp1 to Lp7 are listed in the program)
* and
* | 2 14 | -58.45
* | Lp1*s +...+Lp7*s - R(z) | <= 2
* | |
* Note that 2s = f - s*f = f - hfsq + s*hfsq, where hfsq = f*f/2.
* In order to guarantee error in log below 1ulp, we compute log
* by
* log1p(f) = f - (hfsq - s*(hfsq+R)).
*
* 3. Finally, log1p(x) = k*ln2 + log1p(f).
* = k*ln2_hi+(f-(hfsq-(s*(hfsq+R)+k*ln2_lo)))
* Here ln2 is split into two floating point number:
* ln2_hi + ln2_lo,
* where n*ln2_hi is always exact for |n| < 2000.
*
* Special cases:
* log1p(x) is NaN with signal if x < -1 (including -INF) ;
* log1p(+INF) is +INF; log1p(-1) is -INF with signal;
* log1p(NaN) is that NaN with no signal.
*
* Accuracy:
* according to an error analysis, the error is always less than
* 1 ulp (unit in the last place).
*
* Constants:
* The hexadecimal values are the intended ones for the following
* constants. The decimal values may be used, provided that the
* compiler will convert from decimal to binary accurately enough
* to produce the hexadecimal values shown.
*
* Note: Assuming log() return accurate answer, the following
* algorithm can be used to compute log1p(x) to within a few ULP:
*
* u = 1+x;
* if(u==1.0) return x ; else
* return log(u)*(x/(u-1.0));
*
* See HP-15C Advanced Functions Handbook, p.193.
*/
package kotlin.math.fdlibm
private const val ln2_hi = 6.93147180369123816490e-01 /* 3fe62e42 fee00000 */
private const val ln2_lo = 1.90821492927058770002e-10 /* 3dea39ef 35793c76 */
private const val two54 = 1.80143985094819840000e+16 /* 43500000 00000000 */
private const val Lp1 = 6.666666666666735130e-01 /* 3FE55555 55555593 */
private const val Lp2 = 3.999999999940941908e-01 /* 3FD99999 9997FA04 */
private const val Lp3 = 2.857142874366239149e-01 /* 3FD24924 94229359 */
private const val Lp4 = 2.222219843214978396e-01 /* 3FCC71C5 1D8E78AF */
private const val Lp5 = 1.818357216161805012e-01 /* 3FC74664 96CB03DE */
private const val Lp6 = 1.531383769920937332e-01 /* 3FC39A09 D078C69F */
private const val Lp7 = 1.479819860511658591e-01 /* 3FC2F112 DF3E5244 */
private const val zero = 0.0
internal fun log1p(x: Double): Double {
var hfsq: Double = 0.0
var f: Double = 0.0
var c: Double = 0.0
var s: Double = 0.0
var z: Double = 0.0
var R: Double = 0.0
var u: Double = 0.0
var k: Int = 0
var hx: Int = 0
var hu: Int = 0
var ax: Int = 0
hx = __HI(x) /* high word of x */
ax = hx and 0x7fffffff
k = 1
if (hx < 0x3FDA827A) { /* x < 0.41422 */
if (ax >= 0x3ff00000) { /* x <= -1.0 */
if (x == -1.0) return -two54 / zero /* log1p(-1)=+inf */
else return (x - x) / (x - x) /* log1p(x<-1)=NaN */
}
if (ax < 0x3e200000) { /* |x| < 2**-29 */
if (two54 + x > zero /* raise inexact */
&& ax < 0x3c900000
) /* |x| < 2**-54 */
return x
else
return x - x * x * 0.5
}
if (hx > 0 || hx <= (0xbfd2bec3.toInt())) {
k = 0;f = x;hu = 1
} /* -0.2929<x<0.41422 */
}
if (hx >= 0x7ff00000) return x + x
if (k != 0) {
if (hx < 0x43400000) {
u = 1.0 + x
hu = __HI(u) /* high word of u */
k = (hu shr 20) - 1023
c = if (k > 0) 1.0 - (u - x) else x - (u - 1.0)/* correction term */
c /= u
} else {
u = x
hu = __HI(u) /* high word of u */
k = (hu shr 20) - 1023
c = 0.0
}
hu = hu and 0x000fffff
if (hu < 0x6a09e) {
u = doubleSetWord(d = u, hi = hu or 0x3ff00000) /* normalize u */
} else {
k += 1
u = doubleSetWord(d = u, hi = hu or 0x3fe00000) /* normalize u/2 */
hu = (0x00100000 - hu) shr 2
}
f = u - 1.0
}
hfsq = 0.5 * f * f
if (hu == 0) { /* |f| < 2**-20 */
if (f == zero) if (k == 0) return zero
else {
c += k * ln2_lo; return k * ln2_hi + c
}
R = hfsq * (1.0 - 0.66666666666666666 * f)
if (k == 0) return f - R; else
return k * ln2_hi - ((R - (k * ln2_lo + c)) - f)
}
s = f / (2.0 + f)
z = s * s
R = z * (Lp1 + z * (Lp2 + z * (Lp3 + z * (Lp4 + z * (Lp5 + z * (Lp6 + z * Lp7))))))
if (k == 0) return f - (hfsq - s * (hfsq + R)); else
return k * ln2_hi - ((hfsq - (s * (hfsq + R) + (k * ln2_lo + c))) - f)
}
@@ -0,0 +1,78 @@
/* @(#)s_nextafter.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* IEEE functions
* nextafter(x,y)
* return the next machine floating-point number of x in the
* direction toward y.
* Special cases:
*/
package kotlin.math.fdlibm
internal fun nextafter(x: Double, y: Double): Double {
var x: Double = x
var y: Double = y
var hx: Int = 0
var hy: Int = 0
var ix: Int = 0
var iy: Int = 0
var lx: UInt = 0U
var ly: UInt = 0U
hx = __HI(x) /* high word of x */
lx = __LOu(x) /* low word of x */
hy = __HI(y) /* high word of y */
ly = __LOu(y) /* low word of y */
ix = hx and 0x7fffffff /* |x| */
iy = hy and 0x7fffffff /* |y| */
if (((ix >= 0x7ff00000) && ((ix - 0x7ff00000) or lx.toInt()) != 0) || /* x is nan */
((iy >= 0x7ff00000) && ((iy - 0x7ff00000) or ly.toInt()) != 0)
) /* y is nan */
return x + y
if (x == y) return x /* x=y, return x */
if ((ix or lx.toInt()) == 0) { /* x == 0 */
x = doubleSetWord(d = x, hi = hy and Int.MIN_VALUE) /* return +-minsubnormal */
x = doubleSetWord(d = x, lo = 1)
y = x * x
if (y == x) return y; else return x /* raise underflow flag */
}
if (hx >= 0) { /* x > 0 */
if (hx > hy || ((hx == hy) && (lx > ly))) { /* x > y, x -= ulp */
if (lx == 0U) hx -= 1
lx -= 1U
} else { /* x < y, x += ulp */
lx += 1U
if (lx == 0U) hx += 1
}
} else { /* x < 0 */
if (hy >= 0 || hx > hy || ((hx == hy) && (lx > ly))) {/* x < y, x -= ulp */
if (lx == 0U) hx -= 1
lx -= 1U
} else { /* x > y, x += ulp */
lx += 1U
if (lx == 0U) hx += 1
}
}
hy = hx and 0x7ff00000
if (hy >= 0x7ff00000) return x + x /* overflow */
if (hy < 0x00100000) { /* underflow */
y = x * x
if (y != x) { /* raise underflow flag */
y = doubleSetWord(d = y, hi = hx, lo = lx.toInt())
return y
}
}
x = doubleSetWord(d = x, hi = hx, lo = lx.toInt())
return x
}
@@ -0,0 +1,77 @@
/* @(#)s_rint.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/*
* rint(x)
* Return x rounded to integral value according to the prevailing
* rounding mode.
* Method:
* Using floating addition.
* Exception:
* Inexact flag raised if x not equal to rint(x).
*/
package kotlin.math.fdlibm
private val TWO52 = doubleArrayOf(
4.50359962737049600000e+15, /* 0x43300000, 0x00000000 */
-4.50359962737049600000e+15, /* 0xC3300000, 0x00000000 */
)
internal fun rint(x: Double): Double {
var x: Double = x
var i0: Int = 0
var j0: Int = 0
var sx: Int = 0
var i: UInt = 0U
var i1: UInt = 0U
var w: Double = 0.0
var t: Double = 0.0
i0 = __HI(x)
sx = (i0 shr 31) and 1
i1 = __LOu(x)
j0 = ((i0 shr 20) and 0x7ff) - 0x3ff
if (j0 < 20) {
if (j0 < 0) {
if (((i0 and 0x7fffffff) or i1.toInt()) == 0) return x
i1 = i1 or (i0 and 0x0fffff).toUInt()
i0 = i0 and 0xfffe0000.toInt()
i0 = i0 or (((i1 or i1.negate()) shr 12) and 0x80000.toUInt()).toInt();
x = doubleSetWord(d = x, hi = i0)
w = TWO52[sx] + x
t = w - TWO52[sx]
i0 = __HI(t)
t = doubleSetWord(d = t, hi = (i0 and 0x7fffffff) or (sx shl 31))
return t
} else {
i = ((0x000fffff) shr j0).toUInt()
if (((i0 and i.toInt()) or i1.toInt()) == 0) return x /* x is integral */
i = i shr 1
if (((i0 and i.toInt()) or i1.toInt()) != 0) {
if (j0 == 19) i1 = 0x40000000.toUInt(); else
i0 = (i0 and i.inv().toInt()) or ((0x20000) shr j0)
}
}
} else if (j0 > 51) {
if (j0 == 0x400) return x + x /* inf or NaN */
else return x /* x is integral */
} else {
i = ((0xffffffff.toUInt())) shr (j0 - 20)
if ((i1 and i) == 0U) return x /* x is integral */
i = i shr 1
if ((i1 and i) != 0U) i1 = (i1 and (i.inv())) or ((0x40000000) shr (j0 - 20)).toUInt()
}
x = doubleSetWord(x, hi = i0, lo = i1.toInt())
w = TWO52[sx] + x
return w - TWO52[sx]
}
@@ -0,0 +1,57 @@
/* @(#)s_scalbn.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/*
* scalbn (double x, int n)
* scalbn(x,n) returns x* 2**n computed by exponent
* manipulation rather than by actually performing an
* exponentiation or a multiplication.
*/
package kotlin.math.fdlibm
import kotlin.wasm.internal.wasm_f64_copysign as copysign
private const val two54 = 1.80143985094819840000e+16 /* 0x43500000, 0x00000000 */
private const val twom54 = 5.55111512312578270212e-17 /* 0x3C900000, 0x00000000 */
private const val huge = 1.0e+300
private const val tiny = 1.0e-300
internal fun scalbn(x: Double, n: Int): Double {
var x: Double = x
var k: Int = 0
var hx: Int = 0
var lx: Int = 0
hx = __HI(x)
lx = __LO(x)
k = ((hx and 0x7ff00000) shr 20) /* extract exponent */
if (k == 0) { /* 0 or subnormal x */
if ((lx or (hx and 0x7fffffff)) == 0) return x /* +-0 */
x *= two54
hx = __HI(x)
k = (((hx and 0x7ff00000) shr 20) - 54)
if (n < -50000) return tiny * x /*underflow*/
}
if (k == 0x7ff) return x + x /* NaN or Inf */
k = k + n
if (k > 0x7fe) return huge * copysign(huge, x) /* overflow */
if (k > 0) /* normal result */ {
x = doubleSetWord(d = x, hi = (hx and 0x800fffff.toInt()) or (k shl 20)); return x
}
if (k <= -54)
if (n > 50000) /* in case integer overflow in n+k */
return huge * copysign(huge, x) /*overflow*/
else return tiny * copysign(tiny, x) /*underflow*/
k += 54 /* subnormal result */
x = doubleSetWord(d = x, hi = (hx and 0x800fffff.toInt()) or (k shl 20))
return x * twom54
}
@@ -0,0 +1,72 @@
/* @(#)s_sin.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* sin(x)
* Return sine function of x.
*
* kernel function:
* __kernel_sin ... sine function on [-pi/4,pi/4]
* __kernel_cos ... cose function on [-pi/4,pi/4]
* __ieee754_rem_pio2 ... argument reduction routine
*
* Method.
* Let S,C and T denote the sin, cos and tan respectively on
* [-PI/4, +PI/4]. Reduce the argument x to y1+y2 = x-k*pi/2
* in [-pi/4 , +pi/4], and let n = k mod 4.
* We have
*
* n sin(x) cos(x) tan(x)
* ----------------------------------------------------------
* 0 S C T
* 1 C -S -1/T
* 2 -S -C T
* 3 -C S -1/T
* ----------------------------------------------------------
*
* Special cases:
* Let trig be any of sin, cos, or tan.
* trig(+-INF) is NaN, with signals;
* trig(NaN) is that NaN;
*
* Accuracy:
* TRIG(x) returns trig(x) nearly rounded
*/
package kotlin.math.fdlibm
internal fun sin(x: Double): Double {
val y: DoubleArray = DoubleArray(2)
var z: Double = 0.0
var n: Int = 0
var ix: Int = 0
/* High word of x. */
ix = __HI(x)
/* |x| ~< pi/4 */
ix = ix and 0x7fffffff
if (ix <= 0x3fe921fb) return __kernel_sin(x, z, 0)
/* sin(Inf or NaN) is NaN */
else if (ix >= 0x7ff00000) return x - x
/* argument reduction needed */
else {
n = __ieee754_rem_pio2(x, y)
when (n and 3) {
0 -> return __kernel_sin(y[0], y[1], 1)
1 -> return __kernel_cos(y[0], y[1])
2 -> return -__kernel_sin(y[0], y[1], 1)
else -> return -__kernel_cos(y[0], y[1])
}
}
}
@@ -0,0 +1,67 @@
/* @(#)s_tan.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* tan(x)
* Return tangent function of x.
*
* kernel function:
* __kernel_tan ... tangent function on [-pi/4,pi/4]
* __ieee754_rem_pio2 ... argument reduction routine
*
* Method.
* Let S,C and T denote the sin, cos and tan respectively on
* [-PI/4, +PI/4]. Reduce the argument x to y1+y2 = x-k*pi/2
* in [-pi/4 , +pi/4], and let n = k mod 4.
* We have
*
* n sin(x) cos(x) tan(x)
* ----------------------------------------------------------
* 0 S C T
* 1 C -S -1/T
* 2 -S -C T
* 3 -C S -1/T
* ----------------------------------------------------------
*
* Special cases:
* Let trig be any of sin, cos, or tan.
* trig(+-INF) is NaN, with signals;
* trig(NaN) is that NaN;
*
* Accuracy:
* TRIG(x) returns trig(x) nearly rounded
*/
package kotlin.math.fdlibm
internal fun tan(x: Double): Double {
val y: DoubleArray = DoubleArray(2)
var z: Double = 0.0
var n: Int = 0
var ix: Int = 0
/* High word of x. */
ix = __HI(x)
/* |x| ~< pi/4 */
ix = ix and 0x7fffffff
if (ix <= 0x3fe921fb) return __kernel_tan(x, z, 1)
/* tan(Inf or NaN) is NaN */
else if (ix >= 0x7ff00000) return x - x /* NaN */
/* argument reduction needed */
else {
n = __ieee754_rem_pio2(x, y)
return __kernel_tan(y[0], y[1], 1 - ((n and 1) shl 1)) /* 1 -- n even
-1 -- n odd */
}
}
@@ -0,0 +1,75 @@
/* @(#)s_tanh.c 1.3 95/01/18 */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunSoft, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/* Tanh(x)
* Return the Hyperbolic Tangent of x
*
* Method :
* x -x
* e - e
* 0. tanh(x) is defined to be -----------
* x -x
* e + e
* 1. reduce x to non-negative by tanh(-x) = -tanh(x).
* 2. 0 <= x <= 2**-55 : tanh(x) := x*(one+x)
* -t
* 2**-55 < x <= 1 : tanh(x) := -----; t = expm1(-2x)
* t + 2
* 2
* 1 <= x <= 22.0 : tanh(x) := 1- ----- ; t=expm1(2x)
* t + 2
* 22.0 < x <= INF : tanh(x) := 1.
*
* Special cases:
* tanh(NaN) is NaN;
* only tanh(0)=0 is exact for finite argument.
*/
package kotlin.math.fdlibm
private const val one = 1.0
private const val two = 2.0
private const val tiny = 1.0e-300
internal fun tanh(x: Double): Double {
var t: Double = 0.0
var z: Double = 0.0
var jx: Int = 0
var ix: Int = 0
/* High word of |x|. */
jx = __HI(x)
ix = jx and 0x7fffffff
/* x is INF or NaN */
if (ix >= 0x7ff00000) {
if (jx >= 0) return one / x + one /* tanh(+-inf)=+-1 */
else return one / x - one /* tanh(NaN) = NaN */
}
/* |x| < 22 */
if (ix < 0x40360000) { /* |x|<22 */
if (ix < 0x3c800000) /* |x|<2**-55 */
return x * (one + x) /* tanh(small) = small */
if (ix >= 0x3ff00000) { /* |x|>=1 */
t = expm1(two * fabs(x))
z = one - two / (t + two)
} else {
t = expm1(-two * fabs(x))
z = -t / (t + two)
}
/* |x| > 22, return +-1 */
} else {
z = one - tiny /* raised inexact flag */
}
return if (jx >= 0) z else -z
}
@@ -0,0 +1,19 @@
/*
* Copyright 2010-2020 JetBrains s.r.o. and Kotlin Programming Language contributors.
* Use of this source code is governed by the Apache 2.0 license that can be found in the license/LICENSE.txt file.
*/
package kotlin.math.fdlibm
//#define __HI(x) *(1+(int*)&x)
internal fun __HI(x: Double): Int = (x.toRawBits() ushr 32).toInt()
internal fun __HIu(x: Double): UInt = (x.toRawBits() ushr 32).toUInt()
//#define __LO(x) *(int*)&x
internal fun __LO(x: Double): Int = (x.toRawBits() and 0xFFFFFFFF).toInt()
internal fun __LOu(x: Double): UInt = (x.toRawBits() and 0xFFFFFFFF).toUInt()
internal fun doubleSetWord(d: Double = 0.0, hi: Int = __HI(d), lo: Int = __LO(d)): Double =
Double.fromBits((hi.toLong() shl 32) or (lo.toLong() and 0xFFFFFFFF))
internal fun UInt.negate(): UInt = inv() + 1U
@@ -5,7 +5,6 @@
package kotlin.math
// region ================ Double Math ========================================
/** Computes the sine of the angle [x] given in radians.
@@ -14,7 +13,7 @@ package kotlin.math
* - `sin(NaN|+Inf|-Inf)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun sin(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun sin(x: Double): Double = kotlin.math.fdlibm.sin(x)
/** Computes the cosine of the angle [x] given in radians.
*
@@ -22,7 +21,7 @@ public actual fun sin(x: Double): Double = TODO("Wasm stdlib: Math")
* - `cos(NaN|+Inf|-Inf)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun cos(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun cos(x: Double): Double = kotlin.math.fdlibm.cos(x)
/** Computes the tangent of the angle [x] given in radians.
*
@@ -30,7 +29,7 @@ public actual fun cos(x: Double): Double = TODO("Wasm stdlib: Math")
* - `tan(NaN|+Inf|-Inf)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun tan(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun tan(x: Double): Double = kotlin.math.fdlibm.tan(x)
/**
* Computes the arc sine of the value [x];
@@ -40,17 +39,17 @@ public actual fun tan(x: Double): Double = TODO("Wasm stdlib: Math")
* - `asin(x)` is `NaN`, when `abs(x) > 1` or x is `NaN`
*/
@SinceKotlin("1.2")
public actual fun asin(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun asin(x: Double): Double = kotlin.math.fdlibm.__ieee754_asin(x)
/**
* Computes the arc cosine of the value [x];
* the returned value is an angle in the range from `0.0` to `PI` radians.
*
* Special cases:
* - `acos(x)` is `NaN`, when `abs(x) > 1` or x is `NaN`
* - `acos(x)` is `NaN`, when `abs(x) > 1` or x is `NaN`fasin
*/
@SinceKotlin("1.2")
public actual fun acos(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun acos(x: Double): Double = kotlin.math.fdlibm.__ieee754_acos(x)
/**
* Computes the arc tangent of the value [x];
@@ -60,7 +59,7 @@ public actual fun acos(x: Double): Double = TODO("Wasm stdlib: Math")
* - `atan(NaN)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun atan(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun atan(x: Double): Double = kotlin.math.fdlibm.atan(x)
/**
* Returns the angle `theta` of the polar coordinates `(r, theta)` that correspond
@@ -79,7 +78,7 @@ public actual fun atan(x: Double): Double = TODO("Wasm stdlib: Math")
* - `atan2(NaN, x)` and `atan2(y, NaN)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun atan2(y: Double, x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun atan2(y: Double, x: Double): Double = kotlin.math.fdlibm.__ieee754_atan2(y, x)
/**
* Computes the hyperbolic sine of the value [x].
@@ -90,7 +89,7 @@ public actual fun atan2(y: Double, x: Double): Double = TODO("Wasm stdlib: Math"
* - `sinh(-Inf)` is `-Inf`
*/
@SinceKotlin("1.2")
public actual fun sinh(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun sinh(x: Double): Double = kotlin.math.fdlibm.__ieee754_sinh(x)
/**
* Computes the hyperbolic cosine of the value [x].
@@ -100,7 +99,7 @@ public actual fun sinh(x: Double): Double = TODO("Wasm stdlib: Math")
* - `cosh(+Inf|-Inf)` is `+Inf`
*/
@SinceKotlin("1.2")
public actual fun cosh(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun cosh(x: Double): Double = kotlin.math.fdlibm.__ieee754_cosh(x)
/**
* Computes the hyperbolic tangent of the value [x].
@@ -111,7 +110,7 @@ public actual fun cosh(x: Double): Double = TODO("Wasm stdlib: Math")
* - `tanh(-Inf)` is `-1.0`
*/
@SinceKotlin("1.2")
public actual fun tanh(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun tanh(x: Double): Double = kotlin.math.fdlibm.tanh(x)
/**
* Computes the inverse hyperbolic sine of the value [x].
@@ -124,7 +123,7 @@ public actual fun tanh(x: Double): Double = TODO("Wasm stdlib: Math")
* - `asinh(-Inf)` is `-Inf`
*/
@SinceKotlin("1.2")
public actual fun asinh(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun asinh(x: Double): Double = kotlin.math.fdlibm.asinh(x)
/**
* Computes the inverse hyperbolic cosine of the value [x].
@@ -137,12 +136,12 @@ public actual fun asinh(x: Double): Double = TODO("Wasm stdlib: Math")
* - `acosh(+Inf)` is `+Inf`
*/
@SinceKotlin("1.2")
public actual fun acosh(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun acosh(x: Double): Double = kotlin.math.fdlibm.__ieee754_acosh(x)
/**
* Computes the inverse hyperbolic tangent of the value [x].
*
* The returned value is `y` such that `tanh(y) == x`.
* The returned value is `y` such that `tanh(y) == x `.
*
* Special cases:
* - `tanh(NaN)` is `NaN`
@@ -151,7 +150,7 @@ public actual fun acosh(x: Double): Double = TODO("Wasm stdlib: Math")
* - `tanh(-1.0)` is `-Inf`
*/
@SinceKotlin("1.2")
public actual fun atanh(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun atanh(x: Double): Double = kotlin.math.fdlibm.__ieee754_atanh(x)
/**
* Computes `sqrt(x^2 + y^2)` without intermediate overflow or underflow.
@@ -161,7 +160,7 @@ public actual fun atanh(x: Double): Double = TODO("Wasm stdlib: Math")
* - returns `NaN` if any of arguments is `NaN` and the other is not infinite
*/
@SinceKotlin("1.2")
public actual fun hypot(x: Double, y: Double): Double = TODO("Wasm stdlib: Math")
public actual fun hypot(x: Double, y: Double): Double = kotlin.math.fdlibm.__ieee754_hypot(x, y)
/**
* Computes the positive square root of the value [x].
@@ -170,7 +169,7 @@ public actual fun hypot(x: Double, y: Double): Double = TODO("Wasm stdlib: Math"
* - `sqrt(x)` is `NaN` when `x < 0` or `x` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun sqrt(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun sqrt(x: Double): Double = kotlin.wasm.internal.wasm_f64_sqrt(x)
/**
* Computes Euler's number `e` raised to the power of the value [x].
@@ -181,7 +180,7 @@ public actual fun sqrt(x: Double): Double = TODO("Wasm stdlib: Math")
* - `exp(-Inf)` is `0.0`
*/
@SinceKotlin("1.2")
public actual fun exp(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun exp(x: Double): Double = kotlin.math.fdlibm.__ieee754_exp(x)
/**
* Computes `exp(x) - 1`.
@@ -196,7 +195,7 @@ public actual fun exp(x: Double): Double = TODO("Wasm stdlib: Math")
* @see [exp] function.
*/
@SinceKotlin("1.2")
public actual fun expm1(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun expm1(x: Double): Double = kotlin.math.fdlibm.expm1(x)
/**
* Computes the logarithm of the value [x] to the given [base].
@@ -211,7 +210,13 @@ public actual fun expm1(x: Double): Double = TODO("Wasm stdlib: Math")
* See also logarithm functions for common fixed bases: [ln], [log10] and [log2].
*/
@SinceKotlin("1.2")
public actual fun log(x: Double, base: Double): Double = TODO("Wasm stdlib: Math")
public actual fun log(x: Double, base: Double): Double {
if (x.isNaN() || base.isNaN()) return Double.NaN
if (x < 0.0 || base <= 0.0 || base == 1.0) return Double.NaN
if (x.isInfinite() && base.isInfinite()) return Double.NaN
return kotlin.math.fdlibm.__ieee754_log(x) / kotlin.math.fdlibm.__ieee754_log(base)
}
/**
* Computes the natural logarithm (base `E`) of the value [x].
@@ -223,7 +228,7 @@ public actual fun log(x: Double, base: Double): Double = TODO("Wasm stdlib: Math
* - `ln(0.0)` is `-Inf`
*/
@SinceKotlin("1.2")
public actual fun ln(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun ln(x: Double): Double = kotlin.math.fdlibm.__ieee754_log(x)
/**
* Computes the common logarithm (base 10) of the value [x].
@@ -231,7 +236,7 @@ public actual fun ln(x: Double): Double = TODO("Wasm stdlib: Math")
* @see [ln] function for special cases.
*/
@SinceKotlin("1.2")
public actual fun log10(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun log10(x: Double): Double = kotlin.math.fdlibm.__ieee754_log10(x)
/**
* Computes the binary logarithm (base 2) of the value [x].
@@ -239,7 +244,7 @@ public actual fun log10(x: Double): Double = TODO("Wasm stdlib: Math")
* @see [ln] function for special cases.
*/
@SinceKotlin("1.2")
public actual fun log2(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun log2(x: Double): Double = kotlin.math.fdlibm.__ieee754_log2(x)
/**
* Computes `ln(x + 1)`.
@@ -256,7 +261,7 @@ public actual fun log2(x: Double): Double = TODO("Wasm stdlib: Math")
* @see [expm1] function
*/
@SinceKotlin("1.2")
public actual fun ln1p(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun ln1p(x: Double): Double = kotlin.math.fdlibm.log1p(x)
/**
* Rounds the given value [x] to an integer towards positive infinity.
@@ -267,7 +272,7 @@ public actual fun ln1p(x: Double): Double = TODO("Wasm stdlib: Math")
* - `ceil(x)` is `x` where `x` is `NaN` or `+Inf` or `-Inf` or already a mathematical integer.
*/
@SinceKotlin("1.2")
public actual fun ceil(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun ceil(x: Double): Double = kotlin.wasm.internal.wasm_f64_ceil(x)
/**
* Rounds the given value [x] to an integer towards negative infinity.
@@ -278,7 +283,7 @@ public actual fun ceil(x: Double): Double = TODO("Wasm stdlib: Math")
* - `floor(x)` is `x` where `x` is `NaN` or `+Inf` or `-Inf` or already a mathematical integer.
*/
@SinceKotlin("1.2")
public actual fun floor(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun floor(x: Double): Double = kotlin.wasm.internal.wasm_f64_floor(x)
/**
* Rounds the given value [x] to an integer towards zero.
@@ -289,7 +294,7 @@ public actual fun floor(x: Double): Double = TODO("Wasm stdlib: Math")
* - `truncate(x)` is `x` where `x` is `NaN` or `+Inf` or `-Inf` or already a mathematical integer.
*/
@SinceKotlin("1.2")
public actual fun truncate(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun truncate(x: Double): Double = kotlin.wasm.internal.wasm_f64_truncate(x)
/**
* Rounds the given value [x] towards the closest integer with ties rounded towards even integer.
@@ -298,7 +303,7 @@ public actual fun truncate(x: Double): Double = TODO("Wasm stdlib: Math")
* - `round(x)` is `x` where `x` is `NaN` or `+Inf` or `-Inf` or already a mathematical integer.
*/
@SinceKotlin("1.2")
public actual fun round(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun round(x: Double): Double = kotlin.math.fdlibm.rint(x)
/**
* Returns the absolute value of the given value [x].
@@ -309,7 +314,7 @@ public actual fun round(x: Double): Double = TODO("Wasm stdlib: Math")
* @see absoluteValue extension property for [Double]
*/
@SinceKotlin("1.2")
public actual fun abs(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun abs(x: Double): Double = kotlin.wasm.internal.wasm_f64_abs(x)
/**
* Returns the sign of the given value [x]:
@@ -321,8 +326,12 @@ public actual fun abs(x: Double): Double = TODO("Wasm stdlib: Math")
* - `sign(NaN)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun sign(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun sign(x: Double): Double = when {
x.isNaN() -> Double.NaN
x > 0.0 -> 1.0
x < 0.0 -> -1.0
else -> x
}
/**
* Returns the smaller of two values.
@@ -330,7 +339,7 @@ public actual fun sign(x: Double): Double = TODO("Wasm stdlib: Math")
* If either value is `NaN`, then the result is `NaN`.
*/
@SinceKotlin("1.2")
public actual fun min(a: Double, b: Double): Double = TODO("Wasm stdlib: Math")
public actual fun min(a: Double, b: Double): Double = kotlin.wasm.internal.wasm_f64_min(a, b)
/**
* Returns the greater of two values.
@@ -338,7 +347,7 @@ public actual fun min(a: Double, b: Double): Double = TODO("Wasm stdlib: Math")
* If either value is `NaN`, then the result is `NaN`.
*/
@SinceKotlin("1.2")
public actual fun max(a: Double, b: Double): Double = TODO("Wasm stdlib: Math")
public actual fun max(a: Double, b: Double): Double = kotlin.wasm.internal.wasm_f64_max(a, b)
// extensions
@@ -354,7 +363,7 @@ public actual fun max(a: Double, b: Double): Double = TODO("Wasm stdlib: Math")
* - `b.pow(x)` is `NaN` for `b < 0` and `x` is finite and not an integer
*/
@SinceKotlin("1.2")
public actual fun Double.pow(x: Double): Double = TODO("Wasm stdlib: Math")
public actual fun Double.pow(x: Double): Double = kotlin.math.fdlibm.__ieee754_pow(this, x)
/**
* Raises this value to the integer power [n].
@@ -362,7 +371,7 @@ public actual fun Double.pow(x: Double): Double = TODO("Wasm stdlib: Math")
* See the other overload of [pow] for details.
*/
@SinceKotlin("1.2")
public actual fun Double.pow(n: Int): Double = TODO("Wasm stdlib: Math")
public actual fun Double.pow(n: Int): Double = kotlin.math.fdlibm.__ieee754_pow(this, n.toDouble())
/**
* Returns the absolute value of this value.
@@ -373,7 +382,7 @@ public actual fun Double.pow(n: Int): Double = TODO("Wasm stdlib: Math")
* @see abs function
*/
@SinceKotlin("1.2")
public actual val Double.absoluteValue: Double get() = TODO("Wasm stdlib: Math")
public actual val Double.absoluteValue: Double get() = kotlin.wasm.internal.wasm_f64_abs(this)
/**
* Returns the sign of this value:
@@ -385,7 +394,7 @@ public actual val Double.absoluteValue: Double get() = TODO("Wasm stdlib: Math")
* - `NaN.sign` is `NaN`
*/
@SinceKotlin("1.2")
public actual val Double.sign: Double get() = TODO("Wasm stdlib: Math")
public actual val Double.sign: Double get() = sign(this)
/**
* Returns this value with the sign bit same as of the [sign] value.
@@ -393,13 +402,13 @@ public actual val Double.sign: Double get() = TODO("Wasm stdlib: Math")
* If [sign] is `NaN` the sign of the result is undefined.
*/
@SinceKotlin("1.2")
public actual fun Double.withSign(sign: Double): Double = TODO("Wasm stdlib: Math")
public actual fun Double.withSign(sign: Double): Double = kotlin.wasm.internal.wasm_f64_copysign(this, sign)
/**
* Returns this value with the sign bit same as of the [sign] value.
*/
@SinceKotlin("1.2")
public actual fun Double.withSign(sign: Int): Double = TODO("Wasm stdlib: Math")
public actual fun Double.withSign(sign: Int): Double = kotlin.wasm.internal.wasm_f64_copysign(this, sign.toDouble())
/**
* Returns the ulp (unit in the last place) of this value.
@@ -412,19 +421,32 @@ public actual fun Double.withSign(sign: Int): Double = TODO("Wasm stdlib: Math")
* - `0.0.ulp` is `Double.MIN_VALUE`
*/
@SinceKotlin("1.2")
public actual val Double.ulp: Double get() = TODO("Wasm stdlib: Math")
public actual val Double.ulp: Double get() = when {
this < 0 -> (-this).ulp
this.isNaN() || this == Double.POSITIVE_INFINITY -> this
this == Double.MAX_VALUE -> this - this.nextDown()
else -> this.nextUp() - this
}
/**
* Returns the [Double] value nearest to this value in direction of positive infinity.
*/
@SinceKotlin("1.2")
public actual fun Double.nextUp(): Double = TODO("Wasm stdlib: Math")
public actual fun Double.nextUp(): Double = when {
this.isNaN() || this == Double.POSITIVE_INFINITY -> this
this == 0.0 -> Double.MIN_VALUE
else -> Double.fromBits(this.toRawBits() + if (this > 0) 1 else -1)
}
/**
* Returns the [Double] value nearest to this value in direction of negative infinity.
*/
@SinceKotlin("1.2")
public actual fun Double.nextDown(): Double = TODO("Wasm stdlib: Math")
public actual fun Double.nextDown(): Double = when {
this.isNaN() || this == Double.NEGATIVE_INFINITY -> this
this == 0.0 -> -Double.MIN_VALUE
else -> Double.fromBits(this.toRawBits() + if (this > 0) -1 else 1)
}
/**
* Returns the [Double] value nearest to this value in direction from this value towards the value [to].
@@ -435,7 +457,7 @@ public actual fun Double.nextDown(): Double = TODO("Wasm stdlib: Math")
*
*/
@SinceKotlin("1.2")
public actual fun Double.nextTowards(to: Double): Double = TODO("Wasm stdlib: Math")
public actual fun Double.nextTowards(to: Double): Double = kotlin.math.fdlibm.nextafter(this, to)
/**
* Rounds this [Double] value to the nearest integer and converts the result to [Int].
@@ -448,7 +470,12 @@ public actual fun Double.nextTowards(to: Double): Double = TODO("Wasm stdlib: Ma
* @throws IllegalArgumentException when this value is `NaN`
*/
@SinceKotlin("1.2")
public actual fun Double.roundToInt(): Int = TODO("Wasm stdlib: Math")
public actual fun Double.roundToInt(): Int = when {
isNaN() -> throw IllegalArgumentException("Cannot round NaN value.")
this > Int.MAX_VALUE -> Int.MAX_VALUE
this < Int.MIN_VALUE -> Int.MIN_VALUE
else -> kotlin.math.fdlibm.rint(this).toInt()
}
/**
* Rounds this [Double] value to the nearest integer and converts the result to [Long].
@@ -461,7 +488,12 @@ public actual fun Double.roundToInt(): Int = TODO("Wasm stdlib: Math")
* @throws IllegalArgumentException when this value is `NaN`
*/
@SinceKotlin("1.2")
public actual fun Double.roundToLong(): Long = TODO("Wasm stdlib: Math")
public actual fun Double.roundToLong(): Long = when {
isNaN() -> throw IllegalArgumentException("Cannot round NaN value.")
this > Long.MAX_VALUE -> Long.MAX_VALUE
this < Long.MIN_VALUE -> Long.MIN_VALUE
else -> kotlin.math.fdlibm.rint(this).toLong()
}
// endregion
@@ -475,7 +507,7 @@ public actual fun Double.roundToLong(): Long = TODO("Wasm stdlib: Math")
* - `sin(NaN|+Inf|-Inf)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun sin(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun sin(x: Float): Float = kotlin.math.fdlibm.sin(x.toDouble()).toFloat()
/** Computes the cosine of the angle [x] given in radians.
*
@@ -483,7 +515,7 @@ public actual fun sin(x: Float): Float = TODO("Wasm stdlib: Math")
* - `cos(NaN|+Inf|-Inf)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun cos(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun cos(x: Float): Float = kotlin.math.fdlibm.cos(x.toDouble()).toFloat()
/** Computes the tangent of the angle [x] given in radians.
*
@@ -491,7 +523,7 @@ public actual fun cos(x: Float): Float = TODO("Wasm stdlib: Math")
* - `tan(NaN|+Inf|-Inf)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun tan(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun tan(x: Float): Float = kotlin.math.fdlibm.tan(x.toDouble()).toFloat()
/**
* Computes the arc sine of the value [x];
@@ -501,7 +533,7 @@ public actual fun tan(x: Float): Float = TODO("Wasm stdlib: Math")
* - `asin(x)` is `NaN`, when `abs(x) > 1` or x is `NaN`
*/
@SinceKotlin("1.2")
public actual fun asin(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun asin(x: Float): Float = kotlin.math.fdlibm.__ieee754_asin(x.toDouble()).toFloat()
/**
* Computes the arc cosine of the value [x];
@@ -511,7 +543,7 @@ public actual fun asin(x: Float): Float = TODO("Wasm stdlib: Math")
* - `acos(x)` is `NaN`, when `abs(x) > 1` or x is `NaN`
*/
@SinceKotlin("1.2")
public actual fun acos(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun acos(x: Float): Float = kotlin.math.fdlibm.__ieee754_acos(x.toDouble()).toFloat()
/**
* Computes the arc tangent of the value [x];
@@ -521,7 +553,7 @@ public actual fun acos(x: Float): Float = TODO("Wasm stdlib: Math")
* - `atan(NaN)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun atan(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun atan(x: Float): Float = kotlin.math.fdlibm.atan(x.toDouble()).toFloat()
/**
* Returns the angle `theta` of the polar coordinates `(r, theta)` that correspond
@@ -540,7 +572,7 @@ public actual fun atan(x: Float): Float = TODO("Wasm stdlib: Math")
* - `atan2(NaN, x)` and `atan2(y, NaN)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun atan2(y: Float, x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun atan2(y: Float, x: Float): Float = kotlin.math.fdlibm.__ieee754_atan2(y.toDouble(), x.toDouble()).toFloat()
/**
* Computes the hyperbolic sine of the value [x].
@@ -551,7 +583,7 @@ public actual fun atan2(y: Float, x: Float): Float = TODO("Wasm stdlib: Math")
* - `sinh(-Inf)` is `-Inf`
*/
@SinceKotlin("1.2")
public actual fun sinh(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun sinh(x: Float): Float = kotlin.math.fdlibm.__ieee754_sinh(x.toDouble()).toFloat()
/**
* Computes the hyperbolic cosine of the value [x].
@@ -561,7 +593,7 @@ public actual fun sinh(x: Float): Float = TODO("Wasm stdlib: Math")
* - `cosh(+Inf|-Inf)` is `+Inf`
*/
@SinceKotlin("1.2")
public actual fun cosh(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun cosh(x: Float): Float = kotlin.math.fdlibm.__ieee754_cosh(x.toDouble()).toFloat()
/**
* Computes the hyperbolic tangent of the value [x].
@@ -572,7 +604,7 @@ public actual fun cosh(x: Float): Float = TODO("Wasm stdlib: Math")
* - `tanh(-Inf)` is `-1.0`
*/
@SinceKotlin("1.2")
public actual fun tanh(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun tanh(x: Float): Float = kotlin.math.fdlibm.tanh(x.toDouble()).toFloat()
/**
* Computes the inverse hyperbolic sine of the value [x].
@@ -585,7 +617,7 @@ public actual fun tanh(x: Float): Float = TODO("Wasm stdlib: Math")
* - `asinh(-Inf)` is `-Inf`
*/
@SinceKotlin("1.2")
public actual fun asinh(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun asinh(x: Float): Float = kotlin.math.fdlibm.asinh(x.toDouble()).toFloat()
/**
* Computes the inverse hyperbolic cosine of the value [x].
@@ -598,7 +630,7 @@ public actual fun asinh(x: Float): Float = TODO("Wasm stdlib: Math")
* - `acosh(+Inf)` is `+Inf`
*/
@SinceKotlin("1.2")
public actual fun acosh(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun acosh(x: Float): Float = kotlin.math.fdlibm.__ieee754_acosh(x.toDouble()).toFloat()
/**
* Computes the inverse hyperbolic tangent of the value [x].
@@ -612,7 +644,7 @@ public actual fun acosh(x: Float): Float = TODO("Wasm stdlib: Math")
* - `tanh(-1.0)` is `-Inf`
*/
@SinceKotlin("1.2")
public actual fun atanh(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun atanh(x: Float): Float = kotlin.math.fdlibm.__ieee754_atanh(x.toDouble()).toFloat()
/**
* Computes `sqrt(x^2 + y^2)` without intermediate overflow or underflow.
@@ -622,7 +654,7 @@ public actual fun atanh(x: Float): Float = TODO("Wasm stdlib: Math")
* - returns `NaN` if any of arguments is `NaN` and the other is not infinite
*/
@SinceKotlin("1.2")
public actual fun hypot(x: Float, y: Float): Float = TODO("Wasm stdlib: Math")
public actual fun hypot(x: Float, y: Float): Float = kotlin.math.fdlibm.__ieee754_hypot(x.toDouble(), y.toDouble()).toFloat()
/**
* Computes the positive square root of the value [x].
@@ -631,7 +663,7 @@ public actual fun hypot(x: Float, y: Float): Float = TODO("Wasm stdlib: Math")
* - `sqrt(x)` is `NaN` when `x < 0` or `x` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun sqrt(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun sqrt(x: Float): Float = kotlin.wasm.internal.wasm_f32_sqrt(x)
/**
* Computes Euler's number `e` raised to the power of the value [x].
@@ -642,7 +674,7 @@ public actual fun sqrt(x: Float): Float = TODO("Wasm stdlib: Math")
* - `exp(-Inf)` is `0.0`
*/
@SinceKotlin("1.2")
public actual fun exp(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun exp(x: Float): Float = kotlin.math.fdlibm.__ieee754_exp(x.toDouble()).toFloat()
/**
* Computes `exp(x) - 1`.
@@ -657,7 +689,7 @@ public actual fun exp(x: Float): Float = TODO("Wasm stdlib: Math")
* @see [exp] function.
*/
@SinceKotlin("1.2")
public actual fun expm1(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun expm1(x: Float): Float = kotlin.math.fdlibm.expm1(x.toDouble()).toFloat()
/**
* Computes the logarithm of the value [x] to the given [base].
@@ -672,7 +704,10 @@ public actual fun expm1(x: Float): Float = TODO("Wasm stdlib: Math")
* See also logarithm functions for common fixed bases: [ln], [log10] and [log2].
*/
@SinceKotlin("1.2")
public actual fun log(x: Float, base: Float): Float = TODO("Wasm stdlib: Math")
public actual fun log(x: Float, base: Float): Float {
if (base <= 0.0F || base == 1.0F) return Float.NaN
return log(x.toDouble(), base.toDouble()).toFloat()
}
/**
* Computes the natural logarithm (base `E`) of the value [x].
@@ -684,7 +719,7 @@ public actual fun log(x: Float, base: Float): Float = TODO("Wasm stdlib: Math")
* - `ln(0.0)` is `-Inf`
*/
@SinceKotlin("1.2")
public actual fun ln(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun ln(x: Float): Float = kotlin.math.fdlibm.__ieee754_log(x.toDouble()).toFloat()
/**
* Computes the common logarithm (base 10) of the value [x].
@@ -692,7 +727,7 @@ public actual fun ln(x: Float): Float = TODO("Wasm stdlib: Math")
* @see [ln] function for special cases.
*/
@SinceKotlin("1.2")
public actual fun log10(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun log10(x: Float): Float = kotlin.math.fdlibm.__ieee754_log10(x.toDouble()).toFloat()
/**
* Computes the binary logarithm (base 2) of the value [x].
@@ -700,7 +735,7 @@ public actual fun log10(x: Float): Float = TODO("Wasm stdlib: Math")
* @see [ln] function for special cases.
*/
@SinceKotlin("1.2")
public actual fun log2(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun log2(x: Float): Float = kotlin.math.fdlibm.__ieee754_log2(x.toDouble()).toFloat()
/**
* Computes `ln(a + 1)`.
@@ -717,7 +752,7 @@ public actual fun log2(x: Float): Float = TODO("Wasm stdlib: Math")
* @see [expm1] function
*/
@SinceKotlin("1.2")
public actual fun ln1p(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun ln1p(x: Float): Float = kotlin.math.fdlibm.log1p(x.toDouble()).toFloat()
/**
* Rounds the given value [x] to an integer towards positive infinity.
@@ -728,7 +763,7 @@ public actual fun ln1p(x: Float): Float = TODO("Wasm stdlib: Math")
* - `ceil(x)` is `x` where `x` is `NaN` or `+Inf` or `-Inf` or already a mathematical integer.
*/
@SinceKotlin("1.2")
public actual fun ceil(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun ceil(x: Float): Float = kotlin.wasm.internal.wasm_f32_ceil(x)
/**
* Rounds the given value [x] to an integer towards negative infinity.
@@ -739,7 +774,7 @@ public actual fun ceil(x: Float): Float = TODO("Wasm stdlib: Math")
* - `floor(x)` is `x` where `x` is `NaN` or `+Inf` or `-Inf` or already a mathematical integer.
*/
@SinceKotlin("1.2")
public actual fun floor(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun floor(x: Float): Float = kotlin.wasm.internal.wasm_f32_floor(x)
/**
* Rounds the given value [x] to an integer towards zero.
@@ -750,7 +785,7 @@ public actual fun floor(x: Float): Float = TODO("Wasm stdlib: Math")
* - `truncate(x)` is `x` where `x` is `NaN` or `+Inf` or `-Inf` or already a mathematical integer.
*/
@SinceKotlin("1.2")
public actual fun truncate(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun truncate(x: Float): Float = kotlin.wasm.internal.wasm_f32_truncate(x)
/**
* Rounds the given value [x] towards the closest integer with ties rounded towards even integer.
@@ -759,7 +794,7 @@ public actual fun truncate(x: Float): Float = TODO("Wasm stdlib: Math")
* - `round(x)` is `x` where `x` is `NaN` or `+Inf` or `-Inf` or already a mathematical integer.
*/
@SinceKotlin("1.2")
public actual fun round(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun round(x: Float): Float = round(x.toDouble()).toFloat()
/**
@@ -771,7 +806,7 @@ public actual fun round(x: Float): Float = TODO("Wasm stdlib: Math")
* @see absoluteValue extension property for [Float]
*/
@SinceKotlin("1.2")
public actual fun abs(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun abs(x: Float): Float = kotlin.wasm.internal.wasm_f32_abs(x)
/**
* Returns the sign of the given value [x]:
@@ -783,9 +818,12 @@ public actual fun abs(x: Float): Float = TODO("Wasm stdlib: Math")
* - `sign(NaN)` is `NaN`
*/
@SinceKotlin("1.2")
public actual fun sign(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun sign(x: Float): Float = when {
x.isNaN() -> Float.NaN
x > 0.0f -> 1.0f
x < 0.0f -> -1.0f
else -> x
}
/**
* Returns the smaller of two values.
@@ -793,7 +831,7 @@ public actual fun sign(x: Float): Float = TODO("Wasm stdlib: Math")
* If either value is `NaN`, then the result is `NaN`.
*/
@SinceKotlin("1.2")
public actual fun min(a: Float, b: Float): Float = TODO("Wasm stdlib: Math")
public actual fun min(a: Float, b: Float): Float = kotlin.wasm.internal.wasm_f32_min(a, b)
/**
* Returns the greater of two values.
@@ -801,7 +839,7 @@ public actual fun min(a: Float, b: Float): Float = TODO("Wasm stdlib: Math")
* If either value is `NaN`, then the result is `NaN`.
*/
@SinceKotlin("1.2")
public actual fun max(a: Float, b: Float): Float = TODO("Wasm stdlib: Math")
public actual fun max(a: Float, b: Float): Float = kotlin.wasm.internal.wasm_f32_max(a ,b)
// extensions
@@ -818,7 +856,7 @@ public actual fun max(a: Float, b: Float): Float = TODO("Wasm stdlib: Math")
* - `b.pow(x)` is `NaN` for `b < 0` and `x` is finite and not an integer
*/
@SinceKotlin("1.2")
public actual fun Float.pow(x: Float): Float = TODO("Wasm stdlib: Math")
public actual fun Float.pow(x: Float): Float = kotlin.math.fdlibm.__ieee754_pow(this.toDouble(), x.toDouble()).toFloat()
/**
* Raises this value to the integer power [n].
@@ -826,7 +864,7 @@ public actual fun Float.pow(x: Float): Float = TODO("Wasm stdlib: Math")
* See the other overload of [pow] for details.
*/
@SinceKotlin("1.2")
public actual fun Float.pow(n: Int): Float = TODO("Wasm stdlib: Math")
public actual fun Float.pow(n: Int): Float = kotlin.math.fdlibm.__ieee754_pow(this.toDouble(), n.toDouble()).toFloat()
/**
* Returns the absolute value of this value.
@@ -837,7 +875,7 @@ public actual fun Float.pow(n: Int): Float = TODO("Wasm stdlib: Math")
* @see abs function
*/
@SinceKotlin("1.2")
public actual val Float.absoluteValue: Float get() = TODO("Wasm stdlib: Math")
public actual val Float.absoluteValue: Float get() = kotlin.wasm.internal.wasm_f32_abs(this)
/**
* Returns the sign of this value:
@@ -849,7 +887,7 @@ public actual val Float.absoluteValue: Float get() = TODO("Wasm stdlib: Math")
* - `NaN.sign` is `NaN`
*/
@SinceKotlin("1.2")
public actual val Float.sign: Float get() = TODO("Wasm stdlib: Math")
public actual val Float.sign: Float get() = sign(this)
/**
* Returns this value with the sign bit same as of the [sign] value.
@@ -857,13 +895,13 @@ public actual val Float.sign: Float get() = TODO("Wasm stdlib: Math")
* If [sign] is `NaN` the sign of the result is undefined.
*/
@SinceKotlin("1.2")
public actual fun Float.withSign(sign: Float): Float = TODO("Wasm stdlib: Math")
public actual fun Float.withSign(sign: Float): Float = kotlin.wasm.internal.wasm_f32_copysign(this, sign)
/**
* Returns this value with the sign bit same as of the [sign] value.
*/
@SinceKotlin("1.2")
public actual fun Float.withSign(sign: Int): Float = TODO("Wasm stdlib: Math")
public actual fun Float.withSign(sign: Int): Float = kotlin.wasm.internal.wasm_f32_copysign(this, sign.toFloat())
/**
@@ -877,7 +915,12 @@ public actual fun Float.withSign(sign: Int): Float = TODO("Wasm stdlib: Math")
* @throws IllegalArgumentException when this value is `NaN`
*/
@SinceKotlin("1.2")
public actual fun Float.roundToInt(): Int = TODO("Wasm stdlib: Math")
public actual fun Float.roundToInt(): Int = when {
isNaN() -> throw IllegalArgumentException("Cannot round NaN value.")
this > Int.MAX_VALUE -> Int.MAX_VALUE
this < Int.MIN_VALUE -> Int.MIN_VALUE
else -> floor(this + 0.5f).toInt()
}
/**
* Rounds this [Float] value to the nearest integer and converts the result to [Long].
@@ -890,8 +933,12 @@ public actual fun Float.roundToInt(): Int = TODO("Wasm stdlib: Math")
* @throws IllegalArgumentException when this value is `NaN`
*/
@SinceKotlin("1.2")
public actual fun Float.roundToLong(): Long = TODO("Wasm stdlib: Math")
public actual fun Float.roundToLong(): Long = when {
isNaN() -> throw IllegalArgumentException("Cannot round NaN value.")
this > Long.MAX_VALUE -> Long.MAX_VALUE
this < Long.MIN_VALUE -> Long.MIN_VALUE
else -> floor(this + 0.5f).toLong()
}
// endregion
@@ -919,7 +966,7 @@ public actual fun min(a: Int, b: Int): Int = if (a < b) a else b
* Returns the greater of two values.
*/
@SinceKotlin("1.2")
public actual fun max(a: Int, b: Int): Int = TODO("Wasm stdlib: Math")
public actual fun max(a: Int, b: Int): Int = if (a > b) a else b
/**
* Returns the absolute value of this value.
@@ -930,7 +977,7 @@ public actual fun max(a: Int, b: Int): Int = TODO("Wasm stdlib: Math")
* @see abs function
*/
@SinceKotlin("1.2")
public actual val Int.absoluteValue: Int get() = TODO("Wasm stdlib: Math")
public actual val Int.absoluteValue: Int get() = abs(this)
/**
* Returns the sign of this value:
@@ -939,11 +986,12 @@ public actual val Int.absoluteValue: Int get() = TODO("Wasm stdlib: Math")
* - `1` if the value is positive
*/
@SinceKotlin("1.2")
public actual val Int.sign: Int get() = when {
this < 0 -> -1
this > 0 -> 1
else -> 0
}
public actual val Int.sign: Int
get() = when {
this < 0 -> -1
this > 0 -> 1
else -> 0
}
/**
* Returns the absolute value of the given value [n].
@@ -954,19 +1002,19 @@ public actual val Int.sign: Int get() = when {
* @see absoluteValue extension property for [Long]
*/
@SinceKotlin("1.2")
public actual fun abs(n: Long): Long = if (n >= 0) n else -n
public actual fun abs(n: Long): Long = if (n < 0) -n else n
/**
* Returns the smaller of two values.
*/
@SinceKotlin("1.2")
public actual fun min(a: Long, b: Long): Long = TODO("Wasm stdlib: Math")
public actual fun min(a: Long, b: Long): Long = if (a <= b) a else b
/**
* Returns the greater of two values.
*/
@SinceKotlin("1.2")
public actual fun max(a: Long, b: Long): Long = TODO("Wasm stdlib: Math")
public actual fun max(a: Long, b: Long): Long = if (a >= b) a else b
/**
* Returns the absolute value of this value.
@@ -977,7 +1025,7 @@ public actual fun max(a: Long, b: Long): Long = TODO("Wasm stdlib: Math")
* @see abs function
*/
@SinceKotlin("1.2")
public actual val Long.absoluteValue: Long get() = TODO("Wasm stdlib: Math")
public actual val Long.absoluteValue: Long get() = abs(this)
/**
* Returns the sign of this value:
@@ -987,8 +1035,8 @@ public actual val Long.absoluteValue: Long get() = TODO("Wasm stdlib: Math")
*/
@SinceKotlin("1.2")
public actual val Long.sign: Int get() = when {
this < 0 -> -1
this > 0 -> 1
this < 0L -> -1
this > 0L -> 1
else -> 0
}