1use super::value::{Endpoint, ErrorFlags, Ival, IvalClass, classify};
4use crate::mpfr::{
5 exp2_overflow_threshold, mpfr_get_exp, mpfr_integer, mpfr_log2, mpfr_odd, mpfr_pow, mpfr_pow2,
6 mpfr_sign, zero,
7};
8use rug::{Assign, Float, float::Round};
9
10#[derive(Clone, Copy, Debug, PartialEq, Eq)]
11enum PosIvalClass1 {
12 GreaterOrEqual = 1,
13 Less = -1,
14 Straddles = 0,
15}
16
17impl Ival {
18 pub fn pow_assign(&mut self, x: &Ival, y: &Ival) {
20 let x_lo = x.lo.as_float();
21 let x_hi = x.hi.as_float();
22
23 let hi_is_neg = mpfr_sign(x_hi) == -1 && !x_hi.is_zero();
24 let lo_is_pos = mpfr_sign(x_lo) == 1 || x_lo.is_zero();
25
26 if hi_is_neg {
27 self.pow_neg_assign(x, y);
28 } else if lo_is_pos {
29 self.pow_pos_assign(x, y);
30 } else {
31 let zero_val = zero(x.prec());
32 let (neg, pos) = x.split_at(&zero_val);
33
34 let mut neg_result = Ival::zero(self.prec());
35 let mut pos_result = Ival::zero(self.prec());
36
37 neg_result.pow_neg_assign(&neg, y);
38 pos_result.pow_pos_assign(&pos, y);
39
40 self.assign_from(&neg_result);
41 self.union_assign(pos_result);
42 if !(x.lo.immovable && x.hi.immovable) {
43 self.lo.immovable = false;
44 self.hi.immovable = false;
45 }
46 }
47 }
48
49 pub fn pow2_assign(&mut self, x: &Ival) {
50 let mut abs_x = Ival::zero(x.max_prec());
51 abs_x.pre_fabs_assign(x);
52 self.monotonic_assign(&mpfr_pow2, &abs_x);
53 }
54
55 fn pow_pos_assign(&mut self, x: &Ival, y: &Ival) {
56 let x_class = classify_pos_ival_1(x);
57 let y_class = classify(y, false);
58
59 let mk_pow =
60 |out: &mut Ival, lo_a: &Endpoint, lo_b: &Endpoint, hi_a: &Endpoint, hi_b: &Endpoint| {
61 let prec = out.prec();
62 let mut lo_out = zero(prec);
63 let mut hi_out = zero(prec);
64
65 let lo_imm = eppow(lo_a, lo_b, x_class, y_class, &mut lo_out, Round::Down);
66 let hi_imm = eppow(hi_a, hi_b, x_class, y_class, &mut hi_out, Round::Up);
67
68 let (real_lo_imm, real_hi_imm) = if lo_out.is_zero() || hi_out.is_infinite() {
69 let threshold = exp2_overflow_threshold(prec);
70
71 let log2_sum_exceeds_threshold = |exp: &Float, base: &Float| -> bool {
72 let mut log2_base = zero(prec);
73 mpfr_log2(base, &mut log2_base, Round::Zero);
74 mpfr_get_exp(exp).saturating_add(mpfr_get_exp(&log2_base))
75 > mpfr_get_exp(&threshold)
76 };
77
78 let x_class_i = x_class as i32;
79 let y_class_i = y_class as i32;
80
81 let must_overflow = hi_out.is_infinite()
82 && x_class_i * y_class_i == 1
83 && log2_sum_exceeds_threshold(lo_b.as_float(), lo_a.as_float());
84
85 let must_underflow = lo_out.is_zero()
86 && x_class_i * y_class_i == -1
87 && log2_sum_exceeds_threshold(hi_b.as_float(), hi_a.as_float());
88
89 let new_lo_imm = lo_imm
90 || must_underflow
91 || (lo_out.is_zero() && lo_a.immovable && lo_b.immovable);
92
93 let new_hi_imm = hi_imm
94 || must_underflow
95 || must_overflow
96 || (hi_out.is_infinite() && hi_a.immovable && hi_b.immovable);
97
98 (new_lo_imm, new_hi_imm)
99 } else {
100 (lo_imm, hi_imm)
101 };
102
103 out.lo.as_float_mut().assign(&lo_out);
104 out.lo.immovable = real_lo_imm;
105 out.hi.as_float_mut().assign(&hi_out);
106 out.hi.immovable = real_hi_imm;
107
108 let x_lo_zero = x.lo.as_float().is_zero();
109 out.err = x.err.union(&y.err);
110 if x_lo_zero && !matches!(y_class, IvalClass::Pos) {
111 out.err.partial = true;
112 }
113 if x.hi.as_float().is_zero()
114 && matches!(y_class, IvalClass::Neg)
115 && !y.hi.as_float().is_zero()
116 {
117 out.err.total = true;
118 }
119 };
120
121 let xlo = &x.lo;
122 let xhi = &x.hi;
123 let ylo = &y.lo;
124 let yhi = &y.hi;
125
126 match (x_class, y_class) {
127 (PosIvalClass1::GreaterOrEqual, IvalClass::Pos) => mk_pow(self, xlo, ylo, xhi, yhi),
128 (PosIvalClass1::GreaterOrEqual, IvalClass::Mix) => mk_pow(self, xhi, ylo, xhi, yhi),
129 (PosIvalClass1::GreaterOrEqual, IvalClass::Neg) => mk_pow(self, xhi, ylo, xlo, yhi),
130 (PosIvalClass1::Straddles, IvalClass::Pos) => mk_pow(self, xlo, yhi, xhi, yhi),
131 (PosIvalClass1::Straddles, IvalClass::Neg) => mk_pow(self, xhi, ylo, xlo, ylo),
132 (PosIvalClass1::Less, IvalClass::Pos) => mk_pow(self, xlo, yhi, xhi, ylo),
133 (PosIvalClass1::Less, IvalClass::Mix) => mk_pow(self, xlo, yhi, xlo, ylo),
134 (PosIvalClass1::Less, IvalClass::Neg) => mk_pow(self, xhi, yhi, xlo, ylo),
135 (PosIvalClass1::Straddles, IvalClass::Mix) => {
136 let mut tmp = self.clone();
137 mk_pow(self, xlo, yhi, xhi, yhi);
138 mk_pow(&mut tmp, xhi, ylo, xlo, ylo);
139 self.union_assign(tmp);
140 }
141 }
142 }
143
144 fn pow_neg_assign(&mut self, x: &Ival, y: &Ival) {
145 let y_lo = y.lo.as_float();
146 let y_hi = y.hi.as_float();
147
148 if y_lo == y_hi {
149 if mpfr_integer(y_lo) {
150 let mut abs_x = Ival::zero(x.max_prec());
151 abs_x.exact_fabs_assign(x);
152
153 if mpfr_odd(y_lo) {
154 self.pow_pos_assign(&abs_x, y);
155 self.neg_inplace();
156 } else {
157 self.pow_pos_assign(&abs_x, y);
158 }
159 } else {
160 self.lo.as_float_mut().assign(f64::NAN);
161 self.hi.as_float_mut().assign(f64::NAN);
162 self.lo.immovable = true;
163 self.hi.immovable = true;
164 self.err = ErrorFlags::error();
165 }
166 } else {
167 let mut abs_x = Ival::zero(x.max_prec());
168 abs_x.exact_fabs_assign(x);
169
170 let mut pos_pow = Ival::zero(self.prec());
171 let mut neg_pow = Ival::zero(self.prec());
172
173 pos_pow.pow_pos_assign(&abs_x, y);
174 neg_pow.neg_assign(&pos_pow);
175
176 self.assign_from(&pos_pow);
177 self.union_assign(neg_pow);
178 self.err.partial = true;
179 }
180 }
181}
182
183fn classify_pos_ival_1(x: &Ival) -> PosIvalClass1 {
184 let x_lo = x.lo.as_float();
185 if mpfr_get_exp(x_lo) >= 1 {
186 return PosIvalClass1::GreaterOrEqual;
187 }
188
189 let x_hi = x.hi.as_float();
190 if mpfr_get_exp(x_hi) < 1 && !x_hi.is_infinite() {
191 return PosIvalClass1::Less;
192 }
193
194 PosIvalClass1::Straddles
195}
196
197fn eppow(
198 a: &Endpoint,
199 b: &Endpoint,
200 a_class: PosIvalClass1,
201 b_class: IvalClass,
202 out: &mut Float,
203 rnd: Round,
204) -> bool {
205 let base_val = a.as_float();
206 let exp_val = b.as_float();
207
208 let base = if base_val.is_zero() {
209 &zero(base_val.prec())
210 } else {
211 base_val
212 };
213
214 let exact = mpfr_pow(base, exp_val, out, rnd);
215
216 let a_imm = a.immovable;
217 let b_imm = b.immovable;
218
219 (a_imm && b_imm && exact)
220 || (a_imm && base == &Float::with_val(base.prec(), 1))
221 || (a_imm && base.is_zero() && !matches!(b_class, IvalClass::Mix))
222 || (a_imm && base.is_infinite() && !matches!(b_class, IvalClass::Mix))
223 || (b_imm && exp_val.is_zero())
224 || (b_imm && exp_val.is_infinite() && !matches!(a_class, PosIvalClass1::Straddles))
225}