One of our TAs reported this to me earlier today. (A student reported it to him.)

Sum(sin(1/10000.0*i^2), i = 0 .. 199):
%, value( % ), evalf( % );
            199                                                   
           -----                                                  
            \                                                     
             )      /                 2\                          
            /    sin\0.0001000000000 i /, 80.85387129, 127.2891137
           -----                                                  
           i = 0     

The student used evalf. I suggested that value was a better first option in this case. Since this is just the sum of a fixed set of numbers, an even better alternative is to use add:

add(sin(1/10000.0*i^2), i = 0 .. 199);

                                 80.85387132

The problem disappears when Digits increases to (at least) 15.

Digits := 15:
Sum(sin(1/10000.0*i^2), i = 0 .. 199):
%, value( % ), evalf( % );
    199                                                                  
   -----                                                                 
    \                                                                    
     )      /                      2\                                    
    /    sin\0.000100000000000000 i /, 80.8538712866271, 80.8538712866267
   -----                                                                 
   i = 0                                                                 

add(sin(1/10000.0*i^2), i = 0 .. 199);
                              80.8538712866273 

If you change the 10000.0 to 10000, the same problems are seen in the output.

One the one hand, this has the signs of a floating-point cancellation problem. But, a quick glance at the individual terms shows that there are no really small terms or other indicaton of such a problem.

a := [seq(sin(1/10000.0*i^2), i = 0 .. 199)]:
convert(a,`+`);
                                 80.85387131
for i from 1 to nops(a) do
 convert( a[1..i], `+` ), convert( a[-i..-1], `+`)
end do;

I have the feeling that Maple is able to find a closed-form expression for the sum and its the floating-point evaluation of this expression that is problematic.

Any thoughts?


Please Wait...