Opened 7 years ago

Closed 7 years ago

#11641 closed Bugs (fixed)

boost::multiprecision::fmod does not behave in the same way as std::fmod with neg divisor

Reported by: Rob Stacey <rpstac@…> Owned by: John Maddock
Milestone: To Be Determined Component: multiprecision
Version: Boost 1.57.0 Severity: Problem
Keywords: Cc:

Description

I have discovered that the replacement fmod for use with the high precision types behaves differently to the std::fmod function when using a negative divisor.

My example is boost::multiprecision::fmod(5.5, -3.333...), this gives a remainder of -1.1666... whereas std::fmod gives 2.1666... I'm basing this expectation on the following references.

The description of std::fmod (http://en.cppreference.com/w/cpp/numeric/math/fmod):

"Computes the floating-point remainder of the division operation x/y" "The returned value has the same sign as x and is less or equal to y in magnitude."

And the description of the boost::multiprecision::fmod (http://www.boost.org/doc/libs/1_59_0/libs/multiprecision/doc/html/boost_multiprecision/ref/backendconc.html#boost_multiprecision.ref.backendconc.optional_requirements_on_the_backend_type):

"Performs the equivalent operation to std::fmod on arguments cb and cb2, and store the result in b. Only required when B is an floating-point type. The default version of this function is synthesised from other operations above."

The following code illustrates the problem (with the example case that exposed it for me):

#include <cmath>
#include <iostream>
#include <iomanip>
#include <limits>
#include <boost/multiprecision/cpp_dec_float.hpp>

int main ( int argc , char *argv[] )
{
  // This is the double precision standard impl.                                                                                                                                  
  double d1 = 11.0 / 2.0;
  double d2 = -10.0 / 3.0;
  double d3 = std::fmod(d1, -d2);

  // this is the boost multiprecision impl.                                                                                                                                       
  boost::multiprecision::cpp_dec_float_100 b_11("11.0");
  boost::multiprecision::cpp_dec_float_100 b_12("2.0");
  boost::multiprecision::cpp_dec_float_100 b_21("-10.0");
  boost::multiprecision::cpp_dec_float_100 b_22("3.0");
  boost::multiprecision::cpp_dec_float_100 b1(b_11 / b_12);
  boost::multiprecision::cpp_dec_float_100 b2(b_21 / b_22);
  boost::multiprecision::cpp_dec_float_100 b3 = boost::multiprecision::fmod(b1, -b2);

  std::cout << std::endl << "Boost (Working) =>" << std::endl;
  std::cout << std::setprecision(std::numeric_limits<boost::multiprecision::cpp_dec_float_100>::max_digits10);
  std::cout << " Dividend: " << b1 << std::endl;
  std::cout << " Divisor : " << -b2 << std::endl;
  std::cout << " Result  : " << b3 << std::endl << std::endl;

  std::cout << "Std =>" << std::endl;
  std::cout << std::setprecision(std::numeric_limits<double>::max_digits10);
  std::cout << " Dividend: " << d1 << std::endl;
  std::cout << " Divisor : " << -d2 << std::endl;
  std::cout << " Result  : " << d3 << std::endl << std::endl;

  b3 = boost::multiprecision::fmod(-b1, -b2);
  d3 = std::fmod(-d1, -d2);

  std::cout << "Boost (Working) =>" << std::endl;
  std::cout << std::setprecision(std::numeric_limits<boost::multiprecision::cpp_dec_float_100>::max_digits10);
  std::cout << " Dividend: " << -b1 << std::endl;
  std::cout << " Divisor : " << -b2 << std::endl;
  std::cout << " Result  : " << b3 << std::endl << std::endl;

  std::cout << "Std =>" << std::endl;
  std::cout << std::setprecision(std::numeric_limits<double>::max_digits10);
  std::cout << " Dividend: " << -d1 << std::endl;
  std::cout << " Divisor : " << -d2 << std::endl;
  std::cout << " Result  : " << d3 << std::endl << std::endl;
b3 = boost::multiprecision::fmod(-b1, b2);
  d3 = std::fmod(-d1, d2);

  std::cout << "Boost (Broken) =>" << std::endl;
  std::cout << std::setprecision(std::numeric_limits<boost::multiprecision::cpp_dec_float_100>::max_digits10);
  std::cout << " Dividend: " << -b1 << std::endl;
  std::cout << " Divisor : " << b2 << std::endl;
  std::cout << " Result  : " << b3 << std::endl << std::endl;

  std::cout << "Std =>" << std::endl;
  std::cout << std::setprecision(std::numeric_limits<double>::max_digits10);
  std::cout << " Dividend: " << -d1 << std::endl;
  std::cout << " Divisor : " << d2 << std::endl;
  std::cout << " Result  : " << d3 << std::endl << std::endl;

  b3 = boost::multiprecision::fmod(b1, b2);
  d3 = std::fmod(d1, d2);

  std::cout << "Boost (Broken) =>" << std::endl;
  std::cout << std::setprecision(std::numeric_limits<boost::multiprecision::cpp_dec_float_100>::max_digits10);
  std::cout << " Dividend: " << b1 << std::endl;
  std::cout << " Divisor : " << b2 << std::endl;
  std::cout << " Result  : " << b3 << std::endl << std::endl;

  std::cout << "Std =>" << std::endl;
  std::cout << std::setprecision(std::numeric_limits<double>::max_digits10);
  std::cout << " Dividend: " << d1 << std::endl;
  std::cout << " Divisor : " << d2 << std::endl;
  std::cout << " Result  : " << d3 << std::endl << std::endl;
}

This outputs (on my system):

-bash-4.1$ g++ -std=c++11 calc.cpp -o calc -bash-4.1$ ./calc

Boost (Working) =>

Dividend: 5.5 Divisor : 3.333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333334 Result : 2.166666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666

Std =>

Dividend: 5.5 Divisor : 3.3333333333333335 Result : 2.1666666666666665

Boost (Working) =>

Dividend: -5.5 Divisor : 3.333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333334 Result : -2.166666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666

Std =>

Dividend: -5.5 Divisor : 3.3333333333333335 Result : -2.1666666666666665

Boost (Broken) =>

Dividend: -5.5 Divisor : -3.333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333334 Result : 1.166666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666668

Std =>

Dividend: -5.5 Divisor : -3.3333333333333335 Result : -2.1666666666666665

Boost (Broken) =>

Dividend: 5.5 Divisor : -3.333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333333334 Result : -1.166666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666666668

Std =>

Dividend: 5.5 Divisor : -3.3333333333333335 Result : 2.1666666666666665

I believe that the root cause is line 933 of include/boost/multiprecision/detail/default_ops.hpp (https://github.com/boostorg/multiprecision/blob/master/include/boost/multiprecision/detail/default_ops.hpp)

Changing this from:

   if(eval_get_sign(a) < 0)

to

   if(eval_get_sign(result) < 0)

fixes the issue for this case, I'm currently trying to prepare a merge request but I didn't want to delay someone confirming that this is a issue that needs fixing or whether I have misunderstood something.

Change History (2)

comment:1 by John Maddock, 7 years ago

Thanks, confirmed.

Testing patch and updated tests now, John.

Note: See TracTickets for help on using tickets.