Wilco Dijkstra | 269dc16 | 2018-05-16 15:39:22 +0100 | [diff] [blame] | 1 | /* |
| 2 | * Single-precision cos function. |
| 3 | * |
| 4 | * Copyright (c) 2018, Arm Limited. |
Szabolcs Nagy | 11253b0 | 2018-11-12 11:10:57 +0000 | [diff] [blame] | 5 | * SPDX-License-Identifier: MIT |
Wilco Dijkstra | 269dc16 | 2018-05-16 15:39:22 +0100 | [diff] [blame] | 6 | */ |
| 7 | |
Wilco Dijkstra | 269dc16 | 2018-05-16 15:39:22 +0100 | [diff] [blame] | 8 | #include <stdint.h> |
| 9 | #include <math.h> |
| 10 | #include "math_config.h" |
| 11 | #include "sincosf.h" |
| 12 | |
Wilco Dijkstra | b2fc989 | 2018-08-08 15:03:29 +0100 | [diff] [blame] | 13 | /* Fast cosf implementation. Worst-case ULP is 0.5607, maximum relative |
| 14 | error is 0.5303 * 2^-23. A single-step range reduction is used for |
Wilco Dijkstra | 269dc16 | 2018-05-16 15:39:22 +0100 | [diff] [blame] | 15 | small values. Large inputs have their range reduced using fast integer |
Wilco Dijkstra | b2fc989 | 2018-08-08 15:03:29 +0100 | [diff] [blame] | 16 | arithmetic. */ |
Wilco Dijkstra | 269dc16 | 2018-05-16 15:39:22 +0100 | [diff] [blame] | 17 | float |
| 18 | cosf (float y) |
| 19 | { |
| 20 | double x = y; |
| 21 | double s; |
| 22 | int n; |
Wilco Dijkstra | 3262ef2 | 2018-07-04 17:45:15 +0100 | [diff] [blame] | 23 | const sincos_t *p = &__sincosf_table[0]; |
Wilco Dijkstra | 269dc16 | 2018-05-16 15:39:22 +0100 | [diff] [blame] | 24 | |
| 25 | if (abstop12 (y) < abstop12 (pio4)) |
| 26 | { |
| 27 | double x2 = x * x; |
| 28 | |
| 29 | if (unlikely (abstop12 (y) < abstop12 (0x1p-12f))) |
| 30 | return 1.0f; |
| 31 | |
| 32 | return sinf_poly (x, x2, p, 1); |
| 33 | } |
| 34 | else if (likely (abstop12 (y) < abstop12 (120.0f))) |
| 35 | { |
| 36 | x = reduce_fast (x, p, &n); |
| 37 | |
| 38 | /* Setup the signs for sin and cos. */ |
| 39 | s = p->sign[n & 3]; |
| 40 | |
| 41 | if (n & 2) |
Wilco Dijkstra | 3262ef2 | 2018-07-04 17:45:15 +0100 | [diff] [blame] | 42 | p = &__sincosf_table[1]; |
Wilco Dijkstra | 269dc16 | 2018-05-16 15:39:22 +0100 | [diff] [blame] | 43 | |
| 44 | return sinf_poly (x * s, x * x, p, n ^ 1); |
| 45 | } |
| 46 | else if (abstop12 (y) < abstop12 (INFINITY)) |
| 47 | { |
| 48 | uint32_t xi = asuint (y); |
| 49 | int sign = xi >> 31; |
| 50 | |
| 51 | x = reduce_large (xi, &n); |
| 52 | |
| 53 | /* Setup signs for sin and cos - include original sign. */ |
| 54 | s = p->sign[(n + sign) & 3]; |
| 55 | |
| 56 | if ((n + sign) & 2) |
Wilco Dijkstra | 3262ef2 | 2018-07-04 17:45:15 +0100 | [diff] [blame] | 57 | p = &__sincosf_table[1]; |
Wilco Dijkstra | 269dc16 | 2018-05-16 15:39:22 +0100 | [diff] [blame] | 58 | |
| 59 | return sinf_poly (x * s, x * x, p, n ^ 1); |
| 60 | } |
| 61 | else |
| 62 | return __math_invalidf (y); |
| 63 | } |