| author | |
| committer | |
| log | a6d6af947d09fbd52662c96c6f0a12510ef15f70 |
| tree | 4093cf1120a04a325c8bb5b019b9de9577a83bd8 |
| parent | b5ec3e597e2cd374853dea5244b2bc2bccccd96a |
| signature |
The implementation was ported from `musl`. Unit tests for `f80` and `f128` were also added.
`__cosl` was already implemented in `trig.zig` while working on `sinl`.
The changes were tested by running:
```
$ ./build/stage3/bin/zig build -p stage4 -Denable-llvm -Dno-lib
$ stage4/bin/zig build test-libc -Dlibc-test-path=<LIBC-TEST-PATH> -Dtest-filter='math.cosl' -fqemu -fwasmtime --summary line
Build Summary: 553/553 steps succeeded
```9 files changed, 103 insertions(+), 338 deletions(-)
lib/compiler_rt/cos.zig+103-31| ... | @@ -1,13 +1,23 @@ | ... | @@ -1,13 +1,23 @@ |
| 1 | //! Ported from musl, which is licensed under the MIT license: | ||
| 2 | //! https://git.musl-libc.org/cgit/musl/tree/COPYRIGHT | ||
| 3 | //! | ||
| 4 | //! https://git.musl-libc.org/cgit/musl/tree/src/math/cosf.c | ||
| 5 | //! https://git.musl-libc.org/cgit/musl/tree/src/math/cos.c | ||
| 6 | //! https://git.musl-libc.org/cgit/musl/tree/src/math/cosl.c | ||
| 7 | |||
| 1 | const std = @import("std"); | 8 | const std = @import("std"); |
| 2 | const math = std.math; | 9 | const math = std.math; |
| 3 | const mem = std.mem; | 10 | const mem = std.mem; |
| 4 | const expect = std.testing.expect; | 11 | const expect = std.testing.expect; |
| 12 | const expectApproxEqAbs = std.testing.expectApproxEqAbs; | ||
| 5 | 13 | ||
| 6 | const compiler_rt = @import("../compiler_rt.zig"); | 14 | const compiler_rt = @import("../compiler_rt.zig"); |
| 7 | const symbol = @import("../compiler_rt.zig").symbol; | 15 | const symbol = @import("../compiler_rt.zig").symbol; |
| 8 | const trig = @import("trig.zig"); | 16 | const trig = @import("trig.zig"); |
| 9 | const rem_pio2 = @import("rem_pio2.zig").rem_pio2; | 17 | const rem_pio2 = @import("rem_pio2.zig").rem_pio2; |
| 10 | const rem_pio2f = @import("rem_pio2f.zig").rem_pio2f; | 18 | const rem_pio2f = @import("rem_pio2f.zig").rem_pio2f; |
| 19 | const rem_pio2l = @import("rem_pio2l.zig").rem_pio2l; | ||
| 20 | const utils = @import("math_utils.zig"); | ||
| 11 | 21 | ||
| 12 | comptime { | 22 | comptime { |
| 13 | symbol(&__cosh, "__cosh"); | 23 | symbol(&__cosh, "__cosh"); |
| ... | @@ -112,14 +122,36 @@ pub fn cos(x: f64) callconv(.c) f64 { | ... | @@ -112,14 +122,36 @@ pub fn cos(x: f64) callconv(.c) f64 { |
| 112 | }; | 122 | }; |
| 113 | } | 123 | } |
| 114 | 124 | ||
| 115 | pub fn __cosx(a: f80) callconv(.c) f80 { | 125 | fn coslGeneric(comptime T: type, x: T) T { |
| 116 | // TODO: more efficient implementation | 126 | const se = utils.ldSignExponent(x) & 0x7fff; |
| 117 | return @floatCast(cosq(a)); | 127 | if (se == 0x7fff) { |
| 128 | return x - x; | ||
| 129 | } | ||
| 130 | |||
| 131 | if (@abs(x) < utils.pi_4) { | ||
| 132 | if (se < 0x3fff - math.floatMantissaBits(T)) { | ||
| 133 | // raise inexact if x!=0 | ||
| 134 | return 1.0 + x; | ||
| 135 | } | ||
| 136 | return trig.__cosl(T, x, 0.0); | ||
| 137 | } | ||
| 138 | |||
| 139 | var y: [2]T = undefined; | ||
| 140 | const n = rem_pio2l(T, x, &y); | ||
| 141 | return switch (n & 3) { | ||
| 142 | 0 => trig.__cosl(T, y[0], y[1]), | ||
| 143 | 1 => -trig.__sinl(T, y[0], y[1], 1), | ||
| 144 | 2 => -trig.__cosl(T, y[0], y[1]), | ||
| 145 | else => trig.__sinl(T, y[0], y[1], 1), | ||
| 146 | }; | ||
| 147 | } | ||
| 148 | |||
| 149 | pub fn __cosx(x: f80) callconv(.c) f80 { | ||
| 150 | return coslGeneric(f80, x); | ||
| 118 | } | 151 | } |
| 119 | 152 | ||
| 120 | pub fn cosq(a: f128) callconv(.c) f128 { | 153 | pub fn cosq(x: f128) callconv(.c) f128 { |
| 121 | // TODO: more correct implementation | 154 | return coslGeneric(f128, x); |
| 122 | return cos(@floatCast(a)); | ||
| 123 | } | 155 | } |
| 124 | 156 | ||
| 125 | pub fn cosl(x: c_longdouble) callconv(.c) c_longdouble { | 157 | pub fn cosl(x: c_longdouble) callconv(.c) c_longdouble { |
| ... | @@ -133,38 +165,78 @@ pub fn cosl(x: c_longdouble) callconv(.c) c_longdouble { | ... | @@ -133,38 +165,78 @@ pub fn cosl(x: c_longdouble) callconv(.c) c_longdouble { |
| 133 | } | 165 | } |
| 134 | } | 166 | } |
| 135 | 167 | ||
| 136 | test "cos32" { | 168 | fn testCosSpecial(comptime T: type) !void { |
| 137 | const epsilon = 0.00001; | 169 | const f = switch (T) { |
| 170 | f32 => cosf, | ||
| 171 | f64 => cos, | ||
| 172 | f80 => __cosx, | ||
| 173 | f128 => cosq, | ||
| 174 | else => @compileError("unimplemented"), | ||
| 175 | }; | ||
| 138 | 176 | ||
| 139 | try expect(math.approxEqAbs(f32, cosf(0.0), 1.0, epsilon)); | 177 | try expect(f(0.0) == 1.0); |
| 140 | try expect(math.approxEqAbs(f32, cosf(0.2), 0.980067, epsilon)); | 178 | try expect(f(-0.0) == 1.0); |
| 141 | try expect(math.approxEqAbs(f32, cosf(0.8923), 0.627623, epsilon)); | 179 | try expect(math.isNan(f(math.inf(T)))); |
| 142 | try expect(math.approxEqAbs(f32, cosf(1.5), 0.070737, epsilon)); | 180 | try expect(math.isNan(f(-math.inf(T)))); |
| 143 | try expect(math.approxEqAbs(f32, cosf(-1.5), 0.070737, epsilon)); | 181 | try expect(math.isNan(f(math.nan(T)))); |
| 144 | try expect(math.approxEqAbs(f32, cosf(37.45), 0.969132, epsilon)); | ||
| 145 | try expect(math.approxEqAbs(f32, cosf(89.123), 0.400798, epsilon)); | ||
| 146 | } | 182 | } |
| 147 | 183 | ||
| 148 | test "cos64" { | 184 | test "cos32.normal" { |
| 149 | const epsilon = 0.000001; | 185 | const epsilon = math.floatEps(f32); |
| 150 | 186 | try expectApproxEqAbs(@as(f32, 1.0), cosf(0.0), epsilon); | |
| 151 | try expect(math.approxEqAbs(f64, cos(0.0), 1.0, epsilon)); | 187 | try expectApproxEqAbs(@as(f32, 0.9800666), cosf(0.2), epsilon); |
| 152 | try expect(math.approxEqAbs(f64, cos(0.2), 0.980067, epsilon)); | 188 | try expectApproxEqAbs(@as(f32, 0.6276231), cosf(0.8923), epsilon); |
| 153 | try expect(math.approxEqAbs(f64, cos(0.8923), 0.627623, epsilon)); | 189 | try expectApproxEqAbs(@as(f32, 0.0707372), cosf(1.5), epsilon); |
| 154 | try expect(math.approxEqAbs(f64, cos(1.5), 0.070737, epsilon)); | 190 | try expectApproxEqAbs(@as(f32, 0.0707372), cosf(-1.5), epsilon); |
| 155 | try expect(math.approxEqAbs(f64, cos(-1.5), 0.070737, epsilon)); | 191 | try expectApproxEqAbs(@as(f32, 0.96913195), cosf(37.45), epsilon); |
| 156 | try expect(math.approxEqAbs(f64, cos(37.45), 0.969132, epsilon)); | 192 | try expectApproxEqAbs(@as(f32, 0.40079966), cosf(89.123), epsilon); |
| 157 | try expect(math.approxEqAbs(f64, cos(89.123), 0.40080, epsilon)); | ||
| 158 | } | 193 | } |
| 159 | 194 | ||
| 160 | test "cos32.special" { | 195 | test "cos32.special" { |
| 161 | try expect(math.isNan(cosf(math.inf(f32)))); | 196 | try testCosSpecial(f32); |
| 162 | try expect(math.isNan(cosf(-math.inf(f32)))); | 197 | } |
| 163 | try expect(math.isNan(cosf(math.nan(f32)))); | 198 | |
| 199 | test "cos64.normal" { | ||
| 200 | const epsilon = math.floatEps(f64); | ||
| 201 | try expectApproxEqAbs(@as(f64, 1.0), cos(0.0), epsilon); | ||
| 202 | try expectApproxEqAbs(@as(f64, 0.9800665778412416), cos(0.2), epsilon); | ||
| 203 | try expectApproxEqAbs(@as(f64, 0.6276230983360804), cos(0.8923), epsilon); | ||
| 204 | try expectApproxEqAbs(@as(f64, 0.0707372016677029), cos(1.5), epsilon); | ||
| 205 | try expectApproxEqAbs(@as(f64, 0.0707372016677029), cos(-1.5), epsilon); | ||
| 206 | try expectApproxEqAbs(@as(f64, 0.9691317730707778), cos(37.45), epsilon); | ||
| 207 | try expectApproxEqAbs(@as(f64, 0.4008006809354791), cos(89.123), epsilon); | ||
| 164 | } | 208 | } |
| 165 | 209 | ||
| 166 | test "cos64.special" { | 210 | test "cos64.special" { |
| 167 | try expect(math.isNan(cos(math.inf(f64)))); | 211 | try testCosSpecial(f64); |
| 168 | try expect(math.isNan(cos(-math.inf(f64)))); | 212 | } |
| 169 | try expect(math.isNan(cos(math.nan(f64)))); | 213 | |
| 214 | test "cos80.normal" { | ||
| 215 | const epsilon = math.floatEps(f80); | ||
| 216 | try expectApproxEqAbs(@as(f80, 1.0), __cosx(0.0), epsilon); | ||
| 217 | try expectApproxEqAbs(@as(f80, 0.98006657784124163112419651674816888), __cosx(0.2), epsilon); | ||
| 218 | try expectApproxEqAbs(@as(f80, 0.62762309833608037003563995939286067), __cosx(0.8923), epsilon); | ||
| 219 | try expectApproxEqAbs(@as(f80, 0.070737201667702910088189851434268747), __cosx(1.5), epsilon); | ||
| 220 | try expectApproxEqAbs(@as(f80, 0.070737201667702910088189851434268747), __cosx(-1.5), epsilon); | ||
| 221 | try expectApproxEqAbs(@as(f80, 0.9691317730707771246), __cosx(37.45), epsilon); | ||
| 222 | try expectApproxEqAbs(@as(f80, 0.4008006809354834001), __cosx(89.123), epsilon); | ||
| 223 | } | ||
| 224 | |||
| 225 | test "cos80.special" { | ||
| 226 | try testCosSpecial(f80); | ||
| 227 | } | ||
| 228 | |||
| 229 | test "cos128.normal" { | ||
| 230 | const epsilon = math.floatEps(f128); | ||
| 231 | try expectApproxEqAbs(@as(f128, 1.0), cosq(0.0), epsilon); | ||
| 232 | try expectApproxEqAbs(@as(f128, 0.98006657784124163112419651674816888), cosq(0.2), epsilon); | ||
| 233 | try expectApproxEqAbs(@as(f128, 0.62762309833608037003563995939286067), cosq(0.8923), epsilon); | ||
| 234 | try expectApproxEqAbs(@as(f128, 0.070737201667702910088189851434268747), cosq(1.5), epsilon); | ||
| 235 | try expectApproxEqAbs(@as(f128, 0.070737201667702910088189851434268747), cosq(-1.5), epsilon); | ||
| 236 | try expectApproxEqAbs(@as(f128, 0.96913177307077712443149563847233230), cosq(37.45), epsilon); | ||
| 237 | try expectApproxEqAbs(@as(f128, 0.40080068093548339848199454493704702), cosq(89.123), epsilon); | ||
| 238 | } | ||
| 239 | |||
| 240 | test "cos128.special" { | ||
| 241 | try testCosSpecial(f128); | ||
| 170 | } | 242 | } |
lib/libc/mingw/math/x86/cos.def.h deleted-65| ... | @@ -1,65 +0,0 @@ | ||
| 1 | /* | ||
| 2 | This Software is provided under the Zope Public License (ZPL) Version 2.1. | ||
| 3 | |||
| 4 | Copyright (c) 2009, 2010 by the mingw-w64 project | ||
| 5 | |||
| 6 | See the AUTHORS file for the list of contributors to the mingw-w64 project. | ||
| 7 | |||
| 8 | This license has been certified as open source. It has also been designated | ||
| 9 | as GPL compatible by the Free Software Foundation (FSF). | ||
| 10 | |||
| 11 | Redistribution and use in source and binary forms, with or without | ||
| 12 | modification, are permitted provided that the following conditions are met: | ||
| 13 | |||
| 14 | 1. Redistributions in source code must retain the accompanying copyright | ||
| 15 | notice, this list of conditions, and the following disclaimer. | ||
| 16 | 2. Redistributions in binary form must reproduce the accompanying | ||
| 17 | copyright notice, this list of conditions, and the following disclaimer | ||
| 18 | in the documentation and/or other materials provided with the | ||
| 19 | distribution. | ||
| 20 | 3. Names of the copyright holders must not be used to endorse or promote | ||
| 21 | products derived from this software without prior written permission | ||
| 22 | from the copyright holders. | ||
| 23 | 4. The right to distribute this software or to use it for any purpose does | ||
| 24 | not give you the right to use Servicemarks (sm) or Trademarks (tm) of | ||
| 25 | the copyright holders. Use of them is covered by separate agreement | ||
| 26 | with the copyright holders. | ||
| 27 | 5. If any files are modified, you must cause the modified files to carry | ||
| 28 | prominent notices stating that you changed the files and the date of | ||
| 29 | any change. | ||
| 30 | |||
| 31 | Disclaimer | ||
| 32 | |||
| 33 | THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS ``AS IS'' AND ANY EXPRESSED | ||
| 34 | OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES | ||
| 35 | OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO | ||
| 36 | EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY DIRECT, INDIRECT, | ||
| 37 | INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT | ||
| 38 | LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, | ||
| 39 | OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF | ||
| 40 | LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING | ||
| 41 | NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, | ||
| 42 | EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. | ||
| 43 | */ | ||
| 44 | |||
| 45 | #include "../complex/complex_internal.h" | ||
| 46 | #include <errno.h> | ||
| 47 | |||
| 48 | extern long double __cosl_internal (long double); | ||
| 49 | |||
| 50 | __FLT_TYPE | ||
| 51 | __FLT_ABI(cos) (__FLT_TYPE x) | ||
| 52 | { | ||
| 53 | int x_class = fpclassify (x); | ||
| 54 | if (x_class == FP_NAN) | ||
| 55 | { | ||
| 56 | __FLT_RPT_DOMAIN ("cos", x, 0.0, x); | ||
| 57 | return x; | ||
| 58 | } | ||
| 59 | else if (x_class == FP_INFINITE) | ||
| 60 | { | ||
| 61 | __FLT_RPT_DOMAIN ("cos", x, 0.0, __FLT_NAN); | ||
| 62 | return __FLT_NAN; | ||
| 63 | } | ||
| 64 | return (__FLT_TYPE) __cosl_internal ((long double) x); | ||
| 65 | } | ||
lib/libc/mingw/math/x86/cosl.c deleted-46| ... | @@ -1,46 +0,0 @@ | ||
| 1 | /* | ||
| 2 | This Software is provided under the Zope Public License (ZPL) Version 2.1. | ||
| 3 | |||
| 4 | Copyright (c) 2009, 2010 by the mingw-w64 project | ||
| 5 | |||
| 6 | See the AUTHORS file for the list of contributors to the mingw-w64 project. | ||
| 7 | |||
| 8 | This license has been certified as open source. It has also been designated | ||
| 9 | as GPL compatible by the Free Software Foundation (FSF). | ||
| 10 | |||
| 11 | Redistribution and use in source and binary forms, with or without | ||
| 12 | modification, are permitted provided that the following conditions are met: | ||
| 13 | |||
| 14 | 1. Redistributions in source code must retain the accompanying copyright | ||
| 15 | notice, this list of conditions, and the following disclaimer. | ||
| 16 | 2. Redistributions in binary form must reproduce the accompanying | ||
| 17 | copyright notice, this list of conditions, and the following disclaimer | ||
| 18 | in the documentation and/or other materials provided with the | ||
| 19 | distribution. | ||
| 20 | 3. Names of the copyright holders must not be used to endorse or promote | ||
| 21 | products derived from this software without prior written permission | ||
| 22 | from the copyright holders. | ||
| 23 | 4. The right to distribute this software or to use it for any purpose does | ||
| 24 | not give you the right to use Servicemarks (sm) or Trademarks (tm) of | ||
| 25 | the copyright holders. Use of them is covered by separate agreement | ||
| 26 | with the copyright holders. | ||
| 27 | 5. If any files are modified, you must cause the modified files to carry | ||
| 28 | prominent notices stating that you changed the files and the date of | ||
| 29 | any change. | ||
| 30 | |||
| 31 | Disclaimer | ||
| 32 | |||
| 33 | THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS ``AS IS'' AND ANY EXPRESSED | ||
| 34 | OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES | ||
| 35 | OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO | ||
| 36 | EVENT SHALL THE COPYRIGHT HOLDERS BE LIABLE FOR ANY DIRECT, INDIRECT, | ||
| 37 | INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT | ||
| 38 | LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, | ||
| 39 | OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF | ||
| 40 | LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING | ||
| 41 | NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS SOFTWARE, | ||
| 42 | EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. | ||
| 43 | */ | ||
| 44 | |||
| 45 | #define _NEW_COMPLEX_LDOUBLE 1 | ||
| 46 | #include "cos.def.h" | ||
lib/libc/mingw/math/x86/cosl_internal.S deleted-55| ... | @@ -1,55 +0,0 @@ | ||
| 1 | /** | ||
| 2 | * This file has no copyright assigned and is placed in the Public Domain. | ||
| 3 | * This file is part of the mingw-w64 runtime package. | ||
| 4 | * No warranty is given; refer to the file DISCLAIMER.PD within this package. | ||
| 5 | */ | ||
| 6 | #include <_mingw_mac.h> | ||
| 7 | |||
| 8 | 	.file	"cosl_internal.S" | ||
| 9 | 	.text | ||
| 10 | #ifdef __x86_64__ | ||
| 11 | 	.align 8 | ||
| 12 | #else | ||
| 13 | 	.align 4 | ||
| 14 | #endif | ||
| 15 | .globl __MINGW_USYMBOL(__cosl_internal) | ||
| 16 | 	.def	__MINGW_USYMBOL(__cosl_internal);	.scl	2;	.type	32;	.endef | ||
| 17 | __MINGW_USYMBOL(__cosl_internal): | ||
| 18 | #ifdef __x86_64__ | ||
| 19 | 	fldt	(%rdx) | ||
| 20 | 	fcos | ||
| 21 | 	fnstsw	%ax | ||
| 22 | 	testl	$0x400,%eax | ||
| 23 | 	jz	1f | ||
| 24 | 	fldpi | ||
| 25 | 	fadd	%st(0) | ||
| 26 | 	fxch	%st(1) | ||
| 27 | 2:	fprem1 | ||
| 28 | 	fnstsw	%ax | ||
| 29 | 	testl	$0x400,%eax | ||
| 30 | 	jnz	2b | ||
| 31 | 	fstp	%st(1) | ||
| 32 | 	fcos | ||
| 33 | 1:	movq %rcx,%rax | ||
| 34 | 	movq	$0,8(%rcx) | ||
| 35 | 	fstpt (%rcx) | ||
| 36 | 	ret | ||
| 37 | #else | ||
| 38 | 	fldt	4(%esp) | ||
| 39 | 	fcos | ||
| 40 | 	fnstsw	%ax | ||
| 41 | 	testl	$0x400,%eax | ||
| 42 | 	jnz	1f | ||
| 43 | 	ret | ||
| 44 | 1:	fldpi | ||
| 45 | 	fadd	%st(0) | ||
| 46 | 	fxch	%st(1) | ||
| 47 | 2:	fprem1 | ||
| 48 | 	fnstsw	%ax | ||
| 49 | 	testl	$0x400,%eax | ||
| 50 | 	jnz	2b | ||
| 51 | 	fstp	%st(1) | ||
| 52 | 	fcos | ||
| 53 | 	ret | ||
| 54 | #endif | ||
| 55 | |||
lib/libc/musl/src/math/__cosl.c deleted-96| ... | @@ -1,96 +0,0 @@ | ||
| 1 | /* origin: FreeBSD /usr/src/lib/msun/ld80/k_cosl.c */ | ||
| 2 | /* origin: FreeBSD /usr/src/lib/msun/ld128/k_cosl.c */ | ||
| 3 | /* | ||
| 4 | * ==================================================== | ||
| 5 | * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved. | ||
| 6 | * Copyright (c) 2008 Steven G. Kargl, David Schultz, Bruce D. Evans. | ||
| 7 | * | ||
| 8 | * Developed at SunSoft, a Sun Microsystems, Inc. business. | ||
| 9 | * Permission to use, copy, modify, and distribute this | ||
| 10 | * software is freely granted, provided that this notice | ||
| 11 | * is preserved. | ||
| 12 | * ==================================================== | ||
| 13 | */ | ||
| 14 | |||
| 15 | |||
| 16 | #include "libm.h" | ||
| 17 | |||
| 18 | #if (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384 | ||
| 19 | #if LDBL_MANT_DIG == 64 | ||
| 20 | /* | ||
| 21 | * ld80 version of __cos.c. See __cos.c for most comments. | ||
| 22 | */ | ||
| 23 | /* | ||
| 24 | * Domain [-0.7854, 0.7854], range ~[-2.43e-23, 2.425e-23]: | ||
| 25 | * |cos(x) - c(x)| < 2**-75.1 | ||
| 26 | * | ||
| 27 | * The coefficients of c(x) were generated by a pari-gp script using | ||
| 28 | * a Remez algorithm that searches for the best higher coefficients | ||
| 29 | * after rounding leading coefficients to a specified precision. | ||
| 30 | * | ||
| 31 | * Simpler methods like Chebyshev or basic Remez barely suffice for | ||
| 32 | * cos() in 64-bit precision, because we want the coefficient of x^2 | ||
| 33 | * to be precisely -0.5 so that multiplying by it is exact, and plain | ||
| 34 | * rounding of the coefficients of a good polynomial approximation only | ||
| 35 | * gives this up to about 64-bit precision. Plain rounding also gives | ||
| 36 | * a mediocre approximation for the coefficient of x^4, but a rounding | ||
| 37 | * error of 0.5 ulps for this coefficient would only contribute ~0.01 | ||
| 38 | * ulps to the final error, so this is unimportant. Rounding errors in | ||
| 39 | * higher coefficients are even less important. | ||
| 40 | * | ||
| 41 | * In fact, coefficients above the x^4 one only need to have 53-bit | ||
| 42 | * precision, and this is more efficient. We get this optimization | ||
| 43 | * almost for free from the complications needed to search for the best | ||
| 44 | * higher coefficients. | ||
| 45 | */ | ||
| 46 | static const long double | ||
| 47 | C1 = 0.0416666666666666666136L; /* 0xaaaaaaaaaaaaaa9b.0p-68 */ | ||
| 48 | static const double | ||
| 49 | C2 = -0.0013888888888888874, /* -0x16c16c16c16c10.0p-62 */ | ||
| 50 | C3 = 0.000024801587301571716, /* 0x1a01a01a018e22.0p-68 */ | ||
| 51 | C4 = -0.00000027557319215507120, /* -0x127e4fb7602f22.0p-74 */ | ||
| 52 | C5 = 0.0000000020876754400407278, /* 0x11eed8caaeccf1.0p-81 */ | ||
| 53 | C6 = -1.1470297442401303e-11, /* -0x19393412bd1529.0p-89 */ | ||
| 54 | C7 = 4.7383039476436467e-14; /* 0x1aac9d9af5c43e.0p-97 */ | ||
| 55 | #define POLY(z) (z*(C1+z*(C2+z*(C3+z*(C4+z*(C5+z*(C6+z*C7))))))) | ||
| 56 | #elif LDBL_MANT_DIG == 113 | ||
| 57 | /* | ||
| 58 | * ld128 version of __cos.c. See __cos.c for most comments. | ||
| 59 | */ | ||
| 60 | /* | ||
| 61 | * Domain [-0.7854, 0.7854], range ~[-1.80e-37, 1.79e-37]: | ||
| 62 | * |cos(x) - c(x))| < 2**-122.0 | ||
| 63 | * | ||
| 64 | * 113-bit precision requires more care than 64-bit precision, since | ||
| 65 | * simple methods give a minimax polynomial with coefficient for x^2 | ||
| 66 | * that is 1 ulp below 0.5, but we want it to be precisely 0.5. See | ||
| 67 | * above for more details. | ||
| 68 | */ | ||
| 69 | static const long double | ||
| 70 | C1 = 0.04166666666666666666666666666666658424671L, | ||
| 71 | C2 = -0.001388888888888888888888888888863490893732L, | ||
| 72 | C3 = 0.00002480158730158730158730158600795304914210L, | ||
| 73 | C4 = -0.2755731922398589065255474947078934284324e-6L, | ||
| 74 | C5 = 0.2087675698786809897659225313136400793948e-8L, | ||
| 75 | C6 = -0.1147074559772972315817149986812031204775e-10L, | ||
| 76 | C7 = 0.4779477332386808976875457937252120293400e-13L; | ||
| 77 | static const double | ||
| 78 | C8 = -0.1561920696721507929516718307820958119868e-15, | ||
| 79 | C9 = 0.4110317413744594971475941557607804508039e-18, | ||
| 80 | C10 = -0.8896592467191938803288521958313920156409e-21, | ||
| 81 | C11 = 0.1601061435794535138244346256065192782581e-23; | ||
| 82 | #define POLY(z) (z*(C1+z*(C2+z*(C3+z*(C4+z*(C5+z*(C6+z*(C7+ \ | ||
| 83 | 	z*(C8+z*(C9+z*(C10+z*C11))))))))))) | ||
| 84 | #endif | ||
| 85 | |||
| 86 | long double __cosl(long double x, long double y) | ||
| 87 | { | ||
| 88 | 	long double hz,z,r,w; | ||
| 89 | |||
| 90 | 	z = x*x; | ||
| 91 | 	r = POLY(z); | ||
| 92 | 	hz = 0.5*z; | ||
| 93 | 	w = 1.0-hz; | ||
| 94 | 	return w + (((1.0-w)-hz) + (z*r-x*y)); | ||
| 95 | } | ||
| 96 | #endif | ||
lib/libc/musl/src/math/cosl.c deleted-39| ... | @@ -1,39 +0,0 @@ | ||
| 1 | #include "libm.h" | ||
| 2 | |||
| 3 | #if LDBL_MANT_DIG == 53 && LDBL_MAX_EXP == 1024 | ||
| 4 | long double cosl(long double x) { | ||
| 5 | 	return cos(x); | ||
| 6 | } | ||
| 7 | #elif (LDBL_MANT_DIG == 64 || LDBL_MANT_DIG == 113) && LDBL_MAX_EXP == 16384 | ||
| 8 | long double cosl(long double x) | ||
| 9 | { | ||
| 10 | 	union ldshape u = {x}; | ||
| 11 | 	unsigned n; | ||
| 12 | 	long double y[2], hi, lo; | ||
| 13 | |||
| 14 | 	u.i.se &= 0x7fff; | ||
| 15 | 	if (u.i.se == 0x7fff) | ||
| 16 | 		return x - x; | ||
| 17 | 	x = u.f; | ||
| 18 | 	if (x < M_PI_4) { | ||
| 19 | 		if (u.i.se < 0x3fff - LDBL_MANT_DIG) | ||
| 20 | 			/* raise inexact if x!=0 */ | ||
| 21 | 			return 1.0 + x; | ||
| 22 | 		return __cosl(x, 0); | ||
| 23 | 	} | ||
| 24 | 	n = __rem_pio2l(x, y); | ||
| 25 | 	hi = y[0]; | ||
| 26 | 	lo = y[1]; | ||
| 27 | 	switch (n & 3) { | ||
| 28 | 	case 0: | ||
| 29 | 		return __cosl(hi, lo); | ||
| 30 | 	case 1: | ||
| 31 | 		return -__sinl(hi, lo, 1); | ||
| 32 | 	case 2: | ||
| 33 | 		return -__cosl(hi, lo); | ||
| 34 | 	case 3: | ||
| 35 | 	default: | ||
| 36 | 		return __sinl(hi, lo, 1); | ||
| 37 | 	} | ||
| 38 | } | ||
| 39 | #endif | ||
src/libs/mingw.zig-2| ... | @@ -940,8 +940,6 @@ const mingw32_x86_src = [_][]const u8{ | ... | @@ -940,8 +940,6 @@ const mingw32_x86_src = [_][]const u8{ |
| 940 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "atan2l.c", | 940 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "atan2l.c", |
| 941 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "atanhl.c", | 941 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "atanhl.c", |
| 942 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "atanl.c", | 942 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "atanl.c", |
| 943 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "cosl.c", | ||
| 944 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "cosl_internal.S", | ||
| 945 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "cossinl.c", | 943 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "cossinl.c", |
| 946 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "exp2l.S", | 944 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "exp2l.S", |
| 947 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "expl.c", | 945 | "math" ++ path.sep_str ++ "x86" ++ path.sep_str ++ "expl.c", |
src/libs/musl.zig-2| ... | @@ -811,8 +811,6 @@ const src_files = [_][]const u8{ | ... | @@ -811,8 +811,6 @@ const src_files = [_][]const u8{ |
| 811 | "musl/src/math/__cos.c", | 811 | "musl/src/math/__cos.c", |
| 812 | "musl/src/math/__cosdf.c", | 812 | "musl/src/math/__cosdf.c", |
| 813 | "musl/src/math/coshl.c", | 813 | "musl/src/math/coshl.c", |
| 814 | "musl/src/math/__cosl.c", | ||
| 815 | "musl/src/math/cosl.c", | ||
| 816 | "musl/src/math/erf.c", | 814 | "musl/src/math/erf.c", |
| 817 | "musl/src/math/erff.c", | 815 | "musl/src/math/erff.c", |
| 818 | "musl/src/math/erfl.c", | 816 | "musl/src/math/erfl.c", |
src/libs/wasi_libc.zig-2| ... | @@ -682,8 +682,6 @@ const libc_top_half_src_files = [_][]const u8{ | ... | @@ -682,8 +682,6 @@ const libc_top_half_src_files = [_][]const u8{ |
| 682 | "musl/src/math/__cos.c", | 682 | "musl/src/math/__cos.c", |
| 683 | "musl/src/math/__cosdf.c", | 683 | "musl/src/math/__cosdf.c", |
| 684 | "musl/src/math/coshl.c", | 684 | "musl/src/math/coshl.c", |
| 685 | "musl/src/math/__cosl.c", | ||
| 686 | "musl/src/math/cosl.c", | ||
| 687 | "musl/src/math/erf.c", | 685 | "musl/src/math/erf.c", |
| 688 | "musl/src/math/erff.c", | 686 | "musl/src/math/erff.c", |
| 689 | "musl/src/math/erfl.c", | 687 | "musl/src/math/erfl.c", |