|
From: | Marco Atzeri |
Subject: | Re: odd behaviour with 2x2 matrix floating-point multiplication |
Date: | Tue, 16 Dec 2014 06:12:47 +0100 |
User-agent: | Mozilla/5.0 (Windows NT 6.1; WOW64; rv:31.0) Gecko/20100101 Thunderbird/31.3.0 |
On 12/16/2014 4:47 AM, s.jo wrote:
After Jordi's comment, I updated my cygwin system from 1.7.32 to 1.7.33-2 (ver 2014-11-13). Now I have a perfect computation on zero entry, but a*a+b*b term still shows mismatch with a probability of 0.24 or so: If I try 10000 pairs of a and b, about 2300 cases shows mismatch at the a*a+b*b. I checked that the mismatch occurs at only last bit of hex form, say +/-eps(a*a+b*b), as Dan pointed out. I understand that N pairs term floating-point operation of matrix could be different from the exact math expression. But, in this case, only two terms are involved in machine code that BLAS may generate. I am not sure of how BLAS routine can shuffle the operation order or do anything to show different floating-point results. I suspect that BLAS routines scale the numbers unnecessarily, even trivial cases. Is this a problem related to a release version of BLAS and/or cygwin?
blas dependent surely. cygwin specific I doubt.
-- jo ps. The follow script is to check the probability of a*a+b*b mismatch: %% num_trial=10000; count_err=[0;0]; for ii=1:num_trial a = randn; b = randn; C = [a b; -b a]*[a; b]; D = [a*a+b*b;-b*a+a*b]; err=(C ~= D); count_err = count_err + err; if count_err(1)==100 num2hex(C(1)) num2hex(D(1)) num2hex(D(1)-eps(D(1))) num2hex(D(1)+eps(D(1))) break end end ii count_err ii/count_err
I assume is count_err/ii
ans = 0.240963855421687 0
with blas reference ans = 0.16807 0.99832 with openblas ans = 0.24450 0.00000 so while openblas has ~24% mismatch on C(1)-D(1) has 0% mismatch on C(2)-D(2) Netlib blas has ~16% mismatch on C(1)-D(1) has 100% mismatch on C(2)-D(2) which one do you prefer ? Suggested reading: http://docs.oracle.com/cd/E19957-01/806-3568/ncg_goldberg.html Regards Marco
[Prev in Thread] | Current Thread | [Next in Thread] |