calendrical_calculations/
helpers.rs

1// This file is part of ICU4X.
2//
3// The contents of this file implement algorithms from Calendrical Calculations
4// by Reingold & Dershowitz, Cambridge University Press, 4th edition (2018),
5// which have been released as Lisp code at <https://github.com/EdReingold/calendar-code2/>
6// under the Apache-2.0 license. Accordingly, this file is released under
7// the Apache License, Version 2.0 which can be found at the calendrical_calculations
8// package root or at http://www.apache.org/licenses/LICENSE-2.0.
9
10use crate::astronomy::Location;
11use crate::rata_die::{Moment, RataDie};
12#[allow(unused_imports)]
13use core_maths::*;
14
15pub(crate) trait IntegerRoundings {
16    fn div_ceil(self, rhs: Self) -> Self;
17}
18
19impl IntegerRoundings for i64 {
20    // copied from std
21    fn div_ceil(self, rhs: Self) -> Self {
22        let d = self / rhs;
23        let r = self % rhs;
24        if (r > 0 && rhs > 0) || (r < 0 && rhs < 0) {
25            d + 1
26        } else {
27            d
28        }
29    }
30}
31
32#[test]
33fn test_div_ceil() {
34    assert_eq!(IntegerRoundings::div_ceil(i64::MIN, 1), i64::MIN);
35    assert_eq!(
36        IntegerRoundings::div_ceil(i64::MIN, 2),
37        -4611686018427387904
38    );
39    assert_eq!(
40        IntegerRoundings::div_ceil(i64::MIN, 3),
41        -3074457345618258602
42    );
43
44    assert_eq!(IntegerRoundings::div_ceil(-10, 1), -10);
45    assert_eq!(IntegerRoundings::div_ceil(-10, 2), -5);
46    assert_eq!(IntegerRoundings::div_ceil(-10, 3), -3);
47
48    assert_eq!(IntegerRoundings::div_ceil(-9, 1), -9);
49    assert_eq!(IntegerRoundings::div_ceil(-9, 2), -4);
50    assert_eq!(IntegerRoundings::div_ceil(-9, 3), -3);
51
52    assert_eq!(IntegerRoundings::div_ceil(-8, 1), -8);
53    assert_eq!(IntegerRoundings::div_ceil(-8, 2), -4);
54    assert_eq!(IntegerRoundings::div_ceil(-8, 3), -2);
55
56    assert_eq!(IntegerRoundings::div_ceil(-2, 1), -2);
57    assert_eq!(IntegerRoundings::div_ceil(-2, 2), -1);
58    assert_eq!(IntegerRoundings::div_ceil(-2, 3), 0);
59
60    assert_eq!(IntegerRoundings::div_ceil(-1, 1), -1);
61    assert_eq!(IntegerRoundings::div_ceil(-1, 2), 0);
62    assert_eq!(IntegerRoundings::div_ceil(-1, 3), 0);
63
64    assert_eq!(IntegerRoundings::div_ceil(0, 1), 0);
65    assert_eq!(IntegerRoundings::div_ceil(0, 2), 0);
66    assert_eq!(IntegerRoundings::div_ceil(0, 3), 0);
67
68    assert_eq!(IntegerRoundings::div_ceil(1, 1), 1);
69    assert_eq!(IntegerRoundings::div_ceil(1, 2), 1);
70    assert_eq!(IntegerRoundings::div_ceil(1, 3), 1);
71
72    assert_eq!(IntegerRoundings::div_ceil(2, 1), 2);
73    assert_eq!(IntegerRoundings::div_ceil(2, 2), 1);
74    assert_eq!(IntegerRoundings::div_ceil(2, 3), 1);
75
76    assert_eq!(IntegerRoundings::div_ceil(8, 1), 8);
77    assert_eq!(IntegerRoundings::div_ceil(8, 2), 4);
78    assert_eq!(IntegerRoundings::div_ceil(8, 3), 3);
79
80    assert_eq!(IntegerRoundings::div_ceil(9, 1), 9);
81    assert_eq!(IntegerRoundings::div_ceil(9, 2), 5);
82    assert_eq!(IntegerRoundings::div_ceil(9, 3), 3);
83
84    assert_eq!(IntegerRoundings::div_ceil(10, 1), 10);
85    assert_eq!(IntegerRoundings::div_ceil(10, 2), 5);
86    assert_eq!(IntegerRoundings::div_ceil(10, 3), 4);
87
88    assert_eq!(IntegerRoundings::div_ceil(i64::MAX, 1), 9223372036854775807);
89    assert_eq!(IntegerRoundings::div_ceil(i64::MAX, 2), 4611686018427387904);
90    assert_eq!(IntegerRoundings::div_ceil(i64::MAX, 3), 3074457345618258603);
91
92    for n in -100..100 {
93        for d in 1..5 {
94            let x1 = IntegerRoundings::div_ceil(n, d);
95            let x2 = (n as f64 / d as f64).ceil();
96            assert_eq!(x1, x2 as i64);
97        }
98    }
99}
100
101// Highest power is *last*
102pub(crate) fn poly(x: f64, coeffs: &[f64]) -> f64 {
103    coeffs.iter().rev().fold(0.0, |a, c| a * x + c)
104}
105
106// A generic function that finds a value within an interval
107// where a certain condition is satisfied.
108pub(crate) fn binary_search(
109    mut l: f64,
110    mut h: f64,
111    test: impl Fn(f64) -> bool,
112    epsilon: f64,
113) -> f64 {
114    debug_assert!(l < h);
115
116    loop {
117        let mid = l + (h - l) / 2.0;
118
119        // Determine which direction to go. `test` returns true if we need to go left.
120        (l, h) = if test(mid) { (l, mid) } else { (mid, h) };
121
122        // When the interval width reaches `epsilon`, return the current midpoint.
123        if (h - l) < epsilon {
124            return mid;
125        }
126    }
127}
128
129pub(crate) fn next_moment<F>(mut index: Moment, location: Location, condition: F) -> RataDie
130where
131    F: Fn(Moment, Location) -> bool,
132{
133    loop {
134        if condition(index, location) {
135            return index.as_rata_die();
136        }
137        index += 1.0;
138    }
139}
140
141pub(crate) fn next<F>(mut index: RataDie, condition: F) -> RataDie
142where
143    F: Fn(RataDie) -> bool,
144{
145    loop {
146        if condition(index) {
147            return index;
148        }
149        index += 1;
150    }
151}
152
153pub(crate) fn next_u8<F>(mut index: u8, condition: F) -> u8
154where
155    F: Fn(u8) -> bool,
156{
157    loop {
158        if condition(index) {
159            return index;
160        }
161        index += 1;
162    }
163}
164
165// "Final" is a reserved keyword in rust, which explains the naming convention here.
166pub(crate) fn final_func<F>(mut index: i32, condition: F) -> i32
167where
168    F: Fn(i32) -> bool,
169{
170    while condition(index) {
171        index += 1;
172    }
173    index - 1
174}
175
176#[test]
177fn test_binary_search() {
178    struct TestCase {
179        test_fn: fn(f64) -> bool,
180        range: (f64, f64),
181        expected: f64,
182    }
183
184    let test_cases = [
185        TestCase {
186            test_fn: |x: f64| x >= 4.0,
187            range: (0.0, 10.0),
188            expected: 4.0,
189        },
190        TestCase {
191            test_fn: |x: f64| x * x >= 2.0,
192            range: (0.0, 2.0),
193            expected: 2.0f64.sqrt(),
194        },
195        TestCase {
196            test_fn: |x: f64| x >= -4.0,
197            range: (-10.0, 0.0),
198            expected: -4.0,
199        },
200        TestCase {
201            test_fn: |x: f64| x >= 0.0,
202            range: (0.0, 10.0),
203            expected: 0.0,
204        },
205        TestCase {
206            test_fn: |x: f64| x > 10.0,
207            range: (0.0, 10.0),
208            expected: 10.0,
209        },
210    ];
211
212    for case in test_cases {
213        let result = binary_search(case.range.0, case.range.1, case.test_fn, 1e-4);
214        assert!((result - case.expected).abs() < 0.0001);
215    }
216}
217
218pub(crate) fn invert_angular<F: Fn(f64) -> f64>(f: F, y: f64, r: (f64, f64)) -> f64 {
219    binary_search(r.0, r.1, |x| (f(x) - y).rem_euclid(360.0) < 180.0, 1e-5)
220}
221
222#[test]
223fn test_invert_angular() {
224    struct TestCase {
225        f: fn(f64) -> f64,
226        y: f64,
227        r: (f64, f64),
228        expected: f64,
229    }
230
231    fn f1(x: f64) -> f64 {
232        (2.0 * x).rem_euclid(360.0)
233    }
234
235    fn f2(x: f64) -> f64 {
236        (3.0 * x).rem_euclid(360.0)
237    }
238
239    fn f3(x: f64) -> f64 {
240        (x).rem_euclid(360.0)
241    }
242    // tolerance for comparing floating points.
243    let tolerance = 1e-5;
244
245    let test_cases = [
246        TestCase {
247            f: f1,
248            y: 4.0,
249            r: (0.0, 10.0),
250            expected: 4.0,
251        },
252        TestCase {
253            f: f2,
254            y: 6.0,
255            r: (0.0, 20.0),
256            expected: 6.0,
257        },
258        TestCase {
259            f: f3,
260            y: 400.0,
261            r: (0.0, 10.0),
262            expected: 10.0,
263        },
264        TestCase {
265            f: f3,
266            y: 0.0,
267            r: (0.0, 10.0),
268            expected: 0.0,
269        },
270        TestCase {
271            f: f3,
272            y: 10.0,
273            r: (0.0, 10.0),
274            expected: 10.0,
275        },
276    ];
277
278    for case in test_cases {
279        let x = invert_angular(case.f, case.y, case.r);
280        assert!((((case.f)(x)).rem_euclid(360.0) - case.expected).abs() < tolerance);
281    }
282}
283
284/// Error returned when casting from an i32
285#[derive(Copy, Clone, Debug, displaydoc::Display)]
286#[allow(clippy::exhaustive_enums)] // enum is specific to function and has a closed set of possible values
287pub enum I32CastError {
288    /// Less than i32::MIN
289    BelowMin,
290    /// Greater than i32::MAX
291    AboveMax,
292}
293
294impl core::error::Error for I32CastError {}
295
296impl I32CastError {
297    pub(crate) const fn saturate(self) -> i32 {
298        match self {
299            I32CastError::BelowMin => i32::MIN,
300            I32CastError::AboveMax => i32::MAX,
301        }
302    }
303}
304
305/// Convert an i64 to i32 and with information on which way it was out of bounds if so
306#[inline]
307pub const fn i64_to_i32(input: i64) -> Result<i32, I32CastError> {
308    if input < i32::MIN as i64 {
309        Err(I32CastError::BelowMin)
310    } else if input > i32::MAX as i64 {
311        Err(I32CastError::AboveMax)
312    } else {
313        Ok(input as i32)
314    }
315}
316
317/// Convert an i64 to i32 but saturate at th ebounds
318#[inline]
319pub(crate) fn i64_to_saturated_i32(input: i64) -> i32 {
320    i64_to_i32(input).unwrap_or_else(|i| i.saturate())
321}
322
323#[test]
324fn test_i64_to_saturated_i32() {
325    assert_eq!(i64_to_saturated_i32(i64::MIN), i32::MIN);
326    assert_eq!(i64_to_saturated_i32(-2147483649), -2147483648);
327    assert_eq!(i64_to_saturated_i32(-2147483648), -2147483648);
328    assert_eq!(i64_to_saturated_i32(-2147483647), -2147483647);
329    assert_eq!(i64_to_saturated_i32(-2147483646), -2147483646);
330    assert_eq!(i64_to_saturated_i32(-100), -100);
331    assert_eq!(i64_to_saturated_i32(0), 0);
332    assert_eq!(i64_to_saturated_i32(100), 100);
333    assert_eq!(i64_to_saturated_i32(2147483646), 2147483646);
334    assert_eq!(i64_to_saturated_i32(2147483647), 2147483647);
335    assert_eq!(i64_to_saturated_i32(2147483648), 2147483647);
336    assert_eq!(i64_to_saturated_i32(2147483649), 2147483647);
337    assert_eq!(i64_to_saturated_i32(i64::MAX), i32::MAX);
338}