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}