Azzera filtri
Azzera filtri

inv(A)*B*inv(A) or A\B/A, which is more accurate for a full rank A and half rank B

39 visualizzazioni (ultimi 30 giorni)
My matrices are
A = [0 0 0 0 0 0 0 0 0 1
0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 2 0 0
0 0 0 0 0 0 6 0 0 0
0 0 0 0 0 24 0 0 0 0
1953050.78182336 390611.805629631 78122.6909803108 15624.6041672176 3124.93402773033 624.989444414569 124.998416658843 24.9997888874000 4.99997888869543 1
3515506.25066668 624981.527842487 109372.229170523 18749.6041663820 3124.94722207284 499.993666635371 74.9993666622000 9.99995777739086 1 0
5624833.75058238 874977.833364188 131247.229164674 18749.6833324371 2499.96833317685 299.997466648800 29.9998733321726 2 0 0
7874800.50027769 1049977.83331739 131247.783327059 14999.8099990611 1499.98733324400 119.999493328690 6 0 0 0
9449800.49985652 1049982.26661648 104998.669993428 8999.92399946400 599.997466643452 24 0 0 0 0];
B=[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 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
3515506.25066668 624981.527842487 109372.229170523 18749.6041663820 3124.94722207284 499.993666635371 74.9993666622000 9.99995777739086 1 0
5624833.75058238 874977.833364188 131247.229164674 18749.6833324371 2499.96833317685 299.997466648800 29.9998733321726 2 0 0
7874800.50027769 1049977.83331739 131247.783327059 14999.8099990611 1499.98733324400 119.999493328690 6 0 0 0
9449800.49985652 1049982.26661648 104998.669993428 8999.92399946400 599.997466643452 24 0 0 0 0
9449840.39954828 839989.359947423 62999.4679962480 3599.98479986071 120 0 0 0 0 0];
I tried both inv(A)*B*inv(A) or A\B/A, resulting in very different results for a gradient descent algorithm.
So, which is more accurate?
  2 Commenti
Torsten
Torsten il 11 Giu 2024
Did you think about why you get different results ? Is your matrix A near to singular and badly conditioned ? What do you get with
cond(A)
?
Can you include the matrix here for testing ?
Qiang Li
Qiang Li il 11 Giu 2024
Modificato: Qiang Li il 11 Giu 2024
Hi Torsten, here's cond(A)
K>> cond(A_2stepblock)
ans =
1.439014477720320e+07
I updated the question with the matrices.

Accedi per commentare.

Risposta accettata

