mirrors/cosmopolitan
1
0
Fork 0
mirror of https://github.com/jart/cosmopolitan.git synced 2025-02-07 06:53:33 +00:00
Code Issues Projects Releases Packages Wiki Activity Actions

Get libc/tinymath/ compiling on aarch64

Browse source
This commit is contained in:
Justine Tunney 2023-05-02 18:35:25 -07:00
parent 2b73e72d59
commit 135080fd3e
No known key found for this signature in database
GPG key ID: BE714B4575D6E328
243 changed files with 7773 additions and 4027 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

4
Makefile
View file

@ -288,12 +288,12 @@ o/$(MODE)/hdrs-old.txt: o/$(MODE)/.x $(MAKEFILES) $(call uniq,$(foreach x,$(HDRS
$(file >$@) $(foreach x,$(HDRS) $(INCS),$(file >>$@,$(x)))
TAGS: private .UNSANDBOXED = 1
TAGS: o/$(MODE)/srcs-old.txt $(SRCS) o/$(MODE)/third_party/ctags/ctags.com
TAGS: o/$(MODE)/srcs-old.txt $(SRCS) #o/$(MODE)/third_party/ctags/ctags.com
@$(RM) $@
@o/$(MODE)/third_party/ctags/ctags.com $(TAGSFLAGS) -L $< -o $@
HTAGS: private .UNSANDBOXED = 1
HTAGS: o/$(MODE)/hdrs-old.txt $(HDRS) o/$(MODE)/third_party/ctags/ctags.com
HTAGS: o/$(MODE)/hdrs-old.txt $(HDRS) #o/$(MODE)/third_party/ctags/ctags.com
@$(RM) $@
@build/htags o/$(MODE)/third_party/ctags/ctags.com -L $< -o $@

1
libc/aarch64/aarch64.mk
View file

@ -10,5 +10,4 @@ o/$(MODE)/libc/aarch64/start.o: \
o/$(MODE)/libc/aarch64: \
o/$(MODE)/libc/aarch64/crt.o \
o/$(MODE)/libc/aarch64/fenv.o \
o/$(MODE)/libc/aarch64/start.o

67
libc/aarch64/fenv.S
View file

@ -1,67 +0,0 @@
.global fegetround
.type fegetround,%function
fegetround:
mrs x0,fpcr
and w0,w0,#0xc00000
ret
.global __fesetround
.hidden __fesetround
.type __fesetround,%function
__fesetround:
mrs x1,fpcr
bic w1,w1,#0xc00000
orr w1,w1,w0
msr fpcr,x1
mov w0,#0
ret
.global fetestexcept
.type fetestexcept,%function
fetestexcept:
and w0,w0,#0x1f
mrs x1,fpsr
and w0,w0,w1
ret
.global feclearexcept
.type feclearexcept,%function
feclearexcept:
and w0,w0,#0x1f
mrs x1,fpsr
bic w1,w1,w0
msr fpsr,x1
mov w0,#0
ret
.global feraiseexcept
.type feraiseexcept,%function
feraiseexcept:
and w0,w0,#0x1f
mrs x1,fpsr
orr w1,w1,w0
msr fpsr,x1
mov w0,#0
ret
.global fegetenv
.type fegetenv,%function
fegetenv:
mrs x1,fpcr
mrs x2,fpsr
stp w1,w2,[x0]
mov w0,#0
ret
.global fesetenv
.type fesetenv,%function
fesetenv:
mov x1,#0
mov x2,#0
cmn x0,#1
b.eq 1f
ldp w1,w2,[x0]
1: msr fpcr,x1
msr fpsr,x2
mov w0,#0
ret

60
libc/intrin/fenv.S
View file

@ -28,6 +28,7 @@
#include "libc/macros.internal.h"
feclearexcept:
#ifdef __x86_64__
# maintain exceptions in the sse mxcsr, clear x87 exceptions
mov %edi,%ecx
and $0x3f,%ecx
@ -45,18 +46,36 @@ feclearexcept:
ldmxcsr -8(%rsp)
1: xor %eax,%eax
ret
#elif defined(__aarch64__)
and w0,w0,#0x1f
mrs x1,fpsr
bic w1,w1,w0
msr fpsr,x1
mov w0,#0
ret
#endif
.endfn feclearexcept,globl
feraiseexcept:
#ifdef __x86_64__
and $0x3f,%edi
stmxcsr -8(%rsp)
or %edi,-8(%rsp)
ldmxcsr -8(%rsp)
xor %eax,%eax
ret
#elif defined(__aarch64__)
and w0,w0,#0x1f
mrs x1,fpsr
orr w1,w1,w0
msr fpsr,x1
mov w0,#0
ret
#endif
.endfn feraiseexcept,globl
__fesetround:
#ifdef __x86_64__
push %rax
xor %eax,%eax
mov %edi,%ecx
@ -71,25 +90,48 @@ __fesetround:
ldmxcsr (%rsp)
pop %rcx
ret
#elif defined(__aarch64__)
mrs x1,fpcr
bic w1,w1,#0xc00000
orr w1,w1,w0
msr fpcr,x1
mov w0,#0
ret
#endif
.endfn __fesetround,globl,hidden
fegetround:
#ifdef __x86_64__
push %rax
stmxcsr (%rsp)
pop %rax
shr $3,%eax
and $0xc00,%eax
ret
#elif defined(__aarch64__)
mrs x0,fpcr
and w0,w0,#0xc00000
ret
#endif
.endfn fegetround,globl
fegetenv:
#ifdef __x86_64__
xor %eax,%eax
fnstenv (%rdi)
stmxcsr 28(%rdi)
ret
#elif defined(__aarch64__)
mrs x1,fpcr
mrs x2,fpsr
stp w1,w2,[x0]
mov w0,#0
ret
#endif
.endfn fegetenv,globl
fesetenv:
#ifdef __x86_64__
xor %eax,%eax
inc %rdi
jz 1f
@ -105,9 +147,21 @@ fesetenv:
ldmxcsr (%rsp)
add $40,%rsp
ret
#elif defined(__aarch64__)
mov x1,#0
mov x2,#0
cmn x0,#1
b.eq 1f
ldp w1,w2,[x0]
1: msr fpcr,x1
msr fpsr,x2
mov w0,#0
ret
#endif
.endfn fesetenv,globl
fetestexcept:
#ifdef __x86_64__
and $0x3f,%edi
push %rax
stmxcsr (%rsp)
@ -116,4 +170,10 @@ fetestexcept:
or %esi,%eax
and %edi,%eax
ret
#elif defined(__aarch64__)
and w0,w0,#0x1f
mrs x1,fpsr
and w0,w0,w1
ret
#endif
.endfn fetestexcept,globl

2
libc/intrin/intrin.mk
View file

@ -187,6 +187,8 @@ o/$(MODE)/libc/intrin/memmove.o: private \
-fpie
# these assembly files are safe to build on aarch64
o/$(MODE)/libc/intrin/fenv.o: libc/intrin/fenv.S
@$(COMPILE) -AOBJECTIFY.S $(OBJECTIFY.S) $(OUTPUT_OPTION) -c $<
o/$(MODE)/libc/intrin/kclocknames.o: libc/intrin/kclocknames.S
@$(COMPILE) -AOBJECTIFY.S $(OBJECTIFY.S) $(OUTPUT_OPTION) -c $<
o/$(MODE)/libc/intrin/kdos2errno.o: libc/intrin/kdos2errno.S

1
libc/math.h
View file

@ -255,7 +255,6 @@ long double nearbyintl(long double);
long double nextafterl(long double, long double);
long double nexttowardl(long double, long double);
long double pow10l(long double);
long double powil(long double, int);
long double powl(long double, long double);
long double remainderl(long double, long double);
long double rintl(long double);

2
libc/tinymath/__cexp.c
View file

@ -33,7 +33,7 @@ asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
// clang-format off
/* origin: FreeBSD /usr/src/lib/msun/src/k_exp.c */
/*-

2
libc/tinymath/__cexpf.c
View file

@ -33,7 +33,7 @@ asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
// clang-format off
/* origin: FreeBSD /usr/src/lib/msun/src/k_expf.c */
/*-

26
libc/tinymath/__math_invalidl.c Normal file
View file

@ -0,0 +1,26 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/internal.h"
#if LDBL_MANT_DIG != DBL_MANT_DIG
long double __math_invalidl(long double x) {
return (x - x) / (x - x);
}
#endif

2
libc/tinymath/acos.c
View file

@ -31,7 +31,7 @@ asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
// clang-format off
/* origin: FreeBSD /usr/src/lib/msun/src/e_acos.c */
/*

3
libc/tinymath/acosf.c
View file

@ -27,7 +27,6 @@
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
#ifndef TINY
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
@ -107,5 +106,3 @@ float acosf(float x)
w = R(z)*s+c;
return 2*(df+w);
}
#endif /* TINY */

17
libc/tinymath/acoshl.c
View file

@ -34,9 +34,13 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
/**
* Returns inverse hyperbolic cosine of 𝑥.
*/
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double acoshl(long double x)
{
return acosh(x);
}
#elif LDBL_MANT_DIG == 64 && LDBL_MAX_EXP == 16384
/* acosh(x) = log(x + sqrt(x*x-1)) */
long double acoshl(long double x)
{
union ldshape u = {x};
@ -50,3 +54,10 @@ long double acoshl(long double x)
return logl(2*x - 1/(x+sqrtl(x*x-1)));
return logl(x) + 0.693147180559945309417232121458176568L;
}
#elif LDBL_MANT_DIG == 113 && LDBL_MAX_EXP == 16384
// TODO: broken implementation to make things compile
long double acoshl(long double x)
{
return acosh(x);
}
#endif

43
libc/tinymath/acosl.S
View file

@ -1,43 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 sw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns arc cosine of 𝑥.
//
// @param 𝑥 is an 80-bit long double passed on stack in 16-bytes
// @return result of computation on FPU stack in %st
// @define atan2(fabs(sqrt((1-𝑥)*(1+𝑥))),𝑥)
// @domain -1 ≤ 𝑥 ≤ 1
acosl: pushq %rbp
mov %rsp,%rbp
.profilable
fldt 16(%rbp)
fld1
fld %st
fsub %st(2)
fxch
fadd %st(2)
fmulp
fsqrt
fabs
fxch
fpatan
leave
ret
.endfn acosl,globl

111
libc/tinymath/acosl.c Normal file
View file

@ -0,0 +1,111 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/invtrigl.internal.h"
#include "libc/tinymath/ldshape.internal.h"
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
/* origin: FreeBSD /usr/src/lib/msun/src/e_acosl.c */
/*
* ====================================================
* 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.
* ====================================================
*/
/*
* See comments in acos.c.
* Converted to long double by David Schultz <das@FreeBSD.ORG>.
*/
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double acosl(long double x)
{
return acos(x);
}
#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
#if LDBL_MANT_DIG == 64
#define CLEARBOTTOM(u) (u.i.m &= -1ULL << 32)
#elif LDBL_MANT_DIG == 113
#define CLEARBOTTOM(u) (u.i.lo = 0)
#endif
/**
* Returns arc cosine of 𝑥.
*
* @define atan2(fabs(sqrt((1-𝑥)*(1+𝑥))),𝑥)
* @domain -1 ≤ 𝑥 ≤ 1
*/
long double acosl(long double x)
{
union ldshape u = {x};
long double z, s, c, f;
uint16_t e = u.i.se & 0x7fff;
/* |x| >= 1 or nan */
if (e >= 0x3fff) {
if (x == 1)
return 0;
if (x == -1)
return 2*pio2_hi + 0x1p-120f;
return 0/(x-x);
}
/* |x| < 0.5 */
if (e < 0x3fff - 1) {
if (e < 0x3fff - LDBL_MANT_DIG - 1)
return pio2_hi + 0x1p-120f;
return pio2_hi - (__invtrigl_R(x*x)*x - pio2_lo + x);
}
/* x < -0.5 */
if (u.i.se >> 15) {
z = (1 + x)*0.5;
s = sqrtl(z);
return 2*(pio2_hi - (__invtrigl_R(z)*s - pio2_lo + s));
}
/* x > 0.5 */
z = (1 - x)*0.5;
s = sqrtl(z);
u.f = s;
CLEARBOTTOM(u);
f = u.f;
c = (z - f*f)/(s + f);
return 2*(__invtrigl_R(z)*s + c + f);
}
#endif

6
libc/tinymath/asinf.c
View file

@ -27,8 +27,10 @@
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
#ifndef TINY
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
@ -98,5 +100,3 @@ float asinf(float x)
return -x;
return x;
}
#endif /* TINY */

18
libc/tinymath/asinhl.c
View file

@ -27,6 +27,7 @@
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/feval.internal.h"
#include "libc/tinymath/internal.h"
#include "libc/tinymath/ldshape.internal.h"
asm(".ident\t\"\\n\\n\
@ -35,6 +36,13 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double asinhl(long double x)
{
return asinh(x);
}
#elif LDBL_MANT_DIG == 64 && LDBL_MAX_EXP == 16384
/**
* Returns inverse hyperbolic sine of 𝑥.
* @define asinh(x) = sign(x)*log(|x|+sqrt(x*x+1)) ~= x - x^3/6 + o(x^5)
@ -60,7 +68,15 @@ long double asinhl(long double x)
x = log1pl(x + x*x/(sqrtl(x*x+1)+1));
} else {
/* |x| < 0x1p-32, raise inexact if x!=0 */
fevall(x + 0x1p120f);
FORCE_EVAL(x + 0x1p120f);
}
return s ? -x : x;
}
#elif LDBL_MANT_DIG == 113 && LDBL_MAX_EXP == 16384
// TODO: broken implementation to make things compile
long double asinhl(long double x)
{
return asinh(x);
}
#endif

116
libc/tinymath/asinl.c Normal file
View file

@ -0,0 +1,116 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/internal.h"
#include "libc/tinymath/invtrigl.internal.h"
#include "libc/tinymath/ldshape.internal.h"
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
/* origin: FreeBSD /usr/src/lib/msun/src/e_asinl.c */
/*
* ====================================================
* 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.
* ====================================================
*/
/*
* See comments in asin.c.
* Converted to long double by David Schultz <das@FreeBSD.ORG>.
*/
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double asinl(long double x)
{
return asin(x);
}
#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
#if LDBL_MANT_DIG == 64
#define CLOSETO1(u) (u.i.m>>56 >= 0xf7)
#define CLEARBOTTOM(u) (u.i.m &= -1ULL << 32)
#elif LDBL_MANT_DIG == 113
#define CLOSETO1(u) (u.i.top >= 0xee00)
#define CLEARBOTTOM(u) (u.i.lo = 0)
#endif
/**
* Returns arc sine of 𝑥.
*
* @define atan2(𝑥,sqrt((1-𝑥)*(1+𝑥)))
* @domain -1 ≤ 𝑥 ≤ 1
*/
long double asinl(long double x)
{
union ldshape u = {x};
long double z, r, s;
uint16_t e = u.i.se & 0x7fff;
int sign = u.i.se >> 15;
if (e >= 0x3fff) { /* |x| >= 1 or nan */
/* asin(+-1)=+-pi/2 with inexact */
if (x == 1 || x == -1)
return x*pio2_hi + 0x1p-120f;
return 0/(x-x);
}
if (e < 0x3fff - 1) { /* |x| < 0.5 */
if (e < 0x3fff - (LDBL_MANT_DIG+1)/2) {
/* return x with inexact if x!=0 */
FORCE_EVAL(x + 0x1p120f);
return x;
}
return x + x*__invtrigl_R(x*x);
}
/* 1 > |x| >= 0.5 */
z = (1.0 - fabsl(x))*0.5;
s = sqrtl(z);
r = __invtrigl_R(z);
if (CLOSETO1(u)) {
x = pio2_hi - (2*(s+s*r)-pio2_lo);
} else {
long double f, c;
u.f = s;
CLEARBOTTOM(u);
f = u.f;
c = (z - f*f)/(s + f);
x = 0.5*pio2_hi-(2*s*r - (pio2_lo-2*c) - (0.5*pio2_hi-2*f));
}
return sign ? -x : x;
}
#endif

6
libc/tinymath/atan.c
View file

@ -27,8 +27,10 @@
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
#ifndef TINY
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
@ -152,5 +154,3 @@ double atan(double x)
z = atanhi[id] - (x*(s1+s2) - atanlo[id] - x);
return sign ? -z : z;
}
#endif /* TINY */

