Path: news.mathworks.com!newsfeed-00.mathworks.com!newsfeed2.dallas1.level3.net!news.level3.com!postnews.google.com!v5g2000prm.googlegroups.com!not-for-mail From: "Felipe G. Nievinski" <fgnievinski@gmail.com> Newsgroups: comp.soft-sys.matlab,sci.math.num-analysis Subject: Re: robust truncation at n decimal places Date: Sun, 7 Dec 2008 01:26:53 -0800 (PST) Organization: http://groups.google.com Lines: 109 Message-ID: <fac1f966-62ff-427d-b67c-7cce5297f210@v5g2000prm.googlegroups.com> References: <545015a5-3f51-41e3-a40b-1cf7db531b11@35g2000pry.googlegroups.com> <ghfobu$blk$1@fred.mathworks.com> NNTP-Posting-Host: 128.138.43.113 Mime-Version: 1.0 Content-Type: text/plain; charset=ISO-8859-1 Content-Transfer-Encoding: quoted-printable X-Trace: posting.google.com 1228642013 13697 127.0.0.1 (7 Dec 2008 09:26:53 GMT) X-Complaints-To: groups-abuse@google.com NNTP-Posting-Date: Sun, 7 Dec 2008 09:26:53 +0000 (UTC) Complaints-To: groups-abuse@google.com Injection-Info: v5g2000prm.googlegroups.com; posting-host=128.138.43.113; posting-account=iCBMwwoAAAARtxHOK9kIhRHsTae4FBZX User-Agent: G2/1.0 X-HTTP-UserAgent: Mozilla/5.0 (Windows; U; Windows NT 6.0; en-US) AppleWebKit/525.19 (KHTML, like Gecko) Chrome/0.4.154.29 Safari/525.19,gzip(gfe),gzip(gfe) Xref: news.mathworks.com comp.soft-sys.matlab:505443 sci.math.num-analysis:105694 On Dec 6, 10:53 pm, "Roger Stafford" <ellieandrogerxy...@mindspring.com.invalid> wrote: > "Felipe G. Nievinski" <fgnievin...@gmail.com> wrote in message <545015a5-= 3f51-41e3-a40b-1cf7db531...@35g2000pry.googlegroups.com>... > > > Hi. Given > > x =3D 0.99899998 > > y =3D 0.99900000 > > za =3D sscanf(sprintf('%.7f\n',x), '%f') > > zb =3D fix(x./10^-7).*10^-7 > > za =3D=3D y > > zb =3D=3D y > > then za =3D=3D y is true, but zb =3D=3D y is false. > > As you see, it's possible to perform a robust truncation at n decimal > > places through conversion to/from character string. I was trying to > > achive the same using computationally cheaper double precision > > arithmetic, without success. Is that it or am I missing something? > > > Thanks, > > Felipe. > > Your 'fix' is doing truncation - that is to say, rounding toward zero -= whereas sprintf('%.7f\n',x) is doing a round to nearest. It's not surpris= ing you get different answers. Also for better accuracy you should be mult= iplying by 10^7 instead of dividing by 10^-7 and then dividing by 10^7 inst= ead of multiplying by 10^-7, because the number 10^-7 cannot be represented= precisely with binary floating point numbers. > Thanks for indicating the round-like behavior of sprintf -- I was overlooking that. Now I get the expected result for the particular x value I gave: x =3D 0.99899998 y =3D 0.99900000 za =3D sscanf(sprintf('%.7f\n',x), '%f') zb =3D round(x./10^-7).*10^-7 za =3D=3D y zb =3D=3D y Then I ran the code for another test cases. x =3D 257104079.99899998 y =3D 257104079.9990000 za =3D sscanf(sprintf('%.7f\n',x), '%f') zb =3D round(x.*10^7)./10^7 za =3D=3D y zb =3D=3D y The test above would only pass after I've replaced x*10^-7 by x/(10^7) -- it's interesting that the two operations are not equivalent. But for the test below the algorithm fails: x =3D 257098540.00000003 y =3D 257098540.0000000 za =3D sscanf(sprintf('%.7f\n',x), '%f') zb =3D round(x.*10^7)./10^7 za =3D=3D y zb =3D=3D y :-( > Nevertheless even with these corrections, you probably cannot count on = always getting precisely the same answers clear out to the 53rd binary plac= e unless you know exactly how sprintf's and sscanf's algorithms work. Actually, sprintf can be dismissed in the comparisons above, because I know the answer I expect, e.g., in =3D 257098540.00000003 out =3D 257098540.0000000 I've only included sprintf to demonstrate that what I need is not impossible, even in a finite precision computer. You're right that if sprintf can do it, there must be a way to do it -- I hope without going to/from character strings. Maybe dissecting the floating-point number into exponent and mantissa? I tried this, using [f,e]=3Dlog2(x), where x =3D f * 2^e,: x =3D 257098540.00000003 y =3D 257098540.0000000 za =3D sscanf(sprintf('%.7f\n',x), '%f') zb =3D round(x.*10^7)./10^7 [f,e] =3D log2(x); f =3D f*10^7; temp =3D pow2(f, e); temp =3D round(temp); [f,e] =3D log2(temp); f =3D f/10^7; zc =3D pow2(f, e) za =3D=3D y zb =3D=3D y zc =3D=3D y But it failed the test, too. For the time being I have to keep using the costly sprintf/scanf solution, simply because it's the only one that I know of passing all test cases. > You should be making comparison like abs(za-zb) < tol where tol is someth= ing like 10^-15*abs(x) > Usually I compare two floating point numbers as you prescribe, but for this particular application I need exact equality: I input two vectors containing time tags into intersect() to find the common epochs, and intersect internally relies on ismember(), which checks for exact equality between vector elements. I could try to generalize ismember() -- in the past, I've generalized unique(x, tol) -- but since my input data is given at finite precision, I find it reasonable to expect exact equality, after rounding/truncation. (My time tags are given in the format year month day hour minute seconds -- all fields are integers except seconds, which has seven decimal places. I end up with decimals beyond the seventh place because I convert them through datenum()). Thanks a lot for your ideas. Felipe.