John D'Errico
John D'Errico il 11 Giu 2024
Modificato: John D'Errico il 11 Giu 2024
Why would you want to compute the inverse of A TWICE ANYWAY? Even using both slash and backslash is arguably a bad idea, since again, you are doing way more work than you need.
Instead, decompose A. For example, you might decide to compute a QR factorization for A. Or perhaps an LU. Even then, the condition number of A is large enough that when you do two solves, it still inflates any noise in the system by the SQUARE of the condition number. And since you have given us only numbers to 14 significant digits, that means anything you get will probably be complete crap, no matter what you do.
A = [0 0 0 0 0 0 0 0 0 1
0 0 0 0 0 0 0 0 1 0
0 0 0 0 0 0 0 2 0 0
0 0 0 0 0 0 6 0 0 0
0 0 0 0 0 24 0 0 0 0
1953050.78182336 390611.805629631 78122.6909803108 15624.6041672176 3124.93402773033 624.989444414569 124.998416658843 24.9997888874000 4.99997888869543 1
3515506.25066668 624981.527842487 109372.229170523 18749.6041663820 3124.94722207284 499.993666635371 74.9993666622000 9.99995777739086 1 0
5624833.75058238 874977.833364188 131247.229164674 18749.6833324371 2499.96833317685 299.997466648800 29.9998733321726 2 0 0
7874800.50027769 1049977.83331739 131247.783327059 14999.8099990611 1499.98733324400 119.999493328690 6 0 0 0
9449800.49985652 1049982.26661648 104998.669993428 8999.92399946400 599.997466643452 24 0 0 0 0];
cond(A)
ans = 1.4390e+07
cond(A)^2
ans = 2.0708e+14
Do you see that any noise in the least significant digits of B will potentially be amplified by the SQUARE of the condition number here, and that means, since you have only 14 significant digits in both A and B, the result may be meaningless?
B=[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 0 0
0 0 0 0 0 0 0 0 0 0
0 0 0 0 0 0 0 0 0 0
3515506.25066668 624981.527842487 109372.229170523 18749.6041663820 3124.94722207284 499.993666635371 74.9993666622000 9.99995777739086 1 0
5624833.75058238 874977.833364188 131247.229164674 18749.6833324371 2499.96833317685 299.997466648800 29.9998733321726 2 0 0
7874800.50027769 1049977.83331739 131247.783327059 14999.8099990611 1499.98733324400 119.999493328690 6 0 0 0
9449800.49985652 1049982.26661648 104998.669993428 8999.92399946400 599.997466643452 24 0 0 0 0
9449840.39954828 839989.359947423 62999.4679962480 3599.98479986071 120 0 0 0 0 0];
Even so, we might try a QR or an LU. Simplest would be to do this:
Afac = decomposition(A)
Afac =
decomposition with properties: MatrixSize: [10 10] Type: 'lu' Show all properties
As you can see, it chose an LU by default.
Afac\B/Afac
ans = 10x10
-0.0001 -0.0001 -0.0001 -0.0001 -0.0000 0.0001 -0.0001 0.0001 -0.0001 0.0000 0.0013 0.0029 0.0027 0.0013 0.0003 -0.0013 0.0028 -0.0025 0.0011 -0.0002 -0.0097 -0.0215 -0.0202 -0.0096 -0.0020 0.0097 -0.0200 0.0170 -0.0070 0.0012 0.0323 0.0717 0.0672 0.0320 0.0067 -0.0323 0.0627 -0.0493 0.0184 -0.0027 -0.0403 -0.0896 -0.0840 -0.0400 -0.0083 0.0403 -0.0717 0.0504 -0.0160 0.0017 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 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
That requires MATLAB to perform only one operation in the form of a decomposition. As such, it will be considerably more efficient to perform. It will be a far better solution that using the inverse matrix almost always.
My real question is, can you improve the condition of A? If I look at the last 5 rows of A, it looks like a Vandermonde matrix. That would suggest re-scaling the polynomial variable to be on the order of 1, instead of going out as far as 5. All it is is a variable re-scaling, but it would improve the problem dramatically. For example...
Ahat = [A(1:5,:);A(6:10,:)*diag(5.^(-9:1:0))]
Ahat = 10x10
0 0 0 0 0 0 0 0 0 1.0000 0 0 0 0 0 0 0 0 1.0000 0 0 0 0 0 0 0 0 2.0000 0 0 0 0 0 0 0 0 6.0000 0 0 0 0 0 0 0 0 24.0000 0 0 0 0 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.0000 1.7999 1.6000 1.4000 1.2000 1.0000 0.8000 0.6000 0.4000 0.2000 0 2.8799 2.2399 1.6800 1.2000 0.8000 0.4800 0.2400 0.0800 0 0 4.0319 2.6879 1.6800 0.9600 0.4800 0.1920 0.0480 0 0 0 4.8383 2.6880 1.3440 0.5760 0.1920 0.0384 0 0 0 0
<mw-icon class=""></mw-icon>
<mw-icon class=""></mw-icon>
cond(Ahat)
ans = 7.2583e+04
  9 Commenti
Qiang Li
Qiang Li il 24 Giu 2024 alle 11:51
I think for small matrices as in my case, decomposing the matrix is not worth it, as decomposition costs significantly more than the savings from the inverse. It may work for large matrices though.
Qiang Li
Qiang Li il 26 Giu 2024 alle 0:47
The rescaling suggestion is making a lot of sense now. Definitely need to rescale A to improve the accuracy.

Accedi per commentare.

Più risposte (2)

Matt J
Matt J il 11 Giu 2024
Modificato: Matt J il 11 Giu 2024
I tried both inv(A)*B*inv(A) or A\B/A, resulting in very different results for a gradient descent algorithm.
From this, it would appear that what you are really trying to do is solve A*inv(B)*A*x=y. If so, you should consider lscov,
x=lscov(A,y,B);
  1 Commento
Qiang Li
Qiang Li il 12 Giu 2024
Hi Matt,
I'm not solving for x at the moment. But I may actually need to do it later in my project. So thanks in advance!

Accedi per commentare.


Qiang Li
Qiang Li il 6 Lug 2024 alle 2:13
The short answer for my application: A\B/A.
inv(A)*B*inv(A) introduces inconsistency after iterations.

Categorie

Scopri di più su Linear Algebra in Help Center e File Exchange

Prodotti


Release

R2023b

Community Treasure Hunt

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

Start Hunting!

Translated by