matrix inverse results different between r2023b and 2021b

4 次查看（过去 30 天）
Qiang Li
Qiang Li 2024-6-11
评论： Qiang Li 2024-6-14
I was solving a high-dimensional nonlinear optimization problem with gradient descent method and found out that R2023b and R2021b calculates matrix inverse differently, resulting in completely different gradients and convergence character of the algorithm.
2020b and 2021b gave the following gradient and converged after 18 iterations.
direction = [4.33783e-06 2.33596e-06 -5.33992e-06 -8.36479e-06 -9.19797e-05 1 ]
2023b gave the following gradient and didn't converge after 44 iterations.
direction = [4.17191e-06 2.24661e-06 -5.13617e-06 -7.53861e-06 -5.19154e-05 1 ]
The code is very slow, so I can't let it run forever and find out when it will converge.
I narrowed it down to the inverse of a 10by10 matrix
A_2stepblock =
1.0e+06 *
Columns 1 through 7
0 0 0 0 0 0 0
0 0 0 0 0 0 0
0 0 0 0 0 0 0
0 0 0 0 0 0 0.000006000000000
0 0 0 0 0 0.000024000000000 0
1.953050781823359 0.390611805629631 0.078122690980311 0.015624604167218 0.003124934027730 0.000624989444415 0.000124998416659
3.515506250666680 0.624981527842487 0.109372229170523 0.018749604166382 0.003124947222073 0.000499993666635 0.000074999366662
5.624833750582380 0.874977833364188 0.131247229164674 0.018749683332437 0.002499968333177 0.000299997466649 0.000029999873332
7.874800500277688 1.049977833317391 0.131247783327059 0.014999809999061 0.001499987333244 0.000119999493329 0.000006000000000
9.449800499856519 1.049982266616475 0.104998669993428 0.008999923999464 0.000599997466643 0.000024000000000 0
Columns 8 through 10
0 0 0.000001000000000
0 0.000001000000000 0
0.000002000000000 0 0
0 0 0
0 0 0
0.000024999788887 0.000004999978889 0.000001000000000
0.000009999957777 0.000001000000000 0
0.000002000000000 0 0
0 0 0
0 0 0
r2021b gives
inv(A_2stepblock)
ans =
Columns 1 through 7
-0.000035841361961 -0.000089603026574 -0.000096002837407 -0.000053334684477 -0.000013333614821 0.000035841361961 -0.000089603026574
0.000806427239167 0.002048060531355 0.002240056748025 0.001280027022812 0.000333338963074 -0.000806427239167 0.001984058639750
-0.006912204293321 -0.017920453984200 -0.020160425609289 -0.012000202670662 -0.003333375556299 0.006912204293321 -0.016640421556757
0.026880680976300 0.071681513277471 0.084001418694634 0.053334008900781 0.016666807409588 -0.026880680976300 0.062721324117785
-0.040320851218577 -0.112001891592844 -0.140001773364547 -0.100000844457525 -0.041666842594947 0.040320851218577 -0.089601513274271
0 0 0 0 0.041666666666667 0 0
0 0 0 0.166666666666667 0 0 0
0 0 0.500000000000000 0 0 0 0
0 1.000000000000000 0 0 0 0 0
1.000000000000000 0 0 0 0 0 0
Columns 8 through 10
0.000096002837407 -0.000053334684477 0.000013333614821
-0.002080052694595 0.001120023644960 -0.000266671170459
0.016960358052258 -0.008800148625151 0.002000025333779
-0.061601040376060 0.030667055117945 -0.006666722963834
0.084001064018720 -0.040000337783004 0.008333368518987
0 0 0
0 0 0
0 0 0
0 0 0
0 0 0
and r2023b gives
inv(A_2stepblock)
ans =
Columns 1 through 6
-0.000035841361961 -0.000089603026574 -0.000096002837407 -0.000053334684477 -0.000013333614821 0.000035841361961
0.000806427239167 0.002048060531355 0.002240056748026 0.001280027022812 0.000333338963074 -0.000806427239167
-0.006912204293323 -0.017920453984205 -0.020160425609294 -0.012000202670665 -0.003333375556300 0.006912204293323
0.026880680976306 0.071681513277487 0.084001418694652 0.053334008900791 0.016666807409591 -0.026880680976306
-0.040320851218585 -0.112001891592865 -0.140001773364570 -0.100000844457540 -0.041666842594951 0.040320851218585
0 0 0 0 0.041666666666667 0
0 0 0 0.166666666666667 0 0
0 0 0.500000000000000 0 0 0
0 1.000000000000000 0 0 0 0
1.000000000000000 0 0 0 0 0
Columns 7 through 10
-0.000089603026574 0.000096002837407 -0.000053334684477 0.000013333614821
0.001984058639750 -0.002080052694595 0.001120023644961 -0.000266671170459
-0.016640421556761 0.016960358052262 -0.008800148625154 0.002000025333779
0.062721324117800 -0.061601040376075 0.030667055117952 -0.006666722963836
-0.089601513274291 0.084001064018740 -0.040000337783014 0.008333368518990
0 0 0 0
0 0 0 0
0 0 0 0
0 0 0 0
0 0 0 0
I used
num2str(A,100)
to make sure it's the same matrix to the last digit in both versions.
From the release notes, they did changed quite a few matrix operations from 2022a.
So, which version is more accurate?
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Edit after David Goodmanson's response, the second "r2021b" should be "r2023b", I corrected that.
Running DG's answer directly on my computer to show the difference between the versions, on the r2021b
K>> A = A_2stepblock;
K>> F = inv(A);
K>> max(abs(A*F-eye(10,10)),[],'all')
max(abs(F*A-eye(10,10)),[],'all')
ans =
1.818989403545856e-12
ans =
3.887682601841180e-11
K>> norm(A*F-eye(10,10))
norm(F*A-eye(10,10))
ans =
5.914730386430228e-12
ans =
4.379288551440263e-11
K>> num2hex(A)
ans =
100×16 char array
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'413dcd1ac825935f'
'414ad2392015d885'
'4155750070098aae'
'415e0a3c20048cb5'
'416206290ffed319'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4117d74f38f6f95a'
'412312ab0e415ed4'
'412ab3c3aaaeb5f6'
'41300579d55449dd'
'4130057e4440fa32'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'40f312ab0e415ed4'
'40fab3c3aaaeb5f6'
'41000579d55449dd'
'4100057e4440fa32'
'40f9a26ab84b0751'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'40ce844d5559f487'
'40d24f66aaa978fc'
'40d24f6bbbb7f95e'
'40cd4be7ae0c9aa6'
'40c193f6459d4bb4'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'40a869de38e1f6a6'
'40a869e4fa4aa1d2'
'40a387efc95dbc6f'
'40976ff3077c64f0'
'4082bffacfcaa3cf'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4038000000000000'
'408387ea61d54e42'
'407f3fe60efc60b1'
'4072bff59f96b726'
'405dfff7b2ddd2e6'
'4038000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4018000000000000'
'0000000000000000'
'405f3fe60efc60b1'
'4052bff59f96b726'
'403dfff7b2ddd2e6'
'4018000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4000000000000000'
'0000000000000000'
'0000000000000000'
'4038fff22a1e4988'
'4023fffa773e8c99'
'4000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'3ff0000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4013fffa773e8c99'
'3ff0000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'3ff0000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'3ff0000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
On the r2023b
K>> F = inv(A_2stepblock);
K>> A = A_2stepblock;
K>> max(abs(A*F-eye(10,10)),[],'all')
max(abs(F*A-eye(10,10)),[],'all')
ans =
2.728484105318785e-12
ans =
1.472269235603454e-11
K>> norm(A*F-eye(10,10))
norm(F*A-eye(10,10))
ans =
7.261754408131635e-12
ans =
2.136608381476120e-11
num2hex(A_2stepblock)
ans =
100×16 char array
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'413dcd1ac825935f'
'414ad2392015d885'
'4155750070098aae'
'415e0a3c20048cb5'
'416206290ffed319'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4117d74f38f6f95a'
'412312ab0e415ed4'
'412ab3c3aaaeb5f6'
'41300579d55449dd'
'4130057e4440fa32'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'40f312ab0e415ed4'
'40fab3c3aaaeb5f6'
'41000579d55449dd'
'4100057e4440fa32'
'40f9a26ab84b0751'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'40ce844d5559f487'
'40d24f66aaa978fc'
'40d24f6bbbb7f95e'
'40cd4be7ae0c9aa6'
'40c193f6459d4bb4'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'40a869de38e1f6a6'
'40a869e4fa4aa1d2'
'40a387efc95dbc6f'
'40976ff3077c64f0'
'4082bffacfcaa3cf'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4038000000000000'
'408387ea61d54e42'
'407f3fe60efc60b1'
'4072bff59f96b726'
'405dfff7b2ddd2e6'
'4038000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4018000000000000'
'0000000000000000'
'405f3fe60efc60b1'
'4052bff59f96b726'
'403dfff7b2ddd2e6'
'4018000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4000000000000000'
'0000000000000000'
'0000000000000000'
'4038fff22a1e4988'
'4023fffa773e8c99'
'4000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'3ff0000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'4013fffa773e8c99'
'3ff0000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'3ff0000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'3ff0000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'
'0000000000000000'