34
libc/tinymath/atan2l.S
View file

@ -1,34 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 sw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns arc tangent of 𝑦/𝑥.
//
// @param 𝑥 is an 80-bit long double passed on stack in 16-bytes
// @param 𝑦 is an 80-bit long double passed on stack in 16-bytes
// @return result of computation on FPU stack in %st
atan2l: push %rbp
mov %rsp,%rbp
.profilable
fldt 16(%rbp)
fldt 32(%rbp)
fpatan
pop %rbp
ret
.endfn atan2l,globl

139
libc/tinymath/atan2l.c Normal file
View file

@ -0,0 +1,139 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/internal.h"
#include "libc/tinymath/invtrigl.internal.h"
#include "libc/tinymath/ldshape.internal.h"
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
/* origin: FreeBSD /usr/src/lib/msun/src/e_atan2l.c */
/*
* ====================================================
* 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.
* ====================================================
*
*/
/*
* See comments in atan2.c.
* Converted to long double by David Schultz <das@FreeBSD.ORG>.
*/
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double atan2l(long double y, long double x)
{
return atan2(y, x);
}
#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
/**
* Returns arc tangent of 𝑦/𝑥.
*/
long double atan2l(long double y, long double x)
{
#ifdef __x86__
long double res;
asm("fpatan"
: "=t" (res)
: "0"(x), "u"(y)
: "st(1)");
return res;
#else
union ldshape ux, uy;
long double z;
int m, ex, ey;
if (isnan(x) || isnan(y))
return x+y;
if (x == 1)
return atanl(y);
ux.f = x;
uy.f = y;
ex = ux.i.se & 0x7fff;
ey = uy.i.se & 0x7fff;
m = 2*(ux.i.se>>15) | uy.i.se>>15;
if (y == 0) {
switch(m) {
case 0:
case 1: return y; /* atan(+-0,+anything)=+-0 */
case 2: return 2*pio2_hi; /* atan(+0,-anything) = pi */
case 3: return -2*pio2_hi; /* atan(-0,-anything) =-pi */
}
}
if (x == 0)
return m&1 ? -pio2_hi : pio2_hi;
if (ex == 0x7fff) {
if (ey == 0x7fff) {
switch(m) {
case 0: return pio2_hi/2; /* atan(+INF,+INF) */
case 1: return -pio2_hi/2; /* atan(-INF,+INF) */
case 2: return 1.5*pio2_hi; /* atan(+INF,-INF) */
case 3: return -1.5*pio2_hi; /* atan(-INF,-INF) */
}
} else {
switch(m) {
case 0: return 0.0; /* atan(+...,+INF) */
case 1: return -0.0; /* atan(-...,+INF) */
case 2: return 2*pio2_hi; /* atan(+...,-INF) */
case 3: return -2*pio2_hi; /* atan(-...,-INF) */
}
}
}
if (ex+120 < ey || ey == 0x7fff)
return m&1 ? -pio2_hi : pio2_hi;
/* z = atan(|y/x|) without spurious underflow */
if ((m&2) && ey+120 < ex) /* |y/x| < 0x1p-120, x<0 */
z = 0.0;
else
z = atanl(fabsl(y/x));
switch (m) {
case 0: return z; /* atan(+,+) */
case 1: return -z; /* atan(-,+) */
case 2: return 2*pio2_hi-(z-2*pio2_lo); /* atan(+,-) */
default: /* case 3 */
return (z-2*pio2_lo)-2*pio2_hi; /* atan(-,-) */
}
#endif /* __x86__ */
}
#endif

6
libc/tinymath/atanf.c
View file

@ -27,8 +27,10 @@
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
#ifndef TINY
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
@ -130,5 +132,3 @@ float atanf(float x)
z = atanhi[id] - ((x*(s1+s2) - atanlo[id]) - x);
return sign ? -z : z;
}
#endif /* TINY */

