zenith_float_num/ops/consts/
mod.rs

1mod e;
2mod ln10;
3mod ln2;
4mod pi;
5
6use crate::common::buf::WordBuf;
7use crate::common::util::round_p;
8use crate::mantissa::Mantissa;
9use crate::num::ExactNumNumber;
10use crate::ops::consts::e::ECache;
11use crate::ops::consts::ln10::Ln10Cache;
12use crate::ops::consts::ln2::Ln2Cache;
13use crate::ops::consts::pi::PiCache;
14use crate::Error;
15use crate::ExactNum;
16use crate::RoundingMode;
17use crate::WORD_BIT_SIZE;
18
19#[cfg(not(feature = "std"))]
20use alloc::vec::Vec;
21
22/// Alias for [`Consts`]: a progressive cache of π, e, ln 2, ln 10, √2, φ, and γ.
23pub type ConstCache = Consts;
24
25/// Snapshot of how many mantissa bits of each constant are currently cached.
26#[derive(Clone, Copy, Debug, PartialEq, Eq)]
27pub struct ConstCacheInfo {
28    /// Cached bits of π.
29    pub pi: usize,
30    /// Cached bits of e.
31    pub e: usize,
32    /// Cached bits of ln 2.
33    pub ln2: usize,
34    /// Cached bits of ln 10.
35    pub ln10: usize,
36    /// Cached bits of √2.
37    pub sqrt2: usize,
38    /// Cached bits of φ = (1+√5)/2.
39    pub phi: usize,
40    /// Cached bits of the Euler–Mascheroni constant γ.
41    pub euler: usize,
42}
43
44/// A float stored at extra working precision so later requests at lower (or equal)
45/// precision reuse the cache instead of recomputing.
46#[derive(Clone, Debug)]
47pub struct CachedFBig {
48    inner: ExactNum,
49}
50
51impl CachedFBig {
52    /// Wrap an already-computed value, retaining its current mantissa width.
53    pub fn new(inner: ExactNum) -> Self {
54        CachedFBig { inner }
55    }
56
57    /// Cached mantissa width in bits (`None` for Inf / NaN).
58    pub fn cached_bit_len(&self) -> Option<usize> {
59        self.inner.mantissa_max_bit_len()
60    }
61
62    /// The stored value.
63    pub fn inner(&self) -> &ExactNum {
64        &self.inner
65    }
66
67    /// Round the cached value to `p` bits.
68    pub fn round(&self, p: usize, rm: RoundingMode) -> ExactNum {
69        let mut v = self.inner.clone();
70        let _ = v.set_precision(p, rm);
71        v
72    }
73}
74
75#[derive(Debug)]
76struct ExtraCache {
77    bits: usize,
78    val: Option<ExactNumNumber>,
79}
80
81impl ExtraCache {
82    fn new() -> Self {
83        ExtraCache { bits: 0, val: None }
84    }
85
86    fn cached_bit_len(&self) -> usize {
87        self.bits
88    }
89
90    fn for_prec<F>(
91        &mut self,
92        k: usize,
93        rm: RoundingMode,
94        mut compute: F,
95    ) -> Result<ExactNumNumber, Error>
96    where
97        F: FnMut(usize) -> Result<ExactNumNumber, Error>,
98    {
99        let p = round_p(k);
100        let p_wrk = p.checked_add(WORD_BIT_SIZE).ok_or(Error::InvalidArgument)?;
101        if self.bits >= p {
102            if let Some(v) = &self.val {
103                let mut ret = v.clone()?;
104                ret.set_precision(p, rm)?;
105                return Ok(ret);
106            }
107        }
108        let computed = compute(p_wrk)?;
109        self.bits = computed.mantissa_max_bit_len();
110        self.val = Some(computed.clone()?);
111        let mut ret = computed;
112        ret.set_precision(p, rm)?;
113        Ok(ret)
114    }
115
116    fn install(&mut self, v: ExactNumNumber) {
117        self.bits = v.mantissa_max_bit_len();
118        self.val = Some(v);
119    }
120}
121
122/// Constants cache contains arbitrary-precision mathematical constants.
123#[derive(Debug)]
124pub struct Consts {
125    pi: PiCache,
126    e: ECache,
127    ln2: Ln2Cache,
128    ln10: Ln10Cache,
129    sqrt2: ExtraCache,
130    phi: ExtraCache,
131    euler: ExtraCache,
132    tenpowers: Vec<(WordBuf, WordBuf, usize)>,
133}
134
135/// In an ideal situation, the `Consts` structure is initialized with `Consts::new` only once,
136/// and then used where needed.
137impl Consts {
138    /// Initializes the constants cache.
139    ///
140    /// ## Errors
141    ///
142    ///  - MemoryAllocation: failed to allocate memory for mantissa.
143    pub fn new() -> Result<Self, Error> {
144        Ok(Consts {
145            pi: PiCache::new()?,
146            e: ECache::new()?,
147            ln2: Ln2Cache::new()?,
148            ln10: Ln10Cache::new()?,
149            sqrt2: ExtraCache::new(),
150            phi: ExtraCache::new(),
151            euler: ExtraCache::new(),
152            tenpowers: Vec::new(),
153        })
154    }
155
156    /// Returns the value of the pi number with precision `p` using rounding mode `rm`.
157    /// Precision is rounded upwards to the word size.
158    ///
159    /// ## Errors
160    ///
161    ///  - MemoryAllocation: failed to allocate memory for mantissa.
162    ///  - InvalidArgument: the precision is incorrect.
163    pub(crate) fn pi_num(&mut self, p: usize, rm: RoundingMode) -> Result<ExactNumNumber, Error> {
164        let p = round_p(p);
165        self.pi.for_prec(p, rm)
166    }
167
168    /// Returns the value of the Euler number with precision `p` using rounding mode `rm`.
169    /// Precision is rounded upwards to the word size.
170    ///
171    /// ## Errors
172    ///
173    ///  - MemoryAllocation: failed to allocate memory for mantissa.
174    ///  - InvalidArgument: the precision is incorrect.
175    pub(crate) fn e_num(&mut self, p: usize, rm: RoundingMode) -> Result<ExactNumNumber, Error> {
176        let p = round_p(p);
177        self.e.for_prec(p, rm)
178    }
179
180    /// Returns the value of the natural logarithm of 2 with precision `p` using rounding mode `rm`.
181    /// Precision is rounded upwards to the word size.
182    ///
183    /// ## Errors
184    ///
185    ///  - MemoryAllocation: failed to allocate memory for mantissa.
186    ///  - InvalidArgument: the precision is incorrect.
187    pub(crate) fn ln_2_num(&mut self, p: usize, rm: RoundingMode) -> Result<ExactNumNumber, Error> {
188        let p = round_p(p);
189        self.ln2.for_prec(p, rm)
190    }
191
192    /// Returns the value of the natural logarithm of 10 with precision `p` using rounding mode `rm`.
193    /// Precision is rounded upwards to the word size.
194    ///
195    /// ## Errors
196    ///
197    ///  - MemoryAllocation: failed to allocate memory for mantissa.
198    ///  - InvalidArgument: the precision is incorrect.
199    pub(crate) fn ln_10_num(
200        &mut self,
201        p: usize,
202        rm: RoundingMode,
203    ) -> Result<ExactNumNumber, Error> {
204        let p = round_p(p);
205        self.ln10.for_prec(p, rm)
206    }
207
208    /// Returns the value of the pi number with precision `p` using rounding mode `rm`.
209    /// Precision is rounded upwards to the word size.
210    pub fn pi(&mut self, p: usize, rm: RoundingMode) -> ExactNum {
211        match self.pi_num(p, rm) {
212            Ok(v) => v.into(),
213            Err(e) => ExactNum::nan(Some(e)),
214        }
215    }
216
217    /// Returns the value of the Euler number with precision `p` using rounding mode `rm`.
218    /// Precision is rounded upwards to the word size.
219    pub fn e(&mut self, p: usize, rm: RoundingMode) -> ExactNum {
220        match self.e_num(p, rm) {
221            Ok(v) => v.into(),
222            Err(e) => ExactNum::nan(Some(e)),
223        }
224    }
225
226    /// Returns the value of the natural logarithm of 2 with precision `p` using rounding mode `rm`.
227    /// Precision is rounded upwards to the word size.
228    pub fn ln_2(&mut self, p: usize, rm: RoundingMode) -> ExactNum {
229        match self.ln_2_num(p, rm) {
230            Ok(v) => v.into(),
231            Err(e) => ExactNum::nan(Some(e)),
232        }
233    }
234
235    /// Returns the value of the natural logarithm of 10 with precision `p` using rounding mode `rm`.
236    /// Precision is rounded upwards to the word size.
237    pub fn ln_10(&mut self, p: usize, rm: RoundingMode) -> ExactNum {
238        match self.ln_10_num(p, rm) {
239            Ok(v) => v.into(),
240            Err(e) => ExactNum::nan(Some(e)),
241        }
242    }
243
244    /// Return powers of 10: 100, 10000, 100000000, ...
245    pub(crate) fn tenpowers(&mut self, p: usize) -> Result<&[(WordBuf, WordBuf, usize)], Error> {
246        if p >= self.tenpowers.len() {
247            Mantissa::compute_tenpowers(&mut self.tenpowers, p)?;
248        }
249
250        Ok(&self.tenpowers)
251    }
252
253    /// How many bits of each series / extra constant are currently retained.
254    pub fn cache_info(&self) -> ConstCacheInfo {
255        ConstCacheInfo {
256            pi: self.pi.cached_bit_len(),
257            e: self.e.cached_bit_len(),
258            ln2: self.ln2.cached_bit_len(),
259            ln10: self.ln10.cached_bit_len(),
260            sqrt2: self.sqrt2.cached_bit_len(),
261            phi: self.phi.cached_bit_len(),
262            euler: self.euler.cached_bit_len(),
263        }
264    }
265
266    fn sqrt2_num(&mut self, p: usize, rm: RoundingMode) -> Result<ExactNumNumber, Error> {
267        self.sqrt2.for_prec(p, rm, |p_wrk| {
268            let two = ExactNumNumber::from_word(2, p_wrk)?;
269            two.sqrt(p_wrk, RoundingMode::None)
270        })
271    }
272
273    fn phi_num(&mut self, p: usize, rm: RoundingMode) -> Result<ExactNumNumber, Error> {
274        self.phi.for_prec(p, rm, |p_wrk| {
275            let five = ExactNumNumber::from_word(5, p_wrk)?;
276            let one = ExactNumNumber::from_word(1, p_wrk)?;
277            let two = ExactNumNumber::from_word(2, p_wrk)?;
278            let s = five.sqrt(p_wrk, RoundingMode::None)?;
279            let n = one.add(&s, p_wrk, RoundingMode::None)?;
280            n.div(&two, p_wrk, RoundingMode::None)
281        })
282    }
283
284    /// √2 with precision `p` using rounding mode `rm`.
285    /// Higher requests extend the cache; lower requests reuse it.
286    pub fn sqrt2(&mut self, p: usize, rm: RoundingMode) -> ExactNum {
287        match self.sqrt2_num(p, rm) {
288            Ok(v) => v.into(),
289            Err(e) => ExactNum::nan(Some(e)),
290        }
291    }
292
293    /// Golden ratio φ = (1+√5)/2 with precision `p` using rounding mode `rm`.
294    /// Higher requests extend the cache; lower requests reuse it.
295    pub fn phi(&mut self, p: usize, rm: RoundingMode) -> ExactNum {
296        match self.phi_num(p, rm) {
297            Ok(v) => v.into(),
298            Err(e) => ExactNum::nan(Some(e)),
299        }
300    }
301
302    pub(crate) fn euler_gamma_num(
303        &mut self,
304        p: usize,
305        rm: RoundingMode,
306    ) -> Result<ExactNumNumber, Error> {
307        let p_round = round_p(p);
308        let p_wrk = p_round
309            .checked_add(8 * WORD_BIT_SIZE)
310            .ok_or(Error::InvalidArgument)?;
311        if self.euler.cached_bit_len() < p_round {
312            let ln2 = self.ln_2_num(p_wrk, RoundingMode::None)?;
313            let v = crate::ops::special::euler_mascheroni(p_wrk, &ln2)?;
314            self.euler.install(v);
315        }
316        self.euler.for_prec(p, rm, |_| Err(Error::InvalidArgument))
317    }
318
319    /// Euler–Mascheroni constant γ with precision `p` using rounding mode `rm`.
320    pub fn euler_gamma(&mut self, p: usize, rm: RoundingMode) -> ExactNum {
321        match self.euler_gamma_num(p, rm) {
322            Ok(v) => v.into(),
323            Err(e) => ExactNum::nan(Some(e)),
324        }
325    }
326}
327
328/// Thread-safe [`Consts`] for batch evaluation across threads (`std` only).
329/// Callers still take `&mut Consts` inside [`SharedConsts::with`]; the mutex serializes cache fills.
330#[cfg(feature = "std")]
331#[derive(Debug)]
332pub struct SharedConsts {
333    inner: std::sync::Mutex<Consts>,
334}
335
336#[cfg(feature = "std")]
337impl SharedConsts {
338    /// Allocate an empty progressive constant cache protected by a mutex.
339    pub fn new() -> Result<Self, Error> {
340        Ok(SharedConsts {
341            inner: std::sync::Mutex::new(Consts::new()?),
342        })
343    }
344
345    /// Run `f` with exclusive access to the cache. Recovers from a poisoned mutex by taking the inner value.
346    pub fn with<F, R>(&self, f: F) -> R
347    where
348        F: FnOnce(&mut Consts) -> R,
349    {
350        let mut g = match self.inner.lock() {
351            Ok(g) => g,
352            Err(p) => p.into_inner(),
353        };
354        f(&mut g)
355    }
356}
357
358#[cfg(test)]
359mod tests {
360    use super::*;
361
362    #[test]
363    fn progressive_constant_cache_extends() {
364        let mut cc = Consts::new().expect("constants");
365        let rm = RoundingMode::ToEven;
366        let _ = cc.pi(128, rm);
367        let _ = cc.e(128, rm);
368        let _ = cc.ln_2(128, rm);
369        let _ = cc.ln_10(128, rm);
370        let _ = cc.sqrt2(128, rm);
371        let _ = cc.phi(128, rm);
372        let _ = cc.euler_gamma(128, rm);
373        let before = cc.cache_info();
374        assert!(before.sqrt2 >= 128);
375        assert!(before.phi >= 128);
376        assert!(before.euler >= 128);
377
378        let _ = cc.pi(256, rm);
379        let _ = cc.e(256, rm);
380        let _ = cc.ln_2(256, rm);
381        let _ = cc.ln_10(256, rm);
382        let _ = cc.sqrt2(256, rm);
383        let _ = cc.phi(256, rm);
384        let _ = cc.euler_gamma(256, rm);
385        let after = cc.cache_info();
386        assert!(after.pi >= before.pi);
387        assert!(after.e >= before.e);
388        assert!(after.ln2 >= before.ln2);
389        assert!(after.ln10 >= before.ln10);
390        assert!(after.sqrt2 >= 256);
391        assert!(after.phi >= 256);
392        assert!(after.euler >= 256);
393
394        let a = cc.sqrt2(128, rm);
395        let b = cc.sqrt2(128, rm);
396        assert_eq!(a.cmp(&b), Some(0));
397        let cached = CachedFBig::new(a);
398        assert!(cached.cached_bit_len().unwrap() >= 128);
399        let r = cached.round(64, rm);
400        assert!(!r.is_nan());
401    }
402
403    #[test]
404    fn euler_gamma_matches_known_digits() {
405        let mut cc = Consts::new().expect("constants");
406        let rm = RoundingMode::ToEven;
407        let p = 128;
408        let g = cc.euler_gamma(p, rm);
409        let known = ExactNum::parse(
410            "0.57721566490153286060651209008240243",
411            crate::Radix::Dec,
412            p,
413            rm,
414            &mut cc,
415        );
416        let d = g.sub(&known, p, RoundingMode::None);
417        assert!(d.is_zero() || d.exponent().unwrap() < -((p as i32) / 4));
418    }
419
420    #[cfg(feature = "std")]
421    #[test]
422    fn shared_consts_parallel_pi() {
423        use std::sync::Arc;
424        use std::thread;
425
426        let cc = Arc::new(SharedConsts::new().expect("constants"));
427        let mut hs = Vec::new();
428        for _ in 0..4 {
429            let cc = Arc::clone(&cc);
430            hs.push(thread::spawn(move || {
431                cc.with(|c| c.pi(128, RoundingMode::ToEven))
432            }));
433        }
434        let vals: Vec<_> = hs.into_iter().map(|h| h.join().unwrap()).collect();
435        for v in &vals[1..] {
436            assert_eq!(vals[0].cmp(v), Some(0));
437        }
438    }
439}
zenith_float_num/ops/consts/mod.rs

zenith_float_num/ops/consts/
mod.rs
