LU decomposition algorithm change

9 visualizzazioni (ultimi 30 giorni)
James Tursa
James Tursa il 2 Giu 2023
Commentato: James Tursa il 6 Giu 2023
For this particular example:
A=[ 1 2 3;
4 5 6;
7 8 9;
10 11 12];
I get two different results depending on MATLAB version. E.g.,
Online MATLAB R2023a:
[L,U,P] = lu(A)
L = 4×3
1.0000 0 0 0.1000 1.0000 0 0.7000 0.3333 1.0000 0.4000 0.6667 0
U = 3×3
10.0000 11.0000 12.0000 0 0.9000 1.8000 0 0 0
P = 4×4
0 0 0 1 1 0 0 0 0 0 1 0 0 1 0 0
You can see that the last row of U is exactly 0's, so the L(4,3) spot is arbitrary in that sense but in this case the algorithm produced exactly 0.
Now compare this result with my PCwin MATLAB R2021a:
>> A
A =
1 2 3
4 5 6
7 8 9
10 11 12
>> [L,U,P] = lu(A)
L =
1.0000 0 0
0.1000 1.0000 0
0.7000 0.3333 1.0000
0.4000 0.6667 0.4599
U =
10.0000 11.0000 12.0000
0 0.9000 1.8000
0 0 0.0000
P =
0 0 0 1
1 0 0 0
0 0 1 0
0 1 0 0
>> U(3,3)
ans =
6.9204e-16
>> L(4,3)
ans =
0.4599
You can see that in this case the U(3,3) spot is not exactly 0, and that "arbitrary" L(4,3) spot now has a non-zero value. Any L(4,3) value that is small enough would work equally as well. Slightly different LU algorithms might produce slightly different results due to floating point effects, so both answers are "correct" in that sense and that is not my concern. I am simply wondering when this change took place. I no longer have access to the many MATLAB versions I used to, so I am posing to the group. Which MATLAB version between R2021a and R2023a made the change to report exact 0's for U(3,3) and L(4,3) for this example? Can someone with access to R2021b, R2022a, R2022b, R2023a on desktops run the above code and post what you get? Thanks.

Risposte (2)

chicken vector
chicken vector il 2 Giu 2023
Modificato: chicken vector il 2 Giu 2023
This was introduce in R2022a and might be a promising candidate to answer your question.
Quoting Matlab realease notes:
inv Function: Improved performance when inverting large triangular matrices
The inv function shows improved performance when operating on large triangular matrices.
For example, inverting a 5,000-by-5,000 upper triangular matrix is about 3.7x faster than in the previous release.
function timingInv
rng default
A = randn(5e3);
[~,R] = lu(A);
tic
Y = inv(R);
toc
end
The approximate execution times are:
R2021b: 1.1 s
R2022a: 0.3 s
The code was timed on a Windows 10, Intel Xeon CPU W-2133 @ 3.60 GHz test system by calling the timingInv function.
  7 Commenti
Rik
Rik il 6 Giu 2023
I used the A from your original post. isequal(B,C) returns true on both releases.
R2021b =
10.0000 11.0000 12.0000
0.1000 0.9000 1.8000
0.7000 0.3333 0.0000
0.4000 0.6667 0.4599
R2022a =
10.0000 11.0000 12.0000
0.1000 0.9000 1.8000
0.7000 0.3333 0.0000
0.4000 0.6667 0
James Tursa
James Tursa il 6 Giu 2023
@Rik Thanks. That confirms that the change was in the LAPACK library itself, and not MATLAB's interface code to LAPACK.

Accedi per commentare.


Christine Tobler
Christine Tobler il 2 Giu 2023
Modificato: Christine Tobler il 2 Giu 2023
As you point out, this kind of difference in result is to be expected and is well withing expected round-off behavior.
Another point is that the algorithm change here is most likely not anything deep: The algorithm used for LU has some tuning parameters (the block size) which can change from CPU to CPU for performance reasons. So the reason for the change here is as likely to be that MATLAB Online uses a different type of CPU as it is a change between MATLAB versions.
Round-off can also change based on what registers are available on any CPU, as a performance-oriented implementation will try to use these to best advantage.
If the reason isn't different CPUs, and also isn't different versions of the BLAS or LAPACK library used, the most likely reason is this R2022b release note:
mldivide and pagemldivide Functions: Improved performance with small matrices

Tag

Prodotti


Release

R2021a

Community Treasure Hunt

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

Start Hunting!

Translated by