9
libc/tinymath/atanhl.c
View file

@ -35,6 +35,13 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double atanhl(long double x)
{
return atanh(x);
}
#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
/**
* Returns inverse hyperbolic tangent of 𝑥.
* @define x ? log1p(2 * x / (1 - x)) / 2 : x
@ -64,3 +71,5 @@ long double atanhl(long double x)
}
return s ? -x : x;
}
#endif

38
libc/tinymath/atanl.S
View file

@ -1,38 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 sw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns arc tangent of 𝑥.
//
// 1 3 1 5 1 7 1 9 1 11
// atan(𝑥) = 𝑥 - - 𝑥 + - 𝑥 - - 𝑥 + - 𝑥 - -- 𝑥 ...
// 3 5 7 9 11
//
// @param 𝑥 is an 80-bit long double passed on stack in 16-bytes
// @return result of computation on FPU stack in %st
// @define atan(𝑥) = Σₙ₌₀₋∞ 2²ⁿ(𝑛!)²/(𝟸𝑛+𝟷)!(𝑥²ⁿ⁺¹/(𝑥²+𝟷)ⁿ⁺¹)
atanl: push %rbp
mov %rsp,%rbp
.profilable
fldt 16(%rbp)
fld1
fpatan
pop %rbp
ret
.endfn atanl,globl

234
libc/tinymath/atanl.c Normal file
View file

@ -0,0 +1,234 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/internal.h"
#include "libc/tinymath/ldshape.internal.h"
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
/* origin: FreeBSD /usr/src/lib/msun/src/s_atanl.c */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
/*
* See comments in atan.c.
* Converted to long double by David Schultz <das@FreeBSD.ORG>.
*/
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double atanl(long double x)
{
return atan(x);
}
#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
#if LDBL_MANT_DIG == 64
#define EXPMAN(u) ((u.i.se & 0x7fff)<<8 | (u.i.m>>55 & 0xff))
static const long double atanhi[] = {
4.63647609000806116202e-01L,
7.85398163397448309628e-01L,
9.82793723247329067960e-01L,
1.57079632679489661926e+00L,
};
static const long double atanlo[] = {
1.18469937025062860669e-20L,
-1.25413940316708300586e-20L,
2.55232234165405176172e-20L,
-2.50827880633416601173e-20L,
};
static const long double aT[] = {
3.33333333333333333017e-01L,
-1.99999999999999632011e-01L,
1.42857142857046531280e-01L,
-1.11111111100562372733e-01L,
9.09090902935647302252e-02L,
-7.69230552476207730353e-02L,
6.66661718042406260546e-02L,
-5.88158892835030888692e-02L,
5.25499891539726639379e-02L,
-4.70119845393155721494e-02L,
4.03539201366454414072e-02L,
-2.91303858419364158725e-02L,
1.24822046299269234080e-02L,
};
static long double T_even(long double x)
{
return aT[0] + x * (aT[2] + x * (aT[4] + x * (aT[6] +
x * (aT[8] + x * (aT[10] + x * aT[12])))));
}
static long double T_odd(long double x)
{
return aT[1] + x * (aT[3] + x * (aT[5] + x * (aT[7] +
x * (aT[9] + x * aT[11]))));
}
#elif LDBL_MANT_DIG == 113
#define EXPMAN(u) ((u.i.se & 0x7fff)<<8 | u.i.top>>8)
static const long double atanhi[] = {
4.63647609000806116214256231461214397e-01L,
7.85398163397448309615660845819875699e-01L,
9.82793723247329067985710611014666038e-01L,
1.57079632679489661923132169163975140e+00L,
};
static const long double atanlo[] = {
4.89509642257333492668618435220297706e-36L,
2.16795253253094525619926100651083806e-35L,
-2.31288434538183565909319952098066272e-35L,
4.33590506506189051239852201302167613e-35L,
};
static const long double aT[] = {
3.33333333333333333333333333333333125e-01L,
-1.99999999999999999999999999999180430e-01L,
1.42857142857142857142857142125269827e-01L,
-1.11111111111111111111110834490810169e-01L,
9.09090909090909090908522355708623681e-02L,
-7.69230769230769230696553844935357021e-02L,
6.66666666666666660390096773046256096e-02L,
-5.88235294117646671706582985209643694e-02L,
5.26315789473666478515847092020327506e-02L,
-4.76190476189855517021024424991436144e-02L,
4.34782608678695085948531993458097026e-02L,
-3.99999999632663469330634215991142368e-02L,
3.70370363987423702891250829918659723e-02L,
-3.44827496515048090726669907612335954e-02L,
3.22579620681420149871973710852268528e-02L,
-3.03020767654269261041647570626778067e-02L,
2.85641979882534783223403715930946138e-02L,
-2.69824879726738568189929461383741323e-02L,
2.54194698498808542954187110873675769e-02L,
-2.35083879708189059926183138130183215e-02L,
2.04832358998165364349957325067131428e-02L,
-1.54489555488544397858507248612362957e-02L,
8.64492360989278761493037861575248038e-03L,
-2.58521121597609872727919154569765469e-03L,
};
static long double T_even(long double x)
{
return (aT[0] + x * (aT[2] + x * (aT[4] + x * (aT[6] + x * (aT[8] +
x * (aT[10] + x * (aT[12] + x * (aT[14] + x * (aT[16] +
x * (aT[18] + x * (aT[20] + x * aT[22])))))))))));
}
static long double T_odd(long double x)
{
return (aT[1] + x * (aT[3] + x * (aT[5] + x * (aT[7] + x * (aT[9] +
x * (aT[11] + x * (aT[13] + x * (aT[15] + x * (aT[17] +
x * (aT[19] + x * (aT[21] + x * aT[23])))))))))));
}
#endif
/**
* Returns arc tangent of 𝑥.
*
* 1 3 1 5 1 7 1 9 1 11
* atan(𝑥) = 𝑥 - - 𝑥 + - 𝑥 - - 𝑥 + - 𝑥 - -- 𝑥 ...
* 3 5 7 9 11
*
* @param 𝑥 is an 80-bit long double passed on stack in 16-bytes
* @return result of computation on FPU stack in %st
* @define atan(𝑥) = Σₙ₌₀₋∞ 2²ⁿ(𝑛!)²/(𝟸𝑛+𝟷)!(𝑥²ⁿ⁺¹/(𝑥²+𝟷)ⁿ⁺¹)
*/
long double atanl(long double x)
{
union ldshape u = {x};
long double w, s1, s2, z;
int id;
unsigned e = u.i.se & 0x7fff;
unsigned sign = u.i.se >> 15;
unsigned expman;
if (e >= 0x3fff + LDBL_MANT_DIG + 1) { /* if |x| is large, atan(x)~=pi/2 */
if (isnan(x))
return x;
return sign ? -atanhi[3] : atanhi[3];
}
/* Extract the exponent and the first few bits of the mantissa. */
expman = EXPMAN(u);
if (expman < ((0x3fff - 2) << 8) + 0xc0) { /* |x| < 0.4375 */
if (e < 0x3fff - (LDBL_MANT_DIG+1)/2) { /* if |x| is small, atanl(x)~=x */
/* raise underflow if subnormal */
if (e == 0)
FORCE_EVAL((float)x);
return x;
}
id = -1;
} else {
x = fabsl(x);
if (expman < (0x3fff << 8) + 0x30) { /* |x| < 1.1875 */
if (expman < ((0x3fff - 1) << 8) + 0x60) { /* 7/16 <= |x| < 11/16 */
id = 0;
x = (2.0*x-1.0)/(2.0+x);
} else { /* 11/16 <= |x| < 19/16 */
id = 1;
x = (x-1.0)/(x+1.0);
}
} else {
if (expman < ((0x3fff + 1) << 8) + 0x38) { /* |x| < 2.4375 */
id = 2;
x = (x-1.5)/(1.0+1.5*x);
} else { /* 2.4375 <= |x| */
id = 3;
x = -1.0/x;
}
}
}
/* end of argument reduction */
z = x*x;
w = z*z;
/* break sum aT[i]z**(i+1) into odd and even poly */
s1 = z*T_even(w);
s2 = w*T_odd(w);
if (id < 0)
return x - x*(s1+s2);
z = atanhi[id] - ((x*(s1+s2) - atanlo[id]) - x);
return sign ? -z : z;
}
#endif

47
libc/tinymath/c2rangr.S
View file

@ -1,47 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/runtime/pc.internal.h"
#include "libc/macros.internal.h"
// Computes transcedental trigonometry op w/ reactive scaling.
//
// @param %rdx points to op function
// @param everything else delegates
// @clob %ax
// @see sin,cos,tan
c2rangr:push %rbp
mov %rsp,%rbp
.profilable
call *%rdx
fstsw %ax
test $FPU_C2>>8,%ah
jnz 1f
0: pop %rbp
ret
1: fldpi
fadd %st
fxch
2: fprem1
fnstsw
test $FPU_C2>>8,%ah
jnz 2b
fstp %st(1)
call *%rdx
jmp 0b
.endfn c2rangr,globl,hidden

24
libc/tinymath/cargl.S → libc/tinymath/cabs.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,14 +16,12 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#include "libc/complex.h"
#include "libc/math.h"
cargl: push %rbp
mov %rsp,%rbp
.profilable
fldt 32(%rbp)
fldt 16(%rbp)
fpatan
pop %rbp
ret
.endfn cargl,globl
/**
* Returns absolute value of complex number.
*/
double cabs(double complex z) {
return hypot(creal(z), cimag(z));
}

36
libc/tinymath/cabsf.S
View file

@ -1,36 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
cabsf: push %rbp
mov %rsp,%rbp
.profilable
sub $16,%rsp
movq %xmm0,(%rsp)
movss (%rsp),%xmm0
movss 4(%rsp),%xmm2
movaps %xmm0,%xmm1
mulss %xmm2,%xmm2
mulss %xmm0,%xmm1
movaps %xmm2,%xmm0
addss %xmm1,%xmm0
sqrtss %xmm0,%xmm0
leave
ret
.endfn cabsf,globl

