<html>
<head>
<meta http-equiv="content-type" content="text/html; charset=UTF-8">
</head>
<body text="#000000" bgcolor="#FFFFFF">
<br>
<font face="Courier New">Dear all,<br>
<br>
I wonder if the inversion of a matrix of polynomials is working
properly. I have the following code from the direct resolution of
a two-section constant resistance electrical network::<br>
</font><br>
<font face="Courier New"><span style="color:rgb(0,0,0);"><font
color="#009900">// Data</font><br>
R1</span> <span style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">600</span><span
style="color:rgb(0,0,0);">; </span><span
style="color:rgb(0,0,0);">R2</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">600</span><span
style="color:rgb(0,0,0);">; </span><span
style="color:rgb(0,0,0);">R3</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">600</span><span
style="color:rgb(0,0,0);">;</span> <span
style="color:rgb(0,0,0);">R4</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">600</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);"><br>
R5</span> <span style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">10070</span><span
style="color:rgb(0,0,0);">; </span><span
style="color:rgb(0,0,0);">R6</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">2774</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">R7</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">36</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">R8</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">130</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">R9</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">600</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);"><br>
L1</span> <span style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">2.65</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">L2</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">0.0318</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">L3</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">0.00669</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">L4</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">0.00795</span><span
style="color:rgb(0,0,0);">;</span>
<br>
<span style="color:rgb(0,0,0);">C1</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">18.4e-9</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">C2</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">22.1e-9</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">C3</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">7.3e-6</span><span
style="color:rgb(0,0,0);">;</span>
<span style="color:rgb(0,0,0);">C4</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">88.4e-9</span><span
style="color:rgb(0,0,0);">;<br>
<br>
<span style="color:rgb(0,0,0);"><font color="#009900">// Series
and parallel impedances</font></span><br>
</span><span style="color:rgb(0,0,0);">Z1</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(0,0,0);">R7</span> <span
style="color:rgb(92,92,92);">+</span> <span
style="color:rgb(0,0,0);">L3</span><span
style="color:rgb(92,92,92);">*</span><span
style="color:rgb(218,112,214);">%s</span> <span
style="color:rgb(92,92,92);">+</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(74,85,219);">(</span><span
style="color:rgb(0,0,0);">C3</span><span
style="color:rgb(92,92,92);">*</span><span
style="color:rgb(218,112,214);">%s</span><span
style="color:rgb(74,85,219);">);</span>
<br>
<span style="color:rgb(0,0,0);">Z2</span> <span
style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(74,85,219);">(</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R5</span> <span
style="color:rgb(92,92,92);">+</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(74,85,219);">(</span><span
style="color:rgb(0,0,0);">L1</span><span
style="color:rgb(92,92,92);">*</span><span
style="color:rgb(218,112,214);">%s</span><span
style="color:rgb(74,85,219);">)</span> <span
style="color:rgb(92,92,92);">+</span> <span
style="color:rgb(0,0,0);">C1</span><span
style="color:rgb(92,92,92);">*</span><span
style="color:rgb(218,112,214);">%s</span><span
style="color:rgb(74,85,219);">)</span>;
<span style="color:rgb(0,0,0);"><br>
Z3</span> <span style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(0,0,0);">R8</span> <span
style="color:rgb(92,92,92);">+</span> <span
style="color:rgb(0,0,0);">L4</span><span
style="color:rgb(92,92,92);">*</span><span
style="color:rgb(218,112,214);">%s</span> <span
style="color:rgb(92,92,92);">+</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(74,85,219);">(</span><span
style="color:rgb(0,0,0);">C4</span><span
style="color:rgb(92,92,92);">*</span><span
style="color:rgb(218,112,214);">%s</span><span
style="color:rgb(74,85,219);">);</span>
<span style="color:rgb(0,0,0);"><br>
Z4</span> <span style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(74,85,219);">(</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R6</span> <span
style="color:rgb(92,92,92);">+</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(74,85,219);">(</span><span
style="color:rgb(0,0,0);">L2</span><span
style="color:rgb(92,92,92);">*</span><span
style="color:rgb(218,112,214);">%s</span><span
style="color:rgb(74,85,219);">)</span> <span
style="color:rgb(92,92,92);">+</span> <span
style="color:rgb(0,0,0);">C2</span><span
style="color:rgb(92,92,92);">*</span><span
style="color:rgb(218,112,214);">%s</span><span
style="color:rgb(74,85,219);">);</span><br>
<br>
<span style="color:rgb(0,0,0);"><span style="color:rgb(0,0,0);"><font
color="#009900">// System matrix (node potentials)<br>
</font></span></span><span style="color:rgb(0,0,0);">A</span>
<span style="color:rgb(92,92,92);">=</span> <span
style="color:rgb(74,85,219);">[</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R1</span><span
style="color:rgb(92,92,92);">+</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">Z1</span><span
style="color:rgb(92,92,92);">+</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R2</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(92,92,92);">-</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R2</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(188,106,143);">0</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(188,106,143);">0</span><span
style="color:rgb(0,0,0);">;</span> <span
style="color:rgb(255,170,0);">...</span>
<br>
<span style="color:rgb(92,92,92);">-</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R2</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R2</span><span
style="color:rgb(92,92,92);">+</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">Z2</span><span
style="color:rgb(92,92,92);">+</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R3</span><span
style="color:rgb(92,92,92);">+</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">Z4</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(92,92,92);">-</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R3</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(92,92,92);">-</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">Z4</span><span
style="color:rgb(0,0,0);">;</span> <span
style="color:rgb(255,170,0);">...</span>
<br>
<span style="color:rgb(188,106,143);">0</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(92,92,92);">-</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R3</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R3</span><span
style="color:rgb(92,92,92);">+</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">Z3</span><span
style="color:rgb(92,92,92);">+</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R4</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(92,92,92);">-</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R4</span><span
style="color:rgb(0,0,0);">;</span> <span
style="color:rgb(255,170,0);">...</span> <span
style="color:rgb(188,106,143);"><br>
0</span><span style="color:rgb(0,0,0);">,</span>
<span style="color:rgb(92,92,92);">-</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">Z4</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(92,92,92);">-</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R4</span><span
style="color:rgb(0,0,0);">,</span> <span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R4</span><span
style="color:rgb(92,92,92);">+</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">Z4</span><span
style="color:rgb(92,92,92);">+</span><span
style="color:rgb(188,106,143);">1</span><span
style="color:rgb(92,92,92);">/</span><span
style="color:rgb(0,0,0);">R9</span><span
style="color:rgb(74,85,219);">]</span></font><span
style="color:rgb(0,0,0);"><font face="Courier New">;</font><br>
<br>
The code goes further, but the point is I need to compute inv(A).
When multiplied by A, i.e.<br>
<br>
A*inv(A)<br>
<br>
we should get eye(A). Of course, since both A and inv(A) are
rationals, we don't get exactly eye(A), but we should anyway get
rationals very close to exact num/den cancellation on the
diagonal, and orders of magnitude less outside the diagonal.
However, it may be a bit tricky to recognize these things in a
cluttered output such s the console's, so I attempted replacing
values <br>
<br>
horner(A*inv(A), %i*2*%pi*1000)<br>
<br>
I get the following<br>
<br>
-894.25452 - 52.256219i 22.488572 - 53.503625i 13.243845 -
53.503625i 13.247333 - 53.503625i<br>
55.248428 - 45.072952i -905.16326 + 55.026768i -882.83336 -
2.1937054i -882.84581 - 2.1939413i<br>
22.907079 - 0.7699342i 877.83532 - 1.5398684i 855.90863 +
56.312642i 878.89745 - 1.5491077i<br>
-64.456283 + 100.49087i -877.84725 + 1.5398684i -878.89717 +
1.5491052i -901.87337 + 59.410849i<br>
<br>
It doesn't resemble eye(A) at all.<br>
<br>
Now I tried replacing before inverting:<br>
</span><span style="color:rgb(0,0,0);"><br>
B = horner(inv(A), %i*2*%pi*1000); <br>
C = inv(B);<br>
<br>
Now I compute<br>
<br>
B*C<br>
<br>
The result is much closer to eye(A):<br>
<br>
1. 2.845D-16 0.
0. <br>
-7.216D-16 + 1.332D-15i 1. - 3.886D-16i 1.221D-15i
1.159D-15 + 2.220D-16i<br>
1.249D-16 0. 1.
-1.366D-16i <br>
-2.776D-16 3.775D-15 - 1.665D-16i 0.
1. <br>
<br>
A last experiment: I replace in the inverse:<br>
<br>
D = horner(inv(A), %i*2*%pi*1000); <br>
<br>
Then, computing <br>
<br>
B*D<br>
<br>
yields this completely insane result:<br>
<br>
1.039D+09 + 7.084D+08i 8.026D+09 + 3.899D+09i 8.103D+09 +
3.626D+09i 8.038D+09 + 3.900D+09i<br>
4.565D+09 - 7.139D+09i 1.425D+11 - 2.222D+10i 1.411D+11 -
2.677D+10i 1.426D+11 - 2.231D+10i<br>
-3.751D+08 - 5.719D+08i -3.656D+08 - 8.904D+09i -9.535D+08 -
8.891D+09i -3.712D+08 - 8.916D+09i<br>
5.187D+09 - 6.874D+09i 1.429D+11 - 2.215D+10i 1.415D+11 -
2.671D+10i 1.431D+11 - 2.225D+10i<br>
<br>
<br>
Am I doing anything wrong?<br>
<br>
Regards,<br>
<br>
Federico Miyara <br>
<br>
</span><span style="color:rgb(0,0,0);"><span
style="color:rgb(0,0,0);"></span></span>
</body>
</html>