Check if a double is evenly divisible by another double in C?

How can I check if a double x is evenly divisible by another double y in C? With integers I wou开发者_如何转开发ld just use modulo, but what would be the correct/best way to do it with doubles?

I know floating point numbers carry with them imprecision, but I'm getting the double from standard input. Maybe I should not scan it as a double straight away but as two integers instead, but where would I go from then?

The standard header math.h defines the following functions:

• double fmod(double x, double y);
• float fmodf(float x, float y);
• long double fmodl(long double x, long double y);

These functions return the result of the remainder of x divided by y. The result has the same sign as that of x. You can use r = fmod(x, y); for double numbers x and y, and check if r == 0. If you want to not test for exact divisibility but add some tolerance, then you can check if r is "close enough" to 0 or y (thanks caf).

fmodf() and fmodl() are new in C99.

Edit: C99 also defines a separate remainder(double x, double y) function, that returns the remainder of x/y. From http://docs.sun.com/source/806-3568/ncg_lib.html:

The remainder(x,y) is the operation specified in IEEE Standard 754-1985. The difference between remainder(x,y) and fmod(x,y) is that the sign of the result returned by remainder(x,y) might not agree with the sign of either x or y, whereas fmod(x,y) always returns a result whose sign agrees with x. Both functions return exact results and do not generate inexact exceptions.

...

When y ≠ 0, the remainder r = x REM y is defined regardless of the rounding mode by the mathematical relation r = x - ny, where n is the integer nearest the exact value of x/y; whenever | n - x/y | = 1/2, then n is even. Thus, the remainder is always exact. If r = 0, its sign shall be that of x. This definition is applicable for all implementations.

(Either fmod() or remainder() should work for you.)

The fmod() family of functions give terrible results. Suppose you want to determine if 42 is evenly divisible by 0.4. It is, 105 times. However, fmod does the division and gets a result like 104.99999 which it then rounds down to 104 resulting in a remainder of 0.399999 which gives a false negative result. remainderl(), however, seems to work. Even 0.4 itself is represented inexactly in floating point.

For the folks who don't grok the concept of "evenly divisible", it has nothing to do with the result being an even number - you probably have your etymology backwards. Even numbers are those numbers which are evenly divisible by 2. And the concept of divisibility is entirely valid for non-integers. Evenly divisible means the result of the division is an integer regardless of whether the dividend or divisor are. An example application is if you have a metal lathe with a 3mm pitch leadscrew and are cutting a 0.4mm pitch bolt. 14 threads at 3mm line up with 105 threads at 0.4mm. The divisibility calculation is used to tell where the various moving parts of the lathe sync up again so you can reengage for the next cutting pass. Another example is imperial measurements which have been converted to metric. 50.8mm (2") is evenly divisible by 25.4mm (1"). Even without metric conversions, dimensions are often non-integers yet divisibility is often an issue: 0.5" is evenly divisible by 0.1", 0.125". and 0.250". Converting a floating point number (such as 0.375") to a fractional representation (3/8") is one more application of divisibility to non-integer numbers.

The two alternative calculations in this sample function give the same results for hundreds of different number pairs. However, replacing remainderl() with fmodl() or roundl() with floorl() gives lots of invalid results. I originally used a fuzz of 0.001. Actual calculation error seems to usually be of order 1E-15 so a smaller fuzz can be used. However, comparing the result to 0.0 will give false negative results. You might want to express your fuzz in terms of your denominator in case you are working with very small numbers. divisible(42, 0.4) and divisible(41,0.4) should give the same results as divisible(0.000000042, 0.0000000004) and divisible(0.000000041, 0.0000000004). I.E. are 42nm and 41nm divisible by 0.4nm? With the version of the function given here, they do. With a fixed fuzz, they do not necessarily. However, divisible(42, 0.0000000004) still gives a false negative (error is 1.53003e-15 which is larger than the fuzz of 4E-19) so comparing numbers that differ by 9 orders of magnitude is not reliable. IEEE floating point has its limitations. Notice I used long double calculations to minimize calculation and representation errors. This function was not tested with negative numbers.

int divisible(long double a, long double b) 
{
  int result;
#if 1
   if(fabsl(((roundl(a/b)*b)- a)) <= (1E-9*b) ) {
    result=TRUE;
  } else {
    result=FALSE;
  }
#else
  if( fabsl(remainderl(a,b)) <= (1E-9*b ) ){
    result=TRUE;
  } else {
    result=FALSE;
  }
#endif
  // printf("divisible(%Lg, %Lg): %Lg, %Lg,%d\n", a, b, roundl(a/b), fabsl(((roundl(a/b)*b)-a)), result);
  return(result);
}

If you want to be absolutely precise, you could used fixed-point math. That is, do everything with ints, but ints that are (in your case) some power of 10 of the value they actually represent.

Say the user enters 123.45 and 6789.1. First, you want to make sure you've got the same number of decimal places, so add trailing zeros to the one with fewer decimal places. That gives us 123.45 and 6789.10 (now both with 2 decimal places). Now just remove the decimal point, to get 12345 and 678910. If one divides into the other evenly, that's your answer.

This works because removing the decimal point multiplies both by the same constant (100 in the example above). (x * 100) / (y * 100) == x / y

A few things to be careful about: if you read the integer part and fractional part as ints, be careful that you don't lose leading zeros on the fractional part. (eg: 0.1 and 0.0001 are not the same number!) Also, if there are enough decimal places you can overflow. You probably want to at least use a long.

You could also do the computation with doubles, but it'll be less precise. To do it that way, do the division, and then compare the difference between the result and the rounded result. If within some small tolerance, then it divides evenly.

I am not sure what you are trying to do, but I have used fmod() from math.h in audio synthesis code where I needed my parameters to be floats or doubles and I needed to get a modulo.

1. Scan them in as doubles and call them x1 and x2
2. Find what x1/x2 is using division and call it x3
3. Find x1 - (x2*x3) and see if that number is sufficiently close to zero - if it is then x1 is evenly divisible by x2 - (obviously taking into consideration the the possibility of negative values here)

lol - line 3 fixed :)

How can I check if a double x is evenly dividable by another double y in C? With integers I would just use modulo, but what would be the correct/best way to do it with doubles?

You would include and link to the math library:

#include <math.h>

Then you would call the floating point modulus function fmod:

if (fmod(5.0, 2.5) == 0.0)
  // evenly divisible
else
  // not evenly divisible

You may want to compare the result of fmod with a small value instead of 0.0 depending on your needs.

The concept of "even number" is only defined for integers. You can't apply it to doubles; it does not make mathematical sense. From Wikipedia:

An even number is an integer that is "evenly divisible" by 2, i.e., divisible by 2 without remainder.

I suggest you convert your doubles to ints, by applying whatever method you decide (rounding, truncation) and then use modulo as you suggest.