24
libc/tinymath/cimagf.S → libc/tinymath/cabsf.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,14 +16,12 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#include "libc/complex.h"
#include "libc/math.h"
cimagf: push %rbp
mov %rsp,%rbp
.profilable
sub $16,%rsp
movq %xmm0,(%rsp)
movss 4(%rsp),%xmm0
leave
ret
.endfn cimagf,globl
/**
* Returns absolute value of complex number.
*/
float cabsf(float complex z) {
return hypotf(crealf(z), cimagf(z));
}

29
libc/tinymath/atanf-tiny.S → libc/tinymath/cabsl.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,15 +16,16 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#ifdef TINY
#include "libc/complex.h"
#include "libc/math.h"
// Returns arc tangent of 𝑥.
//
// @param 𝑥 is float scalar in low quarter of %xmm0
// @return float scalar in low quarter of %xmm0
atanf: ezlea atanl,ax
jmp _f2ld2
.endfn atanf,globl
#endif /* TINY */
/**
* Returns absolute value of complex number.
*/
long double cabsl(long double complex z) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
return cabs(z);
#else
return hypotl(creall(z), cimagl(z));
#endif
}

2
libc/tinymath/cacosf.c
View file

@ -35,8 +35,6 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
// FIXME
float complex cacosf(float complex z)

4
libc/tinymath/cacosh.c
View file

@ -35,14 +35,10 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
/* acosh(z) = i acos(z) */
double complex cacosh(double complex z)
{
int zineg = signbit(cimag(z));
z = cacos(z);
if (zineg) return CMPLX(cimag(z), -creal(z));
else return CMPLX(-cimag(z), creal(z));

3
libc/tinymath/cacoshf.c
View file

@ -35,12 +35,9 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
float complex cacoshf(float complex z)
{
int zineg = signbit(cimagf(z));
z = cacosf(z);
if (zineg) return CMPLXF(cimagf(z), -crealf(z));
else return CMPLXF(-cimagf(z), crealf(z));

12
libc/tinymath/cacoshl.c
View file

@ -35,20 +35,14 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
long double complex cacoshl(long double complex z)
{
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double complex cacoshl(long double complex z)
{
return cacosh(z);
}
#else
long double complex cacoshl(long double complex z)
{
int zineg = signbit(cimagl(z));
z = cacosl(z);
if (zineg) return CMPLXL(cimagl(z), -creall(z));
else return CMPLXL(-cimagl(z), creall(z));
}
#endif
}

11
libc/tinymath/cacosl.c
View file

@ -26,7 +26,6 @@
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
@ -35,19 +34,13 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
long double complex cacosl(long double complex z) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double complex cacosl(long double complex z)
{
return cacos(z);
}
#else
// FIXME
#define PI_2 1.57079632679489661923132169163975144L
long double complex cacosl(long double complex z)
{
z = casinl(z);
return CMPLXL(PI_2 - creall(z), -cimagl(z));
}
#endif
}

35
libc/tinymath/carg.S
View file

@ -1,35 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
carg: push %rbp
mov %rsp,%rbp
.profilable
sub $16,%rsp
movsd %xmm0,-16(%rbp)
fldl -16(%rbp)
movsd %xmm1,-16(%rbp)
fldl -16(%rbp)
fxch
fpatan
fstpl -16(%rbp)
movsd -16(%rbp),%xmm0
leave
ret
.endfn carg,globl

15
libc/tinymath/cabs.S → libc/tinymath/carg.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,8 +16,9 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#include "libc/complex.h"
#include "libc/math.h"
// Returns absolute value of complex number.
cabs: jmp hypot
.endfn cabs,globl
double carg(double complex z) {
return atan2(cimag(z), creal(z));
}

24
libc/tinymath/cargf.c Normal file
View file

@ -0,0 +1,24 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
#include "libc/math.h"
float cargf(float complex z) {
return atan2f(cimagf(z), crealf(z));
}

28
libc/tinymath/cargl.c Normal file
View file

@ -0,0 +1,28 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
#include "libc/math.h"
long double cargl(long double complex z) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
return carg(z);
#else
return atan2l(cimagl(z), creall(z));
#endif
}

11
libc/tinymath/casinl.c
View file

@ -35,24 +35,17 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
long double complex casinl(long double complex z) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double complex casinl(long double complex z)
{
return casin(z);
}
#else
// FIXME
long double complex casinl(long double complex z)
{
long double complex w;
long double x, y;
x = creall(z);
y = cimagl(z);
w = CMPLXL(1.0 - (x - y)*(x + y), -2.0*x*y);
long double complex r = clogl(CMPLXL(-y, x) + csqrtl(w));
return CMPLXL(cimagl(r), -creall(r));
}
#endif
}

3
libc/tinymath/catan.c
View file

@ -29,6 +29,9 @@
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
OpenBSD libm (MIT License)\\n\
Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net>\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");

3
libc/tinymath/catanf.c
View file

@ -29,6 +29,9 @@
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
OpenBSD libm (MIT License)\\n\
Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net>\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");

10
libc/tinymath/catanhl.c
View file

@ -35,17 +35,11 @@ Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
long double complex catanhl(long double complex z) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double complex catanhl(long double complex z)
{
return catanh(z);
}
#else
long double complex catanhl(long double complex z)
{
z = catanl(CMPLXL(-cimagl(z), creall(z)));
return CMPLXL(cimagl(z), -creall(z));
}
#endif
}

7
libc/tinymath/catanl.c
View file

@ -29,12 +29,14 @@
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
OpenBSD libm (MIT License)\\n\
Copyright (c) 2008 Stephen L. Moshier <steve@moshier.net>\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
// clang-format off
/* origin: OpenBSD /usr/src/lib/libm/src/s_catanl.c */
/*
@ -146,4 +148,5 @@ long double complex catanl(long double complex z)
w = CMPLXF(w, 0.25L * logl(a));
return w;
}
#endif

138
libc/tinymath/cbrt.c Normal file
View file

@ -0,0 +1,138 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2020 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
// clang-format off
/* origin: FreeBSD /usr/src/lib/msun/src/s_cbrt.c */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
* Optimized by Bruce D. Evans.
*/
static const uint32_t
B1 = 715094163, /* B1 = (1023-1023/3-0.03306235651)*2**20 */
B2 = 696219795; /* B2 = (1023-1023/3-54/3-0.03306235651)*2**20 */
/* |1/cbrt(x) - p(x)| < 2**-23.5 (~[-7.93e-8, 7.929e-8]). */
static const double
P0 = 1.87595182427177009643, /* 0x3ffe03e6, 0x0f61e692 */
P1 = -1.88497979543377169875, /* 0xbffe28e0, 0x92f02420 */
P2 = 1.621429720105354466140, /* 0x3ff9f160, 0x4a49d6c2 */
P3 = -0.758397934778766047437, /* 0xbfe844cb, 0xbee751d9 */
P4 = 0.145996192886612446982; /* 0x3fc2b000, 0xd4e4edd7 */
/**
* Returns cube root of 𝑥.
*/
double cbrt(double x)
{
union {double f; uint64_t i;} u = {x};
double_t r,s,t,w;
uint32_t hx = u.i>>32 & 0x7fffffff;
if (hx >= 0x7ff00000) /* cbrt(NaN,INF) is itself */
return x+x;
/*
* Rough cbrt to 5 bits:
* cbrt(2**e*(1+m) ~= 2**(e/3)*(1+(e%3+m)/3)
* where e is integral and >= 0, m is real and in [0, 1), and "/" and
* "%" are integer division and modulus with rounding towards minus
* infinity. The RHS is always >= the LHS and has a maximum relative
* error of about 1 in 16. Adding a bias of -0.03306235651 to the
* (e%3+m)/3 term reduces the error to about 1 in 32. With the IEEE
* floating point representation, for finite positive normal values,
* ordinary integer divison of the value in bits magically gives
* almost exactly the RHS of the above provided we first subtract the
* exponent bias (1023 for doubles) and later add it back. We do the
* subtraction virtually to keep e >= 0 so that ordinary integer
* division rounds towards minus infinity; this is also efficient.
*/
if (hx < 0x00100000) { /* zero or subnormal? */
u.f = x*0x1p54;
hx = u.i>>32 & 0x7fffffff;
if (hx == 0)
return x; /* cbrt(0) is itself */
hx = hx/3 + B2;
} else
hx = hx/3 + B1;
u.i &= 1ULL<<63;
u.i |= (uint64_t)hx << 32;
t = u.f;
/*
* New cbrt to 23 bits:
* cbrt(x) = t*cbrt(x/t**3) ~= t*P(t**3/x)
* where P(r) is a polynomial of degree 4 that approximates 1/cbrt(r)
* to within 2**-23.5 when |r - 1| < 1/10. The rough approximation
* has produced t such than |t/cbrt(x) - 1| ~< 1/32, and cubing this
* gives us bounds for r = t**3/x.
*
* Try to optimize for parallel evaluation as in __tanf.c.
*/
r = (t*t)*(t/x);
t = t*((P0+r*(P1+r*P2))+((r*r)*r)*(P3+r*P4));
/*
* Round t away from zero to 23 bits (sloppily except for ensuring that
* the result is larger in magnitude than cbrt(x) but not much more than
* 2 23-bit ulps larger). With rounding towards zero, the error bound
* would be ~5/6 instead of ~4/6. With a maximum error of 2 23-bit ulps
* in the rounded t, the infinite-precision error in the Newton
* approximation barely affects third digit in the final error
* 0.667; the error in the rounded t can be up to about 3 23-bit ulps
* before the final error is larger than 0.667 ulps.
*/
u.f = t;
u.i = (u.i + 0x80000000) & 0xffffffffc0000000ULL;
t = u.f;
/* one step Newton iteration to 53 bits with error < 0.667 ulps */
s = t*t; /* t*t is exact */
r = x/s; /* error <= 0.5 ulps; |r| < |t| */
w = t+t; /* t+t is exact */
r = (r-t)/(w+r); /* r-t is exact; w+r ~= 3*t */
t = t+t*r; /* error <= 0.5 + 0.5/3 + epsilon */
return t;
}

106
libc/tinymath/cbrtc.c
View file

