"Richard Crozier" <r.crozier@ed.ac.uk> wrote in message
news:hc3vl6$oag$1@fred.mathworks.com...
> Hello,
>
> I've noticed some undesireable behaviour in some code I've written, which
> may or may not be an error, but for which I would like some advice on how
> to work around.
>
> I have some code which partitions a region of space into 49 sections. On
> the lines running between these sections I am able to determine the exact
> value of a function I am interested in. Then to estimate the value of the
> function between these lines, I then interpolate between the values on
> these lines either side of a point I am interested in. To do this I first
> have to find the correct lines either side of the point I wish to evaluate
> which is where I'm getting an odd error.
>
> By chance one of the points I wish to evaluate has come out as a value
> very close to one of the lines. The nearest line is at the line x=0.4, and
> my actual point has the x position 0.4  5.5511e017. This is just under
> 0.4, which I discovered with the following code:
>
> xyCoords(19,2)  0.4
>
> ans =
>
> 5.5511e017
This is approximately eps(0.4)  about the "spacing" between adjacent
floating point numbers in the vicinity of 0.4.
Read this document for more information about floating point numbers and
EPS:
http://www.mathworks.com/company/newsletters/news_notes/pdf/Fall96Cleve.pdf
> The region is split into 50 lines, and the correct lines either side of
> this point are the lines at x=0.38 and x = 0.4. However, if I do the
> following to find the nearest line below this point:
>
>>> floor(rzCoords(19,2) * 50)/50
> ans =
>
> 0.4000
>
> I get the answer above which is not the line below the point as desired.
> For comparison, if I do:
However, it is the correct answer.
Look at:
format hex
rzCoords(19, 2)
0.4
rzCoords(19, 2)*50
20
The difference between rzCoords(19, 2) and 0.4 is so small that when you
multiply the expression by 50, the difference between 50*rzCoords(19, 2) and
20 is too insignificant (in the floatingpoint arithmetic sense) to be
represented in the result of the calculation.
>>> floor(0.399999999*50)/50
>
> ans =
>
> 0.3800
>
> I get the correct line, the line at 0.38. To check what's going on I can
> also do the following to confirm that the number is definately less than
> 0.4:
>
>>> xyCoords(19,2) == 0.4
>
> ans =
>
> 0
>
>>> xyCoords(19,2) < 0.4
>
> ans =
>
> 1
Suppose I give you a bank account containing $1,000. Not bad, huh? Now I
add $0.001 to your bank account. You don't really care about that, do you?
Now what if I multiplied the contents of your bank account by 1000? Would
you say that you had $1,000,000 or that you had $1,000,001?
> This is a major problem, because if I then use the interp1 function as I
> would like to, as below for instance:
>
> interp1([0.4, 0.42], [0.12154, 0.154894], abs(xyCoords(19,2)))
>
> I get the result NaN as it lies outside the range of the data.
>
> The problem appears to be that xyCoords(19,2) * 50 yields exactly 20, as
> shown below:
>
>>> rzCoords(19,2) * 50
>
> ans =
>
> 20
>
>>> ans == 20
>
> ans =
>
> 1
>
> not a value just under 20 as it should be. What is the best way to fix
> this?
Make your calculations robust to floatingpoint arithmetic.

Steve Lord
slord@mathworks.com
comp.softsys.matlab (CSSM) FAQ: http://matlabwiki.mathworks.com/MATLAB_FAQ
