rational/
extras.rs

1//! Contains some helper functions that can be useful and/or fun.
2
3use crate::Rational;
4
5/// Convenience method for constructing a simple `Rational`.
6///
7/// ## Example
8/// ```rust
9/// # use rational::{Rational, extras::*};
10/// let one_half = r(1, 2);
11/// assert_eq!(one_half, Rational::new(1, 2));
12/// ```
13pub fn r(n: i128, d: i128) -> Rational {
14    Rational::new(n, d)
15}
16
17/// Calculate the greatest common divisor of two numbers using [Stein's algorithm](https://en.wikipedia.org/wiki/Binary_GCD_algorithm).
18///
19/// ## Panics
20/// * If the result does not fit in the `i128` primitive type. This only happens in the cases below.
21/// ```rust,no_run
22/// # use rational::extras::*;
23/// // both of these are equal to `i128::MAX + 1`, which does not fit in an `i128`
24/// gcd(i128::MIN, 0);
25/// gcd(0, i128::MIN);
26/// gcd(i128::MIN, i128::MIN);
27/// ```
28///
29/// ## Example
30/// ```rust
31/// # use rational::extras::*;
32/// assert_eq!(gcd(9, 60), 3);
33/// assert_eq!(gcd(899, 957), 29);
34/// assert_eq!(gcd(-899, 957), 29);
35/// ```
36pub fn gcd(a: i128, b: i128) -> i128 {
37    gcd_checked(a, b)
38        .unwrap_or_else(|| panic!("the gcd of {} and {} is equal to i128::MAX+1, which does not fit in the i128 primitive type", a, b))
39}
40
41/// Calculate the greatest common divisor of two numbers using [Stein's algorithm](https://en.wikipedia.org/wiki/Binary_GCD_algorithm).
42///
43/// ## Returns
44/// * `Some` if the result fits in the `i128` primitive type.
45/// * `None` if the result does not fit in the `i128` primitive type. This only happens in the cases below.
46/// ```rust,no_run
47/// # use rational::extras::*;
48/// gcd(i128::MIN, 0);
49/// gcd(0, i128::MIN);
50/// gcd(i128::MIN, i128::MIN);
51/// ```
52///
53/// ## Example
54/// ```rust
55/// # use rational::extras::*;
56/// assert_eq!(gcd_checked(9, 60), Some(3));
57/// assert_eq!(gcd_checked(899, 957), Some(29));
58/// assert_eq!(gcd_checked(-899, 957), Some(29));
59/// assert_eq!(gcd_checked(i128::MIN, 0), None);
60/// assert_eq!(gcd_checked(0, i128::MIN), None);
61/// assert_eq!(gcd_checked(i128::MIN, i128::MIN), None);
62/// ```
63pub fn gcd_checked(mut a: i128, mut b: i128) -> Option<i128> {
64    if a == 0 || b == 0 {
65        return return_gcd_checked(a | b);
66    }
67
68    let factors_of_two = (a | b).trailing_zeros();
69
70    if a == i128::MIN || b == i128::MIN {
71        return return_gcd_checked(1 << factors_of_two);
72    }
73
74    a = a.abs() >> a.trailing_zeros();
75    b = b.abs() >> b.trailing_zeros();
76
77    while a != b {
78        if a > b {
79            a -= b;
80            a >>= a.trailing_zeros();
81        } else {
82            b -= a;
83            b >>= b.trailing_zeros();
84        }
85    }
86    Some(a << factors_of_two)
87}
88
89/// Internal helper function for returning the correct gcd.
90fn return_gcd_checked(g: i128) -> Option<i128> {
91    if g == i128::MIN {
92        None
93    } else {
94        Some(g.abs())
95    }
96}
97
98/// Calculate the least common multiple of two numbers.
99///
100/// ## Panics
101/// * If overflow occurred
102///
103/// ## Example
104/// ```rust
105/// # use rational::extras::*;
106/// assert_eq!(lcm(6, 4), 12);
107/// assert_eq!(lcm(6, 8), 24);
108/// assert_eq!(lcm(-6, 8), 24);
109/// assert_eq!(lcm(-6, -8), 24);
110/// ```
111pub fn lcm(a: i128, b: i128) -> i128 {
112    let g = gcd(a, b);
113    a.abs() * (b.abs() / g)
114}
115
116/// Calculate the least common multiple of two numbers, returning `None` if overflow occurred.
117///
118/// ## Example
119/// ```rust
120/// # use rational::extras::*;
121/// assert_eq!(lcm_checked(6, 4), Some(12));
122/// assert_eq!(lcm_checked(6, 8), Some(24));
123/// assert_eq!(lcm_checked(-6, 8), Some(24));
124/// assert_eq!(lcm_checked(i128::MAX, i128::MAX - 1), None);
125/// ```
126pub fn lcm_checked(a: i128, b: i128) -> Option<i128> {
127    let g = gcd(a, b);
128    a.abs().checked_mul(b.abs().checked_div(g)?)
129}
130
131/// Checks if `l` and `r` are [coprime](https://en.wikipedia.org/wiki/Coprime_integers). Shorthand for `gcd(l, r) == 1`.
132///
133/// ## Example
134/// ```rust
135/// # use rational::extras::*;
136/// assert!(is_coprime(8, 9));
137/// assert!(is_coprime(7, 9));
138/// assert!(!is_coprime(6, 9));
139/// assert!(is_coprime(-1, 1));
140/// ```
141pub fn is_coprime(l: i128, r: i128) -> bool {
142    gcd(l, r) == 1
143}
144
145/// Create a [continued fraction](https://en.wikipedia.org/wiki/Continued_fraction#Motivation_and_notation).
146///
147/// ## Notes
148/// The size of the numerator and denominator can grow quite quickly with increased lengths of `cont`,
149/// so be careful with integer overflow. If `cont` is empty, then the resulting `Rational` will be equal to `init`.
150///
151/// ## Example
152/// ```rust
153/// # use rational::extras::*;
154/// // to create a rational number that estimates the mathematical constant `e` (~2.7182818...)
155/// // we can use the continued fraction [2;1,2,1,1,4,1,1,6,1,1,8]
156/// let e = continued_fraction(2, &[1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8]);
157/// assert!((e - std::f64::consts::E).abs() < 0.0000001);
158/// ```
159pub fn continued_fraction(init: u64, cont: &[u64]) -> Rational {
160    let mut r = Rational::new(0, 1);
161
162    for &d in cont.iter().rev() {
163        r = Rational::new(1, r + d);
164    }
165
166    r + init
167}
168
169pub fn continued_fraction_iter(init: u64, cont: &[u64]) -> ContinuedFractionIter<'_> {
170    ContinuedFractionIter::new(Rational::new(init, 1), cont)
171}
172
173pub struct ContinuedFractionIter<'a> {
174    init: Rational,
175    fraction: Rational,
176    idx: usize,
177    cont: &'a [u64],
178}
179
180impl<'a> ContinuedFractionIter<'a> {
181    fn new(init: Rational, cont: &'a [u64]) -> Self {
182        Self {
183            init,
184            fraction: Rational::zero(),
185            idx: 0,
186            cont,
187        }
188    }
189
190    pub fn decimals(self) -> impl Iterator<Item = f64> + 'a {
191        self.map(|r| r.decimal_value())
192    }
193}
194
195impl Iterator for ContinuedFractionIter<'_> {
196    type Item = Rational;
197
198    fn next(&mut self) -> Option<Self::Item> {
199        if self.cont.is_empty() {
200            if self.idx == 0 {
201                self.idx += 1;
202                Some(self.init)
203            } else {
204                None
205            }
206        } else if self.idx > self.cont.len() {
207            None
208        } else if self.idx == self.cont.len() {
209            self.idx += 1;
210            Some(self.init + self.fraction)
211        } else {
212            let curr = self.init + self.fraction;
213            self.fraction =
214                Rational::new(1, self.fraction + self.cont[self.cont.len() - self.idx - 1]);
215            self.idx += 1;
216            Some(curr)
217        }
218    }
219}
220
221#[cfg(test)]
222mod tests {
223    use super::*;
224
225    #[test]
226    fn gcd_test() {
227        let eq = |a: i128, b: i128, g: i128| assert_eq!(gcd(a, b), g, "a: {}, b: {}", a, b);
228        eq(1, 2, 1);
229        eq(5, 4, 1);
230        eq(12, 4, 4);
231        eq(-74, 44, 2);
232        eq(-2, -4, 2);
233        eq(i128::MIN, 1, 1);
234    }
235
236    #[test]
237    #[should_panic]
238    fn gcd_should_panic_test_1() {
239        dbg!(gcd(i128::MIN, 0));
240    }
241
242    #[test]
243    #[should_panic]
244    fn gcd_should_panic_test_2() {
245        dbg!(gcd(0, i128::MIN));
246    }
247
248    #[test]
249    #[should_panic]
250    fn gcd_should_panic_test_3() {
251        dbg!(gcd(i128::MIN, i128::MIN));
252    }
253
254    #[test]
255    fn lcm_test() {
256        assert_eq!(lcm(2, 6), 6);
257        assert_eq!(lcm(1, 6), 6);
258    }
259
260    #[test]
261    fn repeated_test() {
262        use std::f64::consts::*;
263        let assert = |init: u64, cont: &[u64], expected: f64| {
264            let actual = continued_fraction(init, cont).decimal_value();
265            assert!(
266                (actual - expected).abs() < 0.000000001,
267                "actual: {}, expected: {}",
268                actual,
269                expected
270            );
271        };
272
273        assert(1, &[2; 15], 2.0_f64.sqrt());
274        assert(1, &[1; 100], 1.6180339887498);
275        assert(4, &[2, 1, 3, 1, 2, 8, 2, 1, 3, 1, 2, 8], 19.0_f64.sqrt());
276        assert(1, &[], 1.0);
277        assert(4, &[], 4.0);
278        assert(2, &[1, 2, 1, 1, 4, 1, 1, 6, 1, 1, 8, 1, 1, 10], E);
279        assert(3, &[7, 15, 1, 292, 1, 1, 1, 2, 1, 3, 1], PI);
280    }
281
282    #[test]
283    fn continued_fraction_iter_test() {
284        let decimal_approx: Vec<_> = continued_fraction_iter(2, &[1, 2, 1, 1, 4, 1, 1, 6])
285            .decimals()
286            .collect();
287        assert_eq!(
288            decimal_approx,
289            vec![
290                2.0,
291                2.1666666666666665,
292                2.857142857142857,
293                2.5384615384615383,
294                2.2203389830508473,
295                2.8194444444444446,
296                2.549618320610687,
297                2.3922155688622753,
298                2.718279569892473,
299            ]
300        );
301    }
302}