On Tue, May 27, 2008 at 12:07:33PM +0200, Friedrich Hagedorn wrote:
>
> On Tue, May 27, 2008 at 11:42:21AM +0200, Ondrej Certik wrote:
> > In [4]: g_dd
> > Out[4]:
> > ⎡ ⎛ 2 2⎞ ⎤
> > ⎢-1 + ω*⎝x + y ⎠ 2*ω*y -2*ω*x 0⎥
> > ⎢ ⎥
> > ⎢ 2*ω*y 1 0 0⎥
> > ⎢ ⎥
> > ⎢ -2*ω*x 0 1 0⎥
> > ⎢ ⎥
> > ⎣ 0 0 0 1⎦
[...]
> > In [3]: g_uu = g_dd.inv()
[...]
> I guess the matrix inversion is done by the gauss algorithm. But my
> experience is that the gauss isnt so good for symbolic calculation.
>
> Try to calc the matrix inverse with the determiants
>
> M^-1 = 1/det(M) * (M_adj).T
>
> with
>
> M_adj := (-1)^(i+j) * det(M_ij)
>
> M_ij := Matrix M without row i and without column j
>
> and the evaluations of the determinat of your matrix is easy. This is
> also (shortly) discribed in
So here is my quick solution with the attached file inversion.py:
In [1]: M = Matrix([(z,y,-x,0),(y, 1, 0, 0), (-x, 0, 1, 0), (0, 0, 0, 1)])
In [2]: M
Out[2]:
⎡ z y -x 0⎤
⎢ y 1 0 0⎥
⎢-x 0 1 0⎥
⎣ 0 0 0 1⎦
In [3]: run inverse.py
In [4]: adj(M).T
Out[4]:
⎡ ⎤
⎢ 1 -y x 0⎥
⎢ 2 ⎥
⎢ -y z - x -x*y 0⎥
⎢ 2 ⎥
⎢ x -x*y z - y 0⎥
⎢ 2 2⎥
⎣ 0 0 0 z - x - y ⎦
In [5]: 1/M.det()
Out[5]:
1
───────────
2 2
z - x - y
If you multiply Out[4] * Out[5] you get the inversion matrix and can test
In [15]: (M*adj(M).T)*(1/M.det())
Out[15]:
⎡1 0 0 0⎤
⎢0 1 0 0⎥
⎢0 0 1 0⎥
⎣0 0 0 1⎦
I hope that helps.
By, Friedrich
--~--~---------~--~----~------------~-------~--~----~
You received this message because you are subscribed to the Google Groups
"sympy" group.
To post to this group, send email to [email protected]
To unsubscribe from this group, send email to [EMAIL PROTECTED]
For more options, visit this group at http://groups.google.com/group/sympy?hl=en
-~----------~----~----~----~------~----~------~--~---
from sympy import *
def adj(M):
rows, cols = M.shape
if rows != cols: return M, 'is not a quadratic matrix'
ret = zero(rows)
for i in range(0, rows):
for j in range(0, cols):
A = M[:,:]
A.row_del(i)
A.col_del(j)
ret[i,j] = (-1)**(i+j)*A.det()
return ret
def inv(M):
return 1/M.det()*adj(M).T
