//===-- Utils which wrap MPC ----------------------------------------------===// // // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions. // See https://llvm.org/LICENSE.txt for license information. // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception // //===----------------------------------------------------------------------===// #include "MPCUtils.h" #include "src/__support/CPP/array.h" #include "src/__support/CPP/stringstream.h" #include "utils/MPFRWrapper/MPCommon.h" #include #include "mpc.h" template using FPBits = LIBC_NAMESPACE::fputil::FPBits; namespace LIBC_NAMESPACE_DECL { namespace testing { namespace mpc { static inline cpp::string str(RoundingMode mode) { switch (mode) { case RoundingMode::Upward: return "MPFR_RNDU"; case RoundingMode::Downward: return "MPFR_RNDD"; case RoundingMode::TowardZero: return "MPFR_RNDZ"; case RoundingMode::Nearest: return "MPFR_RNDN"; } } class MPCNumber { private: unsigned int precision; mpc_t value; mpc_rnd_t mpc_rounding; public: explicit MPCNumber(unsigned int p) : precision(p), mpc_rounding(MPC_RNDNN) { mpc_init2(value, precision); } MPCNumber() : precision(256), mpc_rounding(MPC_RNDNN) { mpc_init2(value, 256); } MPCNumber(unsigned int p, mpc_rnd_t rnd) : precision(p), mpc_rounding(rnd) { mpc_init2(value, precision); } template , bool> = 0> MPCNumber(XType x, unsigned int precision = mpfr::ExtraPrecision::VALUE, RoundingMode rnd = RoundingMode::Nearest) : precision(precision), mpc_rounding(MPC_RND(mpfr::get_mpfr_rounding_mode(rnd), mpfr::get_mpfr_rounding_mode(rnd))) { mpc_init2(value, precision); Complex x_c = cpp::bit_cast>(x); mpfr_t real, imag; mpfr_init2(real, precision); mpfr_init2(imag, precision); mpfr_set_flt(real, x_c.real, mpfr::get_mpfr_rounding_mode(rnd)); mpfr_set_flt(imag, x_c.imag, mpfr::get_mpfr_rounding_mode(rnd)); mpc_set_fr_fr(value, real, imag, mpc_rounding); mpfr_clear(real); mpfr_clear(imag); } template , bool> = 0> MPCNumber(XType x, unsigned int precision = mpfr::ExtraPrecision::VALUE, RoundingMode rnd = RoundingMode::Nearest) : precision(precision), mpc_rounding(MPC_RND(mpfr::get_mpfr_rounding_mode(rnd), mpfr::get_mpfr_rounding_mode(rnd))) { mpc_init2(value, precision); Complex x_c = cpp::bit_cast>(x); mpc_set_d_d(value, x_c.real, x_c.imag, mpc_rounding); } MPCNumber(const MPCNumber &other) : precision(other.precision), mpc_rounding(other.mpc_rounding) { mpc_init2(value, precision); mpc_set(value, other.value, mpc_rounding); } ~MPCNumber() { mpc_clear(value); } MPCNumber &operator=(const MPCNumber &rhs) { precision = rhs.precision; mpc_rounding = rhs.mpc_rounding; mpc_init2(value, precision); mpc_set(value, rhs.value, mpc_rounding); return *this; } void setValue(mpc_t val) const { mpc_set(val, value, mpc_rounding); } mpc_t &getValue() { return value; } MPCNumber carg() const { mpfr_t res; MPCNumber result(precision, mpc_rounding); mpfr_init2(res, precision); mpc_arg(res, value, MPC_RND_RE(mpc_rounding)); mpc_set_fr(result.value, res, mpc_rounding); mpfr_clear(res); return result; } MPCNumber cproj() const { MPCNumber result(precision, mpc_rounding); mpc_proj(result.value, value, mpc_rounding); return result; } }; namespace internal { template cpp::enable_if_t, MPCNumber> unary_operation(Operation op, InputType input, unsigned int precision, RoundingMode rounding) { MPCNumber mpcInput(input, precision, rounding); switch (op) { case Operation::Carg: return mpcInput.carg(); case Operation::Cproj: return mpcInput.cproj(); default: __builtin_unreachable(); } } template bool compare_unary_operation_single_output_same_type(Operation op, InputType input, OutputType libc_result, double ulp_tolerance, RoundingMode rounding) { unsigned int precision = mpfr::get_precision>(ulp_tolerance); MPCNumber mpc_result; mpc_result = unary_operation(op, input, precision, rounding); mpc_t mpc_result_val; mpc_init2(mpc_result_val, precision); mpc_result.setValue(mpc_result_val); mpfr_t real, imag; mpfr_init2(real, precision); mpfr_init2(imag, precision); mpc_real(real, mpc_result_val, mpfr::get_mpfr_rounding_mode(rounding)); mpc_imag(imag, mpc_result_val, mpfr::get_mpfr_rounding_mode(rounding)); mpfr::MPFRNumber mpfr_real(real, precision, rounding); mpfr::MPFRNumber mpfr_imag(imag, precision, rounding); double ulp_real = mpfr_real.ulp( (cpp::bit_cast>>(libc_result)).real); double ulp_imag = mpfr_imag.ulp( (cpp::bit_cast>>(libc_result)).imag); mpc_clear(mpc_result_val); mpfr_clear(real); mpfr_clear(imag); return (ulp_real <= ulp_tolerance) && (ulp_imag <= ulp_tolerance); } template bool compare_unary_operation_single_output_same_type( Operation, _Complex float, _Complex float, double, RoundingMode); template bool compare_unary_operation_single_output_same_type( Operation, _Complex double, _Complex double, double, RoundingMode); template bool compare_unary_operation_single_output_different_type( Operation op, InputType input, OutputType libc_result, double ulp_tolerance, RoundingMode rounding) { unsigned int precision = mpfr::get_precision>(ulp_tolerance); MPCNumber mpc_result; mpc_result = unary_operation(op, input, precision, rounding); mpc_t mpc_result_val; mpc_init2(mpc_result_val, precision); mpc_result.setValue(mpc_result_val); mpfr_t real; mpfr_init2(real, precision); mpc_real(real, mpc_result_val, mpfr::get_mpfr_rounding_mode(rounding)); mpfr::MPFRNumber mpfr_real(real, precision, rounding); double ulp_real = mpfr_real.ulp(libc_result); mpc_clear(mpc_result_val); mpfr_clear(real); return (ulp_real <= ulp_tolerance); } template bool compare_unary_operation_single_output_different_type( Operation, _Complex float, float, double, RoundingMode); template bool compare_unary_operation_single_output_different_type( Operation, _Complex double, double, double, RoundingMode); template void explain_unary_operation_single_output_different_type_error( Operation op, InputType input, OutputType libc_result, double ulp_tolerance, RoundingMode rounding) { unsigned int precision = mpfr::get_precision>(ulp_tolerance); MPCNumber mpc_result; mpc_result = unary_operation(op, input, precision, rounding); mpc_t mpc_result_val; mpc_init2(mpc_result_val, precision); mpc_result.setValue(mpc_result_val); mpfr_t real; mpfr_init2(real, precision); mpc_real(real, mpc_result_val, mpfr::get_mpfr_rounding_mode(rounding)); mpfr::MPFRNumber mpfr_result(real, precision, rounding); mpfr::MPFRNumber mpfrLibcResult(libc_result, precision, rounding); mpfr::MPFRNumber mpfrInputReal( cpp::bit_cast>>(input).real, precision, rounding); mpfr::MPFRNumber mpfrInputImag( cpp::bit_cast>>(input).imag, precision, rounding); cpp::array msg_buf; cpp::StringStream msg(msg_buf); msg << "Match value not within tolerance value of MPFR result:\n" << " Input: " << mpfrInputReal.str() << " + " << mpfrInputImag.str() << "i\n" << " Rounding mode: " << str(rounding) << '\n' << " Libc: " << mpfrLibcResult.str() << '\n' << " MPC: " << mpfr_result.str() << '\n' << '\n' << " ULP error: " << mpfr_result.ulp_as_mpfr_number(libc_result).str() << '\n'; tlog << msg.str(); mpc_clear(mpc_result_val); mpfr_clear(real); } template void explain_unary_operation_single_output_different_type_error( Operation, _Complex float, float, double, RoundingMode); template void explain_unary_operation_single_output_different_type_error( Operation, _Complex double, double, double, RoundingMode); template void explain_unary_operation_single_output_same_type_error( Operation op, InputType input, OutputType libc_result, double ulp_tolerance, RoundingMode rounding) { unsigned int precision = mpfr::get_precision>(ulp_tolerance); MPCNumber mpc_result; mpc_result = unary_operation(op, input, precision, rounding); mpc_t mpc_result_val; mpc_init2(mpc_result_val, precision); mpc_result.setValue(mpc_result_val); mpfr_t real, imag; mpfr_init2(real, precision); mpfr_init2(imag, precision); mpc_real(real, mpc_result_val, mpfr::get_mpfr_rounding_mode(rounding)); mpc_imag(imag, mpc_result_val, mpfr::get_mpfr_rounding_mode(rounding)); mpfr::MPFRNumber mpfr_real(real, precision, rounding); mpfr::MPFRNumber mpfr_imag(imag, precision, rounding); mpfr::MPFRNumber mpfrLibcResultReal( cpp::bit_cast>>(libc_result).real, precision, rounding); mpfr::MPFRNumber mpfrLibcResultImag( cpp::bit_cast>>(libc_result).imag, precision, rounding); mpfr::MPFRNumber mpfrInputReal( cpp::bit_cast>>(input).real, precision, rounding); mpfr::MPFRNumber mpfrInputImag( cpp::bit_cast>>(input).imag, precision, rounding); cpp::array msg_buf; cpp::StringStream msg(msg_buf); msg << "Match value not within tolerance value of MPFR result:\n" << " Input: " << mpfrInputReal.str() << " + " << mpfrInputImag.str() << "i\n" << " Rounding mode: " << str(rounding) << " , " << str(rounding) << '\n' << " Libc: " << mpfrLibcResultReal.str() << " + " << mpfrLibcResultImag.str() << "i\n" << " MPC: " << mpfr_real.str() << " + " << mpfr_imag.str() << "i\n" << '\n' << " ULP error: " << mpfr_real .ulp_as_mpfr_number( cpp::bit_cast>>(libc_result) .real) .str() << " , " << mpfr_imag .ulp_as_mpfr_number( cpp::bit_cast>>(libc_result) .imag) .str() << '\n'; tlog << msg.str(); mpc_clear(mpc_result_val); mpfr_clear(real); mpfr_clear(imag); } template void explain_unary_operation_single_output_same_type_error( Operation, _Complex float, _Complex float, double, RoundingMode); template void explain_unary_operation_single_output_same_type_error( Operation, _Complex double, _Complex double, double, RoundingMode); } // namespace internal } // namespace mpc } // namespace testing } // namespace LIBC_NAMESPACE_DECL