@ -1,106 +0,0 @@
/* origin: FreeBSD /usr/src/lib/msun/src/s_cbrt.c */
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
* Optimized by Bruce D. Evans.
*/
#include "libc/math.h"
asm(".ident\t\"\\n\\n\
fdlibm\\n\
Copyright 1993 Sun Microsystems, Inc.\\n\
Developed at SunPro, a Sun Microsystems, Inc. business.\"");
static const uint32_t
B1 = 715094163, /* B1 = (1023-1023/3-0.03306235651)*2**20 */
B2 = 696219795; /* B2 = (1023-1023/3-54/3-0.03306235651)*2**20 */
/* |1/cbrt(x) - p(x)| < 2**-23.5 (~[-7.93e-8, 7.929e-8]). */
static const double P0 = 1.87595182427177009643, /* 0x3ffe03e6, 0x0f61e692 */
P1 = -1.88497979543377169875, /* 0xbffe28e0, 0x92f02420 */
P2 = 1.621429720105354466140, /* 0x3ff9f160, 0x4a49d6c2 */
P3 = -0.758397934778766047437, /* 0xbfe844cb, 0xbee751d9 */
P4 = 0.145996192886612446982; /* 0x3fc2b000, 0xd4e4edd7 */
/**
* Returns cube root of 𝑥.
*/
double __cbrt(double x) {
union {
double f;
uint64_t i;
} u = {x};
double_t r, s, t, w;
uint32_t hx = u.i >> 32 & 0x7fffffff;
if (hx >= 0x7ff00000) /* cbrt(NaN,INF) is itself */
return x + x;
/*
* Rough cbrt to 5 bits:
* cbrt(2**e*(1+m) ~= 2**(e/3)*(1+(e%3+m)/3)
* where e is integral and >= 0, m is real and in [0, 1), and "/" and
* "%" are integer division and modulus with rounding towards minus
* infinity. The RHS is always >= the LHS and has a maximum relative
* error of about 1 in 16. Adding a bias of -0.03306235651 to the
* (e%3+m)/3 term reduces the error to about 1 in 32. With the IEEE
* floating point representation, for finite positive normal values,
* ordinary integer division of the value in bits magically gives
* almost exactly the RHS of the above provided we first subtract the
* exponent bias (1023 for doubles) and later add it back. We do the
* subtraction virtually to keep e >= 0 so that ordinary integer
* division rounds towards minus infinity; this is also efficient.
*/
if (hx < 0x00100000) { /* zero or subnormal? */
u.f = x * 0x1p54;
hx = u.i >> 32 & 0x7fffffff;
if (hx == 0) return x; /* cbrt(0) is itself */
hx = hx / 3 + B2;
} else
hx = hx / 3 + B1;
u.i &= 1ULL << 63;
u.i |= (uint64_t)hx << 32;
t = u.f;
/*
* New cbrt to 23 bits:
* cbrt(x) = t*cbrt(x/t**3) ~= t*P(t**3/x)
* where P(r) is a polynomial of degree 4 that approximates 1/cbrt(r)
* to within 2**-23.5 when |r - 1| < 1/10. The rough approximation
* has produced t such than |t/cbrt(x) - 1| ~< 1/32, and cubing this
* gives us bounds for r = t**3/x.
*
* Try to optimize for parallel evaluation as in __tanf.c.
*/
r = (t * t) * (t / x);
t = t * ((P0 + r * (P1 + r * P2)) + ((r * r) * r) * (P3 + r * P4));
/*
* Round t away from zero to 23 bits (sloppily except for ensuring that
* the result is larger in magnitude than cbrt(x) but not much more than
* 2 23-bit ulps larger). With rounding towards zero, the error bound
* would be ~5/6 instead of ~4/6. With a maximum error of 2 23-bit ulps
* in the rounded t, the infinite-precision error in the Newton
* approximation barely affects third digit in the final error
* 0.667; the error in the rounded t can be up to about 3 23-bit ulps
* before the final error is larger than 0.667 ulps.
*/
u.f = t;
u.i = (u.i + 0x80000000) & 0xffffffffc0000000ULL;
t = u.f;
/* one step Newton iteration to 53 bits with error < 0.667 ulps */
s = t * t; /* t*t is exact */
r = x / s; /* error <= 0.5 ulps; |r| < |t| */
w = t + t; /* t+t is exact */
r = (r - t) / (w + r); /* r-t is exact; w+r ~= 3*t */
t = t + t * r; /* error <= 0.5 + 0.5/3 + epsilon */
return t;
}

32
libc/tinymath/cbrtf.S
View file

@ -1,32 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns cube root of 𝑥.
//
// @param %xmm0 holds binary32 number
// @return %xmm0 holds binary32 result
cbrtf: pushq %rbp
mov %rsp,%rbp
cvtss2sd %xmm0,%xmm0
call __cbrt
cvtsd2ss %xmm0,%xmm0
popq %rbp
ret
.endfn cbrtf,globl

101
libc/tinymath/cbrtf.c Normal file
View file

@ -0,0 +1,101 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2020 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
// clang-format off
/* origin: FreeBSD /usr/src/lib/msun/src/s_cbrtf.c */
/*
* Conversion to float by Ian Lance Taylor, Cygnus Support, ian@cygnus.com.
* Debugged and optimized by Bruce D. Evans.
*/
/*
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*/
static const unsigned
B1 = 709958130, /* B1 = (127-127.0/3-0.03306235651)*2**23 */
B2 = 642849266; /* B2 = (127-127.0/3-24/3-0.03306235651)*2**23 */
/**
* Returns cube root of 𝑥.
*/
float cbrtf(float x)
{
double_t r,T;
union {float f; uint32_t i;} u = {x};
uint32_t hx = u.i & 0x7fffffff;
if (hx >= 0x7f800000) /* cbrt(NaN,INF) is itself */
return x + x;
/* rough cbrt to 5 bits */
if (hx < 0x00800000) { /* zero or subnormal? */
if (hx == 0)
return x; /* cbrt(+-0) is itself */
u.f = x*0x1p24f;
hx = u.i & 0x7fffffff;
hx = hx/3 + B2;
} else
hx = hx/3 + B1;
u.i &= 0x80000000;
u.i |= hx;
/*
* First step Newton iteration (solving t*t-x/t == 0) to 16 bits. In
* double precision so that its terms can be arranged for efficiency
* without causing overflow or underflow.
*/
T = u.f;
r = T*T*T;
T = T*((double_t)x+x+r)/(x+r+r);
/*
* Second step Newton iteration to 47 bits. In double precision for
* efficiency and accuracy.
*/
r = T*T*T;
T = T*((double_t)x+x+r)/(x+r+r);
/* rounding to 24 bits is perfect in round-to-nearest mode */
return T;
}

36
libc/tinymath/cbrtl.S
View file

@ -1,36 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns cube root of 𝑥.
//
// @param 𝑥 is an 80-bit long double passed on stack in 16-bytes
// @return result of computation on FPU stack in %st
cbrtl: pushq %rbp
mov %rsp,%rbp
sub $16,%rsp
fldt 16(%rbp)
fstpl -8(%rbp)
movsd -8(%rbp),%xmm0
call __cbrt
movsd %xmm0,-8(%rbp)
fldl -8(%rbp)
leave
ret
.endfn cbrtl,globl

168
libc/tinymath/cbrtl.c Normal file
View file

