Re: A problem with range objects and floating point numbers

[Daniel J Sebald]

On 09/29/2014 11:36 AM, s.jo wrote:
> Well, I don't think that way on the printf format:
>
> octave:18>  sprintf('%.43g',0.1)
> ans = 0.1000000000000000055511151231257827021181583
>
> Do you think that the above ans is rounded around 12th position?
> Conversely, the decimation at 12th position is the way how Matlab shows the
> number in display
> only when the explicit precision number is not specified such as '%12f' or
> '%.43g'.
> Even though Matlab shows that way on display, any raw bits are NOT
> discarded.
> You can check it by using num2hex command in matlab.

I understand what you are looking at.


> Anyhow, it is fun to check the detail behaviors:
> Here are the test result up-to test 8:
>
> % generate ranges :
> %  1. range object, 2. matrix object, 3. linspace, 4. perfect decimation
> N=1000;
> range_test1=(-N:0.1:N);
> range_test2=[-N:0.1:N];
> range_test3=linspace(-N,N,20*N+1);
> range_test4=(str2num(sprintf('%.12f ',-N:0.1:N)));
> range_test5=(-N/0.1:1:N/0.1)*0.1;
> range_test6=[-N/0.1:1:N/0.1]*0.1;
> % same as range_test6
> NN=20*N+1;
> range_test7=zeros(1,NN);
> for i=(1:1:NN)
>     range_test7(i)=(i-1-N/0.1)*0.1;
> endfor
> % the linear spacing (interpolation) algorithm
> range_test8=zeros(1,NN); A=-N; B=N;
> for i=(1:1:NN)
>     range_test8(i)=(A*(NN-i) + B*(i-1))/(NN-1);
> endfor
>
> + sum (abs (sin (range_test1 * pi))<  eps (range_test1 * pi))
> ans =       1055
> + sum (abs (sin (range_test2 * pi))<  eps (range_test2 * pi))
> ans =       1168
> + sum (abs (sin (range_test3 * pi))<  eps (range_test3 * pi))
> ans =       2001
> + sum (abs (sin (range_test4 * pi))<  eps (range_test4 * pi))
> ans =       2001
> + sum (abs (sin (range_test5 * pi))<  eps (range_test5 * pi))
> ans =       1055
> + sum (abs (sin (range_test6 * pi))<  eps (range_test6 * pi))
> ans =       2001
> + sum (abs (sin (range_test7 * pi))<  eps (range_test7 * pi))
> ans =       2001
> + sum (abs (sin (range_test8 * pi))<  eps (range_test8 * pi))
> ans =       2001
>
> The results are as you expected.

Great!  So this may be the way to go with linspace():

result(i) = (A*(N-i) + B*i)/N;

And if so, perhaps the range operator can be changed to utilize the 
linear-spacing algorithm.  That is

Step 1) First compute what the final B value is for the range, i.e.,

  N = floor((limit-base)/increment);
  B = base + N * increment;

Step 2) Apply the linear spacing routine, i.e., the equivalent of 
linspace(base, B, N).

Now that first step for computing B could be fraught with it's own 
numerical peculiarity, so would need some thought in that as well.

Dan