CODEWITHSUNDEEP
pythonsympy
QuantumFool
QuantumFool

Reputation: 609

Why does sympy think that a function of only real variables is complex?

I've been trying to get sympy to simplify this complicated algebraic expression for me:
Expression
Where Delta, Gamma, t, and hbar are real. The code I used to generate this expression is:

from __future__ import division
from pylab import *
from sympy import *

def main():
    d, g, t = symbols("Delta Gamma t", real=True)
    hbar = symbols("hbar", positive=True, real=True)
    dg = sqrt(d**2 + g**2)
    S = simplify(Matrix([[(dg - g) / d, -(g + dg)/d], [1, 1]]))
    D = simplify(Matrix([[exp(I * t * dg / hbar), 0], [0, exp(-I * t * dg / hbar)]]))
    Sinverse = simplify(Matrix([[d / (2*dg), (g / dg + 1)/2], [-d / (2 * dg), 1/2 - g / (2*dg)]]))
    U = simplify(S * D * Sinverse)
    initial = Matrix([[1], [0]])
    later = simplify(U * initial)
    P1 = simplify(abs(later[0])**2)
    preview(P1)

if __name__=="__main__":
    main()

For reasons I don't understand, sympy refuses to recognize that the magnitude of a phase is 1 (i.e. that we can eliminate the multiplicative exponential at the far right of the expression inside absolute value signs). I tested whether sympy would simplify this scenario if we had a very simple exponential in absolute value signs. It seems not to:

>>> from pylab import *
>>> from sympy import *
>>> t = symbols("t", real=True)
>>> z = exp(t*I)
>>> abs(z).simplify()
Abs(exp(I*t))

I've tried two remedies for this problem, neither one of which has worked:

(1) Replace the absolute value signs squared with the argument times the argument's complex conjugate. When I alter the end of the main() function to read:

    ...
    P1 = simplify(later[0] * later[0].conjugate())
    preview(P1)

I get the even uglier expression:

Expression2

This did fix the above made-up scenaro:

>>> from pylab import *
>>> from sympy import *
>>> t = symbols("t", real=True)
>>> z = exp(t*I)
>>> abs(z).simplify()
Abs(exp(I*t))
>>> z * z.conjugate()
1

(2) Expand the magnitude with the keyword complex=True, and then simplify. This method also solves the made-up scenario:

>>> from pylab import *
>>> from sympy import *
>>> t = symbols("t", real=True)
>>> z = exp(t*I)
>>> abs(z).simplify()
Abs(exp(I*t))
>>> abs(z).expand(complex=True).simplify()
1

However, it fails comically for my actual expression. When I tweak my main() function to use this method:

    ...
    P1 = (abs(later[0])**2).expand(complex=True).simplify()
    preview(P1)

The program crashes. When I take off the .simplify() and change the preview to a print, I get 613.8 kB worth of text output! Looking at the first and last page of the output file, it actually seems like one really gigantic expression (i.e. I don't think that it's some stupidly long error message or something like that). No wonder the program crashed when I tried to simplify it :)

I have no idea what went wrong with the second approach, but looking at the output of the first, it seems that sympy doesn't realize that the square root of the sum of squares of real variables is also real. What should I do to fix that? Is there some argument I need to pass to some function to tell it that the square root of the sum of squares of Gamma and Delta is a real number?
Any help will be appreciated!

Upvotes: 4

Views: 1704

Answers (1)

asmeurer
asmeurer

Reputation: 91500

It looks like the development version of SymPy (soon to be released as 1.0) can do this kind of simplification:

In [1]: t = Symbol('t', real=True)

In [2]: abs(exp(-t*I))
Out[2]: 1

For your example I get

In [13]: P1
Out[13]:
                                                            2
 │                                                _________│
 │                                               ╱  2    2 │
 │                                       2⋅ⅈ⋅t⋅╲╱  Δ  + Γ  │
 │       _________   ⎛       _________⎞  ──────────────────│
 │      ╱  2    2    ⎜      ╱  2    2 ⎟          h̅        │
-│Γ + ╲╱  Δ  + Γ   - ⎝Γ - ╲╱  Δ  + Γ  ⎠⋅ℯ                  │
──────────────────────────────────────────────────────────────
                            2      2
                         4⋅Δ  + 4⋅Γ

However I checked and the latest stable version (0.7.6.1) doesn't do it, so you'll either need to wait for 1.0 or use the git version.

Upvotes: 3

Related Questions