matrix inverse results different between r2023b and 2021b

Question

Qiang Li 2024-6-11

1
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2127476-matrix-inverse-results-different-between-r2023b-and-2021b

评论： 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'

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

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

Answer 1

David Goodmanson 2024-6-11

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/2127476-matrix-inverse-results-different-between-r2023b-and-2021b#answer_1470591

编辑：David Goodmanson 2024-6-13

在 MATLAB Online 中打开

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更早的评论隐藏 1更早的评论

David Goodmanson 2024-6-13

在 MATLAB Online 中打开

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 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.

请先登录，再进行评论。