xref: /llvm-project/libc/utils/MPFRWrapper/MPFRUtils.cpp (revision 7f37b34d31914120a5bb6bd341e7616773df7613)
1 //===-- Utils which wrap MPFR ---------------------------------------------===//
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 "MPFRUtils.h"
10 #include "MPCommon.h"
11 
12 #include "src/__support/CPP/array.h"
13 #include "src/__support/CPP/stringstream.h"
14 #include "src/__support/FPUtil/fpbits_str.h"
15 #include "src/__support/macros/config.h"
16 #include "src/__support/macros/properties/types.h"
17 
18 namespace LIBC_NAMESPACE_DECL {
19 namespace testing {
20 namespace mpfr {
21 namespace internal {
22 
23 template <typename InputType>
24 cpp::enable_if_t<cpp::is_floating_point_v<InputType>, MPFRNumber>
25 unary_operation(Operation op, InputType input, unsigned int precision,
26                 RoundingMode rounding) {
27   MPFRNumber mpfrInput(input, precision, rounding);
28   switch (op) {
29   case Operation::Abs:
30     return mpfrInput.abs();
31   case Operation::Acos:
32     return mpfrInput.acos();
33   case Operation::Acosh:
34     return mpfrInput.acosh();
35   case Operation::Asin:
36     return mpfrInput.asin();
37   case Operation::Asinh:
38     return mpfrInput.asinh();
39   case Operation::Atan:
40     return mpfrInput.atan();
41   case Operation::Atanh:
42     return mpfrInput.atanh();
43   case Operation::Cbrt:
44     return mpfrInput.cbrt();
45   case Operation::Ceil:
46     return mpfrInput.ceil();
47   case Operation::Cos:
48     return mpfrInput.cos();
49   case Operation::Cosh:
50     return mpfrInput.cosh();
51   case Operation::Cospi:
52     return mpfrInput.cospi();
53   case Operation::Erf:
54     return mpfrInput.erf();
55   case Operation::Exp:
56     return mpfrInput.exp();
57   case Operation::Exp2:
58     return mpfrInput.exp2();
59   case Operation::Exp2m1:
60     return mpfrInput.exp2m1();
61   case Operation::Exp10:
62     return mpfrInput.exp10();
63   case Operation::Exp10m1:
64     return mpfrInput.exp10m1();
65   case Operation::Expm1:
66     return mpfrInput.expm1();
67   case Operation::Floor:
68     return mpfrInput.floor();
69   case Operation::Log:
70     return mpfrInput.log();
71   case Operation::Log2:
72     return mpfrInput.log2();
73   case Operation::Log10:
74     return mpfrInput.log10();
75   case Operation::Log1p:
76     return mpfrInput.log1p();
77   case Operation::Mod2PI:
78     return mpfrInput.mod_2pi();
79   case Operation::ModPIOver2:
80     return mpfrInput.mod_pi_over_2();
81   case Operation::ModPIOver4:
82     return mpfrInput.mod_pi_over_4();
83   case Operation::Round:
84     return mpfrInput.round();
85   case Operation::RoundEven:
86     return mpfrInput.roundeven();
87   case Operation::Sin:
88     return mpfrInput.sin();
89   case Operation::Sinpi:
90     return mpfrInput.sinpi();
91   case Operation::Sinh:
92     return mpfrInput.sinh();
93   case Operation::Sqrt:
94     return mpfrInput.sqrt();
95   case Operation::Tan:
96     return mpfrInput.tan();
97   case Operation::Tanh:
98     return mpfrInput.tanh();
99   case Operation::Tanpi:
100     return mpfrInput.tanpi();
101   case Operation::Trunc:
102     return mpfrInput.trunc();
103   default:
104     __builtin_unreachable();
105   }
106 }
107 
108 template <typename InputType>
109 cpp::enable_if_t<cpp::is_floating_point_v<InputType>, MPFRNumber>
110 unary_operation_two_outputs(Operation op, InputType input, int &output,
111                             unsigned int precision, RoundingMode rounding) {
112   MPFRNumber mpfrInput(input, precision, rounding);
113   switch (op) {
114   case Operation::Frexp:
115     return mpfrInput.frexp(output);
116   default:
117     __builtin_unreachable();
118   }
119 }
120 
121 template <typename InputType>
122 cpp::enable_if_t<cpp::is_floating_point_v<InputType>, MPFRNumber>
123 binary_operation_one_output(Operation op, InputType x, InputType y,
124                             unsigned int precision, RoundingMode rounding) {
125   MPFRNumber inputX(x, precision, rounding);
126   MPFRNumber inputY(y, precision, rounding);
127   switch (op) {
128   case Operation::Add:
129     return inputX.add(inputY);
130   case Operation::Atan2:
131     return inputX.atan2(inputY);
132   case Operation::Div:
133     return inputX.div(inputY);
134   case Operation::Fmod:
135     return inputX.fmod(inputY);
136   case Operation::Hypot:
137     return inputX.hypot(inputY);
138   case Operation::Mul:
139     return inputX.mul(inputY);
140   case Operation::Pow:
141     return inputX.pow(inputY);
142   case Operation::Sub:
143     return inputX.sub(inputY);
144   default:
145     __builtin_unreachable();
146   }
147 }
148 
149 template <typename InputType>
150 cpp::enable_if_t<cpp::is_floating_point_v<InputType>, MPFRNumber>
151 binary_operation_two_outputs(Operation op, InputType x, InputType y,
152                              int &output, unsigned int precision,
153                              RoundingMode rounding) {
154   MPFRNumber inputX(x, precision, rounding);
155   MPFRNumber inputY(y, precision, rounding);
156   switch (op) {
157   case Operation::RemQuo:
158     return inputX.remquo(inputY, output);
159   default:
160     __builtin_unreachable();
161   }
162 }
163 
164 template <typename InputType>
165 cpp::enable_if_t<cpp::is_floating_point_v<InputType>, MPFRNumber>
166 ternary_operation_one_output(Operation op, InputType x, InputType y,
167                              InputType z, unsigned int precision,
168                              RoundingMode rounding) {
169   // For FMA function, we just need to compare with the mpfr_fma with the same
170   // precision as InputType.  Using higher precision as the intermediate results
171   // to compare might incorrectly fail due to double-rounding errors.
172   MPFRNumber inputX(x, precision, rounding);
173   MPFRNumber inputY(y, precision, rounding);
174   MPFRNumber inputZ(z, precision, rounding);
175   switch (op) {
176   case Operation::Fma:
177     return inputX.fma(inputY, inputZ);
178   default:
179     __builtin_unreachable();
180   }
181 }
182 
183 // Remark: For all the explain_*_error functions, we will use std::stringstream
184 // to build the complete error messages before sending it to the outstream `OS`
185 // once at the end.  This will stop the error messages from interleaving when
186 // the tests are running concurrently.
187 template <typename InputType, typename OutputType>
188 void explain_unary_operation_single_output_error(Operation op, InputType input,
189                                                  OutputType matchValue,
190                                                  double ulp_tolerance,
191                                                  RoundingMode rounding) {
192   unsigned int precision = get_precision<InputType>(ulp_tolerance);
193   MPFRNumber mpfrInput(input, precision);
194   MPFRNumber mpfr_result;
195   mpfr_result = unary_operation(op, input, precision, rounding);
196   MPFRNumber mpfrMatchValue(matchValue);
197   cpp::array<char, 1024> msg_buf;
198   cpp::StringStream msg(msg_buf);
199   msg << "Match value not within tolerance value of MPFR result:\n"
200       << "  Input decimal: " << mpfrInput.str() << '\n';
201   msg << "     Input bits: " << str(FPBits<InputType>(input)) << '\n';
202   msg << '\n' << "  Match decimal: " << mpfrMatchValue.str() << '\n';
203   msg << "     Match bits: " << str(FPBits<OutputType>(matchValue)) << '\n';
204   msg << '\n' << "    MPFR result: " << mpfr_result.str() << '\n';
205   msg << "   MPFR rounded: "
206       << str(FPBits<OutputType>(mpfr_result.as<OutputType>())) << '\n';
207   msg << '\n';
208   msg << "      ULP error: " << mpfr_result.ulp_as_mpfr_number(matchValue).str()
209       << '\n';
210   if (msg.overflow())
211     __builtin_unreachable();
212   tlog << msg.str();
213 }
214 
215 template void explain_unary_operation_single_output_error(Operation op, float,
216                                                           float, double,
217                                                           RoundingMode);
218 template void explain_unary_operation_single_output_error(Operation op, double,
219                                                           double, double,
220                                                           RoundingMode);
221 template void explain_unary_operation_single_output_error(Operation op,
222                                                           long double,
223                                                           long double, double,
224                                                           RoundingMode);
225 template void explain_unary_operation_single_output_error(Operation op, double,
226                                                           float, double,
227                                                           RoundingMode);
228 template void explain_unary_operation_single_output_error(Operation op,
229                                                           long double, float,
230                                                           double, RoundingMode);
231 template void explain_unary_operation_single_output_error(Operation op,
232                                                           long double, double,
233                                                           double, RoundingMode);
234 
235 #ifdef LIBC_TYPES_HAS_FLOAT16
236 template void explain_unary_operation_single_output_error(Operation op, float16,
237                                                           float16, double,
238                                                           RoundingMode);
239 template void explain_unary_operation_single_output_error(Operation op, float,
240                                                           float16, double,
241                                                           RoundingMode);
242 template void explain_unary_operation_single_output_error(Operation op, double,
243                                                           float16, double,
244                                                           RoundingMode);
245 template void explain_unary_operation_single_output_error(Operation op,
246                                                           long double, float16,
247                                                           double, RoundingMode);
248 #ifdef LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
249 template void explain_unary_operation_single_output_error(Operation op,
250                                                           float128, float16,
251                                                           double, RoundingMode);
252 #endif // LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
253 #endif // LIBC_TYPES_HAS_FLOAT16
254 
255 #ifdef LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
256 template void explain_unary_operation_single_output_error(Operation op,
257                                                           float128, float128,
258                                                           double, RoundingMode);
259 template void explain_unary_operation_single_output_error(Operation op,
260                                                           float128, float,
261                                                           double, RoundingMode);
262 template void explain_unary_operation_single_output_error(Operation op,
263                                                           float128, double,
264                                                           double, RoundingMode);
265 template void explain_unary_operation_single_output_error(Operation op,
266                                                           float128, long double,
267                                                           double, RoundingMode);
268 #endif // LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
269 
270 template <typename T>
271 void explain_unary_operation_two_outputs_error(
272     Operation op, T input, const BinaryOutput<T> &libc_result,
273     double ulp_tolerance, RoundingMode rounding) {
274   unsigned int precision = get_precision<T>(ulp_tolerance);
275   MPFRNumber mpfrInput(input, precision);
276   int mpfrIntResult;
277   MPFRNumber mpfr_result = unary_operation_two_outputs(op, input, mpfrIntResult,
278                                                        precision, rounding);
279 
280   if (mpfrIntResult != libc_result.i) {
281     tlog << "MPFR integral result: " << mpfrIntResult << '\n'
282          << "Libc integral result: " << libc_result.i << '\n';
283   } else {
284     tlog << "Integral result from libc matches integral result from MPFR.\n";
285   }
286 
287   MPFRNumber mpfrMatchValue(libc_result.f);
288   tlog
289       << "Libc floating point result is not within tolerance value of the MPFR "
290       << "result.\n\n";
291 
292   tlog << "            Input decimal: " << mpfrInput.str() << "\n\n";
293 
294   tlog << "Libc floating point value: " << mpfrMatchValue.str() << '\n';
295   tlog << " Libc floating point bits: " << str(FPBits<T>(libc_result.f))
296        << '\n';
297   tlog << "\n\n";
298 
299   tlog << "              MPFR result: " << mpfr_result.str() << '\n';
300   tlog << "             MPFR rounded: " << str(FPBits<T>(mpfr_result.as<T>()))
301        << '\n';
302   tlog << '\n'
303        << "                ULP error: "
304        << mpfr_result.ulp_as_mpfr_number(libc_result.f).str() << '\n';
305 }
306 
307 template void explain_unary_operation_two_outputs_error<float>(
308     Operation, float, const BinaryOutput<float> &, double, RoundingMode);
309 template void explain_unary_operation_two_outputs_error<double>(
310     Operation, double, const BinaryOutput<double> &, double, RoundingMode);
311 template void explain_unary_operation_two_outputs_error<long double>(
312     Operation, long double, const BinaryOutput<long double> &, double,
313     RoundingMode);
314 
315 template <typename T>
316 void explain_binary_operation_two_outputs_error(
317     Operation op, const BinaryInput<T> &input,
318     const BinaryOutput<T> &libc_result, double ulp_tolerance,
319     RoundingMode rounding) {
320   unsigned int precision = get_precision<T>(ulp_tolerance);
321   MPFRNumber mpfrX(input.x, precision);
322   MPFRNumber mpfrY(input.y, precision);
323   int mpfrIntResult;
324   MPFRNumber mpfr_result = binary_operation_two_outputs(
325       op, input.x, input.y, mpfrIntResult, precision, rounding);
326   MPFRNumber mpfrMatchValue(libc_result.f);
327 
328   tlog << "Input decimal: x: " << mpfrX.str() << " y: " << mpfrY.str() << '\n'
329        << "MPFR integral result: " << mpfrIntResult << '\n'
330        << "Libc integral result: " << libc_result.i << '\n'
331        << "Libc floating point result: " << mpfrMatchValue.str() << '\n'
332        << "               MPFR result: " << mpfr_result.str() << '\n';
333   tlog << "Libc floating point result bits: " << str(FPBits<T>(libc_result.f))
334        << '\n';
335   tlog << "              MPFR rounded bits: "
336        << str(FPBits<T>(mpfr_result.as<T>())) << '\n';
337   tlog << "ULP error: " << mpfr_result.ulp_as_mpfr_number(libc_result.f).str()
338        << '\n';
339 }
340 
341 template void explain_binary_operation_two_outputs_error<float>(
342     Operation, const BinaryInput<float> &, const BinaryOutput<float> &, double,
343     RoundingMode);
344 template void explain_binary_operation_two_outputs_error<double>(
345     Operation, const BinaryInput<double> &, const BinaryOutput<double> &,
346     double, RoundingMode);
347 template void explain_binary_operation_two_outputs_error<long double>(
348     Operation, const BinaryInput<long double> &,
349     const BinaryOutput<long double> &, double, RoundingMode);
350 
351 template <typename InputType, typename OutputType>
352 void explain_binary_operation_one_output_error(
353     Operation op, const BinaryInput<InputType> &input, OutputType libc_result,
354     double ulp_tolerance, RoundingMode rounding) {
355   unsigned int precision = get_precision<InputType>(ulp_tolerance);
356   MPFRNumber mpfrX(input.x, precision);
357   MPFRNumber mpfrY(input.y, precision);
358   FPBits<InputType> xbits(input.x);
359   FPBits<InputType> ybits(input.y);
360   MPFRNumber mpfr_result =
361       binary_operation_one_output(op, input.x, input.y, precision, rounding);
362   MPFRNumber mpfrMatchValue(libc_result);
363 
364   tlog << "Input decimal: x: " << mpfrX.str() << " y: " << mpfrY.str() << '\n';
365   tlog << "First input bits: " << str(FPBits<InputType>(input.x)) << '\n';
366   tlog << "Second input bits: " << str(FPBits<InputType>(input.y)) << '\n';
367 
368   tlog << "Libc result: " << mpfrMatchValue.str() << '\n'
369        << "MPFR result: " << mpfr_result.str() << '\n';
370   tlog << "Libc floating point result bits: "
371        << str(FPBits<OutputType>(libc_result)) << '\n';
372   tlog << "              MPFR rounded bits: "
373        << str(FPBits<OutputType>(mpfr_result.as<OutputType>())) << '\n';
374   tlog << "ULP error: " << mpfr_result.ulp_as_mpfr_number(libc_result).str()
375        << '\n';
376 }
377 
378 template void
379 explain_binary_operation_one_output_error(Operation, const BinaryInput<float> &,
380                                           float, double, RoundingMode);
381 template void explain_binary_operation_one_output_error(
382     Operation, const BinaryInput<double> &, float, double, RoundingMode);
383 template void explain_binary_operation_one_output_error(
384     Operation, const BinaryInput<double> &, double, double, RoundingMode);
385 template void explain_binary_operation_one_output_error(
386     Operation, const BinaryInput<long double> &, float, double, RoundingMode);
387 template void explain_binary_operation_one_output_error(
388     Operation, const BinaryInput<long double> &, double, double, RoundingMode);
389 template void
390 explain_binary_operation_one_output_error(Operation,
391                                           const BinaryInput<long double> &,
392                                           long double, double, RoundingMode);
393 #ifdef LIBC_TYPES_HAS_FLOAT16
394 template void explain_binary_operation_one_output_error(
395     Operation, const BinaryInput<float16> &, float16, double, RoundingMode);
396 template void
397 explain_binary_operation_one_output_error(Operation, const BinaryInput<float> &,
398                                           float16, double, RoundingMode);
399 template void explain_binary_operation_one_output_error(
400     Operation, const BinaryInput<double> &, float16, double, RoundingMode);
401 template void explain_binary_operation_one_output_error(
402     Operation, const BinaryInput<long double> &, float16, double, RoundingMode);
403 #endif
404 
405 template <typename InputType, typename OutputType>
406 void explain_ternary_operation_one_output_error(
407     Operation op, const TernaryInput<InputType> &input, OutputType libc_result,
408     double ulp_tolerance, RoundingMode rounding) {
409   unsigned int precision = get_precision<InputType>(ulp_tolerance);
410   MPFRNumber mpfrX(input.x, precision);
411   MPFRNumber mpfrY(input.y, precision);
412   MPFRNumber mpfrZ(input.z, precision);
413   FPBits<InputType> xbits(input.x);
414   FPBits<InputType> ybits(input.y);
415   FPBits<InputType> zbits(input.z);
416   MPFRNumber mpfr_result = ternary_operation_one_output(
417       op, input.x, input.y, input.z, precision, rounding);
418   MPFRNumber mpfrMatchValue(libc_result);
419 
420   tlog << "Input decimal: x: " << mpfrX.str() << " y: " << mpfrY.str()
421        << " z: " << mpfrZ.str() << '\n';
422   tlog << " First input bits: " << str(FPBits<InputType>(input.x)) << '\n';
423   tlog << "Second input bits: " << str(FPBits<InputType>(input.y)) << '\n';
424   tlog << " Third input bits: " << str(FPBits<InputType>(input.z)) << '\n';
425 
426   tlog << "Libc result: " << mpfrMatchValue.str() << '\n'
427        << "MPFR result: " << mpfr_result.str() << '\n';
428   tlog << "Libc floating point result bits: "
429        << str(FPBits<OutputType>(libc_result)) << '\n';
430   tlog << "              MPFR rounded bits: "
431        << str(FPBits<OutputType>(mpfr_result.as<OutputType>())) << '\n';
432   tlog << "ULP error: " << mpfr_result.ulp_as_mpfr_number(libc_result).str()
433        << '\n';
434 }
435 
436 template void explain_ternary_operation_one_output_error(
437     Operation, const TernaryInput<float> &, float, double, RoundingMode);
438 template void explain_ternary_operation_one_output_error(
439     Operation, const TernaryInput<double> &, float, double, RoundingMode);
440 template void explain_ternary_operation_one_output_error(
441     Operation, const TernaryInput<double> &, double, double, RoundingMode);
442 template void explain_ternary_operation_one_output_error(
443     Operation, const TernaryInput<long double> &, float, double, RoundingMode);
444 template void explain_ternary_operation_one_output_error(
445     Operation, const TernaryInput<long double> &, double, double, RoundingMode);
446 template void
447 explain_ternary_operation_one_output_error(Operation,
448                                            const TernaryInput<long double> &,
449                                            long double, double, RoundingMode);
450 
451 #ifdef LIBC_TYPES_HAS_FLOAT16
452 template void explain_ternary_operation_one_output_error(
453     Operation, const TernaryInput<float> &, float16, double, RoundingMode);
454 template void explain_ternary_operation_one_output_error(
455     Operation, const TernaryInput<double> &, float16, double, RoundingMode);
456 template void
457 explain_ternary_operation_one_output_error(Operation,
458                                            const TernaryInput<long double> &,
459                                            float16, double, RoundingMode);
460 #endif
461 
462 template <typename InputType, typename OutputType>
463 bool compare_unary_operation_single_output(Operation op, InputType input,
464                                            OutputType libc_result,
465                                            double ulp_tolerance,
466                                            RoundingMode rounding) {
467   unsigned int precision = get_precision<InputType>(ulp_tolerance);
468   MPFRNumber mpfr_result;
469   mpfr_result = unary_operation(op, input, precision, rounding);
470   double ulp = mpfr_result.ulp(libc_result);
471   return (ulp <= ulp_tolerance);
472 }
473 
474 template bool compare_unary_operation_single_output(Operation, float, float,
475                                                     double, RoundingMode);
476 template bool compare_unary_operation_single_output(Operation, double, double,
477                                                     double, RoundingMode);
478 template bool compare_unary_operation_single_output(Operation, long double,
479                                                     long double, double,
480                                                     RoundingMode);
481 template bool compare_unary_operation_single_output(Operation, double, float,
482                                                     double, RoundingMode);
483 template bool compare_unary_operation_single_output(Operation, long double,
484                                                     float, double,
485                                                     RoundingMode);
486 template bool compare_unary_operation_single_output(Operation, long double,
487                                                     double, double,
488                                                     RoundingMode);
489 #ifdef LIBC_TYPES_HAS_FLOAT16
490 template bool compare_unary_operation_single_output(Operation, float16, float16,
491                                                     double, RoundingMode);
492 template bool compare_unary_operation_single_output(Operation, float, float16,
493                                                     double, RoundingMode);
494 template bool compare_unary_operation_single_output(Operation, double, float16,
495                                                     double, RoundingMode);
496 template bool compare_unary_operation_single_output(Operation, long double,
497                                                     float16, double,
498                                                     RoundingMode);
499 #ifdef LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
500 template bool compare_unary_operation_single_output(Operation, float128,
501                                                     float16, double,
502                                                     RoundingMode);
503 #endif // LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
504 #endif // LIBC_TYPES_HAS_FLOAT16
505 
506 #ifdef LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
507 template bool compare_unary_operation_single_output(Operation, float128,
508                                                     float128, double,
509                                                     RoundingMode);
510 template bool compare_unary_operation_single_output(Operation, float128, float,
511                                                     double, RoundingMode);
512 template bool compare_unary_operation_single_output(Operation, float128, double,
513                                                     double, RoundingMode);
514 template bool compare_unary_operation_single_output(Operation, float128,
515                                                     long double, double,
516                                                     RoundingMode);
517 #endif // LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
518 
519 template <typename T>
520 bool compare_unary_operation_two_outputs(Operation op, T input,
521                                          const BinaryOutput<T> &libc_result,
522                                          double ulp_tolerance,
523                                          RoundingMode rounding) {
524   int mpfrIntResult;
525   unsigned int precision = get_precision<T>(ulp_tolerance);
526   MPFRNumber mpfr_result = unary_operation_two_outputs(op, input, mpfrIntResult,
527                                                        precision, rounding);
528   double ulp = mpfr_result.ulp(libc_result.f);
529 
530   if (mpfrIntResult != libc_result.i)
531     return false;
532 
533   return (ulp <= ulp_tolerance);
534 }
535 
536 template bool compare_unary_operation_two_outputs<float>(
537     Operation, float, const BinaryOutput<float> &, double, RoundingMode);
538 template bool compare_unary_operation_two_outputs<double>(
539     Operation, double, const BinaryOutput<double> &, double, RoundingMode);
540 template bool compare_unary_operation_two_outputs<long double>(
541     Operation, long double, const BinaryOutput<long double> &, double,
542     RoundingMode);
543 
544 template <typename T>
545 bool compare_binary_operation_two_outputs(Operation op,
546                                           const BinaryInput<T> &input,
547                                           const BinaryOutput<T> &libc_result,
548                                           double ulp_tolerance,
549                                           RoundingMode rounding) {
550   int mpfrIntResult;
551   unsigned int precision = get_precision<T>(ulp_tolerance);
552   MPFRNumber mpfr_result = binary_operation_two_outputs(
553       op, input.x, input.y, mpfrIntResult, precision, rounding);
554   double ulp = mpfr_result.ulp(libc_result.f);
555 
556   if (mpfrIntResult != libc_result.i) {
557     if (op == Operation::RemQuo) {
558       if ((0x7 & mpfrIntResult) != (0x7 & libc_result.i))
559         return false;
560     } else {
561       return false;
562     }
563   }
564 
565   return (ulp <= ulp_tolerance);
566 }
567 
568 template bool compare_binary_operation_two_outputs<float>(
569     Operation, const BinaryInput<float> &, const BinaryOutput<float> &, double,
570     RoundingMode);
571 template bool compare_binary_operation_two_outputs<double>(
572     Operation, const BinaryInput<double> &, const BinaryOutput<double> &,
573     double, RoundingMode);
574 template bool compare_binary_operation_two_outputs<long double>(
575     Operation, const BinaryInput<long double> &,
576     const BinaryOutput<long double> &, double, RoundingMode);
577 
578 template <typename InputType, typename OutputType>
579 bool compare_binary_operation_one_output(Operation op,
580                                          const BinaryInput<InputType> &input,
581                                          OutputType libc_result,
582                                          double ulp_tolerance,
583                                          RoundingMode rounding) {
584   unsigned int precision = get_precision<InputType>(ulp_tolerance);
585   MPFRNumber mpfr_result =
586       binary_operation_one_output(op, input.x, input.y, precision, rounding);
587   double ulp = mpfr_result.ulp(libc_result);
588 
589   return (ulp <= ulp_tolerance);
590 }
591 
592 template bool compare_binary_operation_one_output(Operation,
593                                                   const BinaryInput<float> &,
594                                                   float, double, RoundingMode);
595 template bool compare_binary_operation_one_output(Operation,
596                                                   const BinaryInput<double> &,
597                                                   double, double, RoundingMode);
598 template bool
599 compare_binary_operation_one_output(Operation, const BinaryInput<long double> &,
600                                     float, double, RoundingMode);
601 template bool
602 compare_binary_operation_one_output(Operation, const BinaryInput<long double> &,
603                                     double, double, RoundingMode);
604 template bool
605 compare_binary_operation_one_output(Operation, const BinaryInput<long double> &,
606                                     long double, double, RoundingMode);
607 
608 template bool compare_binary_operation_one_output(Operation,
609                                                   const BinaryInput<double> &,
610                                                   float, double, RoundingMode);
611 #ifdef LIBC_TYPES_HAS_FLOAT16
612 template bool compare_binary_operation_one_output(Operation,
613                                                   const BinaryInput<float16> &,
614                                                   float16, double,
615                                                   RoundingMode);
616 template bool compare_binary_operation_one_output(Operation,
617                                                   const BinaryInput<float> &,
618                                                   float16, double,
619                                                   RoundingMode);
620 template bool compare_binary_operation_one_output(Operation,
621                                                   const BinaryInput<double> &,
622                                                   float16, double,
623                                                   RoundingMode);
624 template bool
625 compare_binary_operation_one_output(Operation, const BinaryInput<long double> &,
626                                     float16, double, RoundingMode);
627 #endif
628 
629 template <typename InputType, typename OutputType>
630 bool compare_ternary_operation_one_output(Operation op,
631                                           const TernaryInput<InputType> &input,
632                                           OutputType libc_result,
633                                           double ulp_tolerance,
634                                           RoundingMode rounding) {
635   unsigned int precision = get_precision<InputType>(ulp_tolerance);
636   MPFRNumber mpfr_result = ternary_operation_one_output(
637       op, input.x, input.y, input.z, precision, rounding);
638   double ulp = mpfr_result.ulp(libc_result);
639 
640   return (ulp <= ulp_tolerance);
641 }
642 
643 template bool compare_ternary_operation_one_output(Operation,
644                                                    const TernaryInput<float> &,
645                                                    float, double, RoundingMode);
646 template bool compare_ternary_operation_one_output(Operation,
647                                                    const TernaryInput<double> &,
648                                                    float, double, RoundingMode);
649 template bool compare_ternary_operation_one_output(Operation,
650                                                    const TernaryInput<double> &,
651                                                    double, double,
652                                                    RoundingMode);
653 template bool compare_ternary_operation_one_output(
654     Operation, const TernaryInput<long double> &, float, double, RoundingMode);
655 template bool compare_ternary_operation_one_output(
656     Operation, const TernaryInput<long double> &, double, double, RoundingMode);
657 template bool
658 compare_ternary_operation_one_output(Operation,
659                                      const TernaryInput<long double> &,
660                                      long double, double, RoundingMode);
661 
662 #ifdef LIBC_TYPES_HAS_FLOAT16
663 template bool compare_ternary_operation_one_output(Operation,
664                                                    const TernaryInput<float> &,
665                                                    float16, double,
666                                                    RoundingMode);
667 template bool compare_ternary_operation_one_output(Operation,
668                                                    const TernaryInput<double> &,
669                                                    float16, double,
670                                                    RoundingMode);
671 template bool
672 compare_ternary_operation_one_output(Operation,
673                                      const TernaryInput<long double> &, float16,
674                                      double, RoundingMode);
675 #endif
676 
677 } // namespace internal
678 
679 template <typename T> bool round_to_long(T x, long &result) {
680   MPFRNumber mpfr(x);
681   return mpfr.round_to_long(result);
682 }
683 
684 template bool round_to_long<float>(float, long &);
685 template bool round_to_long<double>(double, long &);
686 template bool round_to_long<long double>(long double, long &);
687 
688 #ifdef LIBC_TYPES_HAS_FLOAT16
689 template bool round_to_long<float16>(float16, long &);
690 #endif // LIBC_TYPES_HAS_FLOAT16
691 
692 #ifdef LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
693 template bool round_to_long<float128>(float128, long &);
694 #endif // LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
695 
696 template <typename T> bool round_to_long(T x, RoundingMode mode, long &result) {
697   MPFRNumber mpfr(x);
698   return mpfr.round_to_long(get_mpfr_rounding_mode(mode), result);
699 }
700 
701 template bool round_to_long<float>(float, RoundingMode, long &);
702 template bool round_to_long<double>(double, RoundingMode, long &);
703 template bool round_to_long<long double>(long double, RoundingMode, long &);
704 
705 #ifdef LIBC_TYPES_HAS_FLOAT16
706 template bool round_to_long<float16>(float16, RoundingMode, long &);
707 #endif // LIBC_TYPES_HAS_FLOAT16
708 
709 #ifdef LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
710 template bool round_to_long<float128>(float128, RoundingMode, long &);
711 #endif // LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
712 
713 template <typename T> T round(T x, RoundingMode mode) {
714   MPFRNumber mpfr(x);
715   MPFRNumber result = mpfr.rint(get_mpfr_rounding_mode(mode));
716   return result.as<T>();
717 }
718 
719 template float round<float>(float, RoundingMode);
720 template double round<double>(double, RoundingMode);
721 template long double round<long double>(long double, RoundingMode);
722 
723 #ifdef LIBC_TYPES_HAS_FLOAT16
724 template float16 round<float16>(float16, RoundingMode);
725 #endif // LIBC_TYPES_HAS_FLOAT16
726 
727 #ifdef LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
728 template float128 round<float128>(float128, RoundingMode);
729 #endif // LIBC_TYPES_FLOAT128_IS_NOT_LONG_DOUBLE
730 
731 } // namespace mpfr
732 } // namespace testing
733 } // namespace LIBC_NAMESPACE_DECL
734