1use crate::{IBig, RBig, UBig};
6use oxinum_core::{OxiNumError, OxiNumResult};
7
8pub fn continued_fraction(x: &RBig) -> Vec<IBig> {
34 let mut coeffs = Vec::new();
35 let mut num = x.numerator().clone();
36 let mut den: IBig = x.denominator().clone().into();
37
38 if den == IBig::ZERO {
39 return coeffs;
40 }
41
42 loop {
43 let (q, r) = div_floor(&num, &den);
45 coeffs.push(q);
46 if r == IBig::ZERO {
47 break;
48 }
49 num = den;
50 den = r;
51 }
52 coeffs
53}
54
55pub fn from_continued_fraction(coeffs: &[IBig]) -> OxiNumResult<RBig> {
71 if coeffs.is_empty() {
72 return Err(OxiNumError::Parse("empty continued fraction".into()));
73 }
74 let mut result = RBig::from(coeffs[coeffs.len() - 1].clone());
76 for coeff in coeffs[..coeffs.len() - 1].iter().rev() {
77 if result == RBig::ZERO {
79 return Err(OxiNumError::DivByZero);
80 }
81 let reciprocal = rational_reciprocal(&result)?;
82 result = RBig::from(coeff.clone()) + reciprocal;
83 }
84 Ok(result)
85}
86
87pub fn best_rational_approximation(x: &RBig, max_denom: &UBig) -> RBig {
104 if *max_denom == UBig::ZERO {
105 return RBig::from(x.floor());
106 }
107 let coeffs = continued_fraction(x);
108 if coeffs.is_empty() {
109 return RBig::ZERO;
110 }
111
112 let mut h_prev2 = IBig::ZERO; let mut h_prev1 = IBig::ONE; let mut k_prev2 = IBig::ONE; let mut k_prev1 = IBig::ZERO; let max_denom_i: IBig = max_denom.clone().into();
124 let mut best = RBig::from(coeffs[0].clone());
125
126 for coeff in &coeffs {
127 let h_n = coeff * &h_prev1 + &h_prev2;
128 let k_n = coeff * &k_prev1 + &k_prev2;
129
130 if k_n > max_denom_i {
131 break;
134 }
135
136 best = rbig_from_signed(&h_n, &k_n);
138
139 h_prev2 = h_prev1;
140 h_prev1 = h_n;
141 k_prev2 = k_prev1;
142 k_prev1 = k_n;
143 }
144
145 best
146}
147
148pub fn to_decimal_string(x: &RBig, decimal_places: usize) -> String {
169 let num = x.numerator().clone();
170 let den: IBig = x.denominator().clone().into();
171
172 let is_negative = num < IBig::ZERO;
173 let abs_num = if is_negative { -num } else { num };
174
175 let (int_part, mut remainder) = div_floor(&abs_num, &den);
177
178 let mut result = String::new();
179 if is_negative && (int_part != IBig::ZERO || decimal_places > 0) {
180 result.push('-');
181 }
182 result.push_str(&int_part.to_string());
183
184 if decimal_places > 0 {
185 result.push('.');
186 let ten = IBig::from(10);
187 for _ in 0..decimal_places {
188 remainder *= &ten;
189 let (digit, new_rem) = div_floor(&remainder, &den);
190 result.push_str(&digit.to_string());
191 remainder = new_rem;
192 }
193 }
194 result
195}
196
197pub fn mediant(a: &RBig, b: &RBig) -> RBig {
218 let a_num = a.numerator().clone();
219 let a_den: IBig = a.denominator().clone().into();
220 let b_num = b.numerator().clone();
221 let b_den: IBig = b.denominator().clone().into();
222
223 let num = a_num + b_num;
224 let den = a_den + b_den;
225 rbig_from_signed(&num, &den)
226}
227
228pub fn mixed_number(x: &RBig) -> (IBig, RBig) {
249 let whole = x.trunc();
250 let frac = x.fract();
251 (whole, frac)
252}
253
254pub fn rational_floor(x: &RBig) -> IBig {
260 x.floor()
261}
262
263pub fn rational_ceil(x: &RBig) -> IBig {
265 x.ceil()
266}
267
268pub fn rational_round(x: &RBig) -> IBig {
270 x.round()
271}
272
273pub fn rational_truncate(x: &RBig) -> IBig {
275 x.trunc()
276}
277
278pub fn rational_from_integer(n: &IBig) -> RBig {
296 RBig::from_parts(n.clone(), UBig::ONE)
297}
298
299pub fn rational_is_integer(x: &RBig) -> bool {
316 *x.denominator() == UBig::ONE
317}
318
319pub fn rational_to_integer(x: &RBig) -> Option<IBig> {
334 if rational_is_integer(x) {
335 Some(x.numerator().clone())
336 } else {
337 None
338 }
339}
340
341pub fn rational_abs(x: &RBig) -> RBig {
357 if x.numerator() < &IBig::ZERO {
358 -x.clone()
359 } else {
360 x.clone()
361 }
362}
363
364pub fn rational_signum(x: &RBig) -> IBig {
376 if x.is_zero() {
377 IBig::ZERO
378 } else if x.numerator() < &IBig::ZERO {
379 IBig::from(-1)
380 } else {
381 IBig::ONE
382 }
383}
384
385pub fn rational_reciprocal(x: &RBig) -> OxiNumResult<RBig> {
402 if x.is_zero() {
403 return Err(OxiNumError::DivByZero);
404 }
405 let num = x.numerator().clone();
408 let den: IBig = x.denominator().clone().into();
409 let sign_negative = num < IBig::ZERO;
410 let abs_num = if sign_negative { -&num } else { num };
412 let new_num = if sign_negative { -den } else { den };
413 Ok(rbig_from_signed(&new_num, &abs_num))
414}
415
416pub fn rational_pow(x: &RBig, n: i32) -> OxiNumResult<RBig> {
439 if n == 0 {
440 return Ok(RBig::ONE);
441 }
442 if n > 0 {
443 Ok(x.pow(n as usize))
444 } else {
445 let reciprocal = rational_reciprocal(x)?;
446 let abs_n = n.unsigned_abs() as usize;
447 Ok(reciprocal.pow(abs_n))
448 }
449}
450
451pub(crate) fn div_floor(num: &IBig, den: &IBig) -> (IBig, IBig) {
458 use dashu_base::DivRem;
459 let (mut q, mut r) = num.clone().div_rem(den.clone());
460 if (r < IBig::ZERO && *den > IBig::ZERO) || (r > IBig::ZERO && *den < IBig::ZERO) {
462 q -= IBig::ONE;
463 r += den;
464 }
465 (q, r)
466}
467
468pub(crate) fn rbig_from_signed(num: &IBig, den: &IBig) -> RBig {
470 RBig::from_parts_signed(num.clone(), den.clone())
471}
472
473#[cfg(test)]
478mod tests {
479 use super::*;
480
481 #[test]
482 fn cf_355_over_113() {
483 let r = RBig::from_parts(IBig::from(355), UBig::from(113u32));
485 let cf = continued_fraction(&r);
486 assert_eq!(cf, vec![IBig::from(3), IBig::from(7), IBig::from(16)]);
487 }
488
489 #[test]
490 fn cf_integer() {
491 let r = RBig::from(5u32);
493 let cf = continued_fraction(&r);
494 assert_eq!(cf, vec![IBig::from(5)]);
495 }
496
497 #[test]
498 fn cf_half() {
499 let r = RBig::from_parts(IBig::from(1), UBig::from(2u32));
501 let cf = continued_fraction(&r);
502 assert_eq!(cf, vec![IBig::ZERO, IBig::from(2)]);
503 }
504
505 #[test]
506 fn cf_roundtrip() {
507 let r = RBig::from_parts(IBig::from(355), UBig::from(113u32));
508 let cf = continued_fraction(&r);
509 let reconstructed = from_continued_fraction(&cf).expect("ok");
510 assert_eq!(reconstructed, r);
511 }
512
513 #[test]
514 fn cf_roundtrip_various() {
515 let cases = [(22, 7u32), (1, 7u32), (100, 3u32), (17, 12u32), (1, 1u32)];
516 for (n, d) in cases {
517 let r = RBig::from_parts(IBig::from(n), UBig::from(d));
518 let cf = continued_fraction(&r);
519 let back = from_continued_fraction(&cf).expect("ok");
520 assert_eq!(back, r, "roundtrip failed for {n}/{d}");
521 }
522 }
523
524 #[test]
525 fn from_cf_pi_approx() {
526 let cf = vec![IBig::from(3), IBig::from(7), IBig::from(16)];
528 let r = from_continued_fraction(&cf).expect("ok");
529 assert_eq!(r, RBig::from_parts(IBig::from(355), UBig::from(113u32)));
530 }
531
532 #[test]
533 fn from_cf_empty_errors() {
534 assert!(from_continued_fraction(&[]).is_err());
535 }
536
537 #[test]
538 fn best_approx_pi() {
539 let pi = RBig::from_parts(IBig::from(355), UBig::from(113u32));
541 let approx = best_rational_approximation(&pi, &UBig::from(100u32));
542 assert_eq!(approx, RBig::from_parts(IBig::from(22), UBig::from(7u32)));
543 }
544
545 #[test]
546 fn best_approx_small_denom() {
547 let pi = RBig::from_parts(IBig::from(355), UBig::from(113u32));
549 let approx = best_rational_approximation(&pi, &UBig::from(10u32));
550 assert_eq!(approx, RBig::from_parts(IBig::from(22), UBig::from(7u32)));
551 }
552
553 #[test]
554 fn best_approx_denom_one() {
555 let r = RBig::from_parts(IBig::from(7), UBig::from(2u32));
557 let approx = best_rational_approximation(&r, &UBig::ONE);
558 assert_eq!(approx, RBig::from(3u32));
559 }
560
561 #[test]
562 fn decimal_third() {
563 let third = RBig::from_parts(IBig::from(1), UBig::from(3u32));
564 assert_eq!(to_decimal_string(&third, 6), "0.333333");
565 }
566
567 #[test]
568 fn decimal_quarter() {
569 let quarter = RBig::from_parts(IBig::from(1), UBig::from(4u32));
570 assert_eq!(to_decimal_string(&quarter, 4), "0.2500");
571 }
572
573 #[test]
574 fn decimal_negative() {
575 let neg = RBig::from_parts(IBig::from(-1), UBig::from(2u32));
576 assert_eq!(to_decimal_string(&neg, 2), "-0.50");
577 }
578
579 #[test]
580 fn decimal_improper() {
581 let r = RBig::from_parts(IBig::from(7), UBig::from(4u32));
583 assert_eq!(to_decimal_string(&r, 2), "1.75");
584 }
585
586 #[test]
587 fn decimal_zero_places() {
588 let r = RBig::from_parts(IBig::from(7), UBig::from(4u32));
589 assert_eq!(to_decimal_string(&r, 0), "1");
590 }
591
592 #[test]
593 fn mediant_basic() {
594 let a = RBig::from_parts(IBig::from(1), UBig::from(2u32));
596 let b = RBig::from_parts(IBig::from(1), UBig::from(3u32));
597 let m = mediant(&a, &b);
598 assert_eq!(m, RBig::from_parts(IBig::from(2), UBig::from(5u32)));
599 }
600
601 #[test]
602 fn mediant_between() {
603 let a = RBig::ZERO;
605 let b = RBig::ONE;
606 let m = mediant(&a, &b);
607 assert_eq!(m, RBig::from_parts(IBig::from(1), UBig::from(2u32)));
608 }
609
610 #[test]
611 fn mixed_number_basic() {
612 let r = RBig::from_parts(IBig::from(7), UBig::from(3u32));
614 let (whole, frac) = mixed_number(&r);
615 assert_eq!(whole, IBig::from(2));
616 assert_eq!(frac, RBig::from_parts(IBig::from(1), UBig::from(3u32)));
617 }
618
619 #[test]
620 fn mixed_number_negative() {
621 let r = RBig::from_parts(IBig::from(-7), UBig::from(3u32));
623 let (whole, frac) = mixed_number(&r);
624 assert_eq!(whole, IBig::from(-2));
625 assert_eq!(frac, RBig::from_parts(IBig::from(-1), UBig::from(3u32)));
626 }
627
628 #[test]
629 fn floor_ceil_round_trunc() {
630 let r = RBig::from_parts(IBig::from(7), UBig::from(3u32));
632 assert_eq!(rational_floor(&r), IBig::from(2));
633 assert_eq!(rational_ceil(&r), IBig::from(3));
634 assert_eq!(rational_round(&r), IBig::from(2));
635 assert_eq!(rational_truncate(&r), IBig::from(2));
636 }
637
638 #[test]
639 fn floor_ceil_negative() {
640 let r = RBig::from_parts(IBig::from(-7), UBig::from(3u32));
642 assert_eq!(rational_floor(&r), IBig::from(-3));
643 assert_eq!(rational_ceil(&r), IBig::from(-2));
644 assert_eq!(rational_truncate(&r), IBig::from(-2));
645 }
646
647 #[test]
648 fn abs_basic() {
649 let r = RBig::from_parts(IBig::from(-3), UBig::from(4u32));
650 assert_eq!(
651 rational_abs(&r),
652 RBig::from_parts(IBig::from(3), UBig::from(4u32))
653 );
654 }
655
656 #[test]
657 fn signum_basic() {
658 let pos = RBig::from_parts(IBig::from(3), UBig::from(4u32));
659 let neg = RBig::from_parts(IBig::from(-3), UBig::from(4u32));
660 assert_eq!(rational_signum(&pos), IBig::ONE);
661 assert_eq!(rational_signum(&neg), IBig::from(-1));
662 assert_eq!(rational_signum(&RBig::ZERO), IBig::ZERO);
663 }
664
665 #[test]
666 fn reciprocal_basic() {
667 let r = RBig::from_parts(IBig::from(3), UBig::from(4u32));
668 let recip = rational_reciprocal(&r).expect("ok");
669 assert_eq!(recip, RBig::from_parts(IBig::from(4), UBig::from(3u32)));
670 }
671
672 #[test]
673 fn reciprocal_negative() {
674 let r = RBig::from_parts(IBig::from(-3), UBig::from(4u32));
675 let recip = rational_reciprocal(&r).expect("ok");
676 assert_eq!(recip, RBig::from_parts(IBig::from(-4), UBig::from(3u32)));
677 }
678
679 #[test]
680 fn reciprocal_zero_errors() {
681 assert!(rational_reciprocal(&RBig::ZERO).is_err());
682 }
683
684 #[test]
685 fn reciprocal_roundtrip() {
686 let r = RBig::from_parts(IBig::from(7), UBig::from(13u32));
688 let recip = rational_reciprocal(&r).expect("ok");
689 assert_eq!(r * recip, RBig::ONE);
690 }
691
692 #[test]
693 fn pow_positive() {
694 let r = RBig::from_parts(IBig::from(2), UBig::from(3u32));
695 let result = rational_pow(&r, 2).expect("ok");
696 assert_eq!(result, RBig::from_parts(IBig::from(4), UBig::from(9u32)));
697 }
698
699 #[test]
700 fn pow_negative() {
701 let r = RBig::from_parts(IBig::from(2), UBig::from(3u32));
702 let result = rational_pow(&r, -1).expect("ok");
703 assert_eq!(result, RBig::from_parts(IBig::from(3), UBig::from(2u32)));
704 }
705
706 #[test]
707 fn pow_zero() {
708 let r = RBig::from_parts(IBig::from(5), UBig::from(7u32));
709 assert_eq!(rational_pow(&r, 0).expect("ok"), RBig::ONE);
710 }
711
712 #[test]
713 fn pow_negative_squared() {
714 let r = RBig::from_parts(IBig::from(2), UBig::from(3u32));
716 let result = rational_pow(&r, -2).expect("ok");
717 assert_eq!(result, RBig::from_parts(IBig::from(9), UBig::from(4u32)));
718 }
719
720 #[test]
721 fn div_floor_positive() {
722 let (q, r) = div_floor(&IBig::from(7), &IBig::from(3));
723 assert_eq!(q, IBig::from(2));
724 assert_eq!(r, IBig::from(1));
725 }
726
727 #[test]
728 fn div_floor_negative_dividend() {
729 let (q, r) = div_floor(&IBig::from(-7), &IBig::from(3));
731 assert_eq!(q, IBig::from(-3));
732 assert_eq!(r, IBig::from(2));
733 }
734
735 #[test]
738 fn from_integer_basic() {
739 let r = rational_from_integer(&IBig::from(42));
740 assert_eq!(r.numerator(), &IBig::from(42));
741 assert_eq!(r.denominator(), &UBig::ONE);
742 }
743
744 #[test]
745 fn from_integer_negative() {
746 let r = rational_from_integer(&IBig::from(-7));
747 assert_eq!(r.numerator(), &IBig::from(-7));
748 assert_eq!(r.denominator(), &UBig::ONE);
749 }
750
751 #[test]
752 fn from_integer_zero() {
753 let r = rational_from_integer(&IBig::ZERO);
754 assert_eq!(r, RBig::ZERO);
755 }
756
757 #[test]
758 fn is_integer_true_for_whole() {
759 let r = rbig_from_signed(&IBig::from(3), &IBig::from(1));
760 assert!(rational_is_integer(&r));
761 }
762
763 #[test]
764 fn is_integer_false_for_fraction() {
765 let r = rbig_from_signed(&IBig::from(3), &IBig::from(2));
766 assert!(!rational_is_integer(&r));
767 }
768
769 #[test]
770 fn is_integer_after_simplification() {
771 let r = RBig::from_parts(IBig::from(10), UBig::from(5u32));
773 assert!(rational_is_integer(&r));
774 }
775
776 #[test]
777 fn to_integer_round_trip() {
778 let n = IBig::from(42);
779 let r = rational_from_integer(&n);
780 assert_eq!(rational_to_integer(&r), Some(n));
781 }
782
783 #[test]
784 fn to_integer_round_trip_negative() {
785 let n = IBig::from(-12345);
786 let r = rational_from_integer(&n);
787 assert_eq!(rational_to_integer(&r), Some(n));
788 }
789
790 #[test]
791 fn to_integer_none_for_fraction() {
792 let r = RBig::from_parts(IBig::from(3), UBig::from(2u32));
793 assert_eq!(rational_to_integer(&r), None);
794 }
795}