@ -0,0 +1,168 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2020 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/ldshape.internal.h"
asm(".ident\t\"\\n\\n\
FreeBSD libm (BSD-2 License)\\n\
Copyright (c) 2005-2011, Bruce D. Evans, Steven G. Kargl, David Schultz.\"");
asm(".ident\t\"\\n\\n\
fdlibm (fdlibm license)\\n\
Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
// clang-format off
/* origin: FreeBSD /usr/src/lib/msun/src/s_cbrtl.c */
/*-
* ====================================================
* Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
* Copyright (c) 2009-2011, Bruce D. Evans, Steven G. Kargl, David Schultz.
*
* Developed at SunPro, a Sun Microsystems, Inc. business.
* Permission to use, copy, modify, and distribute this
* software is freely granted, provided that this notice
* is preserved.
* ====================================================
*
* The argument reduction and testing for exceptional cases was
* written by Steven G. Kargl with input from Bruce D. Evans
* and David A. Schultz.
*/
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double cbrtl(long double x)
{
return cbrt(x);
}
#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
static const unsigned B1 = 709958130; /* B1 = (127-127.0/3-0.03306235651)*2**23 */
/**
* Returns cube root of 𝑥.
*/
long double cbrtl(long double x)
{
union ldshape u = {x}, v;
union {float f; uint32_t i;} uft;
long double r, s, t, w;
double_t dr, dt, dx;
float_t ft;
int e = u.i.se & 0x7fff;
int sign = u.i.se & 0x8000;
/*
* If x = +-Inf, then cbrt(x) = +-Inf.
* If x = NaN, then cbrt(x) = NaN.
*/
if (e == 0x7fff)
return x + x;
if (e == 0) {
/* Adjust subnormal numbers. */
u.f *= 0x1p120;
e = u.i.se & 0x7fff;
/* If x = +-0, then cbrt(x) = +-0. */
if (e == 0)
return x;
e -= 120;
}
e -= 0x3fff;
u.i.se = 0x3fff;
x = u.f;
switch (e % 3) {
case 1:
case -2:
x *= 2;
e--;
break;
case 2:
case -1:
x *= 4;
e -= 2;
break;
}
v.f = 1.0;
v.i.se = sign | (0x3fff + e/3);
/*
* The following is the guts of s_cbrtf, with the handling of
* special values removed and extra care for accuracy not taken,
* but with most of the extra accuracy not discarded.
*/
/* ~5-bit estimate: */
uft.f = x;
uft.i = (uft.i & 0x7fffffff)/3 + B1;
ft = uft.f;
/* ~16-bit estimate: */
dx = x;
dt = ft;
dr = dt * dt * dt;
dt = dt * (dx + dx + dr) / (dx + dr + dr);
/* ~47-bit estimate: */
dr = dt * dt * dt;
dt = dt * (dx + dx + dr) / (dx + dr + dr);
#if LDBL_MANT_DIG == 64
/*
* dt is cbrtl(x) to ~47 bits (after x has been reduced to 1 <= x < 8).
* Round it away from zero to 32 bits (32 so that t*t is exact, and
* away from zero for technical reasons).
*/
t = dt + (0x1.0p32L + 0x1.0p-31L) - 0x1.0p32;
#elif LDBL_MANT_DIG == 113
/*
* Round dt away from zero to 47 bits. Since we don't trust the 47,
* add 2 47-bit ulps instead of 1 to round up. Rounding is slow and
* might be avoidable in this case, since on most machines dt will
* have been evaluated in 53-bit precision and the technical reasons
* for rounding up might not apply to either case in cbrtl() since
* dt is much more accurate than needed.
*/
t = dt + 0x2.0p-46 + 0x1.0p60L - 0x1.0p60;
#endif
/*
* Final step Newton iteration to 64 or 113 bits with
* error < 0.667 ulps
*/
s = t*t; /* t*t is exact */
r = x/s; /* error <= 0.5 ulps; |r| < |t| */
w = t+t; /* t+t is exact */
r = (r-t)/(w+r); /* r-t is exact; w+r ~= 3*t */
t = t+t*r; /* error <= 0.5 + 0.5/3 + epsilon */
t *= v.f;
return t;
}
#endif

3
libc/tinymath/ccosh.c
View file

@ -29,6 +29,9 @@
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
FreeBSD libm (BSD-2 License)\\n\
Copyright (c) 2005-2011, Bruce D. Evans, Steven G. Kargl, David Schultz.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");

3
libc/tinymath/ccoshf.c
View file

@ -29,6 +29,9 @@
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
FreeBSD libm (BSD-2 License)\\n\
Copyright (c) 2005-2011, Bruce D. Evans, Steven G. Kargl, David Schultz.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");

9
libc/tinymath/ccosl.c
View file

@ -37,14 +37,11 @@ asm(".include \"libc/disclaimer.inc\"");
long double complex ccosl(long double complex z)
{
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
long double complex ccosl(long double complex z)
{
return ccos(z);
}
#else
long double complex ccosl(long double complex z)
{
return ccoshl(CMPLXL(-cimagl(z), creall(z)));
}
#endif
}

57
libc/tinymath/ceil.S
View file

@ -1,57 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns smallest integral not less than 𝑥.
//
// @param 𝑥 is double scalar in low half of %xmm0
// @return double scalar in low half of %xmm0
// @see round(),rint(),nearbyint()
// @see vroundsd $_MM_FROUND_TO_POS_INF|_MM_FROUND_NO_EXC,%xmm0,%xmm0,%xmm0
ceil: .leafprologue
.profilable
movsd 4f(%rip),%xmm3
movsd 2f(%rip),%xmm4
movapd %xmm0,%xmm2
movapd %xmm0,%xmm1
andpd %xmm3,%xmm2
ucomisd %xmm2,%xmm4
jbe 1f
cvttsd2siq %xmm0,%rax
pxor %xmm2,%xmm2
movsd 3f(%rip),%xmm4
andnpd %xmm1,%xmm3
cvtsi2sdq %rax,%xmm2
cmpnlesd %xmm2,%xmm0
andpd %xmm4,%xmm0
addsd %xmm2,%xmm0
orpd %xmm3,%xmm0
1: .leafepilogue
.endfn ceil,globl
.rodata.cst8
2: .long 0x00000000
.long 0x43300000
3: .long 0x00000000
.long 0x3ff00000
.rodata.cst16
4: .long 0xffffffff
.long 0x7fffffff
.long 0x00000000
.long 0x00000000

65
libc/tinymath/ceil.c Normal file
View file

@ -0,0 +1,65 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/internal.h"
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
// clang-format off
#if FLT_EVAL_METHOD==0 || FLT_EVAL_METHOD==1
#define EPS DBL_EPSILON
#elif FLT_EVAL_METHOD==2
#define EPS LDBL_EPSILON
#endif
static const double_t toint = 1/EPS;
double ceil(double x)
{
union {double f; uint64_t i;} u = {x};
int e = u.i >> 52 & 0x7ff;
double_t y;
if (e >= 0x3ff+52 || x == 0)
return x;
/* y = int(x) - x, where int(x) is an integer neighbor of x */
if (u.i >> 63)
y = x - toint + toint - x;
else
y = x + toint - toint - x;
/* special case because of non-nearest rounding modes */
if (e <= 0x3ff-1) {
FORCE_EVAL(y);
return u.i >> 63 ? -0.0 : 1;
}
if (y < 0)
return x + y + 1;
return x + y;
}

55
libc/tinymath/ceilf.S
View file

@ -1,55 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns smallest integral not less than 𝑥.
//
// @param 𝑥 is float scalar in low quarter of %xmm0
// @return float scalar in low quarter of %xmm0
// @see round(),rint(),nearbyint()
// @see vroundss $_MM_FROUND_TO_POS_INF|_MM_FROUND_NO_EXC,%xmm0,%xmm0,%xmm0
ceilf: .leafprologue
.profilable
movss 4f(%rip),%xmm3
movss 2f(%rip),%xmm4
movaps %xmm0,%xmm2
movaps %xmm0,%xmm1
andps %xmm3,%xmm2
ucomiss %xmm2,%xmm4
jbe 1f
cvttss2sil %xmm0,%eax
pxor %xmm2,%xmm2
movss 3f(%rip),%xmm4
andnps %xmm1,%xmm3
cvtsi2ssl %eax,%xmm2
cmpnless %xmm2,%xmm0
andps %xmm4,%xmm0
addss %xmm2,%xmm0
orps %xmm3,%xmm0
1: .leafepilogue
.endfn ceilf,globl
.rodata.cst4
2: .long 0x4b000000
3: .long 0x3f800000
.rodata.cst16
4: .long 0x7fffffff
.long 0x00000000
.long 0x00000000
.long 0x00000000

61
libc/tinymath/ceilf.c Normal file
View file

@ -0,0 +1,61 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/internal.h"
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
// clang-format off
float ceilf(float x)
{
union {float f; uint32_t i;} u = {x};
int e = (int)(u.i >> 23 & 0xff) - 0x7f;
uint32_t m;
if (e >= 23)
return x;
if (e >= 0) {
m = 0x007fffff >> e;
if ((u.i & m) == 0)
return x;
FORCE_EVAL(x + 0x1p120f);
if (u.i >> 31 == 0)
u.i += m;
u.i &= ~m;
} else {
FORCE_EVAL(x + 0x1p120f);
if (u.i >> 31)
u.f = -0.0;
else if (u.i << 1)
u.f = 1.0;
}
return u.f;
}

40
libc/tinymath/ceill.S
View file

@ -1,40 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns smallest integral not less than 𝑥.
//
// @param 𝑥 is long double passed on stack
// @return long double in %st
ceill: pushq %rbp
mov %rsp,%rbp
.profilable
sub $16,%rsp
fnstcw -2(%rbp)
fldt 16(%rbp)
movzwl -2(%rbp),%eax
and $-13,%ah
or $8,%ah
mov %ax,-4(%rbp)
fldcw -4(%rbp)
frndint
fldcw -2(%rbp)
leave
ret
.endfn ceill,globl

63
libc/tinymath/ceill.c Normal file
View file

@ -0,0 +1,63 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
#include "libc/tinymath/internal.h"
#include "libc/tinymath/ldshape.internal.h"
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
// clang-format off
long double ceill(long double x) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
return ceil(x);
#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
static const long double toint = 1/LDBL_EPSILON;
union ldshape u = {x};
int e = u.i.se & 0x7fff;
long double y;
if (e >= 0x3fff+LDBL_MANT_DIG-1 || x == 0)
return x;
/* y = int(x) - x, where int(x) is an integer neighbor of x */
if (u.i.se >> 15)
y = x - toint + toint - x;
else
y = x + toint - toint - x;
/* special case because of non-nearest rounding modes */
if (e <= 0x3fff-1) {
FORCE_EVAL(y);
return u.i.se >> 15 ? -0.0 : 1;
}
if (y < 0)
return x + y + 1;
return x + y;
#endif
}

14
libc/tinymath/fmodf.S → libc/tinymath/cimag.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,8 +16,8 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#include "libc/complex.h"
fmodf: ezlea fmodl,ax
jmp _f2ld2
.endfn fmodf,globl
double(cimag)(double complex z) {
return cimag(z);
}

15
libc/tinymath/lrint.S → libc/tinymath/cimagf.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,9 +16,8 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#include "libc/complex.h"
lrint: cvtsd2siq %xmm0,%rax
ret
.endfn lrint,globl
.alias lrint,llrint
float(cimagf)(float complex z) {
return cimagf(z);
}

15
libc/tinymath/lrintf.S → libc/tinymath/cimagl.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,9 +16,8 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#include "libc/complex.h"
lrintf: cvtss2siq %xmm0,%rax
ret
.endfn lrintf,globl
.alias lrintf,llrintf
long double(cimagl)(long double complex z) {
return cimagl(z);
}

31
libc/tinymath/conj.S
View file

@ -1,31 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
conj: .leafprologue
.profilable
xorpd .L1(%rip),%xmm1
.leafepilogue
.endfn conj,globl
.rodata.cst16
.L1: .long 0
.long -2147483648
.long 0
.long 0

23
libc/tinymath/conj.c Normal file
View file

@ -0,0 +1,23 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
double complex conj(double complex z) {
return CMPLX(creal(z), -cimag(z));
}

38
libc/tinymath/conjf.S
View file

@ -1,38 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
conjf: .profilable
sub $16,%rsp
movq %xmm0,8(%rsp)
movss 12(%rsp),%xmm0
xorps .L1(%rip),%xmm0
movss 8(%rsp),%xmm1
movss %xmm1,(%rsp)
movss %xmm0,4(%rsp)
movq (%rsp),%xmm0
add $16,%rsp
ret
.endfn conjf,globl
.rodata.cst16
.L1: .long 2147483648
.long 0
.long 0
.long 0

23
libc/tinymath/conjf.c Normal file
View file

@ -0,0 +1,23 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
float complex conjf(float complex z) {
return CMPLXF(crealf(z), -cimagf(z));
}

42
libc/tinymath/conjl.S
View file

