Dr. StrangeNumbers or: How I Learned to Stop Worrying and Love Floating Point Arithmetic

Question

The following is the program.

test[t_, dt_] := 
 Module[{}, For[ti = dt, ti <= t, ti = ti + dt, Print[ti];]; 
  Print[MemberQ[{0.01, 0.02}, ti]]; Return[0];]    
test[0.01, 0.001]

(* 0.001
...
False
0 *)

Obviously the result is wrong. Copy the above result, you will be surprised to see:

(* ...
0.009000000000000001`
0.010000000000000002`
... *)

Why is this so?

Answer 1

While it is a bad idea to compare floating point numbers, in this case I think something simpler is going on: you have an off-by-1 problem in your loop. See what your code does:

test[t_, dt_] := Module[{},
  For[ti = dt, ti <= t, ti = ti + dt, Print[ti];];
  Print[MemberQ[{0.01, 0.02}, ti]];
  Return[0];
  ]

then, after

test[0.01, .001];

look at ti:

ti
(*0.011*)

If you use < rather than <= in the condition, the final value will indeed be what you expect it to be:

test2[t_, dt_] := Module[{},
  For[ti = dt, ti < t, ti = ti + dt, Print[ti];];
  Print[MemberQ[{0.01, 0.02}, ti]];
  Return[0];
  ]

test2[0.01, .001]
ti

prints a list up to 0.009 and returns True.

In general, though, don't do things like that with floating-point numbers. Even if it works here, it will eventually fail. Observe:

(1 + $MaxMachineNumber) == $MaxMachineNumber
(*True*)

or, as JM points out,

1 + $MachineEpsilon == 1
(*True*)

Answer 2

This seems to work fine and is more in the spirit of the Mathematica way and perhaps illustrates why a Functional approach tends to lead to fewer bugs than the procedural approach:

test[t_, dt_] := {#, MemberQ[{0.01, 0.02}, #]} & /@ Range[dt, t, dt];

Which when run,

test[0.01, 0.001] // TableForm

gives: