zig

fork of https://codeberg.org/ziglang/zig
Log | Files | Refs | README | LICENSE

tan.zig (6175B) - Raw


      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/tanf.c
      5 //! https://git.musl-libc.org/cgit/musl/tree/src/math/tan.c
      6 //! https://golang.org/src/math/tan.go
      7 
      8 const std = @import("std");
      9 const builtin = @import("builtin");
     10 const math = std.math;
     11 const mem = std.mem;
     12 const expect = std.testing.expect;
     13 
     14 const kernel = @import("trig.zig");
     15 const rem_pio2 = @import("rem_pio2.zig").rem_pio2;
     16 const rem_pio2f = @import("rem_pio2f.zig").rem_pio2f;
     17 
     18 const arch = builtin.cpu.arch;
     19 const common = @import("common.zig");
     20 
     21 pub const panic = common.panic;
     22 
     23 comptime {
     24     @export(&__tanh, .{ .name = "__tanh", .linkage = common.linkage, .visibility = common.visibility });
     25     @export(&tanf, .{ .name = "tanf", .linkage = common.linkage, .visibility = common.visibility });
     26     @export(&tan, .{ .name = "tan", .linkage = common.linkage, .visibility = common.visibility });
     27     @export(&__tanx, .{ .name = "__tanx", .linkage = common.linkage, .visibility = common.visibility });
     28     if (common.want_ppc_abi) {
     29         @export(&tanq, .{ .name = "tanf128", .linkage = common.linkage, .visibility = common.visibility });
     30     }
     31     @export(&tanq, .{ .name = "tanq", .linkage = common.linkage, .visibility = common.visibility });
     32     @export(&tanl, .{ .name = "tanl", .linkage = common.linkage, .visibility = common.visibility });
     33 }
     34 
     35 pub fn __tanh(x: f16) callconv(.c) f16 {
     36     // TODO: more efficient implementation
     37     return @floatCast(tanf(x));
     38 }
     39 
     40 pub fn tanf(x: f32) callconv(.c) f32 {
     41     // Small multiples of pi/2 rounded to double precision.
     42     const t1pio2: f64 = 1.0 * math.pi / 2.0; // 0x3FF921FB, 0x54442D18
     43     const t2pio2: f64 = 2.0 * math.pi / 2.0; // 0x400921FB, 0x54442D18
     44     const t3pio2: f64 = 3.0 * math.pi / 2.0; // 0x4012D97C, 0x7F3321D2
     45     const t4pio2: f64 = 4.0 * math.pi / 2.0; // 0x401921FB, 0x54442D18
     46 
     47     var ix: u32 = @bitCast(x);
     48     const sign = ix >> 31 != 0;
     49     ix &= 0x7fffffff;
     50 
     51     if (ix <= 0x3f490fda) { // |x| ~<= pi/4
     52         if (ix < 0x39800000) { // |x| < 2**-12
     53             // raise inexact if x!=0 and underflow if subnormal
     54             if (common.want_float_exceptions) {
     55                 if (ix < 0x00800000) {
     56                     mem.doNotOptimizeAway(x / 0x1p120);
     57                 } else {
     58                     mem.doNotOptimizeAway(x + 0x1p120);
     59                 }
     60             }
     61             return x;
     62         }
     63         return kernel.__tandf(x, false);
     64     }
     65     if (ix <= 0x407b53d1) { // |x| ~<= 5*pi/4
     66         if (ix <= 0x4016cbe3) { // |x| ~<= 3pi/4
     67             return kernel.__tandf((if (sign) x + t1pio2 else x - t1pio2), true);
     68         } else {
     69             return kernel.__tandf((if (sign) x + t2pio2 else x - t2pio2), false);
     70         }
     71     }
     72     if (ix <= 0x40e231d5) { // |x| ~<= 9*pi/4
     73         if (ix <= 0x40afeddf) { // |x| ~<= 7*pi/4
     74             return kernel.__tandf((if (sign) x + t3pio2 else x - t3pio2), true);
     75         } else {
     76             return kernel.__tandf((if (sign) x + t4pio2 else x - t4pio2), false);
     77         }
     78     }
     79 
     80     // tan(Inf or NaN) is NaN
     81     if (ix >= 0x7f800000) {
     82         return x - x;
     83     }
     84 
     85     var y: f64 = undefined;
     86     const n = rem_pio2f(x, &y);
     87     return kernel.__tandf(y, n & 1 != 0);
     88 }
     89 
     90 pub fn tan(x: f64) callconv(.c) f64 {
     91     var ix = @as(u64, @bitCast(x)) >> 32;
     92     ix &= 0x7fffffff;
     93 
     94     // |x| ~< pi/4
     95     if (ix <= 0x3fe921fb) {
     96         if (ix < 0x3e400000) { // |x| < 2**-27
     97             // raise inexact if x!=0 and underflow if subnormal
     98             if (common.want_float_exceptions) {
     99                 if (ix < 0x00100000) {
    100                     mem.doNotOptimizeAway(x / 0x1p120);
    101                 } else {
    102                     mem.doNotOptimizeAway(x + 0x1p120);
    103                 }
    104             }
    105             return x;
    106         }
    107         return kernel.__tan(x, 0.0, false);
    108     }
    109 
    110     // tan(Inf or NaN) is NaN
    111     if (ix >= 0x7ff00000) {
    112         return x - x;
    113     }
    114 
    115     var y: [2]f64 = undefined;
    116     const n = rem_pio2(x, &y);
    117     return kernel.__tan(y[0], y[1], n & 1 != 0);
    118 }
    119 
    120 pub fn __tanx(x: f80) callconv(.c) f80 {
    121     // TODO: more efficient implementation
    122     return @floatCast(tanq(x));
    123 }
    124 
    125 pub fn tanq(x: f128) callconv(.c) f128 {
    126     // TODO: more correct implementation
    127     return tan(@floatCast(x));
    128 }
    129 
    130 pub fn tanl(x: c_longdouble) callconv(.c) c_longdouble {
    131     switch (@typeInfo(c_longdouble).float.bits) {
    132         16 => return __tanh(x),
    133         32 => return tanf(x),
    134         64 => return tan(x),
    135         80 => return __tanx(x),
    136         128 => return tanq(x),
    137         else => @compileError("unreachable"),
    138     }
    139 }
    140 
    141 test "tan" {
    142     try expect(tan(@as(f32, 0.0)) == tanf(0.0));
    143     try expect(tan(@as(f64, 0.0)) == tan(0.0));
    144 }
    145 
    146 test "tan32" {
    147     const epsilon = 0.00001;
    148 
    149     try expect(math.approxEqAbs(f32, tanf(0.0), 0.0, epsilon));
    150     try expect(math.approxEqAbs(f32, tanf(0.2), 0.202710, epsilon));
    151     try expect(math.approxEqAbs(f32, tanf(0.8923), 1.240422, epsilon));
    152     try expect(math.approxEqAbs(f32, tanf(1.5), 14.101420, epsilon));
    153     try expect(math.approxEqAbs(f32, tanf(37.45), -0.254397, epsilon));
    154     try expect(math.approxEqAbs(f32, tanf(89.123), 2.285852, epsilon));
    155 }
    156 
    157 test "tan64" {
    158     const epsilon = 0.000001;
    159 
    160     try expect(math.approxEqAbs(f64, tan(0.0), 0.0, epsilon));
    161     try expect(math.approxEqAbs(f64, tan(0.2), 0.202710, epsilon));
    162     try expect(math.approxEqAbs(f64, tan(0.8923), 1.240422, epsilon));
    163     try expect(math.approxEqAbs(f64, tan(1.5), 14.101420, epsilon));
    164     try expect(math.approxEqAbs(f64, tan(37.45), -0.254397, epsilon));
    165     try expect(math.approxEqAbs(f64, tan(89.123), 2.2858376, epsilon));
    166 }
    167 
    168 test "tan32.special" {
    169     try expect(tanf(0.0) == 0.0);
    170     try expect(tanf(-0.0) == -0.0);
    171     try expect(math.isNan(tanf(math.inf(f32))));
    172     try expect(math.isNan(tanf(-math.inf(f32))));
    173     try expect(math.isNan(tanf(math.nan(f32))));
    174 }
    175 
    176 test "tan64.special" {
    177     try expect(tan(0.0) == 0.0);
    178     try expect(tan(-0.0) == -0.0);
    179     try expect(math.isNan(tan(math.inf(f64))));
    180     try expect(math.isNan(tan(-math.inf(f64))));
    181     try expect(math.isNan(tan(math.nan(f64))));
    182 }