num_irrational/cont_frac/
general.rs

1//! Implementation of general continued fractions
2
3use super::block::Block;
4use super::infinite::InfiniteContinuedFraction;
5use crate::traits::{Approximation, Computable};
6use num_integer::Integer;
7use num_rational::Ratio;
8use num_traits::{CheckedAdd, CheckedMul, FromPrimitive, Num, NumRef, RefNum, Signed, ToPrimitive};
9
10/// This struct defines utility functions for generalized continued fraction number
11/// `b_1 + a_2 / (b_2 + a_3 / (b_3 + a_4 / .. ))`.
12/// 
13/// It can be converted from any iterator that returns a pair of number (by using [from()][GeneralContinuedFraction::from]).
14/// The first value will be regarded as a_k while the second value as b_k. You need to make
15/// sure that a_1 = 1.
16/// 
17#[derive(Clone, Copy)]
18pub struct GeneralContinuedFraction<I>(I);
19
20impl<I: Iterator> From<I> for GeneralContinuedFraction<I> {
21    fn from(iter: I) -> Self {
22        Self(iter)
23    }
24}
25impl<I: Iterator<Item = (T, T)>, T> IntoIterator for GeneralContinuedFraction<I> {
26    type Item = (T, T);
27    type IntoIter = I;    
28    fn into_iter(self) -> Self::IntoIter {
29        self.0
30    }
31}
32
33impl<I: Iterator<Item = (T, T)>, T: Integer + NumRef> GeneralContinuedFraction<I>
34where
35    for<'r> &'r T: RefNum<T>,
36{
37    /// Compute the convergents of the generalized continued fraction
38    pub fn convergents(self) -> Convergents<I, T> {
39        Convergents {
40            block: Block::identity(),
41            g_coeffs: self.0,
42        }
43    }
44
45    /// Simplify the generalized continued fraction to an `InfiniteContinuedFraction`
46    /// Note that usually this function will generate a simple continued fraction with most
47    /// coefficients being positive. If you want to ensure the positiveness, you can either
48    /// call collect() on `InfiniteContinuedFraction`, or call simplify() again.
49    pub fn simplify(self) -> InfiniteContinuedFraction<Simplified<I, T>> {
50        Simplified {
51            block: Block::identity(),
52            g_coeffs: self.0,
53        }.into()
54    }
55
56    /// Retrieve the decimal representation of the number, as an iterator of digits.
57    /// The iterator will stop if the capacity of T is reached
58    pub fn decimals(self) -> DecimalDigits<I, T> {
59        DecimalDigits {
60            block: Block::identity(),
61            g_coeffs: self.0,
62        }
63    }
64
65    // XXX: we can also implement homo and bihomo function on general continued fraction
66    //      however the result will still be an InfiniteContinuedFraction
67}
68
69// pub trait GeneralContinuedFraction<T: Integer + NumRef>: Iterator<Item = (T, T)>
70// where
71//     for<'r> &'r T: RefNum<T>,
72// {
73//     /// Compute the convergents of the generalized continued fraction
74//     fn convergents(self) -> Convergents<Self, T>;
75
76//     /// Simplify the generalized continued fraction to an `InfiniteContinuedFraction`
77//     /// Note that usually this function will generate a simple continued fraction with most
78//     /// coefficients being positive. If you want to ensure the positiveness, you can either
79//     /// call collect() on `InfiniteContinuedFraction`, or call simplify() again.
80//     fn simplify(self) -> InfiniteContinuedFraction<Simplified<Self, T>>
81//     where
82//         Self: Sized;
83
84//     /// Retrieve the decimal representation of the number, as an iterator of digits.
85//     /// The iterator will stop if the capacity of T is reached
86//     fn decimals(self) -> DecimalDigits<Self, T>;
87
88//     // XXX: we can also implement homo and bihomo function on general continued fraction
89//     //      however the result will still be an InfiniteContinuedFraction
90// }
91
92/// Iterator of convergents of [GeneralContinuedFraction]
93#[derive(Debug, Clone)]
94pub struct Convergents<I: Iterator<Item = (T, T)> + ?Sized, T> {
95    block: Block<T>,
96    g_coeffs: I,
97}
98
99/// Iterator of [GeneralContinuedFraction::simplify()] result
100#[derive(Debug, Clone)]
101pub struct Simplified<I: Iterator<Item = (T, T)> + ?Sized, T> {
102    block: Block<T>,
103    g_coeffs: I,
104}
105
106/// Iterator of [GeneralContinuedFraction::decimals()] result
107#[derive(Debug, Clone)]
108pub struct DecimalDigits<I: Iterator<Item = (T, T)> + ?Sized, T> {
109    block: Block<T>,
110    g_coeffs: I,
111}
112
113impl<I: Iterator<Item = (T, T)>, T: Integer + NumRef + CheckedAdd + CheckedMul + Clone> Iterator
114    for Convergents<I, T>
115where
116    for<'r> &'r T: RefNum<T>,
117{
118    type Item = Ratio<T>;
119
120    fn next(&mut self) -> Option<Ratio<T>> {
121        let (a, b) = self.g_coeffs.next()?;
122        let (p, q) = self.block.checked_gmove(a, b)?;
123        self.block.update(p.clone(), q.clone());
124
125        Some(Ratio::new(p, q))
126    }
127}
128
129impl<I: Iterator<Item = (T, T)>, T: Integer + NumRef> Iterator for Simplified<I, T>
130where
131    for<'r> &'r T: RefNum<T>,
132{
133    type Item = T;
134
135    fn next(&mut self) -> Option<T> {
136        loop {
137            match self.block.reduce_recip() {
138                Some(i) => break Some(i),
139                None => match self.g_coeffs.next() {
140                    Some((a, b)) => self.block.gmove(a, b),
141                    None => break None,
142                },
143            }
144        }
145    }
146}
147
148impl<
149        I: Iterator<Item = (T, T)>,
150        T: Integer + NumRef + CheckedAdd + CheckedMul + FromPrimitive + ToPrimitive,
151    > Iterator for DecimalDigits<I, T>
152where
153    for<'r> &'r T: RefNum<T>,
154{
155    type Item = u8;
156
157    fn next(&mut self) -> Option<u8> {
158        // TODO: print point and handle signed & integer > 10 case
159        loop {
160            match self.block.reduce_mul(T::from_u8(10).unwrap()) {
161                Some(i) => break Some(i.to_u8().unwrap() + b'0'),
162                None => match self.g_coeffs.next() {
163                    Some((a, b)) => {
164                        let (p, q) = self.block.checked_gmove(a, b)?;
165                        self.block.update(p, q);
166                    }
167                    None => break None,
168                },
169            }
170        }
171    }
172}
173
174impl<I: Iterator<Item = (T, T)> + Clone, T: Integer + NumRef + Clone + CheckedAdd + CheckedMul>
175    Computable<T> for GeneralContinuedFraction<I>
176where
177    for<'r> &'r T: RefNum<T>,
178{
179    fn approximated(&self, limit: &T) -> Approximation<Ratio<T>> {
180        let mut convergents = self.clone().convergents();
181        let mut last_conv = convergents.next().unwrap();
182        if last_conv.denom() > limit {
183            return Approximation::Approximated(last_conv);
184        }
185        loop {
186            last_conv = match convergents.next() {
187                Some(v) => {
188                    if v.denom() < limit {
189                        v
190                    } else {
191                        return Approximation::Approximated(last_conv);
192                    }
193                }
194                None => return Approximation::Exact(last_conv),
195            }
196        }
197    }
198}
199
200// TODO: implement continued fraction for various functions
201// REF: https://crypto.stanford.edu/pbc/notes/contfrac/cheat.html
202// List: coth(1/n), tan(1/n), arctan(n), tanh(n), tan(n), log(1+x), exp(2m/n), exp(1/n)
203
204/// Iterator of simple continued fraction coefficients generated by [exp()]
205#[derive(Debug, Clone, Copy)]
206pub struct ExpCoefficients<T> {
207    exponent: T,
208    i: T,
209}
210
211impl<T: Num + NumRef + Clone + Signed> Iterator for ExpCoefficients<T>
212where
213    for<'r> &'r T: RefNum<T>,
214{
215    type Item = (T, T);
216
217    fn next(&mut self) -> Option<Self::Item> {
218        let result = if self.i.is_zero() {
219            Some((T::one(), T::one()))
220        } else if self.i.is_one() {
221            Some((self.exponent.clone(), T::one()))
222        } else {
223            Some((
224                -((&self.i - T::one()) * &self.exponent),
225                &self.i + &self.exponent,
226            ))
227        };
228
229        self.i = &self.i + T::one();
230        result
231    }
232}
233
234pub fn exp<T: Num + NumRef + Clone + Signed>(target: T) -> GeneralContinuedFraction<ExpCoefficients<T>>
235where
236    for<'r> &'r T: RefNum<T>,
237{
238    ExpCoefficients {
239        exponent: target,
240        i: T::zero(),
241    }.into()
242}
243
244// TODO: implement operators to caculate HAKMEM Constant (and add as a example)
245// https://crypto.stanford.edu/pbc/notes/contfrac/hakmem.html
246
247#[cfg(test)]
248mod tests {
249    use super::*;
250    use crate::symbols::{Pi, E};
251
252    #[test]
253    fn general_conf_frac_test() {
254        let e = E {};
255        assert_eq!(e.cfrac::<i32>().take(10), exp(1i32).simplify().take(10));
256
257        let pi = Pi {};
258        assert_eq!(
259            pi.gfrac()
260                .convergents()
261                .skip(1)
262                .take(5)
263                .map(|c| c.into())
264                .collect::<Vec<(i32, i32)>>(),
265            vec![(4, 1), (3, 1), (19, 6), (160, 51), (1744, 555)]
266        );
267        assert_eq!(
268            pi.gfrac().simplify().into_iter().take(6).collect::<Vec<u64>>(),
269            vec![3, 7, 15, 1, 292, 1]
270        );
271        assert_eq!(
272            pi.gfrac::<u32>().decimals().take(20).collect::<Vec<_>>(),
273            b"31415926"
274        );
275    }
276
277    #[test]
278    fn functions_test() {
279        assert_eq!(
280            exp(1i32).into_iter().take(5).collect::<Vec<_>>(),
281            vec![(1i32, 1), (1, 1), (-1, 3), (-2, 4), (-3, 5)]
282        )
283    }
284}
num_irrational/cont_frac/general.rs

num_irrational/cont_frac/
general.rs
