+- HP Forums (https://www.hpmuseum.org/forum)
+-- Forum: HP Calculators (and very old HP Computers) (/forum-3.html)
+--- Forum: HP Prime (/forum-5.html)
+--- Thread: Problem inverting complex matrix (/thread-14204.html)
Problem inverting complex matrix - jesuscf - 12-21-2019 09:48 AM
Consider this complex matrix:
M1 =
\begin{bmatrix}
0.0005 & 0 & 0 & -0.0005 & 0 & 0 & 0 & 1 \\
0 & 9.69696969697E-4+0.00125i & 0 & -0.00125i & -6.66666666667E-4 & -3.0303030303E-4 & 0 & 0 \\
0 & 0 & 3.57142857143E-4-0.0025i & 0 & 0.0025i & 0 & -3.57142857143E-4 & 0 \\
-0.0005 & -0.00125i & 0 & 0.0005+0.00125i & 0 & 0 & 0 & 0 \\
0 & -6.66666666667E-4 & 0.0025i & 0 & 6.66666666667E-4-0.0025i & 0 & 0 & 0 \\
0 & -3.0303030303E-4 & 0 & 0 & 0 & 3.0303030303E-4-1.66666666667E-3i & 0 & 0 \\
0 & 0 & -3.57142857143E-4 & 0 & 0 & 0 & 3.57142857143E-4+0.005i & 0 \\
1 & 0 & 0 & 0 & 0 & 0 & 0 & 0 \\
\end{bmatrix}
Inverting with the HP-Prime using 1/M1 works:
[attachment=7938]
Inverting it using inv(M1) fails.
[attachment=7939]
Is this a bug or there is something I failed to consider?
If you want to try, M1.hpmat is included in this zip file:
[attachment=7940]
RE: Problem inverting complex matrix - Werner - 12-21-2019 12:17 PM
On a 48GX, matrix inversion works ok, and 1/M is not allowed.
Strangely, COND causes a warmstart, but SVD reveals that the matrix is of rank 8.
The smallest singular value is in the order of 2e-4, so the matrix is far from being singular.
Cheers, Werner
RE: Problem inverting complex matrix - parisse - 12-21-2019 01:02 PM
inv calls the CAS, which works better with exact data. You can run approx(inv(exact(...))). If you want me to make a more complete diagnostic, please paste your input as plain text (other formats are more complicated to handle), using the following syntax (replace ... by appropriate data):
m1:=[[0.0005,0,...],[...]]; inv(m1);
RE: Problem inverting complex matrix - Joe Horn - 12-21-2019 01:15 PM
(12-21-2019 01:02 PM)parisse Wrote: inv calls the CAS, which works better with exact data. You can run approx(inv(exact(...))). If you want me to make a more complete diagnostic, please paste your input as plain text (other formats are more complicated to handle), using the following syntax (replace ... by appropriate data):
m1:=[[0.0005,0,...],[...]]; inv(m1);
Here's his matrix in plain-text format, for easy copy & paste:
M1:=[[0.0005,0,0,-0.0005,0,0,0,1],[0,9.69696969697E-4+0.00125*i,0,-0.00125*i,-6.66666666667E-4,-3.0303030303E-4,0,0],[0,0,3.57142857143E-4-0.0025*i,0,0.0025*i,0,-3.57142857143E-4,0],[-0.0005,-0.00125*i,0,0.0005+0.00125*i,0,0,0,0],[0,-6.66666666667E-4,0.0025*i,0,6.66666666667E-4-0.0025*i,0,0,0],[0,-3.0303030303E-4,0,0,0,3.0303030303E-4-1.66666666667E-3*i,0,0],[0,0,-3.57142857143E-4,0,0,0,3.57142857143E-4+0.005*i,0],[1,0,0,0,0,0,0,0]]
RE: Problem inverting complex matrix - parisse - 12-21-2019 04:52 PM
det(exact(m)) returns:
(-1245591211435133830166343338512439683+51639354804136441655581854007061939*i)/44328554791945219310918948796875252340096000000000000, evalf() to -2.80990710679e-17+1.16492303993e-18*i.
The smallest singular values of m is not too small, but there are 6 very small singular values, this is the reason why the determinant is small, and that raises an error at the end of the Gauss algorithm.
As I said, approx(inv(exact(...))) does the job.
RE: Problem inverting complex matrix - jesuscf - 12-21-2019 08:43 PM
(12-21-2019 04:52 PM)parisse Wrote: As I said, approx(inv(exact(...))) does the job.
approx(inv(exact(...))) takes over 15 seconds to produce the result.
RE: Problem inverting complex matrix - Jacob Wall - 12-21-2019 10:17 PM
Slightly off-topic, but anyone know of a simpler way to start editing a matrix than to type EDITMAT(M1)? I always loved the 50g where a cursor down with matrix on the stack just opened the matrix editor. Even if when a matrix is selected on the screen if there was something like an "EDIT" option on the menu that would launch the matrix editor. As is I find it cumbersome to enter and edit matrices, unless I'm missing something obvious?
RE: Problem inverting complex matrix - jesuscf - 12-21-2019 10:42 PM
(12-21-2019 10:17 PM)Jacob Wall Wrote: Slightly off-topic, but anyone know of a simpler way to start editing a matrix than to type EDITMAT(M1)?
Shift, Matrix, M1.
Shift, Mem, Matrices, M1.
RE: Problem inverting complex matrix - Jacob Wall - 12-21-2019 10:57 PM
Quote:Shift, Matrix, M1.
Right, and if I don't want/need a named matrix?
Perhaps even more so for RPN mode, it would be very handy to have:
- Shift Matrix - opens matrix editor
- OK or ENTER pops the matrix on the stack
- From the stack: a) select to reveal EDIT option and b) cursor DOWN to edit
This way to enter two matrices and perform some basic calculations would be IMHO way faster and more intuitive.
Shift Mem Matrices M0-M9 to edit those matrices is fine.
RE: Problem inverting complex matrix - parisse - 12-22-2019 08:54 AM
(12-21-2019 08:43 PM)jesuscf Wrote:(12-21-2019 04:52 PM)parisse Wrote: As I said, approx(inv(exact(...))) does the job.
approx(inv(exact(...))) takes over 15 seconds to produce the result.
This is expected since computations are done with rationals.
After more tracing with the numeric algorithm, I discovered that it failed because checking of non zero pivot failed: inv(M1) fails while inv(1000*M1) works. This should now be fixed in the source code.
RE: Problem inverting complex matrix - jesuscf - 12-22-2019 05:35 PM
(12-22-2019 08:54 AM)parisse Wrote:(12-21-2019 08:43 PM)jesuscf Wrote: approx(inv(exact(...))) takes over 15 seconds to produce the result.
This is expected since computations are done with rationals.
After more tracing with the numeric algorithm, I discovered that it failed because checking of non zero pivot failed: inv(M1) fails while inv(1000*M1) works. This should now be fixed in the source code.
Thank you Parisse! I also want to clarify that to get the correct result, inv(1000*M1) needs to be multiplied by 1000: inv(1000*M1)*1000.