fpmath 0.1.1

A floating point math library
Documentation 
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83

use super::log::log_hi_lo_inner;
use super::Log;
use crate::double::SemiDouble;
use crate::traits::Int as _;

pub(crate) fn atanh<F: Log>(x: F) -> F {
    let e = x.raw_exp();
    let absx = x.abs();
    if (e == F::MAX_RAW_EXP && x.raw_mant() != F::Raw::ZERO)
        || e <= (F::EXP_OFFSET - F::RawExp::from(F::MANT_BITS))
    {
        // propagate NaN
        // or
        // very small, includes subnormal and zero
        // atanh(x) ~= x
        // also handles atanh(-0) = -0
        x
    } else if absx == F::one() {
        // atanh(±1) = ±inf
        F::INFINITY.copysign(x)
    } else if x.abs() > F::one() {
        // |x| > 1, return NaN
        F::NAN
    } else {
        atanh_inner(x)
    }
}

fn atanh_inner<F: Log>(x: F) -> F {
    // t1 = 2 * x / (1 - x)
    let t1 = SemiDouble::new(F::two() * x) / SemiDouble::new_qsub11(F::one(), x);

    // t2 = (1 + x) / (1 - x) = t1 + 1
    let t2 = t1 + F::one();
    let t2 = t2.to_norm();

    // atanh(x) = 0.5 * log((1 + x) / (1 - x))
    F::half() * log_hi_lo_inner(t2.hi(), t2.lo())
}

#[cfg(test)]
mod tests {
    use crate::traits::Float;
    use crate::FloatMath;

    fn test<F: Float + FloatMath>() {
        use crate::atanh;

        let f = F::parse;

        assert_is_nan!(atanh(F::NAN));
        assert_is_nan!(atanh(f("1.5")));
        assert_is_nan!(atanh(f("-1.5")));
        assert_is_nan!(atanh(F::INFINITY));
        assert_is_nan!(atanh(F::neg_infinity()));
        assert_total_eq!(atanh(F::ZERO), F::ZERO);
        assert_total_eq!(atanh(-F::ZERO), -F::ZERO);
        assert_total_eq!(atanh(F::one()), F::INFINITY);
        assert_total_eq!(atanh(-F::one()), F::neg_infinity());
    }

    #[test]
    fn test_f32() {
        test::<f32>();
    }

    #[cfg(feature = "soft-float")]
    #[test]
    fn test_soft_f32() {
        test::<crate::SoftF32>();
    }

    #[test]
    fn test_f64() {
        test::<f64>();
    }

    #[cfg(feature = "soft-float")]
    #[test]
    fn test_soft_f64() {
        test::<crate::SoftF64>();
    }
}