matrix inverse results different between r2023b and 2021b
조회 수: 5 (최근 30일)
이전 댓글 표시
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?
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
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'
댓글 수: 0
채택된 답변
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
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
추가 답변 (0개)
참고 항목
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!