Problem inverting complex matrix
12-21-2019, 09:48 AM
Post: #1
 jesuscf Junior Member Posts: 9 Joined: Oct 2017
Problem inverting complex matrix
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:

Inverting it using inv(M1) fails.

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:

12-21-2019, 12:17 PM
Post: #2
 Werner Senior Member Posts: 891 Joined: Dec 2013
RE: Problem inverting complex matrix
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

41CV†,42S,48GX,49G,DM42,DM41X,17BII,15CE,DM15L,12C,16CE
12-21-2019, 01:02 PM
Post: #3
 parisse Senior Member Posts: 1,327 Joined: Dec 2013
RE: Problem inverting complex matrix
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);
12-21-2019, 01:15 PM
Post: #4
 Joe Horn Senior Member Posts: 2,011 Joined: Dec 2013
RE: Problem inverting complex matrix
(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]]

<0|ɸ|0>
-Joe-
12-21-2019, 04:52 PM
Post: #5
 parisse Senior Member Posts: 1,327 Joined: Dec 2013
RE: Problem inverting complex matrix
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.
12-21-2019, 08:43 PM
Post: #6
 jesuscf Junior Member Posts: 9 Joined: Oct 2017
RE: Problem inverting complex matrix
(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.
12-21-2019, 10:17 PM
Post: #7
 Jacob Wall Member Posts: 164 Joined: Dec 2013
RE: Problem inverting complex matrix
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?
12-21-2019, 10:42 PM (This post was last modified: 12-21-2019 10:44 PM by jesuscf.)
Post: #8
 jesuscf Junior Member Posts: 9 Joined: Oct 2017
RE: Problem inverting complex matrix
(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.
12-21-2019, 10:57 PM
Post: #9
 Jacob Wall Member Posts: 164 Joined: Dec 2013
RE: Problem inverting complex matrix
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.
12-22-2019, 08:54 AM
Post: #10
 parisse Senior Member Posts: 1,327 Joined: Dec 2013
RE: Problem inverting complex matrix
(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.
12-22-2019, 05:35 PM
Post: #11
 jesuscf Junior Member Posts: 9 Joined: Oct 2017
RE: Problem inverting complex matrix
(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.
 « Next Oldest | Next Newest »

User(s) browsing this thread: 1 Guest(s)