<div dir="ltr">The matrix has three rows columns equal to zero, which makes it singular. Several of the eigenvalues are zero and that's probably why inverse_matrix doesn't work and returns all missing values. You will have the same problem if you just replace all zeros with the same non-zero number. In that case you have several eigenvalues that are the same and the solution to the linear system that Will suggested does not have a unique solution.<div><br></div><div>I attached a matlab code that computes the three examples Dennis sent and I added one more that just replaces the zeros with 7 and you still have a singular matrix. So non-invertable.<br></div><div><br></div><div>I don't know what you are trying to do, but you could just delete those three column/rows that have only zeros, since they don't provide any additional information about the matrix, and invert the rest of the matrix. So, shrinking you matrix and inverting only on the subspace where the matrix is non-singular. </div><div><br></div><div>Maria</div></div><div class="gmail_extra"><br><div class="gmail_quote">On Mon, Nov 10, 2014 at 4:25 PM, Walter Hannah <span dir="ltr"><<a href="mailto:whannah@atmos.colostate.edu" target="_blank">whannah@atmos.colostate.edu</a>></span> wrote:<br><blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div bgcolor="#FFFFFF" text="#000000">
Will,<br>
<br>
I tried your suggestion, but all I got back was missing values.<br>
<br>
After Dennis's suggestion, I tried replacing the zeros with tiny
random numbers, using this code:<br>
<blockquote><tt>R = random_normal( 0.,1e-12, (/nlev*2,nlev*2/))</tt><tt><br>
</tt><tt>M = where(M.eq.0,R,M)</tt><tt><br>
</tt><tt>iM = inverse_matrix(M)</tt><br>
</blockquote>
Interestingly, this worked perfectly, and I got what I expected.<br>
I also tried using this approach to trying Will's suggestion, and
that also worked with the random numbers!<br>
<br>
I can't figure out a mathematical explanation, but this seems like a
pretty nice fix. Might be good to add this to the online
documentation.<br>
<br>
Thanks,<br>
Walter<div><div class="h5"><br>
<br>
<div>On 11/10/14 6:17 PM, Will Hobbs wrote:<br>
</div>
<blockquote type="cite">
<div>Walter</div>
<div><br>
</div>
<div>Instead of using the inverse matrix() function, try using
solve_linsys(A, B) instead, for A = your matrix and B = I (the
identity matrix). Since A * inv(A) = I, the output will be the
inverse.</div>
<div><br>
</div>
<div>Mathematically this is no different to inverse_matrix(), but
the code is a little different. In particular, if the matrix is
singular (i.e. the determinant = 0, and it cannot be inverted)
then it will return an error code.</div>
<div><a href="https://www.ncl.ucar.edu/Document/Functions/Built-in/solve_linsys.shtml" target="_blank">https://www.ncl.ucar.edu/Document/Functions/Built-in/solve_linsys.shtml</a></div>
<div><br>
</div>
<div>This might help you figure out if the matrix itself is the
problem (I suspect, as Dennis has suggested, that may be the
issue).</div>
<div><br>
</div>
<div>Will</div>
<div><br>
</div>
<div><br>
</div>
<span>
<div style="font-family:Calibri;font-size:11pt;text-align:left;color:black;BORDER-BOTTOM:medium none;BORDER-LEFT:medium none;PADDING-BOTTOM:0in;PADDING-LEFT:0in;PADDING-RIGHT:0in;BORDER-TOP:#b5c4df 1pt solid;BORDER-RIGHT:medium none;PADDING-TOP:3pt">
<span style="font-weight:bold">From: </span>Walter Hannah
<<a href="mailto:whannah@atmos.colostate.edu" target="_blank">whannah@atmos.colostate.edu</a>><br>
<span style="font-weight:bold">Date: </span>Tuesday, 11
November 2014 9:52 AM<br>
<span style="font-weight:bold">To: </span>Dennis Shea <<a href="mailto:shea@ucar.edu" target="_blank">shea@ucar.edu</a>><br>
<span style="font-weight:bold">Cc: </span>"<a href="mailto:ncl-talk@ucar.edu" target="_blank">ncl-talk@ucar.edu</a>"
<<a href="mailto:ncl-talk@ucar.edu" target="_blank">ncl-talk@ucar.edu</a>><br>
<span style="font-weight:bold">Subject: </span>Re: [ncl-talk]
inverse_matrix problem<br>
</div>
<div><br>
</div>
<div>
<div bgcolor="#FFFFFF" text="#000000">Dennis,<br>
<br>
Wow! Thanks for going above and beyond! This helps a lot.
I'll try to look for a way to get a modified matrix that
doesn't have all those zeros.<br>
<br>
Walter<br>
<br>
<div>On 11/10/14 5:50 PM, Dennis
Shea wrote:<br>
</div>
<blockquote type="cite">
<div dir="ltr">
<div>
<div>
<div>You have rows of 0.0 ... I think that is
causing issues. I'm sure some 'math-whiz' could
answer but I attached a simple test code that
reads your data (thx for providing). It calculates
the inverse and the identity matrix. Note that NCL
is using a direct call to the documented LAPACK
routines.<br>
<br>
===<br>
</div>
TEST1: duplicates your findings<br>
<br>
</div>
TEST 2: create a bogus 20x20 array of random numbers.
No row consists of all constants. Identity matrix is
successfully calculated.<br>
<br>
</div>
TEST 3: wherever the original array had a row of all
0.0, substitute random numbers. The identity matrix is
successfully calculated.</div>
<div class="gmail_extra"><br>
<div class="gmail_quote">On Mon, Nov 10, 2014 at 1:57
PM, Walter Hannah <span dir="ltr">
<<a href="mailto:whannah@atmos.colostate.edu" target="_blank">whannah@atmos.colostate.edu</a>></span>
wrote:<br>
<blockquote class="gmail_quote" style="margin:0 0 0 .8ex;border-left:1px #ccc solid;padding-left:1ex">
<div bgcolor="#FFFFFF" text="#000000">Hi,<br>
<br>
I'm trying to use the "inverse_matrix" function.
The matrix doesn't have any missing values, but I
still get this error message:<br>
<blockquote><tt>Variable: M</tt><tt><br>
</tt><tt>Type: float</tt><tt><br>
</tt><tt>Total Size: 1600 bytes</tt><tt><br>
</tt><tt> 400 values</tt><tt><br>
</tt><tt>Number of Dimensions: 2</tt><tt><br>
</tt><tt>Dimensions and sizes: [pertlev | 20]
x [tendlev | 20]</tt><tt><br>
</tt><tt>Coordinates:</tt><tt><br>
</tt><tt> pertlev: [50..1950]</tt><tt><br>
</tt><tt> tendlev: [50..1950]</tt><tt><br>
</tt><tt>Number Of Attributes: 2</tt><tt><br>
</tt><tt> missing_value : 9.96921e+36</tt><tt><br>
</tt><tt> _FillValue : 9.96921e+36</tt><tt><br>
</tt><tt>(0)
----------------------------------------------------------------------------</tt><tt><br>
</tt><tt><br>
</tt><tt> 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00<br>
0.00 -0.31 -0.06 0.15 -0.27 -0.56
-1.21 -1.38 -1.36 1.49 0.00 0.00
4.65 0.37 0.46 0.47 0.52 0.55
0.77 2.49<br>
0.00 0.26 -0.35 -0.29 -0.65 -1.13
-2.24 -2.45 -2.47 2.91 0.00 0.00
4.07 2.45 0.95 0.87 0.90 1.05
1.55 4.81<br>
0.00 -0.11 1.43 -0.91 -1.69 -1.52
-2.63 -2.81 -2.69 3.17 0.00 0.00
-2.35 2.61 1.49 1.01 1.00 1.18
1.78 5.32<br>
0.00 -0.02 0.55 1.08 -1.74 -2.04
-2.59 -2.71 -2.39 2.63 0.00 0.00
-0.75 -1.10 1.35 1.02 0.89 1.00
1.57 4.87<br>
0.00 -0.06 0.16 0.52 0.72 -1.86
-2.88 -2.74 -2.29 2.31 0.00 0.00
0.31 -1.30 -0.26 0.87 0.75 0.82
1.27 4.50<br>
0.00 -0.13 -0.07 0.20 0.54 0.27
-2.74 -3.10 -2.27 2.04 0.00 0.00
0.08 -0.49 -0.47 0.14 0.72 0.69
1.13 4.24<br>
0.00 -0.15 -0.06 0.20 -0.17 0.02
-0.77 -3.87 -3.13 2.39 0.00 0.00
0.69 -0.08 0.09 0.03 0.36 0.90
1.29 4.59<br>
0.00 -0.14 0.05 0.09 -0.17 -0.29
-0.14 -0.66 -5.36 2.00 0.00 0.00
0.39 0.17 0.14 0.07 -0.01 0.22
1.00 4.05<br>
0.00 0.02 -0.07 -0.01 0.00 0.02
0.43 0.94 0.77 -2.90 0.00 0.00
-0.09 0.09 -0.04 -0.07 -0.17 -0.26
-0.49 0.54<br>
0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00<br>
0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00 0.00 0.00 0.00 0.00
0.00 0.00<br>
0.00 -0.01 -0.04 0.12 -0.02 -0.01
-0.03 -0.05 -0.05 0.04 0.00 0.00
-2.17 -0.25 -0.14 0.01 0.02 0.02
0.02 0.12<br>
0.00 0.00 -0.19 0.11 0.28 0.06
0.11 0.06 0.07 -0.13 0.00 0.00
2.75 -1.97 -0.04 -0.12 -0.05 -0.05
-0.08 -0.10<br>
0.00 -0.02 -0.14 -0.08 0.31 0.21
0.13 0.08 0.06 -0.08 0.00 0.00
0.48 1.11 -1.09 -0.13 -0.10 -0.06
-0.09 -0.12<br>
0.00 -0.01 -0.03 -0.11 -0.07 0.32
0.27 0.12 0.18 -0.17 0.00 0.00
0.02 0.63 0.52 -0.94 -0.02 -0.15
-0.08 -0.25<br>
0.00 0.02 0.08 -0.06 -0.24 0.10
0.43 0.26 0.16 -0.17 0.00 0.00
0.07 0.31 0.35 0.40 -1.06 -0.06
-0.20 -0.26<br>
0.00 0.06 0.03 -0.10 0.07 0.11
1.03 1.38 0.53 -1.10 0.00 0.00
-0.44 0.14 -0.03 0.08 0.19 -1.61
-0.76 -1.16<br>
0.00 0.11 -0.03 -0.08 0.07 0.16
0.46 1.28 1.65 -1.21 0.00 0.00
-0.25 0.07 -0.06 0.00 0.14 0.41
-3.04 -0.66<br>
0.00 0.12 -0.30 -0.15 0.80 1.58
2.86 3.62 5.02 -3.15 0.00 0.00
-1.23 -0.91 -0.80 -0.94 -0.85 -0.66
0.21 -10.64<br>
<br>
(0)
----------------------------------------------------------------------------<br>
(0) # missing values: 0<br>
(0) # values = 0.0 : 111</tt><tt><br>
</tt><tt>warning:inverse_matrix: info = 1;
missing values not allowed</tt><br>
</blockquote>
<br>
I thought it might have something to do with those
zero values, but if I replace them with non-zero
values I get a similar, but not identical, result:<br>
<blockquote><tt>...<br>
(0)
----------------------------------------------------------------------------</tt><tt><br>
</tt><tt>(0) # missing values: 0</tt><tt><br>
</tt><tt>(0) # values = 0.0 : 0</tt><tt><br>
</tt><tt>warning:inverse_matrix: info = 20;
missing values not allowed</tt><br>
</blockquote>
<br>
<br>
Thanks,<br>
Walter </div>
<br>
_______________________________________________<br>
ncl-talk mailing list<br>
List instructions, subscriber options, unsubscribe:<br>
<a href="http://mailman.ucar.edu/mailman/listinfo/ncl-talk" target="_blank">http://mailman.ucar.edu/mailman/listinfo/ncl-talk</a><br>
<br>
</blockquote>
</div>
<br>
</div>
</blockquote>
<br>
</div>
</div>
</span>
</blockquote>
<br>
</div></div></div>
<br>_______________________________________________<br>
ncl-talk mailing list<br>
List instructions, subscriber options, unsubscribe:<br>
<a href="http://mailman.ucar.edu/mailman/listinfo/ncl-talk" target="_blank">http://mailman.ucar.edu/mailman/listinfo/ncl-talk</a><br>
<br></blockquote></div><br></div>