1//===-- Utils used by both MPCWrapper and MPFRWrapper----------------------===//
2//
3// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
4// See https://llvm.org/LICENSE.txt for license information.
5// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
6//
7//===----------------------------------------------------------------------===//
8
9#include "MPCommon.h"
10
11#include "src/__support/CPP/string_view.h"
12#include "src/__support/FPUtil/bfloat16.h"
13#include "src/__support/FPUtil/cast.h"
14#include "src/__support/macros/config.h"
15#include "src/__support/macros/properties/types.h"
16
17namespace LIBC_NAMESPACE_DECL {
18namespace testing {
19namespace mpfr {
20
21MPFRNumber::MPFRNumber() : mpfr_precision(256), mpfr_rounding(MPFR_RNDN) {
22 mpfr_init2(value, mpfr_precision);
23}
24
25MPFRNumber::MPFRNumber(const MPFRNumber &other)
26 : mpfr_precision(other.mpfr_precision), mpfr_rounding(other.mpfr_rounding) {
27 mpfr_init2(value, mpfr_precision);
28 mpfr_set(value, other.value, mpfr_rounding);
29}
30
31MPFRNumber::MPFRNumber(const MPFRNumber &other, unsigned int precision)
32 : mpfr_precision(precision), mpfr_rounding(other.mpfr_rounding) {
33 mpfr_init2(value, mpfr_precision);
34 mpfr_set(value, other.value, mpfr_rounding);
35}
36
37MPFRNumber::MPFRNumber(const mpfr_t x, unsigned int precision,
38 RoundingMode rounding)
39 : mpfr_precision(precision),
40 mpfr_rounding(get_mpfr_rounding_mode(rounding)) {
41 mpfr_init2(value, mpfr_precision);
42 mpfr_set(value, x, mpfr_rounding);
43}
44
45MPFRNumber::~MPFRNumber() { mpfr_clear(value); }
46
47MPFRNumber &MPFRNumber::operator=(const MPFRNumber &rhs) {
48 mpfr_precision = rhs.mpfr_precision;
49 mpfr_rounding = rhs.mpfr_rounding;
50 mpfr_set(value, rhs.value, mpfr_rounding);
51 return *this;
52}
53
54bool MPFRNumber::is_nan() const { return mpfr_nan_p(value); }
55
56MPFRNumber MPFRNumber::abs() const {
57 MPFRNumber result(*this);
58 mpfr_abs(result.value, value, mpfr_rounding);
59 return result;
60}
61
62MPFRNumber MPFRNumber::acos() const {
63 MPFRNumber result(*this);
64 mpfr_acos(result.value, value, mpfr_rounding);
65 return result;
66}
67
68MPFRNumber MPFRNumber::acosh() const {
69 MPFRNumber result(*this);
70 mpfr_acosh(result.value, value, mpfr_rounding);
71 return result;
72}
73
74MPFRNumber MPFRNumber::acospi() const {
75 MPFRNumber result(*this);
76
77#if MPFR_VERSION >= MPFR_VERSION_NUM(4, 2, 0)
78 mpfr_acospi(result.value, value, mpfr_rounding);
79 return result;
80#else
81 MPFRNumber value_acos(0.0, 1280);
82 mpfr_acos(value_acos.value, value, MPFR_RNDN);
83 MPFRNumber value_pi(0.0, 1280);
84 mpfr_const_pi(value_pi.value, MPFR_RNDN);
85 mpfr_div(result.value, value_acos.value, value_pi.value, mpfr_rounding);
86 return result;
87#endif
88}
89
90MPFRNumber MPFRNumber::add(const MPFRNumber &b) const {
91 MPFRNumber result(*this);
92 mpfr_add(result.value, value, b.value, mpfr_rounding);
93 return result;
94}
95
96MPFRNumber MPFRNumber::asin() const {
97 MPFRNumber result(*this);
98 mpfr_asin(result.value, value, mpfr_rounding);
99 return result;
100}
101
102MPFRNumber MPFRNumber::asinh() const {
103 MPFRNumber result(*this);
104 mpfr_asinh(result.value, value, mpfr_rounding);
105 return result;
106}
107
108MPFRNumber MPFRNumber::atan() const {
109 MPFRNumber result(*this);
110 mpfr_atan(result.value, value, mpfr_rounding);
111 return result;
112}
113
114MPFRNumber MPFRNumber::atan2(const MPFRNumber &b) {
115 MPFRNumber result(*this);
116 mpfr_atan2(result.value, value, b.value, mpfr_rounding);
117 return result;
118}
119
120MPFRNumber MPFRNumber::atanh() const {
121 MPFRNumber result(*this);
122 mpfr_atanh(result.value, value, mpfr_rounding);
123 return result;
124}
125
126MPFRNumber MPFRNumber::cbrt() const {
127 MPFRNumber result(*this);
128 mpfr_cbrt(result.value, value, mpfr_rounding);
129 return result;
130}
131
132MPFRNumber MPFRNumber::ceil() const {
133 MPFRNumber result(*this);
134 mpfr_ceil(result.value, value);
135 return result;
136}
137
138MPFRNumber MPFRNumber::cos() const {
139 MPFRNumber result(*this);
140 mpfr_cos(result.value, value, mpfr_rounding);
141 return result;
142}
143
144MPFRNumber MPFRNumber::cosh() const {
145 MPFRNumber result(*this);
146 mpfr_cosh(result.value, value, mpfr_rounding);
147 return result;
148}
149
150MPFRNumber MPFRNumber::cospi() const {
151 MPFRNumber result(*this);
152
153#if MPFR_VERSION >= MPFR_VERSION_NUM(4, 2, 0)
154 mpfr_cospi(result.value, value, mpfr_rounding);
155 return result;
156#else
157 if (mpfr_integer_p(value)) {
158 mpz_t integer;
159 mpz_init(integer);
160 mpfr_get_z(integer, value, mpfr_rounding);
161
162 int d = mpz_tstbit(integer, 0);
163 mpfr_set_si(result.value, d ? -1 : 1, mpfr_rounding);
164 mpz_clear(integer);
165 return result;
166 }
167
168 MPFRNumber value_pi(0.0, 1280);
169 mpfr_const_pi(value_pi.value, MPFR_RNDN);
170 mpfr_mul(value_pi.value, value_pi.value, value, MPFR_RNDN);
171 mpfr_cos(result.value, value_pi.value, mpfr_rounding);
172
173 return result;
174#endif
175}
176
177MPFRNumber MPFRNumber::erf() const {
178 MPFRNumber result(*this);
179 mpfr_erf(result.value, value, mpfr_rounding);
180 return result;
181}
182
183MPFRNumber MPFRNumber::exp() const {
184 MPFRNumber result(*this);
185 mpfr_exp(result.value, value, mpfr_rounding);
186 return result;
187}
188
189MPFRNumber MPFRNumber::exp2() const {
190 MPFRNumber result(*this);
191 mpfr_exp2(result.value, value, mpfr_rounding);
192 return result;
193}
194
195MPFRNumber MPFRNumber::exp2m1() const {
196 // TODO: Only use mpfr_exp2m1 once CI and buildbots get MPFR >= 4.2.0.
197#if MPFR_VERSION >= MPFR_VERSION_NUM(4, 2, 0)
198 MPFRNumber result(*this);
199 mpfr_exp2m1(result.value, value, mpfr_rounding);
200 return result;
201#else
202 unsigned int prec = mpfr_precision * 3;
203 MPFRNumber result(*this, prec);
204
205 float f = mpfr_get_flt(abs().value, mpfr_rounding);
206 if (f > 0.5f && f < 0x1.0p30f) {
207 mpfr_exp2(result.value, value, mpfr_rounding);
208 mpfr_sub_ui(result.value, result.value, 1, mpfr_rounding);
209 return result;
210 }
211
212 MPFRNumber ln2(2.0f, prec);
213 // log(2)
214 mpfr_log(ln2.value, ln2.value, mpfr_rounding);
215 // x * log(2)
216 mpfr_mul(result.value, value, ln2.value, mpfr_rounding);
217 // e^(x * log(2)) - 1
218 int ex = mpfr_expm1(result.value, result.value, mpfr_rounding);
219 mpfr_subnormalize(result.value, ex, mpfr_rounding);
220 return result;
221#endif
222}
223
224MPFRNumber MPFRNumber::exp10() const {
225 MPFRNumber result(*this);
226 mpfr_exp10(result.value, value, mpfr_rounding);
227 return result;
228}
229
230MPFRNumber MPFRNumber::exp10m1() const {
231 // TODO: Only use mpfr_exp10m1 once CI and buildbots get MPFR >= 4.2.0.
232#if MPFR_VERSION >= MPFR_VERSION_NUM(4, 2, 0)
233 MPFRNumber result(*this);
234 mpfr_exp10m1(result.value, value, mpfr_rounding);
235 return result;
236#else
237 unsigned int prec = mpfr_precision * 3;
238 MPFRNumber result(*this, prec);
239
240 MPFRNumber ln10(10.0f, prec);
241 // log(10)
242 mpfr_log(ln10.value, ln10.value, mpfr_rounding);
243 // x * log(10)
244 mpfr_mul(result.value, value, ln10.value, mpfr_rounding);
245 // e^(x * log(10)) - 1
246 int ex = mpfr_expm1(result.value, result.value, mpfr_rounding);
247 mpfr_subnormalize(result.value, ex, mpfr_rounding);
248 return result;
249#endif
250}
251
252MPFRNumber MPFRNumber::expm1() const {
253 MPFRNumber result(*this);
254 mpfr_expm1(result.value, value, mpfr_rounding);
255 return result;
256}
257
258MPFRNumber MPFRNumber::div(const MPFRNumber &b) const {
259 MPFRNumber result(*this);
260 mpfr_div(result.value, value, b.value, mpfr_rounding);
261 return result;
262}
263
264MPFRNumber MPFRNumber::floor() const {
265 MPFRNumber result(*this);
266 mpfr_floor(result.value, value);
267 return result;
268}
269
270MPFRNumber MPFRNumber::fmod(const MPFRNumber &b) {
271 MPFRNumber result(*this);
272 mpfr_fmod(result.value, value, b.value, mpfr_rounding);
273 return result;
274}
275
276MPFRNumber MPFRNumber::frexp(int &exp) {
277 MPFRNumber result(*this);
278 mpfr_exp_t resultExp;
279 mpfr_frexp(&resultExp, result.value, value, mpfr_rounding);
280 exp = static_cast<int>(resultExp);
281 return result;
282}
283
284MPFRNumber MPFRNumber::hypot(const MPFRNumber &b) {
285 MPFRNumber result(*this);
286 mpfr_hypot(result.value, value, b.value, mpfr_rounding);
287 return result;
288}
289
290MPFRNumber MPFRNumber::log() const {
291 MPFRNumber result(*this);
292 mpfr_log(result.value, value, mpfr_rounding);
293 return result;
294}
295
296MPFRNumber MPFRNumber::log2() const {
297 MPFRNumber result(*this);
298 mpfr_log2(result.value, value, mpfr_rounding);
299 return result;
300}
301
302MPFRNumber MPFRNumber::log10() const {
303 MPFRNumber result(*this);
304 mpfr_log10(result.value, value, mpfr_rounding);
305 return result;
306}
307
308MPFRNumber MPFRNumber::log1p() const {
309 MPFRNumber result(*this);
310 mpfr_log1p(result.value, value, mpfr_rounding);
311 return result;
312}
313
314MPFRNumber MPFRNumber::pow(const MPFRNumber &b) {
315 MPFRNumber result(*this);
316 mpfr_pow(result.value, value, b.value, mpfr_rounding);
317 return result;
318}
319
320MPFRNumber MPFRNumber::remquo(const MPFRNumber &divisor, int &quotient) {
321 MPFRNumber remainder(*this);
322 long q;
323 mpfr_remquo(remainder.value, &q, value, divisor.value, mpfr_rounding);
324 quotient = static_cast<int>(q);
325 return remainder;
326}
327
328MPFRNumber MPFRNumber::round() const {
329 MPFRNumber result(*this);
330 mpfr_round(result.value, value);
331 return result;
332}
333
334MPFRNumber MPFRNumber::roundeven() const {
335 MPFRNumber result(*this);
336#if MPFR_VERSION_MAJOR >= 4
337 mpfr_roundeven(result.value, value);
338#else
339 mpfr_rint(result.value, value, MPFR_RNDN);
340#endif
341 return result;
342}
343
344bool MPFRNumber::round_to_long(long &result) const {
345 // We first calculate the rounded value. This way, when converting
346 // to long using mpfr_get_si, the rounding direction of MPFR_RNDN
347 // (or any other rounding mode), does not have an influence.
348 MPFRNumber roundedValue = round();
349 mpfr_clear_erangeflag();
350 result = mpfr_get_si(roundedValue.value, MPFR_RNDN);
351 return mpfr_erangeflag_p();
352}
353
354bool MPFRNumber::round_to_long(mpfr_rnd_t rnd, long &result) const {
355 MPFRNumber rint_result(*this);
356 mpfr_rint(rint_result.value, value, rnd);
357 return rint_result.round_to_long(result);
358}
359
360MPFRNumber MPFRNumber::rint(mpfr_rnd_t rnd) const {
361 MPFRNumber result(*this);
362 mpfr_rint(result.value, value, rnd);
363 return result;
364}
365
366MPFRNumber MPFRNumber::mod_2pi() const {
367 MPFRNumber result(0.0, 1280);
368 MPFRNumber _2pi(0.0, 1280);
369 mpfr_const_pi(_2pi.value, MPFR_RNDN);
370 mpfr_mul_si(_2pi.value, _2pi.value, 2, MPFR_RNDN);
371 mpfr_fmod(result.value, value, _2pi.value, MPFR_RNDN);
372 return result;
373}
374
375MPFRNumber MPFRNumber::mod_pi_over_2() const {
376 MPFRNumber result(0.0, 1280);
377 MPFRNumber pi_over_2(0.0, 1280);
378 mpfr_const_pi(pi_over_2.value, MPFR_RNDN);
379 mpfr_mul_d(pi_over_2.value, pi_over_2.value, 0.5, MPFR_RNDN);
380 mpfr_fmod(result.value, value, pi_over_2.value, MPFR_RNDN);
381 return result;
382}
383
384MPFRNumber MPFRNumber::mod_pi_over_4() const {
385 MPFRNumber result(0.0, 1280);
386 MPFRNumber pi_over_4(0.0, 1280);
387 mpfr_const_pi(pi_over_4.value, MPFR_RNDN);
388 mpfr_mul_d(pi_over_4.value, pi_over_4.value, 0.25, MPFR_RNDN);
389 mpfr_fmod(result.value, value, pi_over_4.value, MPFR_RNDN);
390 return result;
391}
392
393MPFRNumber MPFRNumber::sin() const {
394 MPFRNumber result(*this);
395 mpfr_sin(result.value, value, mpfr_rounding);
396 return result;
397}
398
399MPFRNumber MPFRNumber::sinpi() const {
400 MPFRNumber result(*this);
401
402#if MPFR_VERSION >= MPFR_VERSION_NUM(4, 2, 0)
403 mpfr_sinpi(result.value, value, mpfr_rounding);
404 return result;
405#else
406 if (mpfr_integer_p(value)) {
407 mpfr_set_si(result.value, 0, mpfr_rounding);
408 return result;
409 }
410
411 MPFRNumber value_mul_two(*this);
412 mpfr_mul_si(value_mul_two.value, value, 2, MPFR_RNDN);
413
414 if (mpfr_integer_p(value_mul_two.value)) {
415 auto d = mpfr_get_si(value, MPFR_RNDD);
416 mpfr_set_si(result.value, (d & 1) ? -1 : 1, mpfr_rounding);
417 return result;
418 }
419
420 MPFRNumber value_pi(0.0, 1280);
421 mpfr_const_pi(value_pi.value, MPFR_RNDN);
422 mpfr_mul(value_pi.value, value_pi.value, value, MPFR_RNDN);
423 mpfr_sin(result.value, value_pi.value, mpfr_rounding);
424 return result;
425#endif
426}
427
428MPFRNumber MPFRNumber::sinh() const {
429 MPFRNumber result(*this);
430 mpfr_sinh(result.value, value, mpfr_rounding);
431 return result;
432}
433
434MPFRNumber MPFRNumber::sqrt() const {
435 MPFRNumber result(*this);
436 mpfr_sqrt(result.value, value, mpfr_rounding);
437 return result;
438}
439
440MPFRNumber MPFRNumber::sub(const MPFRNumber &b) const {
441 MPFRNumber result(*this);
442 mpfr_sub(result.value, value, b.value, mpfr_rounding);
443 return result;
444}
445
446MPFRNumber MPFRNumber::tan() const {
447 MPFRNumber result(*this);
448 mpfr_tan(result.value, value, mpfr_rounding);
449 return result;
450}
451
452MPFRNumber MPFRNumber::tanh() const {
453 MPFRNumber result(*this);
454 mpfr_tanh(result.value, value, mpfr_rounding);
455 return result;
456}
457
458MPFRNumber MPFRNumber::tanpi() const {
459 MPFRNumber result(*this);
460
461#if MPFR_VERSION >= MPFR_VERSION_NUM(4, 2, 0)
462 mpfr_tanpi(result.value, value, mpfr_rounding);
463 return result;
464#else
465 MPFRNumber value_ret_exact(*this);
466 MPFRNumber value_one(*this);
467 mpfr_set_si(value_one.value, 1, MPFR_RNDN);
468 mpfr_fmod(value_ret_exact.value, value, value_one.value, mpfr_rounding);
469 mpfr_mul_si(value_ret_exact.value, value_ret_exact.value, 4, MPFR_RNDN);
470
471 if (mpfr_integer_p(value_ret_exact.value)) {
472 int mod = mpfr_get_si(value_ret_exact.value, MPFR_RNDN);
473 mod = (mod < 0 ? -1 * mod : mod);
474
475 switch (mod) {
476 case 0:
477 mpfr_set_si(result.value, 0, mpfr_rounding);
478 break;
479 case 1:
480 mpfr_set_si(result.value, (mpfr_signbit(value) ? -1 : 1), mpfr_rounding);
481 break;
482 case 2: {
483 auto d = mpfr_get_si(value, MPFR_RNDZ);
484 d += mpfr_sgn(value) > 0 ? 0 : 1;
485 mpfr_set_inf(result.value, (d & 1) ? -1 : 1);
486 break;
487 }
488 case 3:
489 mpfr_set_si(result.value, (mpfr_signbit(value) ? 1 : -1), mpfr_rounding);
490 break;
491 }
492
493 return result;
494 }
495
496 MPFRNumber value_pi(0.0, 1280);
497 mpfr_const_pi(value_pi.value, MPFR_RNDN);
498 mpfr_mul(value_pi.value, value_pi.value, value, MPFR_RNDN);
499 mpfr_tan(result.value, value_pi.value, mpfr_rounding);
500 return result;
501#endif
502}
503
504MPFRNumber MPFRNumber::trunc() const {
505 MPFRNumber result(*this);
506 mpfr_trunc(result.value, value);
507 return result;
508}
509
510MPFRNumber MPFRNumber::fma(const MPFRNumber &b, const MPFRNumber &c) {
511 MPFRNumber result(*this);
512 mpfr_fma(result.value, value, b.value, c.value, mpfr_rounding);
513 return result;
514}
515
516MPFRNumber MPFRNumber::mul(const MPFRNumber &b) {
517 MPFRNumber result(*this);
518 mpfr_mul(result.value, value, b.value, mpfr_rounding);
519 return result;
520}
521
522cpp::string MPFRNumber::str() const {
523 // 200 bytes should be more than sufficient to hold a 100-digit number
524 // plus additional bytes for the decimal point, '-' sign etc.
525 constexpr size_t printBufSize = 200;
526 char buffer[printBufSize];
527 mpfr_snprintf(buffer, printBufSize, "%100.50Rf", value);
528 cpp::string_view view(buffer);
529 // Trim whitespaces
530 const char whitespace = ' ';
531 while (!view.empty() && view.front() == whitespace)
532 view.remove_prefix(1);
533 while (!view.empty() && view.back() == whitespace)
534 view.remove_suffix(1);
535 return cpp::string(view.data());
536}
537
538void MPFRNumber::dump(const char *msg) const {
539 mpfr_printf("%s%.128g\n", msg, value);
540}
541
542template <> float MPFRNumber::as<float>() const {
543 return mpfr_get_flt(value, mpfr_rounding);
544}
545
546template <> double MPFRNumber::as<double>() const {
547 return mpfr_get_d(value, mpfr_rounding);
548}
549
550template <> long double MPFRNumber::as<long double>() const {
551 return mpfr_get_ld(value, mpfr_rounding);
552}
553
554#ifdef LIBC_TYPES_HAS_FLOAT16
555template <> float16 MPFRNumber::as<float16>() const {
556 // TODO: Either prove that this cast won't cause double-rounding errors, or
557 // find a better way to get a float16.
558 return fputil::cast<float16>(mpfr_get_d(value, mpfr_rounding));
559}
560#endif
561
562#ifdef LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
563template <> float128 MPFRNumber::as<float128>() const {
564 return mpfr_get_float128(value, mpfr_rounding);
565}
566#endif // LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
567
568template <> bfloat16 MPFRNumber::as<bfloat16>() const {
569 return fputil::cast<bfloat16>(mpfr_get_flt(value, mpfr_rounding));
570}
571
572} // namespace mpfr
573} // namespace testing
574} // namespace LIBC_NAMESPACE_DECL
575

source code of libc/utils/MPFRWrapper/MPCommon.cpp