Hi,
to improve computations speed in case of symbolic matrices (as g_dd
is) it would be convenient to use M.berkowitz_det() rather that
M.det():
[EMAIL PROTECTED] ~/hg/sympy $ SYMPY_USE_CACHE=no ./bin/isympy -q
Python 2.5.2 console for SymPy 0.5.15-hg (cache: off)
In [1]: var("t omega")
Out[1]: (t, ω)
In [2]: g_dd = Matrix([(-1+omega*(x**2+y**2),2*omega*y,-2*omega*x,0),
...: (2*omega*y, 1, 0, 0), (-2*omega*x, 0, 1, 0), (0, 0, 0, 1)])
In [3]: time a = g_dd.berkowitz_det()
CPU times: user 0.26 s, sys: 0.00 s, total: 0.27 s
Wall time: 0.28 s
In [5]: time b = g_dd.det()
CPU times: user 6.99 s, sys: 0.06 s, total: 7.04 s
Wall time: 7.19 s
In [7]: a == b
Out[7]: True
Mateusz
2008/5/27 Ondrej Certik <[EMAIL PROTECTED]>:
> 2008/5/27 Friedrich Hagedorn <[EMAIL PROTECTED]>:
>> 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.
>
> Wow, indeed, that helps a lot!
>
> With my original problem, I am getting this:
>
>
> In [18]: e = adj(M).T * (1/M.det())
>
> In [19]: one = M*e
>
> In [20]: one[0,0]
> Out[20]:
> 2
> 1 ω*x
> - ──────────────────────────────────── + ────────────────────────────────────
> 2 2 2 2 2 2 2 2 2 2 2 2
> -1 + ω*x + ω*y - 4*ω *x - 4*ω *y -1 + ω*x + ω*y - 4*ω *x - 4*ω *y
>
> 2 2 2
> ω*y 4*ω *x
> + ──────────────────────────────────── - ────────────────────────────────────
> 2 2 2 2 2 2 2 2 2 2 2 2
> -1 + ω*x + ω*y - 4*ω *x - 4*ω *y -1 + ω*x + ω*y - 4*ω *x - 4*ω *y
>
> 2 2
> 4*ω *y
> - ────────────────────────────────────
> 2 2 2 2 2 2
> -1 + ω*x + ω*y - 4*ω *x - 4*ω *y
>
> In [21]: simplify(one[0,0])
> Out[21]: 1
>
>
> So first the result looks much nicer and second simplify can handle
> that. Very cool, I created a new issue for your code to get in:
>
> http://code.google.com/p/sympy/issues/detail?id=866
>
> Ondrej
>
> >
>
--~--~---------~--~----~------------~-------~--~----~
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
-~----------~----~----~----~------~----~------~--~---