请先登录，再回答此问题。

采纳的回答

David Goodmanson
David Goodmanson 2024-6-11
编辑：David Goodmanson 2024-6-13
Hi QL.
you called both the inverse matrices '2021b' so I called the first one F and the second one S (and the orginal matrix is A). So why not find out. Sounds like you were careful to get all the digits in place but there will still be a bit of uncertainty in copying numbers off of the display as opposed to having them in memory. Num2hex is the real solution because it shows exactly what is stored in memory**. But with that caveat and on the basis of how large [A*inverse(A) - I] (or the other way round) is, F is slightly better in this particular case. In both measures below, S*A is the outlier.
cond(A)
ans =
1.4390e+07 % pretty close to singular
format compact
% largest nonzero element of the matrix that should be zero
max(abs(A*F-eye(10,10)),[],'all')
max(abs(F*A-eye(10,10)),[],'all')
max(abs(A*S-eye(10,10)),[],'all')
max(abs(S*A-eye(10,10)),[],'all')
ans =
2.7399e-09
ans =
5.1108e-09
ans =
2.9945e-09
ans =
9.9112e-09
% norm of the matrix that should be zero
norm(A*F-eye(10,10))
norm(F*A-eye(10,10))
norm(A*S-eye(10,10))
norm(S*A-eye(10,10))
ans =
9.2492e-09
ans =
6.7615e-09
ans =
8.6579e-09
ans =
1.3678e-08
** for example
num2hex(1/100) ans = '3f847ae147ae147b'
num2hex(1/128) ans = '3f80000000000000'
illustrates the difference between storing inverse powers of 10 and inverse powers of 2.
  3 个评论
显示 1更早的评论
David Goodmanson
David Goodmanson 2024-6-13
HI Ql, It doesn't seem to. If you run the following example a few times there are cases where b*a is bigger in both measures, smaller in both measures, or one of each.
a = 2*rand(10,10) -1;
b = inv(a);
maxx(abs(a*b-eye(size(a))))
maxx(abs(b*a-eye(size(a))))
norm(a*b-eye(size(a)))
norm(b*a-eye(size(a)))
function b = maxx(a)
% maximum value of the elements of an n-dimensional array
%
% b = maxx(a)
b = max(a(:));
end
Qiang Li
Qiang Li 2024-6-14
I've improved my code over the past couple of days. The difference is still there, but the code is more robust now and I don't feel like that it matters that much any more. Thanks for all the help.

请先登录，再进行评论。

更多回答（0 个）

类别

在 Help Center 和 File Exchange 中查找有关 Logical 的更多信息

标签

产品


版本

R2023b

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by