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);
}