@ -1,42 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
conjl: .profilable
sub $24,%rsp
fldt 32(%rsp)
fnstcw 14(%rsp)
movzwl 14(%rsp),%eax
orb $12,%ah
movw %ax,12(%rsp)
fldcw 12(%rsp)
fistpq (%rsp)
fldcw 14(%rsp)
movq (%rsp),%rsi
fldt 48(%rsp)
mov %rsi,%rax
fchs
fldcw 12(%rsp)
fistpq (%rsp)
fldcw 14(%rsp)
movq (%rsp),%rcx
add $24,%rsp
mov %rcx,%rdx
ret
.endfn conjl,globl

23
libc/tinymath/conjl.c Normal file
View file

@ -0,0 +1,23 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
long double complex conjl(long double complex z) {
return CMPLXL(creall(z), -cimagl(z));
}

44
libc/tinymath/copysign.S
View file

@ -1,44 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns 𝑥 with same sign as 𝑦.
//
// @param 𝑥 is double scalar in low half of %xmm0
// @param 𝑦 is double scalar in low half of %xmm1
// @return double scalar in low half of %xmm0
copysign:
.leafprologue
.profilable
movapd %xmm1,%xmm2
andpd .Lnan(%rip),%xmm0
andpd .Lneg0(%rip),%xmm2
orpd %xmm2,%xmm0
.leafepilogue
.endfn copysign,globl
.rodata.cst16
.Lnan: .long 0xffffffff
.long 0x7fffffff
.long 0x00000000
.long 0x00000000
.Lneg0: .long 0x00000000
.long 0x80000000
.long 0x00000000
.long 0x00000000

32
libc/tinymath/copysign.c Normal file
View file

@ -0,0 +1,32 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
/**
* Returns 𝑥 with same sign as 𝑦.
*/
double copysign(double x, double y) {
union {
double f;
uint64_t i;
} ux = {x}, uy = {y};
ux.i &= -1ULL / 2;
ux.i |= uy.i & 1ULL << 63;
return ux.f;
}

44
libc/tinymath/copysignf.S
View file

@ -1,44 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns 𝑥 with same sign as 𝑦.
//
// @param 𝑦 is float scalar in low quarter of %xmm0
// @param 𝑥 is float scalar in low quarter of %xmm1
// @return float scalar in low quarter of %xmm0
copysignf:
.leafprologue
.profilable
movaps %xmm1,%xmm2
andps .LC8(%rip),%xmm0
andps .LC10(%rip),%xmm2
orps %xmm2,%xmm0
.leafepilogue
.endfn copysignf,globl
.rodata.cst16
.LC8: .long 2147483647
.long 0
.long 0
.long 0
.LC10: .long 2147483648
.long 0
.long 0
.long 0

33
libc/tinymath/cargf.S → libc/tinymath/copysignf.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,18 +16,17 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#include "libc/math.h"
cargf: push %rbp
mov %rsp,%rbp
.profilable
sub $16,%rsp
movq %xmm0,8(%rsp)
flds 12(%rsp)
flds 8(%rsp)
fpatan
fstps 4(%rsp)
movss 4(%rsp),%xmm0
leave
ret
.endfn cargf,globl
/**
* Returns 𝑥 with same sign as 𝑦.
*/
float copysignf(float x, float y) {
union {
float f;
uint32_t i;
} ux = {x}, uy = {y};
ux.i &= 0x7fffffff;
ux.i |= uy.i & 0x80000000;
return ux.f;
}

42
libc/tinymath/copysignl.S
View file

@ -1,42 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/runtime/pc.internal.h"
#include "libc/macros.internal.h"
// Returns 𝑥 with same sign as 𝑦.
//
// @param 𝑥 is an 80-bit long double passed on stack in 16-bytes
// @param 𝑦 is the power, also pushed on stack, in reverse order
// @return result on FPU stack in %st
copysignl:
push %rbp
mov %rsp,%rbp
.profilable
fldt 32(%rbp)
fxam
fnstsw
fstp %st
fldt 16(%rbp)
test $FPU_C1>>8,%ah
fabs
je 1f
fchs
1: pop %rbp
ret
.endfn copysignl,globl

43
libc/tinymath/asinl.S → libc/tinymath/copysignl.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 sw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,26 +16,19 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#include "libc/math.h"
#include "libc/tinymath/ldshape.internal.h"
// Returns arc sine of 𝑥.
//
// @param 𝑥 is an 80-bit long double passed on stack in 16-bytes
// @return result of computation on FPU stack in %st
// @define atan2(𝑥,sqrt((1-𝑥)*(1+𝑥)))
// @domain -1 ≤ 𝑥 ≤ 1
asinl: pushq %rbp
mov %rsp,%rbp
.profilable
fldt 16(%rbp)
fld1
fld %st
fsub %st(2)
fxch
fadd %st(2)
fmulp
fsqrt
fpatan
leave
ret
.endfn asinl,globl
/**
* Returns 𝑥 with same sign as 𝑦.
*/
long double copysignl(long double x, long double y) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
return copysign(x, y);
#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
union ldshape ux = {x}, uy = {y};
ux.i.se &= 0x7fff;
ux.i.se |= uy.i.se & 0x8000;
return ux.f;
#endif
}

10
libc/tinymath/coshl.c
View file

@ -43,8 +43,10 @@ asm(".include \"libc/disclaimer.inc\"");
* = 1 + 0.5*(exp(x)-1)*(exp(x)-1)/exp(x)
* = 1 + x*x/2 + o(x^4)
*/
long double coshl(long double x)
{
long double coshl(long double x) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
return cosh(x);
#elif LDBL_MANT_DIG == 64 && LDBL_MAX_EXP == 16384
union ldshape u = {x};
unsigned ex = u.i.se & 0x7fff;
uint32_t w;
@ -74,4 +76,8 @@ long double coshl(long double x)
/* |x| > log(LDBL_MAX) or nan */
t = expl(0.5*x);
return 0.5*t*t;
#elif LDBL_MANT_DIG == 113 && LDBL_MAX_EXP == 16384
// TODO: broken implementation to make things compile
return cosh(x);
#endif
}

40
libc/tinymath/cosl.S
View file

@ -1,40 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 sw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns cosine of 𝑥.
//
// @param 𝑥 is an 80-bit long double passed on stack in 16-bytes
// @domain -(3π/8) < 𝑥 < 3π/8 for best accuracy
// @return %st stores result
cosl: push %rbp
mov %rsp,%rbp
.profilable
fldt 16(%rbp)
ezlea _cos,dx
call c2rangr
pop %rbp
ret
.endfn cosl,globl
_cos: .leafprologue
.profilable
fcos
.leafepilogue
.endfn _cos

72
libc/tinymath/cosl.c Normal file
View file

@ -0,0 +1,72 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
#include "libc/math.h"
#include "libc/tinymath/kernel.internal.h"
#include "libc/tinymath/ldshape.internal.h"
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
// clang-format off
long double cosl(long double x) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
return cos(x);
#elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384
union ldshape u = {x};
unsigned n;
long double y[2], hi, lo;
u.i.se &= 0x7fff;
if (u.i.se == 0x7fff)
return x - x;
x = u.f;
if (x < M_PI_4) {
if (u.i.se < 0x3fff - LDBL_MANT_DIG)
/* raise inexact if x!=0 */
return 1.0 + x;
return __cosl(x, 0);
}
n = __rem_pio2l(x, y);
hi = y[0];
lo = y[1];
switch (n & 3) {
case 0:
return __cosl(hi, lo);
case 1:
return -__sinl(hi, lo, 1);
case 2:
return -__cosl(hi, lo);
case 3:
default:
return __sinl(hi, lo, 1);
}
#endif
}

58
libc/tinymath/cproj.c
View file

