Adding upstream version 0.0~git20221206.8b7148f.upstream/0.0_git20221206.8b7148f

Signed-off-by: Daniel Baumann <daniel.baumann@progress-linux.org>

73 files changed, 1024 insertions, 0 deletions

diff --git a/libc-top-half/musl/src/math/i386/__invtrigl.s b/libc-top-half/musl/src/math/i386/__invtrigl.s
new file mode 100644
index 0000000..e69de29
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/__invtrigl.s
diff --git a/libc-top-half/musl/src/math/i386/acos.s b/libc-top-half/musl/src/math/i386/acos.s
new file mode 100644
index 0000000..af423a2
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/acos.s
@@ -0,0 +1,18 @@
+# use acos(x) = atan2(fabs(sqrt((1-x)*(1+x))), x)
+
+.global acos
+.type acos,@function
+acos:
+	fldl 4(%esp)
+	fld %st(0)
+	fld1
+	fsub %st(0),%st(1)
+	fadd %st(2)
+	fmulp
+	fsqrt
+	fabs         # fix sign of zero (matters in downward rounding mode)
+	fxch %st(1)
+	fpatan
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/acosf.s b/libc-top-half/musl/src/math/i386/acosf.s
new file mode 100644
index 0000000..d2cdfdb
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/acosf.s
@@ -0,0 +1,16 @@
+.global acosf
+.type acosf,@function
+acosf:
+	flds 4(%esp)
+	fld %st(0)
+	fld1
+	fsub %st(0),%st(1)
+	fadd %st(2)
+	fmulp
+	fsqrt
+	fabs         # fix sign of zero (matters in downward rounding mode)
+	fxch %st(1)
+	fpatan
+	fstps 4(%esp)
+	flds 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/acosl.s b/libc-top-half/musl/src/math/i386/acosl.s
new file mode 100644
index 0000000..599c823
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/acosl.s
@@ -0,0 +1,14 @@
+.global acosl
+.type acosl,@function
+acosl:
+	fldt 4(%esp)
+	fld %st(0)
+	fld1
+	fsub %st(0),%st(1)
+	fadd %st(2)
+	fmulp
+	fsqrt
+	fabs         # fix sign of zero (matters in downward rounding mode)
+	fxch %st(1)
+	fpatan
+	ret
diff --git a/libc-top-half/musl/src/math/i386/asin.s b/libc-top-half/musl/src/math/i386/asin.s
new file mode 100644
index 0000000..2bc8356
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/asin.s
@@ -0,0 +1,21 @@
+.global asin
+.type asin,@function
+asin:
+	fldl 4(%esp)
+	mov 8(%esp),%eax
+	add %eax,%eax
+	cmp $0x00200000,%eax
+	jb 1f
+	fld %st(0)
+	fld1
+	fsub %st(0),%st(1)
+	fadd %st(2)
+	fmulp
+	fsqrt
+	fpatan
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	ret
+		# subnormal x, return x with underflow
+1:	fsts 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/asinf.s b/libc-top-half/musl/src/math/i386/asinf.s
new file mode 100644
index 0000000..0590975
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/asinf.s
@@ -0,0 +1,23 @@
+.global asinf
+.type asinf,@function
+asinf:
+	flds 4(%esp)
+	mov 4(%esp),%eax
+	add %eax,%eax
+	cmp $0x01000000,%eax
+	jb 1f
+	fld %st(0)
+	fld1
+	fsub %st(0),%st(1)
+	fadd %st(2)
+	fmulp
+	fsqrt
+	fpatan
+	fstps 4(%esp)
+	flds 4(%esp)
+	ret
+		# subnormal x, return x with underflow
+1:	fld %st(0)
+	fmul %st(1)
+	fstps 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/asinl.s b/libc-top-half/musl/src/math/i386/asinl.s
new file mode 100644
index 0000000..e973fc8
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/asinl.s
@@ -0,0 +1,12 @@
+.global asinl
+.type asinl,@function
+asinl:
+	fldt 4(%esp)
+	fld %st(0)
+	fld1
+	fsub %st(0),%st(1)
+	fadd %st(2)
+	fmulp
+	fsqrt
+	fpatan
+	ret
diff --git a/libc-top-half/musl/src/math/i386/atan.s b/libc-top-half/musl/src/math/i386/atan.s
new file mode 100644
index 0000000..2c57f6b
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/atan.s
@@ -0,0 +1,16 @@
+.global atan
+.type atan,@function
+atan:
+	fldl 4(%esp)
+	mov 8(%esp),%eax
+	add %eax,%eax
+	cmp $0x00200000,%eax
+	jb 1f
+	fld1
+	fpatan
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	ret
+		# subnormal x, return x with underflow
+1:	fsts 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/atan2.s b/libc-top-half/musl/src/math/i386/atan2.s
new file mode 100644
index 0000000..8bc441b
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/atan2.s
@@ -0,0 +1,15 @@
+.global atan2
+.type atan2,@function
+atan2:
+	fldl 4(%esp)
+	fldl 12(%esp)
+	fpatan
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	mov 8(%esp),%eax
+	add %eax,%eax
+	cmp $0x00200000,%eax
+	jae 1f
+		# subnormal x, return x with underflow
+	fsts 4(%esp)
+1:	ret
diff --git a/libc-top-half/musl/src/math/i386/atan2f.s b/libc-top-half/musl/src/math/i386/atan2f.s
new file mode 100644
index 0000000..3908c86
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/atan2f.s
@@ -0,0 +1,17 @@
+.global atan2f
+.type atan2f,@function
+atan2f:
+	flds 4(%esp)
+	flds 8(%esp)
+	fpatan
+	fstps 4(%esp)
+	flds 4(%esp)
+	mov 4(%esp),%eax
+	add %eax,%eax
+	cmp $0x01000000,%eax
+	jae 1f
+		# subnormal x, return x with underflow
+	fld %st(0)
+	fmul %st(1)
+	fstps 4(%esp)
+1:	ret
diff --git a/libc-top-half/musl/src/math/i386/atan2l.s b/libc-top-half/musl/src/math/i386/atan2l.s
new file mode 100644
index 0000000..adf6e10
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/atan2l.s
@@ -0,0 +1,7 @@
+.global atan2l
+.type atan2l,@function
+atan2l:
+	fldt 4(%esp)
+	fldt 16(%esp)
+	fpatan
+	ret
diff --git a/libc-top-half/musl/src/math/i386/atanf.s b/libc-top-half/musl/src/math/i386/atanf.s
new file mode 100644
index 0000000..c2cbe2e
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/atanf.s
@@ -0,0 +1,18 @@
+.global atanf
+.type atanf,@function
+atanf:
+	flds 4(%esp)
+	mov 4(%esp),%eax
+	add %eax,%eax
+	cmp $0x01000000,%eax
+	jb 1f
+	fld1
+	fpatan
+	fstps 4(%esp)
+	flds 4(%esp)
+	ret
+		# subnormal x, return x with underflow
+1:	fld %st(0)
+	fmul %st(1)
+	fstps 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/atanl.s b/libc-top-half/musl/src/math/i386/atanl.s
new file mode 100644
index 0000000..c508bc4
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/atanl.s
@@ -0,0 +1,7 @@
+.global atanl
+.type atanl,@function
+atanl:
+	fldt 4(%esp)
+	fld1
+	fpatan
+	ret
diff --git a/libc-top-half/musl/src/math/i386/ceil.s b/libc-top-half/musl/src/math/i386/ceil.s
new file mode 100644
index 0000000..bc29f15
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/ceil.s
@@ -0,0 +1 @@
+# see floor.s
diff --git a/libc-top-half/musl/src/math/i386/ceilf.s b/libc-top-half/musl/src/math/i386/ceilf.s
new file mode 100644
index 0000000..bc29f15
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/ceilf.s
@@ -0,0 +1 @@
+# see floor.s
diff --git a/libc-top-half/musl/src/math/i386/ceill.s b/libc-top-half/musl/src/math/i386/ceill.s
new file mode 100644
index 0000000..bc29f15
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/ceill.s
@@ -0,0 +1 @@
+# see floor.s
diff --git a/libc-top-half/musl/src/math/i386/exp2l.s b/libc-top-half/musl/src/math/i386/exp2l.s
new file mode 100644
index 0000000..8125761
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/exp2l.s
@@ -0,0 +1 @@
+# see exp_ld.s
diff --git a/libc-top-half/musl/src/math/i386/exp_ld.s b/libc-top-half/musl/src/math/i386/exp_ld.s
new file mode 100644
index 0000000..99cba01
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/exp_ld.s
@@ -0,0 +1,93 @@
+.global expm1l
+.type expm1l,@function
+expm1l:
+	fldt 4(%esp)
+	fldl2e
+	fmulp
+	mov $0xc2820000,%eax
+	push %eax
+	flds (%esp)
+	pop %eax
+	fucomp %st(1)
+	fnstsw %ax
+	sahf
+	fld1
+	jb 1f
+		# x*log2e < -65, return -1 without underflow
+	fstp %st(1)
+	fchs
+	ret
+1:	fld %st(1)
+	fabs
+	fucom %st(1)
+	fnstsw %ax
+	fstp %st(0)
+	fstp %st(0)
+	sahf
+	ja 1f
+	f2xm1
+	ret
+1:	call 1f
+	fld1
+	fsubrp
+	ret
+
+.global exp2l
+.global __exp2l
+.hidden __exp2l
+.type exp2l,@function
+exp2l:
+__exp2l:
+	fldt 4(%esp)
+1:	sub $12,%esp
+	fld %st(0)
+	fstpt (%esp)
+	mov 8(%esp),%ax
+	and $0x7fff,%ax
+	cmp $0x3fff+13,%ax
+	jb 4f             # |x| < 8192
+	cmp $0x3fff+15,%ax
+	jae 3f            # |x| >= 32768
+	fsts (%esp)
+	cmpl $0xc67ff800,(%esp)
+	jb 2f             # x > -16382
+	movl $0x5f000000,(%esp)
+	flds (%esp)       # 0x1p63
+	fld %st(1)
+	fsub %st(1)
+	faddp
+	fucomp %st(1)
+	fnstsw
+	sahf
+	je 2f             # x - 0x1p63 + 0x1p63 == x
+	movl $1,(%esp)
+	flds (%esp)       # 0x1p-149
+	fdiv %st(1)
+	fstps (%esp)      # raise underflow
+2:	fld1
+	fld %st(1)
+	frndint
+	fxch %st(2)
+	fsub %st(2)       # st(0)=x-rint(x), st(1)=1, st(2)=rint(x)
+	f2xm1
+	faddp             # 2^(x-rint(x))
+1:	fscale
+	fstp %st(1)
+	add $12,%esp
+	ret
+3:	xor %eax,%eax
+4:	cmp $0x3fff-64,%ax
+	fld1
+	jb 1b             # |x| < 0x1p-64
+	fstpt (%esp)
+	fistl 8(%esp)
+	fildl 8(%esp)
+	fsubrp %st(1)
+	addl $0x3fff,8(%esp)
+	f2xm1
+	fld1
+	faddp             # 2^(x-rint(x))
+	fldt (%esp)       # 2^rint(x)
+	fmulp
+	add $12,%esp
+	ret
diff --git a/libc-top-half/musl/src/math/i386/expl.s b/libc-top-half/musl/src/math/i386/expl.s
new file mode 100644
index 0000000..b5124e8
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/expl.s
@@ -0,0 +1,101 @@
+# exp(x) = 2^hi + 2^hi (2^lo - 1)
+# where hi+lo = log2e*x with 128bit precision
+# exact log2e*x calculation depends on nearest rounding mode
+# using the exact multiplication method of Dekker and Veltkamp
+
+.global expl
+.type expl,@function
+expl:
+	fldt 4(%esp)
+
+		# interesting case: 0x1p-32 <= |x| < 16384
+		# check if (exponent|0x8000) is in [0xbfff-32, 0xbfff+13]
+	mov 12(%esp), %ax
+	or $0x8000, %ax
+	sub $0xbfdf, %ax
+	cmp $45, %ax
+	jbe 2f
+	test %ax, %ax
+	fld1
+	js 1f
+		# if |x|>=0x1p14 or nan return 2^trunc(x)
+	fscale
+	fstp %st(1)
+	ret
+		# if |x|<0x1p-32 return 1+x
+1:	faddp
+	ret
+
+		# should be 0x1.71547652b82fe178p0L == 0x3fff b8aa3b29 5c17f0bc
+		# it will be wrong on non-nearest rounding mode
+2:	fldl2e
+	subl $44, %esp
+		# hi = log2e_hi*x
+		# 2^hi = exp2l(hi)
+	fmul %st(1),%st
+	fld %st(0)
+	fstpt (%esp)
+	fstpt 16(%esp)
+	fstpt 32(%esp)
+.hidden __exp2l
+	call __exp2l
+		# if 2^hi == inf return 2^hi
+	fld %st(0)
+	fstpt (%esp)
+	cmpw $0x7fff, 8(%esp)
+	je 1f
+	fldt 32(%esp)
+	fldt 16(%esp)
+		# fpu stack: 2^hi x hi
+		# exact mult: x*log2e
+	fld %st(1)
+		# c = 0x1p32+1
+	pushl $0x41f00000
+	pushl $0x00100000
+	fldl (%esp)
+		# xh = x - c*x + c*x
+		# xl = x - xh
+	fmulp
+	fld %st(2)
+	fsub %st(1), %st
+	faddp
+	fld %st(2)
+	fsub %st(1), %st
+		# yh = log2e_hi - c*log2e_hi + c*log2e_hi
+	pushl $0x3ff71547
+	pushl $0x65200000
+	fldl (%esp)
+		# fpu stack: 2^hi x hi xh xl yh
+		# lo = hi - xh*yh + xl*yh
+	fld %st(2)
+	fmul %st(1), %st
+	fsubp %st, %st(4)
+	fmul %st(1), %st
+	faddp %st, %st(3)
+		# yl = log2e_hi - yh
+	pushl $0x3de705fc
+	pushl $0x2f000000
+	fldl (%esp)
+		# fpu stack: 2^hi x lo xh xl yl
+		# lo += xh*yl + xl*yl
+	fmul %st, %st(2)
+	fmulp %st, %st(1)
+	fxch %st(2)
+	faddp
+	faddp
+		# log2e_lo
+	pushl $0xbfbe
+	pushl $0x82f0025f
+	pushl $0x2dc582ee
+	fldt (%esp)
+	addl $36,%esp
+		# fpu stack: 2^hi x lo log2e_lo
+		# lo += log2e_lo*x
+		# return 2^hi + 2^hi (2^lo - 1)
+	fmulp %st, %st(2)
+	faddp
+	f2xm1
+	fmul %st(1), %st
+	faddp
+1:	addl $44, %esp
+	ret
diff --git a/libc-top-half/musl/src/math/i386/expm1l.s b/libc-top-half/musl/src/math/i386/expm1l.s
new file mode 100644
index 0000000..8125761
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/expm1l.s
@@ -0,0 +1 @@
+# see exp_ld.s
diff --git a/libc-top-half/musl/src/math/i386/fabs.c b/libc-top-half/musl/src/math/i386/fabs.c
new file mode 100644
index 0000000..3967278
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/fabs.c
@@ -0,0 +1,7 @@
+#include <math.h>
+
+double fabs(double x)
+{
+	__asm__ ("fabs" : "+t"(x));
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/fabsf.c b/libc-top-half/musl/src/math/i386/fabsf.c
new file mode 100644
index 0000000..d882eee
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/fabsf.c
@@ -0,0 +1,7 @@
+#include <math.h>
+
+float fabsf(float x)
+{
+	__asm__ ("fabs" : "+t"(x));
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/fabsl.c b/libc-top-half/musl/src/math/i386/fabsl.c
new file mode 100644
index 0000000..cc1c9ed
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/fabsl.c
@@ -0,0 +1,7 @@
+#include <math.h>
+
+long double fabsl(long double x)
+{
+	__asm__ ("fabs" : "+t"(x));
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/floor.s b/libc-top-half/musl/src/math/i386/floor.s
new file mode 100644
index 0000000..46ba88d
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/floor.s
@@ -0,0 +1,67 @@
+.global floorf
+.type floorf,@function
+floorf:
+	flds 4(%esp)
+	jmp 1f
+
+.global floorl
+.type floorl,@function
+floorl:
+	fldt 4(%esp)
+	jmp 1f
+
+.global floor
+.type floor,@function
+floor:
+	fldl 4(%esp)
+1:	mov $0x7,%al
+1:	fstcw 4(%esp)
+	mov 5(%esp),%ah
+	mov %al,5(%esp)
+	fldcw 4(%esp)
+	frndint
+	mov %ah,5(%esp)
+	fldcw 4(%esp)
+	ret
+
+.global ceil
+.type ceil,@function
+ceil:
+	fldl 4(%esp)
+	mov $0xb,%al
+	jmp 1b
+
+.global ceilf
+.type ceilf,@function
+ceilf:
+	flds 4(%esp)
+	mov $0xb,%al
+	jmp 1b
+
+.global ceill
+.type ceill,@function
+ceill:
+	fldt 4(%esp)
+	mov $0xb,%al
+	jmp 1b
+
+.global trunc
+.type trunc,@function
+trunc:
+	fldl 4(%esp)
+	mov $0xf,%al
+	jmp 1b
+
+.global truncf
+.type truncf,@function
+truncf:
+	flds 4(%esp)
+	mov $0xf,%al
+	jmp 1b
+
+.global truncl
+.type truncl,@function
+truncl:
+	fldt 4(%esp)
+	mov $0xf,%al
+	jmp 1b
diff --git a/libc-top-half/musl/src/math/i386/floorf.s b/libc-top-half/musl/src/math/i386/floorf.s
new file mode 100644
index 0000000..bc29f15
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/floorf.s
@@ -0,0 +1 @@
+# see floor.s
diff --git a/libc-top-half/musl/src/math/i386/floorl.s b/libc-top-half/musl/src/math/i386/floorl.s
new file mode 100644
index 0000000..bc29f15
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/floorl.s
@@ -0,0 +1 @@
+# see floor.s
diff --git a/libc-top-half/musl/src/math/i386/fmod.c b/libc-top-half/musl/src/math/i386/fmod.c
new file mode 100644
index 0000000..ea0c58d
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/fmod.c
@@ -0,0 +1,10 @@
+#include <math.h>
+
+double fmod(double x, double y)
+{
+	unsigned short fpsr;
+	// fprem does not introduce excess precision into x
+	do __asm__ ("fprem; fnstsw %%ax" : "+t"(x), "=a"(fpsr) : "u"(y));
+	while (fpsr & 0x400);
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/fmodf.c b/libc-top-half/musl/src/math/i386/fmodf.c
new file mode 100644
index 0000000..90b56ab
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/fmodf.c
@@ -0,0 +1,10 @@
+#include <math.h>
+
+float fmodf(float x, float y)
+{
+	unsigned short fpsr;
+	// fprem does not introduce excess precision into x
+	do __asm__ ("fprem; fnstsw %%ax" : "+t"(x), "=a"(fpsr) : "u"(y));
+	while (fpsr & 0x400);
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/fmodl.c b/libc-top-half/musl/src/math/i386/fmodl.c
new file mode 100644
index 0000000..3daeab0
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/fmodl.c
@@ -0,0 +1,9 @@
+#include <math.h>
+
+long double fmodl(long double x, long double y)
+{
+	unsigned short fpsr;
+	do __asm__ ("fprem; fnstsw %%ax" : "+t"(x), "=a"(fpsr) : "u"(y));
+	while (fpsr & 0x400);
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/hypot.s b/libc-top-half/musl/src/math/i386/hypot.s
new file mode 100644
index 0000000..299c2e1
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/hypot.s
@@ -0,0 +1,45 @@
+.global hypot
+.type hypot,@function
+hypot:
+	mov 8(%esp),%eax
+	mov 16(%esp),%ecx
+	add %eax,%eax
+	add %ecx,%ecx
+	and %eax,%ecx
+	cmp $0xffe00000,%ecx
+	jae 2f
+	or 4(%esp),%eax
+	jnz 1f
+	fldl 12(%esp)
+	fabs
+	ret
+1:	mov 16(%esp),%eax
+	add %eax,%eax
+	or 12(%esp),%eax
+	jnz 1f
+	fldl 4(%esp)
+	fabs
+	ret
+1:	fldl 4(%esp)
+	fld %st(0)
+	fmulp
+	fldl 12(%esp)
+	fld %st(0)
+	fmulp
+	faddp
+	fsqrt
+	ret
+2:	sub $0xffe00000,%eax
+	or 4(%esp),%eax
+	jnz 1f
+	fldl 4(%esp)
+	fabs
+	ret
+1:	mov 16(%esp),%eax
+	add %eax,%eax
+	sub $0xffe00000,%eax
+	or 12(%esp),%eax
+	fldl 12(%esp)
+	jnz 1f
+	fabs
+1:	ret
diff --git a/libc-top-half/musl/src/math/i386/hypotf.s b/libc-top-half/musl/src/math/i386/hypotf.s
new file mode 100644
index 0000000..068935e
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/hypotf.s
@@ -0,0 +1,42 @@
+.global hypotf
+.type hypotf,@function
+hypotf:
+	mov 4(%esp),%eax
+	mov 8(%esp),%ecx
+	add %eax,%eax
+	add %ecx,%ecx
+	and %eax,%ecx
+	cmp $0xff000000,%ecx
+	jae 2f
+	test %eax,%eax
+	jnz 1f
+	flds 8(%esp)
+	fabs
+	ret
+1:	mov 8(%esp),%eax
+	add %eax,%eax
+	jnz 1f
+	flds 4(%esp)
+	fabs
+	ret
+1:	flds 4(%esp)
+	fld %st(0)
+	fmulp
+	flds 8(%esp)
+	fld %st(0)
+	fmulp
+	faddp
+	fsqrt
+	ret
+2:	cmp $0xff000000,%eax
+	jnz 1f
+	flds 4(%esp)
+	fabs
+	ret
+1:	mov 8(%esp),%eax
+	add %eax,%eax
+	cmp $0xff000000,%eax
+	flds 8(%esp)
+	jnz 1f
+	fabs
+1:	ret
diff --git a/libc-top-half/musl/src/math/i386/ldexp.s b/libc-top-half/musl/src/math/i386/ldexp.s
new file mode 100644
index 0000000..c430f78
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/ldexp.s
@@ -0,0 +1 @@
+# see scalbn.s
diff --git a/libc-top-half/musl/src/math/i386/ldexpf.s b/libc-top-half/musl/src/math/i386/ldexpf.s
new file mode 100644
index 0000000..3f8e4b9
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/ldexpf.s
@@ -0,0 +1 @@
+# see scalbnf.s
diff --git a/libc-top-half/musl/src/math/i386/ldexpl.s b/libc-top-half/musl/src/math/i386/ldexpl.s
new file mode 100644
index 0000000..86fe562
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/ldexpl.s
@@ -0,0 +1 @@
+# see scalbnl.s
diff --git a/libc-top-half/musl/src/math/i386/llrint.c b/libc-top-half/musl/src/math/i386/llrint.c
new file mode 100644
index 0000000..aa40081
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/llrint.c
@@ -0,0 +1,8 @@
+#include <math.h>
+
+long long llrint(double x)
+{
+	long long r;
+	__asm__ ("fistpll %0" : "=m"(r) : "t"(x) : "st");
+	return r;
+}
diff --git a/libc-top-half/musl/src/math/i386/llrintf.c b/libc-top-half/musl/src/math/i386/llrintf.c
new file mode 100644
index 0000000..c41a317
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/llrintf.c
@@ -0,0 +1,8 @@
+#include <math.h>
+
+long long llrintf(float x)
+{
+	long long r;
+	__asm__ ("fistpll %0" : "=m"(r) : "t"(x) : "st");
+	return r;
+}
diff --git a/libc-top-half/musl/src/math/i386/llrintl.c b/libc-top-half/musl/src/math/i386/llrintl.c
new file mode 100644
index 0000000..c439ef2
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/llrintl.c
@@ -0,0 +1,8 @@
+#include <math.h>
+
+long long llrintl(long double x)
+{
+	long long r;
+	__asm__ ("fistpll %0" : "=m"(r) : "t"(x) : "st");
+	return r;
+}
diff --git a/libc-top-half/musl/src/math/i386/log.s b/libc-top-half/musl/src/math/i386/log.s
new file mode 100644
index 0000000..08c5992
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log.s
@@ -0,0 +1,9 @@
+.global log
+.type log,@function
+log:
+	fldln2
+	fldl 4(%esp)
+	fyl2x
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/log10.s b/libc-top-half/musl/src/math/i386/log10.s
new file mode 100644
index 0000000..120e91e
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log10.s
@@ -0,0 +1,9 @@
+.global log10
+.type log10,@function
+log10:
+	fldlg2
+	fldl 4(%esp)
+	fyl2x
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/log10f.s b/libc-top-half/musl/src/math/i386/log10f.s
new file mode 100644
index 0000000..b055493
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log10f.s
@@ -0,0 +1,9 @@
+.global log10f
+.type log10f,@function
+log10f:
+	fldlg2
+	flds 4(%esp)
+	fyl2x
+	fstps 4(%esp)
+	flds 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/log10l.s b/libc-top-half/musl/src/math/i386/log10l.s
new file mode 100644
index 0000000..aaa44f2
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log10l.s
@@ -0,0 +1,7 @@
+.global log10l
+.type log10l,@function
+log10l:
+	fldlg2
+	fldt 4(%esp)
+	fyl2x
+	ret
diff --git a/libc-top-half/musl/src/math/i386/log1p.s b/libc-top-half/musl/src/math/i386/log1p.s
new file mode 100644
index 0000000..f3c95f8
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log1p.s
@@ -0,0 +1,25 @@
+.global log1p
+.type log1p,@function
+log1p:
+	mov 8(%esp),%eax
+	fldln2
+	and $0x7fffffff,%eax
+	fldl 4(%esp)
+	cmp $0x3fd28f00,%eax
+	ja 1f
+	cmp $0x00100000,%eax
+	jb 2f
+	fyl2xp1
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	ret
+1:	fld1
+	faddp
+	fyl2x
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	ret
+		# subnormal x, return x with underflow
+2:	fsts 4(%esp)
+	fstp %st(1)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/log1pf.s b/libc-top-half/musl/src/math/i386/log1pf.s
new file mode 100644
index 0000000..9f13d95
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log1pf.s
@@ -0,0 +1,26 @@
+.global log1pf
+.type log1pf,@function
+log1pf:
+	mov 4(%esp),%eax
+	fldln2
+	and $0x7fffffff,%eax
+	flds 4(%esp)
+	cmp $0x3e940000,%eax
+	ja 1f
+	cmp $0x00800000,%eax
+	jb 2f
+	fyl2xp1
+	fstps 4(%esp)
+	flds 4(%esp)
+	ret
+1:	fld1
+	faddp
+	fyl2x
+	fstps 4(%esp)
+	flds 4(%esp)
+	ret
+		# subnormal x, return x with underflow
+2:	fxch
+	fmul %st(1)
+	fstps 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/log1pl.s b/libc-top-half/musl/src/math/i386/log1pl.s
new file mode 100644
index 0000000..a048ab6
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log1pl.s
@@ -0,0 +1,15 @@
+.global log1pl
+.type log1pl,@function
+log1pl:
+	mov 10(%esp),%eax
+	fldln2
+	and $0x7fffffff,%eax
+	fldt 4(%esp)
+	cmp $0x3ffd9400,%eax
+	ja 1f
+	fyl2xp1
+	ret
+1:	fld1
+	faddp
+	fyl2x
+	ret
diff --git a/libc-top-half/musl/src/math/i386/log2.s b/libc-top-half/musl/src/math/i386/log2.s
new file mode 100644
index 0000000..7eff0b6
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log2.s
@@ -0,0 +1,9 @@
+.global log2
+.type log2,@function
+log2:
+	fld1
+	fldl 4(%esp)
+	fyl2x
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/log2f.s b/libc-top-half/musl/src/math/i386/log2f.s
new file mode 100644
index 0000000..b32fa2f
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log2f.s
@@ -0,0 +1,9 @@
+.global log2f
+.type log2f,@function
+log2f:
+	fld1
+	flds 4(%esp)
+	fyl2x
+	fstps 4(%esp)
+	flds 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/log2l.s b/libc-top-half/musl/src/math/i386/log2l.s
new file mode 100644
index 0000000..c58f56f
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/log2l.s
@@ -0,0 +1,7 @@
+.global log2l
+.type log2l,@function
+log2l:
+	fld1
+	fldt 4(%esp)
+	fyl2x
+	ret
diff --git a/libc-top-half/musl/src/math/i386/logf.s b/libc-top-half/musl/src/math/i386/logf.s
new file mode 100644
index 0000000..4d0346a
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/logf.s
@@ -0,0 +1,9 @@
+.global logf
+.type logf,@function
+logf:
+	fldln2
+	flds 4(%esp)
+	fyl2x
+	fstps 4(%esp)
+	flds 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/logl.s b/libc-top-half/musl/src/math/i386/logl.s
new file mode 100644
index 0000000..d4e3339
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/logl.s
@@ -0,0 +1,7 @@
+.global logl
+.type logl,@function
+logl:
+	fldln2
+	fldt 4(%esp)
+	fyl2x
+	ret
diff --git a/libc-top-half/musl/src/math/i386/lrint.c b/libc-top-half/musl/src/math/i386/lrint.c
new file mode 100644
index 0000000..89563ab
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/lrint.c
@@ -0,0 +1,8 @@
+#include <math.h>
+
+long lrint(double x)
+{
+	long r;
+	__asm__ ("fistpl %0" : "=m"(r) : "t"(x) : "st");
+	return r;
+}
diff --git a/libc-top-half/musl/src/math/i386/lrintf.c b/libc-top-half/musl/src/math/i386/lrintf.c
new file mode 100644
index 0000000..0bbf29d
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/lrintf.c
@@ -0,0 +1,8 @@
+#include <math.h>
+
+long lrintf(float x)
+{
+	long r;
+	__asm__ ("fistpl %0" : "=m"(r) : "t"(x) : "st");
+	return r;
+}
diff --git a/libc-top-half/musl/src/math/i386/lrintl.c b/libc-top-half/musl/src/math/i386/lrintl.c
new file mode 100644
index 0000000..eb8c090
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/lrintl.c
@@ -0,0 +1,8 @@
+#include <math.h>
+
+long lrintl(long double x)
+{
+	long r;
+	__asm__ ("fistpl %0" : "=m"(r) : "t"(x) : "st");
+	return r;
+}
diff --git a/libc-top-half/musl/src/math/i386/remainder.c b/libc-top-half/musl/src/math/i386/remainder.c
new file mode 100644
index 0000000..c083df9
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/remainder.c
@@ -0,0 +1,12 @@
+#include <math.h>
+
+double remainder(double x, double y)
+{
+	unsigned short fpsr;
+	// fprem1 does not introduce excess precision into x
+	do __asm__ ("fprem1; fnstsw %%ax" : "+t"(x), "=a"(fpsr) : "u"(y));
+	while (fpsr & 0x400);
+	return x;
+}
+
+weak_alias(remainder, drem);
diff --git a/libc-top-half/musl/src/math/i386/remainderf.c b/libc-top-half/musl/src/math/i386/remainderf.c
new file mode 100644
index 0000000..280207d
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/remainderf.c
@@ -0,0 +1,12 @@
+#include <math.h>
+
+float remainderf(float x, float y)
+{
+	unsigned short fpsr;
+	// fprem1 does not introduce excess precision into x
+	do __asm__ ("fprem1; fnstsw %%ax" : "+t"(x), "=a"(fpsr) : "u"(y));
+	while (fpsr & 0x400);
+	return x;
+}
+
+weak_alias(remainderf, dremf);
diff --git a/libc-top-half/musl/src/math/i386/remainderl.c b/libc-top-half/musl/src/math/i386/remainderl.c
new file mode 100644
index 0000000..8cf7507
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/remainderl.c
@@ -0,0 +1,9 @@
+#include <math.h>
+
+long double remainderl(long double x, long double y)
+{
+	unsigned short fpsr;
+	do __asm__ ("fprem1; fnstsw %%ax" : "+t"(x), "=a"(fpsr) : "u"(y));
+	while (fpsr & 0x400);
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/remquo.s b/libc-top-half/musl/src/math/i386/remquo.s
new file mode 100644
index 0000000..598e754
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/remquo.s
@@ -0,0 +1,50 @@
+.global remquof
+.type remquof,@function
+remquof:
+	mov 12(%esp),%ecx
+	flds 8(%esp)
+	flds 4(%esp)
+	mov 11(%esp),%dh
+	xor 7(%esp),%dh
+	jmp 1f
+
+.global remquol
+.type remquol,@function
+remquol:
+	mov 28(%esp),%ecx
+	fldt 16(%esp)
+	fldt 4(%esp)
+	mov 25(%esp),%dh
+	xor 13(%esp),%dh
+	jmp 1f
+
+.global remquo
+.type remquo,@function
+remquo:
+	mov 20(%esp),%ecx
+	fldl 12(%esp)
+	fldl 4(%esp)
+	mov 19(%esp),%dh
+	xor 11(%esp),%dh
+1:	fprem1
+	fnstsw %ax
+	sahf
+	jp 1b
+	fstp %st(1)
+	mov %ah,%dl
+	shr %dl
+	and $1,%dl
+	mov %ah,%al
+	shr $5,%al
+	and $2,%al
+	or %al,%dl
+	mov %ah,%al
+	shl $2,%al
+	and $4,%al
+	or %al,%dl
+	test %dh,%dh
+	jns 1f
+	neg %dl
+1:	movsbl %dl,%edx
+	mov %edx,(%ecx)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/remquof.s b/libc-top-half/musl/src/math/i386/remquof.s
new file mode 100644
index 0000000..511a6bc
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/remquof.s
@@ -0,0 +1 @@
+# see remquo.s
diff --git a/libc-top-half/musl/src/math/i386/remquol.s b/libc-top-half/musl/src/math/i386/remquol.s
new file mode 100644
index 0000000..511a6bc
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/remquol.s
@@ -0,0 +1 @@
+# see remquo.s
diff --git a/libc-top-half/musl/src/math/i386/rint.c b/libc-top-half/musl/src/math/i386/rint.c
new file mode 100644
index 0000000..a5276a6
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/rint.c
@@ -0,0 +1,7 @@
+#include <math.h>
+
+double rint(double x)
+{
+	__asm__ ("frndint" : "+t"(x));
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/rintf.c b/libc-top-half/musl/src/math/i386/rintf.c
new file mode 100644
index 0000000..bb4121a
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/rintf.c
@@ -0,0 +1,7 @@
+#include <math.h>
+
+float rintf(float x)
+{
+	__asm__ ("frndint" : "+t"(x));
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/rintl.c b/libc-top-half/musl/src/math/i386/rintl.c
new file mode 100644
index 0000000..e1a9207
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/rintl.c
@@ -0,0 +1,7 @@
+#include <math.h>
+
+long double rintl(long double x)
+{
+	__asm__ ("frndint" : "+t"(x));
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/scalbln.s b/libc-top-half/musl/src/math/i386/scalbln.s
new file mode 100644
index 0000000..c430f78
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/scalbln.s
@@ -0,0 +1 @@
+# see scalbn.s
diff --git a/libc-top-half/musl/src/math/i386/scalblnf.s b/libc-top-half/musl/src/math/i386/scalblnf.s
new file mode 100644
index 0000000..3f8e4b9
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/scalblnf.s
@@ -0,0 +1 @@
+# see scalbnf.s
diff --git a/libc-top-half/musl/src/math/i386/scalblnl.s b/libc-top-half/musl/src/math/i386/scalblnl.s
new file mode 100644
index 0000000..86fe562
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/scalblnl.s
@@ -0,0 +1 @@
+# see scalbnl.s
diff --git a/libc-top-half/musl/src/math/i386/scalbn.s b/libc-top-half/musl/src/math/i386/scalbn.s
new file mode 100644
index 0000000..8bf302f
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/scalbn.s
@@ -0,0 +1,33 @@
+.global ldexp
+.type ldexp,@function
+ldexp:
+	nop
+
+.global scalbln
+.type scalbln,@function
+scalbln:
+	nop
+
+.global scalbn
+.type scalbn,@function
+scalbn:
+	mov 12(%esp),%eax
+	add $0x3ffe,%eax
+	cmp $0x7ffd,%eax
+	jb 1f
+	sub $0x3ffe,%eax
+	sar $31,%eax
+	xor $0xfff,%eax
+	add $0x3ffe,%eax
+1:	inc %eax
+	fldl 4(%esp)
+	mov %eax,12(%esp)
+	mov $0x80000000,%eax
+	mov %eax,8(%esp)
+	xor %eax,%eax
+	mov %eax,4(%esp)
+	fldt 4(%esp)
+	fmulp
+	fstpl 4(%esp)
+	fldl 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/scalbnf.s b/libc-top-half/musl/src/math/i386/scalbnf.s
new file mode 100644
index 0000000..9cb9ef5
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/scalbnf.s
@@ -0,0 +1,32 @@
+.global ldexpf
+.type ldexpf,@function
+ldexpf:
+	nop
+
+.global scalblnf
+.type scalblnf,@function
+scalblnf:
+	nop
+
+.global scalbnf
+.type scalbnf,@function
+scalbnf:
+	mov 8(%esp),%eax
+	add $0x3fe,%eax
+	cmp $0x7fd,%eax
+	jb 1f
+	sub $0x3fe,%eax
+	sar $31,%eax
+	xor $0x1ff,%eax
+	add $0x3fe,%eax
+1:	inc %eax
+	shl $20,%eax
+	flds 4(%esp)
+	mov %eax,8(%esp)
+	xor %eax,%eax
+	mov %eax,4(%esp)
+	fldl 4(%esp)
+	fmulp
+	fstps 4(%esp)
+	flds 4(%esp)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/scalbnl.s b/libc-top-half/musl/src/math/i386/scalbnl.s
new file mode 100644
index 0000000..54414c2
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/scalbnl.s
@@ -0,0 +1,32 @@
+.global ldexpl
+.type ldexpl,@function
+ldexpl:
+	nop
+
+.global scalblnl
+.type scalblnl,@function
+scalblnl:
+	nop
+
+.global scalbnl
+.type scalbnl,@function
+scalbnl:
+	mov 16(%esp),%eax
+	add $0x3ffe,%eax
+	cmp $0x7ffd,%eax
+	jae 1f
+	inc %eax
+	fldt 4(%esp)
+	mov %eax,12(%esp)
+	mov $0x80000000,%eax
+	mov %eax,8(%esp)
+	xor %eax,%eax
+	mov %eax,4(%esp)
+	fldt 4(%esp)
+	fmulp
+	ret
+1:	fildl 16(%esp)
+	fldt 4(%esp)
+	fscale
+	fstp %st(1)
+	ret
diff --git a/libc-top-half/musl/src/math/i386/sqrt.c b/libc-top-half/musl/src/math/i386/sqrt.c
new file mode 100644
index 0000000..934fbcc
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/sqrt.c
@@ -0,0 +1,15 @@
+#include "libm.h"
+
+double sqrt(double x)
+{
+	union ldshape ux;
+	unsigned fpsr;
+	__asm__ ("fsqrt; fnstsw %%ax": "=t"(ux.f), "=a"(fpsr) : "0"(x));
+	if ((ux.i.m & 0x7ff) != 0x400)
+		return (double)ux.f;
+	/* Rounding to double would have encountered an exact halfway case.
+	   Adjust mantissa downwards if fsqrt rounded up, else upwards.
+	   (result of fsqrt could not have been exact) */
+	ux.i.m ^= (fpsr & 0x200) + 0x300;
+	return (double)ux.f;
+}
diff --git a/libc-top-half/musl/src/math/i386/sqrtf.c b/libc-top-half/musl/src/math/i386/sqrtf.c
new file mode 100644
index 0000000..41c65c2
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/sqrtf.c
@@ -0,0 +1,12 @@
+#include <math.h>
+
+float sqrtf(float x)
+{
+	long double t;
+	/* The long double result has sufficient precision so that
+	 * second rounding to float still keeps the returned value
+	 * correctly rounded, see Pierre Roux, "Innocuous Double
+	 * Rounding of Basic Arithmetic Operations". */
+	__asm__ ("fsqrt" : "=t"(t) : "0"(x));
+	return (float)t;
+}
diff --git a/libc-top-half/musl/src/math/i386/sqrtl.c b/libc-top-half/musl/src/math/i386/sqrtl.c
new file mode 100644
index 0000000..864cfcc
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/sqrtl.c
@@ -0,0 +1,7 @@
+#include <math.h>
+
+long double sqrtl(long double x)
+{
+	__asm__ ("fsqrt" : "+t"(x));
+	return x;
+}
diff --git a/libc-top-half/musl/src/math/i386/trunc.s b/libc-top-half/musl/src/math/i386/trunc.s
new file mode 100644
index 0000000..bc29f15
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/trunc.s
@@ -0,0 +1 @@
+# see floor.s
diff --git a/libc-top-half/musl/src/math/i386/truncf.s b/libc-top-half/musl/src/math/i386/truncf.s
new file mode 100644
index 0000000..bc29f15
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/truncf.s
@@ -0,0 +1 @@
+# see floor.s
diff --git a/libc-top-half/musl/src/math/i386/truncl.s b/libc-top-half/musl/src/math/i386/truncl.s
new file mode 100644
index 0000000..bc29f15
--- /dev/null
+++ b/libc-top-half/musl/src/math/i386/truncl.s
@@ -0,0 +1 @@
+# see floor.s