1use num_bigint::BigInt;
10
11use crate::algebraic::AlgebraicNumber;
12use crate::basis::Basis;
13use crate::computable::ComputableReal;
14use crate::primes::factorize;
15use crate::rational::RnsRational;
16use crate::rns::RnsInt;
17use crate::symbolic::{IdentityGraph, SymbolicExpr};
18
19#[derive(Debug, Clone, Copy, PartialEq, Eq, PartialOrd, Ord)]
21pub enum TowerLevel {
22 Integer = 0,
23 Rational = 1,
24 Algebraic = 2,
25 Computable = 3,
26 Symbolic = 4,
27}
28
29#[derive(Clone)]
31pub enum TowerValue {
32 Integer(RnsInt),
33 Rational(RnsRational),
34 Algebraic(AlgebraicNumber),
35 Computable(ComputableReal),
36 Symbolic(SymbolicExpr),
37}
38
39impl TowerValue {
40 pub fn level(&self) -> TowerLevel {
42 match self {
43 TowerValue::Integer(_) => TowerLevel::Integer,
44 TowerValue::Rational(_) => TowerLevel::Rational,
45 TowerValue::Algebraic(_) => TowerLevel::Algebraic,
46 TowerValue::Computable(_) => TowerLevel::Computable,
47 TowerValue::Symbolic(_) => TowerLevel::Symbolic,
48 }
49 }
50
51 pub fn channels(&self) -> Basis {
53 match self {
54 TowerValue::Integer(i) => i.basis.clone(),
55 TowerValue::Rational(r) => r.channels.clone(),
56 TowerValue::Algebraic(a) => a.channels.clone(),
57 TowerValue::Computable(c) => c.channels(),
58 TowerValue::Symbolic(_) => Basis::standard(),
59 }
60 }
61
62 pub fn reduce(&self) -> TowerValue {
64 let mut current = self.clone();
65 loop {
66 let next = current.reduce_once();
67 if next.level() == current.level() {
68 return next;
69 }
70 current = next;
71 }
72 }
73
74 fn reduce_once(&self) -> TowerValue {
75 match self {
76 TowerValue::Algebraic(a) if a.degree() == 1 => {
77 TowerValue::Rational(a.to_rational().unwrap())
78 }
79 TowerValue::Rational(r) if r.is_integer() => {
80 let p = r.to_pair().0;
81 TowerValue::Integer(RnsInt::from_bigint(&p, r.channels.clone()))
82 }
83 TowerValue::Symbolic(e) => {
84 let simplified = IdentityGraph::standard().simplify(e.clone());
85 match symbolic_to_value(&simplified, &self.channels()) {
86 Some(v) => v,
87 None => TowerValue::Symbolic(simplified),
88 }
89 }
90 other => other.clone(),
91 }
92 }
93
94 pub fn elevate_to(&self, target: TowerLevel) -> TowerValue {
97 let mut current = self.clone();
98 while current.level() < target {
99 current = current.elevate_once();
100 }
101 current
102 }
103
104 fn elevate_once(&self) -> TowerValue {
105 match self {
106 TowerValue::Integer(i) => {
107 TowerValue::Rational(RnsRational::new(i.to_bigint(), BigInt::from(1), i.basis.clone()))
108 }
109 TowerValue::Rational(r) => TowerValue::Algebraic(AlgebraicNumber::from_rational(r.clone())),
110 TowerValue::Algebraic(a) => TowerValue::Computable(a.to_computable()),
111 TowerValue::Computable(_) => self.clone(), TowerValue::Symbolic(_) => self.clone(),
113 }
114 }
115
116 pub fn to_f64(&self) -> Option<f64> {
118 match self {
119 TowerValue::Integer(i) => {
120 Some(RnsRational::new(i.to_bigint(), BigInt::from(1), i.basis.clone()).to_f64())
121 }
122 TowerValue::Rational(r) => Some(r.to_f64()),
123 TowerValue::Algebraic(a) => Some(a.to_f64()),
124 TowerValue::Computable(c) => Some(c.evaluate_f64()),
125 TowerValue::Symbolic(_) => self.reduce().non_symbolic_to_f64(),
126 }
127 }
128
129 fn non_symbolic_to_f64(&self) -> Option<f64> {
130 match self {
131 TowerValue::Symbolic(e) => {
132 let ch = self.channels();
134 symbolic_to_computable(e, &ch).map(|c| c.evaluate_f64())
135 }
136 other => other.to_f64(),
137 }
138 }
139
140 pub fn digits(&self, precision: u64) -> RnsRational {
142 match self {
143 TowerValue::Integer(i) => {
144 RnsRational::new(i.to_bigint(), BigInt::from(1), i.basis.clone())
145 }
146 TowerValue::Rational(r) => r.clone(),
147 TowerValue::Algebraic(a) => {
148 let mut clone = a.clone();
149 let target = RnsRational::new(
150 BigInt::from(1),
151 BigInt::from(10u8).pow((precision + 1) as u32),
152 a.channels.clone(),
153 );
154 clone.refine_interval(&target);
155 clone.interval.0.midpoint(&clone.interval.1)
156 }
157 TowerValue::Computable(c) => c.evaluate(precision),
158 TowerValue::Symbolic(e) => {
159 let ch = self.channels();
160 let simplified = IdentityGraph::standard().simplify(e.clone());
161 if let Some(v) = symbolic_to_value(&simplified, &ch) {
162 return v.digits(precision);
163 }
164 match symbolic_to_computable(&simplified, &ch) {
165 Some(c) => c.evaluate(precision),
166 None => panic!("cannot produce digits for irreducible symbolic expression"),
167 }
168 }
169 }
170 }
171
172 pub fn add(&self, other: &TowerValue) -> TowerValue {
176 self.binop(other, Op::Add).reduce()
177 }
178
179 pub fn mul(&self, other: &TowerValue) -> TowerValue {
181 self.binop(other, Op::Mul).reduce()
182 }
183
184 fn binop(&self, other: &TowerValue, op: Op) -> TowerValue {
185 let lvl = self.level().max(other.level());
186 if lvl == TowerLevel::Symbolic {
187 let a = self.to_symbolic();
188 let b = other.to_symbolic();
189 let expr = match op {
190 Op::Add => SymbolicExpr::Add(vec![a, b]),
191 Op::Mul => SymbolicExpr::Mul(vec![a, b]),
192 };
193 return TowerValue::Symbolic(IdentityGraph::standard().simplify(expr));
194 }
195 let a = self.elevate_to(lvl);
196 let b = other.elevate_to(lvl);
197 match (a, b) {
198 (TowerValue::Integer(x), TowerValue::Integer(y)) => TowerValue::Integer(match op {
199 Op::Add => x.add(&y),
200 Op::Mul => x.mul(&y),
201 }),
202 (TowerValue::Rational(x), TowerValue::Rational(y)) => TowerValue::Rational(match op {
203 Op::Add => x.add(&y),
204 Op::Mul => x.mul(&y),
205 }),
206 (TowerValue::Algebraic(x), TowerValue::Algebraic(y)) => TowerValue::Algebraic(match op {
207 Op::Add => x.add(&y),
208 Op::Mul => x.mul(&y),
209 }),
210 (TowerValue::Computable(x), TowerValue::Computable(y)) => TowerValue::Computable(match op {
211 Op::Add => x.add(&y),
212 Op::Mul => x.mul(&y),
213 }),
214 _ => unreachable!("levels were equalized before the operation"),
215 }
216 }
217
218 pub fn sqrt(&self) -> TowerValue {
221 let ch = self.channels();
222 if let TowerValue::Integer(i) = self {
223 let n = i.to_bigint();
224 if let Some(nu) = bigint_to_u64(&n) {
225 let root = (nu as f64).sqrt().round() as u64;
226 if root * root == nu {
227 return TowerValue::Integer(RnsInt::from_bigint(&BigInt::from(root), ch));
228 }
229 return TowerValue::Algebraic(AlgebraicNumber::sqrt(nu, ch)).reduce();
230 }
231 }
232 let v = self.to_f64().unwrap_or(f64::NAN);
234 if v >= 0.0 && v.fract() == 0.0 {
235 return TowerValue::Algebraic(AlgebraicNumber::sqrt(v as u64, ch)).reduce();
236 }
237 panic!("sqrt of non-integer values is not supported at the tower level yet");
238 }
239
240 pub fn sin(&self) -> TowerValue {
242 let expr = SymbolicExpr::Sin(Box::new(self.to_symbolic()));
243 TowerValue::Symbolic(IdentityGraph::standard().simplify(expr)).reduce()
244 }
245
246 fn to_symbolic(&self) -> SymbolicExpr {
247 match self {
248 TowerValue::Integer(i) => SymbolicExpr::Integer(i.to_bigint()),
249 TowerValue::Rational(r) => {
250 let (p, q) = r.to_pair();
251 SymbolicExpr::Rational(p, q)
252 }
253 TowerValue::Symbolic(e) => e.clone(),
254 TowerValue::Algebraic(a) => {
256 if let Some(r) = a.to_rational() {
257 let (p, q) = r.to_pair();
258 SymbolicExpr::Rational(p, q)
259 } else {
260 panic!("cannot lift this algebraic number to symbolic form")
261 }
262 }
263 TowerValue::Computable(_) => panic!("cannot lift a computable real to symbolic form"),
264 }
265 }
266}
267
268#[derive(Clone, Copy)]
269enum Op {
270 Add,
271 Mul,
272}
273
274fn bigint_to_u64(n: &BigInt) -> Option<u64> {
275 use num_traits::ToPrimitive;
276 n.to_u64()
277}
278
279fn symbolic_to_value(e: &SymbolicExpr, ch: &Basis) -> Option<TowerValue> {
281 match e {
282 SymbolicExpr::Integer(n) => Some(TowerValue::Integer(RnsInt::from_bigint(n, ch.clone()))),
283 SymbolicExpr::Rational(p, q) => {
284 Some(TowerValue::Rational(RnsRational::new(p.clone(), q.clone(), ch.clone())))
285 }
286 SymbolicExpr::Sqrt { radicand } => {
287 let nu = bigint_to_u64(radicand)?;
288 Some(TowerValue::Algebraic(AlgebraicNumber::sqrt(nu, ch.clone())))
289 }
290 SymbolicExpr::ScaledSqrt { coeff: (a, b), rad } => {
291 let nu = bigint_to_u64(rad)?;
292 let s = AlgebraicNumber::sqrt(nu, ch.clone());
293 let coeff = RnsRational::new(a.clone(), b.clone(), ch.clone());
294 Some(TowerValue::Algebraic(s).mul(&TowerValue::Rational(coeff)))
295 }
296 _ => None,
297 }
298}
299
300fn symbolic_to_computable(e: &SymbolicExpr, ch: &Basis) -> Option<ComputableReal> {
303 match e {
304 SymbolicExpr::Integer(n) => Some(ComputableReal::from_rational(RnsRational::new(
305 n.clone(),
306 BigInt::from(1),
307 ch.clone(),
308 ))),
309 SymbolicExpr::Rational(p, q) => Some(ComputableReal::from_rational(RnsRational::new(
310 p.clone(),
311 q.clone(),
312 ch.clone(),
313 ))),
314 SymbolicExpr::Pi => Some(ComputableReal::pi(ch.clone())),
315 SymbolicExpr::E => Some(ComputableReal::e(ch.clone())),
316 SymbolicExpr::Sqrt { radicand } => {
317 let r = RnsRational::new(radicand.clone(), BigInt::from(1), ch.clone());
318 Some(ComputableReal::sqrt(r))
319 }
320 SymbolicExpr::Add(terms) => {
321 let mut acc: Option<ComputableReal> = None;
322 for t in terms {
323 let c = symbolic_to_computable(t, ch)?;
324 acc = Some(match acc {
325 Some(a) => a.add(&c),
326 None => c,
327 });
328 }
329 acc
330 }
331 SymbolicExpr::Mul(factors) => {
332 let mut acc: Option<ComputableReal> = None;
333 for f in factors {
334 let c = symbolic_to_computable(f, ch)?;
335 acc = Some(match acc {
336 Some(a) => a.mul(&c),
337 None => c,
338 });
339 }
340 acc
341 }
342 _ => None,
343 }
344}
345
346pub fn is_perfect_square(n: u64) -> bool {
348 let r = (n as f64).sqrt().round() as u64;
349 r * r == n
350}
351
352#[allow(dead_code)]
353fn _uses_factorize() {
354 let _ = factorize(12);
355}
356
357#[cfg(test)]
358mod tests {
359 use super::*;
360
361 fn ch() -> Basis {
362 Basis::standard()
363 }
364
365 #[test]
366 fn rational_level() {
367 let v = TowerValue::Rational(RnsRational::from_fraction(2, 3, ch()));
368 assert_eq!(v.level(), TowerLevel::Rational);
369 }
370
371 #[test]
372 fn rational_with_denom_one_reduces_to_integer() {
373 let v = TowerValue::Rational(RnsRational::from_fraction(6, 2, ch()));
374 assert_eq!(v.reduce().level(), TowerLevel::Integer);
375 }
376
377 #[test]
378 fn sqrt2_times_sqrt2_drops_to_integer() {
379 let s = TowerValue::Algebraic(AlgebraicNumber::sqrt(2, ch()));
380 let prod = s.mul(&s);
381 assert_eq!(prod.level(), TowerLevel::Integer);
382 assert_eq!(prod.to_f64().unwrap().round(), 2.0);
383 }
384
385 #[test]
386 fn mixed_level_add() {
387 let a = TowerValue::Integer(RnsInt::from_i64(1, ch()));
389 let b = TowerValue::Rational(RnsRational::from_fraction(1, 2, ch()));
390 let sum = a.add(&b);
391 assert_eq!(sum.level(), TowerLevel::Rational);
392 assert!((sum.to_f64().unwrap() - 1.5).abs() < 1e-12);
393 }
394
395 #[test]
396 fn pi_plus_one_digits() {
397 let pi = TowerValue::Symbolic(SymbolicExpr::Pi);
399 let one = TowerValue::Integer(RnsInt::from_i64(1, ch()));
400 let sum = pi.add(&one);
401 let d = sum.digits(20);
402 assert!((d.to_f64() - (std::f64::consts::PI + 1.0)).abs() < 1e-12);
403 }
404
405 #[test]
406 fn sin_pi_is_zero() {
407 let pi = TowerValue::Symbolic(SymbolicExpr::Pi);
408 let s = pi.sin();
409 assert_eq!(s.to_f64().unwrap(), 0.0);
410 }
411}