@ -1,45 +1,27 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
/* clang-format off */
double complex cproj(double complex z)
{
if (isinf(creal(z)) || isinf(cimag(z)))
return CMPLX(INFINITY, copysign(0.0, creal(z)));
return z;
double complex cproj(double complex z) {
if (isinf(creal(z)) || isinf(cimag(z))) {
return CMPLX(INFINITY, copysign(0.0, creal(z)));
}
return z;
}

33
libc/tinymath/cprojf.S
View file

@ -1,33 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
cprojf: push %rbp
mov %rsp,%rbp
.profilable
sub $16,%rsp
movq %xmm0,8(%rsp)
movss 8(%rsp),%xmm0
movss %xmm0,(%rsp)
movss 12(%rsp),%xmm0
movss %xmm0,4(%rsp)
movq (%rsp),%xmm0
leave
ret
.endfn cprojf,globl

21
libc/tinymath/cbrt.S → libc/tinymath/cprojf.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,11 +16,12 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#include "libc/complex.h"
#include "libc/math.h"
// Returns cube root of 𝑥.
//
// @param %xmm0 holds binary64 number
// @return %xmm0 holds binary64 result
cbrt: jmp __cbrt
.endfn cbrt,globl
float complex cprojf(float complex z) {
if (isinf(crealf(z)) || isinf(cimagf(z))) {
return CMPLXF(INFINITY, copysignf(0.0, cimagf(z)));
}
return z;
}

37
libc/tinymath/cprojl.S
View file

@ -1,37 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Projects into Rienmann sphere.
//
// @param z is complex long double passed on stack
// @note needs sse3
cprojl: .profilable
sub $24,%rsp
fldt 32(%rsp)
fisttpq 8(%rsp)
fldt 48(%rsp)
movq 8(%rsp),%rsi
mov %rsi,%rax
fisttpq 8(%rsp)
mov 8(%rsp),%rcx
add $24,%rsp
mov %rcx,%rdx
ret
.endfn cprojl,globl

29
libc/tinymath/asinf-tiny.S → libc/tinymath/cprojl.c
View file

@ -1,7 +1,7 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
@ -16,15 +16,16 @@
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#ifdef TINY
#include "libc/complex.h"
#include "libc/math.h"
// Returns arc sine of 𝑥.
//
// @param 𝑥 is float scalar in low quarter of %xmm0
// @return float scalar in low quarter of %xmm0
asinf: ezlea asinl,ax
jmp _f2ld2
.endfn asinf,globl
#endif /* TINY */
long double complex cprojl(long double complex z) {
#if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024
return cproj(z);
#else
if (isinf(creall(z)) || isinf(cimagl(z))) {
return CMPLXL(INFINITY, copysignl(0.0, cimagl(z)));
}
return z;
#endif
}

23
libc/tinymath/creal.c Normal file
View file

@ -0,0 +1,23 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
double(creal)(double complex z) {
return creal(z);
}

29
libc/tinymath/crealf.S
View file

@ -1,29 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
crealf: push %rbp
mov %rsp,%rbp
.profilable
push %rax
movq %xmm0,(%rsp)
movss (%rsp),%xmm0
leave
ret
.endfn crealf,globl

23
libc/tinymath/crealf.c Normal file
View file

@ -0,0 +1,23 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
float(crealf)(float complex z) {
return crealf(z);
}

27
libc/tinymath/creall.S
View file

@ -1,27 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
creall: push %rbp
mov %rsp,%rbp
.profilable
fldt 16(%rbp)
pop %rbp
ret
.endfn creall,globl

23
libc/tinymath/creall.c Normal file
View file

@ -0,0 +1,23 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/complex.h"
long double(creall)(long double complex z) {
return creall(z);
}

3
libc/tinymath/csinh.c
View file

@ -29,6 +29,9 @@
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
FreeBSD libm (BSD-2 License)\\n\
Copyright (c) 2005-2011, Bruce D. Evans, Steven G. Kargl, David Schultz.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");

3
libc/tinymath/csinhf.c
View file

@ -29,6 +29,9 @@
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
FreeBSD libm (BSD-2 License)\\n\
Copyright (c) 2005-2011, Bruce D. Evans, Steven G. Kargl, David Schultz.\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");

3
libc/tinymath/ctanh.c
View file

@ -29,6 +29,9 @@
#include "libc/math.h"
#include "libc/tinymath/complex.internal.h"
asm(".ident\t\"\\n\\n\
ctahnh (BSD-2 License)\\n\
Copyright (c) 2011 David Schultz <das@FreeBSD.ORG>\"");
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");

41
libc/tinymath/d2ld2.S
View file

@ -1,41 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Thunks double(*fn)(double,double) -> long double fn.
//
// @param %xmm0[0] contains double param
// @return %xmm0[0] contains double result
// @note 100% negligible overhead
_d2ld2: push %rbp
mov %rsp,%rbp
.profilable
sub $32,%rsp
movsd %xmm0,-32(%rbp)
fldl -32(%rbp)
fstpt -32(%rbp)
movsd %xmm1,-16(%rbp)
fldl -16(%rbp)
fstpt -16(%rbp)
call *%rax
fstpl -16(%rbp)
movsd -16(%rbp),%xmm0
leave
ret
.endfn _d2ld2,globl,hidden

23
libc/tinymath/drem.c Normal file
View file

@ -0,0 +1,23 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
double drem(double x, double y) {
return remainder(x, y);
}

23
libc/tinymath/dremf.c Normal file
View file

@ -0,0 +1,23 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
float dremf(float x, float y) {
return remainderf(x, y);
}

23
libc/tinymath/dreml.c Normal file
View file

@ -0,0 +1,23 @@
/*-*- mode:c;indent-tabs-mode:nil;c-basic-offset:2;tab-width:8;coding:utf-8 -*-│
│vi: set net ft=c ts=2 sts=2 sw=2 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2023 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
long double dreml(long double x, long double y) {
return remainderl(x, y);
}

30
libc/tinymath/exp-tiny.S
View file

@ -1,30 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
#ifdef TINY
// Returns 𝑒^x.
//
// @param 𝑥 is float scalar in low quarter of %xmm0
// @return float scalar in low quarter of %xmm0
exp: ezlea expl,ax
jmp _d2ld2
.endfn exp,globl
#endif /* TINY */

7
libc/tinymath/exp.c
View file

@ -2,8 +2,8 @@
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ Optimized Routines │
│ Copyright (c) 1999-2022, Arm Limited. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
@ -29,7 +29,6 @@
#include "libc/math.h"
#include "libc/tinymath/exp_data.internal.h"
#include "libc/tinymath/internal.h"
#ifndef TINY
asm(".ident\t\"\\n\\n\
Double-precision math functions (MIT License)\\n\
@ -169,5 +168,3 @@ double exp(double x)
is no spurious underflow here even without fma. */
return eval_as_double(scale + scale * tmp);
}
#endif /* !TINY */

29
libc/tinymath/exp10.S
View file

@ -1,29 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns 10^x.
//
// @param 𝑥 is double scalar in low half of %xmm0
// @return double scalar in low half of %xmm0
// @see pow(), exp()
exp10: ezlea exp10l,ax
jmp _d2ld2
.endfn exp10,globl
.alias exp10,pow10

53
libc/tinymath/exp10.c Normal file
View file

@ -0,0 +1,53 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
// clang-format off
double exp10(double x)
{
static const double p10[] = {
1e-15, 1e-14, 1e-13, 1e-12, 1e-11, 1e-10,
1e-9, 1e-8, 1e-7, 1e-6, 1e-5, 1e-4, 1e-3, 1e-2, 1e-1,
1, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7, 1e8, 1e9,
1e10, 1e11, 1e12, 1e13, 1e14, 1e15
};
double n, y = modf(x, &n);
union {double f; uint64_t i;} u = {n};
/* fabs(n) < 16 without raising invalid on nan */
if ((u.i>>52 & 0x7ff) < 0x3ff+4) {
if (!y) return p10[(int)n+15];
y = exp2(3.32192809488736234787031942948939 * y);
return y * p10[(int)n+15];
}
return pow(10.0, x);
}

29
libc/tinymath/exp10f.S
View file

@ -1,29 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns 10^x.
//
// @param 𝑥 is float scalar in low quarter of %xmm0
// @return float scalar in low quarter of %xmm0
exp10f:
ezlea exp10l,ax
jmp _f2ld2
.endfn exp10f,globl
.alias exp10f,pow10f

51
libc/tinymath/exp10f.c Normal file
View file

@ -0,0 +1,51 @@
/*-*- mode:c;indent-tabs-mode:t;c-basic-offset:8;tab-width:8;coding:utf-8 -*-│
│vi: set et ft=c ts=8 tw=8 fenc=utf-8 :vi│
╚──────────────────────────────────────────────────────────────────────────────╝
│ │
│ Musl Libc │
│ Copyright © 2005-2014 Rich Felker, et al. │
│ │
│ Permission is hereby granted, free of charge, to any person obtaining │
│ a copy of this software and associated documentation files (the │
│ "Software"), to deal in the Software without restriction, including │
│ without limitation the rights to use, copy, modify, merge, publish, │
│ distribute, sublicense, and/or sell copies of the Software, and to │
│ permit persons to whom the Software is furnished to do so, subject to │
│ the following conditions: │
│ │
│ The above copyright notice and this permission notice shall be │
│ included in all copies or substantial portions of the Software. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, │
│ EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF │
│ MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. │
│ IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY │
│ CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT, │
│ TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE │
│ SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. │
│ │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/math.h"
asm(".ident\t\"\\n\\n\
Musl libc (MIT License)\\n\
Copyright 2005-2014 Rich Felker, et. al.\"");
asm(".include \"libc/disclaimer.inc\"");
// clang-format off
float exp10f(float x)
{
static const float p10[] = {
1e-7f, 1e-6f, 1e-5f, 1e-4f, 1e-3f, 1e-2f, 1e-1f,
1, 1e1, 1e2, 1e3, 1e4, 1e5, 1e6, 1e7
};
float n, y = modff(x, &n);
union {float f; uint32_t i;} u = {n};
/* fabsf(n) < 8 without raising invalid on nan */
if ((u.i>>23 & 0xff) < 0x7f+3) {
if (!y) return p10[(int)n+7];
y = exp2f(3.32192809488736234787031942948939f * y);
return y * p10[(int)n+7];
}
return exp2(3.32192809488736234787031942948939 * x);
}

54
libc/tinymath/exp10l.S
View file

@ -1,54 +0,0 @@
/*-*- mode:unix-assembly; indent-tabs-mode:t; tab-width:8; coding:utf-8 -*-│
│vi: set et ft=asm ts=8 tw=8 fenc=utf-8 :vi│
╞══════════════════════════════════════════════════════════════════════════════╡
│ Copyright 2020 Justine Alexandra Roberts Tunney │
│ │
│ Permission to use, copy, modify, and/or distribute this software for │
│ any purpose with or without fee is hereby granted, provided that the │
│ above copyright notice and this permission notice appear in all copies. │
│ │
│ THE SOFTWARE IS PROVIDED "AS IS" AND THE AUTHOR DISCLAIMS ALL │
│ WARRANTIES WITH REGARD TO THIS SOFTWARE INCLUDING ALL IMPLIED │
│ WARRANTIES OF MERCHANTABILITY AND FITNESS. IN NO EVENT SHALL THE │
│ AUTHOR BE LIABLE FOR ANY SPECIAL, DIRECT, INDIRECT, OR CONSEQUENTIAL │
│ DAMAGES OR ANY DAMAGES WHATSOEVER RESULTING FROM LOSS OF USE, DATA OR │
│ PROFITS, WHETHER IN AN ACTION OF CONTRACT, NEGLIGENCE OR OTHER │
│ TORTIOUS ACTION, ARISING OUT OF OR IN CONNECTION WITH THE USE OR │
│ PERFORMANCE OF THIS SOFTWARE. │
╚─────────────────────────────────────────────────────────────────────────────*/
#include "libc/macros.internal.h"
// Returns 10^x.
//
// @param 𝑥 is an 80-bit long double passed on stack in 16-bytes
// @return result of exponentiation on FPU stack in %st
exp10l: push %rbp
mov %rsp,%rbp
.profilable
fldt 16(%rbp)
fxam # isinf(x)
fstsw %ax
mov %ah,%al
and $0x45,%ah
cmp $5,%ah
je 1f
fldl2t
fmulp %st,%st(1)
fld %st
frndint
fsubr %st,%st(1)
fxch
f2xm1
fld1
faddp
fscale
0: fstp %st(1)
pop %rbp
ret
1: test $2,%al # signbit(x)
jz 0b
fstp %st
fldz
jmp 0b
.endfn exp10l,globl
.alias exp10l,pow10l

Some files were not shown because too many files have changed in this diff Show more