On Sun, Oct 5, 2008 at 1:54 PM, glt <t.gu...@gmail.com> wrote:
>
> Hello all,
>
> I am looking for a quick and short way to implement computations with
> the Virasoro algebra in sympy. What i want to achieve is pretty basic:
> i want to define symbolic noncommutative objects L(n), where n is an
> integer, and have sympy handle products of L(n)'s as follows:
>
> * On two-terms products L(n)L(m):
> * if n<=m, do nothing
> * otherwise, replace L(n)L(m) by L(m)L(n) + [L(n),L(m)] where [ ]
> is the commutator of the Virasoro operators (expressed as another L
> operator plus possibly another term)
>
> * On products L(n)L(m)...L(p): recursively apply the two-terms rule,
> so that n, m, .. p are in increasing order
Thanks for doing it, that's very useful. What other operations would
you like to do with it?
Looking here, it's just a Lie algebra, right?
http://en.wikipedia.org/wiki/Virasoro_algebra
>
> So i started by subclassing Basic and ArithMeths in an L class, with
> is_commutative=False. Now i'm stuck the the definition of the product
> rule. I guess this should happen somewhere in Mul.flatten, but what is
> the proper way of doing it? Subclassing Mul to define my own
> multiplication? How can then i couple it to the L objects?
>
> Another thing i have thought about would be representing the product
> of L objects as a Basic object itself, with its own rules, but i'm not
> too sure about it.
>
> How would you guys do it?
> Thanks in advance
I would suggest you either subclass Mul, or just "fork it", e.g.
create VirasoroMul by copying the current Mul and then customizing it
to do what you want. If you have problems, just create a function that
takes some set of L(n)L(m)...L(p) and returns the simplifications
using the algorithm you described. After you have some working
prototype, we'll help you integrate it well in sympy, with the patch
and getting it in.
In the past, we played with having a NCMul like ginac, then we
switched to just having one Mul and the commutative=False assumption,
currently attached to the symbols directly, in the future this will be
handled by the assumptions system, just like in Mathematica. If you
have a better idea how this could/should be handled, please share it.
Looking forward,
Ondrej
But note that Mathematica handles noncommutative multiplication with a
separate "class", not with assumptions. I think a NCMul class would
be a better solution for SymPy too.
Fredrik
Right, thanks for the info.
>
>>
>> Looking here, it's just a Lie algebra, right?
>>
>> http://en.wikipedia.org/wiki/Virasoro_algebra
>>
>>
>
> Yes it is. The idea of sorting the product terms by increasing order
> is that it makes it easy to see which terms vanish on a lowest weight
> representation vector of the algebra. You just have to discard all
> products ending with an L_i with i>0.
I see. I had a similar problem with calculating cross-sections in QED
manually using creation and anihilation operators, one also needs to
employ commutation/anticommutation relations to simplify it. So SymPy
should be able to do this kind of manipulations.
>
>>
>> I would suggest you either subclass Mul, or just "fork it", e.g.
>> create VirasoroMul by copying the current Mul and then customizing it
>> to do what you want. If you have problems, just create a function that
>> takes some set of L(n)L(m)...L(p) and returns the simplifications
>> using the algorithm you described. After you have some working
>> prototype, we'll help you integrate it well in sympy, with the patch
>> and getting it in.
>
> What i have done right now is very basic but works and does the job
> (sort of). I've defined an L class which basically defines a product
> L(a)L(b)...L(p)L(q), constructed in python as L(a,b,...,p,q).
>
> Here are a few lines demonstrating how it works:
>
> ############
>
> from sympy import *
>
> class L(Basic):
> is_commutative = False
>
> def sort_virasoro(self, c):
> """Returns the product of Virasoro operators as a sum of
> products of L operators in order of increasing indices"""
> idx = self.args
>
> # Find index of first pair of operators to swap
> try:
> iswap = [p>q for p,q in zip(idx[:-1],idx[1:])].index(True)
> except:
> return self
^^^ Don't use empty "except:", because this catches all exceptions,
e.g. it hides problems. I guess you are catching the IndexError from
.index()?
>
> # Split L(a)L(b)...L(m)L(n)...L(p)L(q) into 4 parts:
> # pre = L(a)L(b)..., L(m), L(n), and post = ...L(p)L(q)
> pre = idx[:iswap]
> post = idx[iswap+2:]
> m = idx[iswap]
> n = idx[iswap+1]
>
> # This is the swapped term L(a)L(b)...L(n)L(m)...L(p)L(q)
> out = L(*(pre + (n,m) + post)).sort_virasoro(c)
>
> # This is the commutator L(a)L(b)...L(m+n)...L(p)L(q)*(m-n)
> if m-n!=0:
> out += (m-n)*L(*(pre+(m+n,)+post)).sort_virasoro(c)
>
> # This is the central term
> if m+n==0:
> if len(pre+post)==0:
> Lterm = 1
> else:
> Lterm = L(*(pre + post)).sort_virasoro(c)
>
> out += Lterm * c/12*(m**3-m)
>
> return out
>
> if __name__=="__main__":
> c = Symbol("c") # The central charge
> foo = L(-7,3,-3,-4,-5).sort_virasoro(c)
> print foo.expand()
Otherwise looks good.
>
> ############
>
>
> Note that it's not optimal in many aspects (e.g. the sorting is
> totally recursive and suboptimal) but it works.
Excellent, thanks for the code. I created an issue to get this into sympy:
http://code.google.com/p/sympy/issues/detail?id=1138
>
> I don't know how this kind of code could be integrated tighter with,
> say Mul, so that you could just define L(3)*L(2)*L(1), and have sympy
> automatically put it into this kind of "canonical" ordering for
> example. Do you have insights on this? I've tried to look into the Mul
> code but have not made it very far yet.
Couple questions first:
1) should this happen automatically, or maybe only after calling some
function on the product?
If the latter, then all that is needed is to parse the noncummative
Mul expression, plug it into your function and construct the result
back. And integrate it nicely in sympy somewhere.
2) How is this implemented in other symbolic software, e.g.
Mathematica or Maple? (If you have experience with them.)
I learned it pays off to investigate other systems first, before
trying to come up with something ourselves.
> My education is much more in physics than in programming, but I'm
> quite interested in sympy and its flexible use for this kind of
> symbolic computation.
The same here, I also studied physics and not programming. What is
your application of this? E.g. what kind of stuff do you want to
calculate? My own little aim with SymPy is to be able to do basically
any calculation that one does in theoretical physics, but we still
have a long way to go. So I am glad you need this noncommutative
stuff, because it's just waiting to get polished and more widely used,
so if we can get SymPy do what you want, it will help us a lot.
Ondrej
Not yet, but I want to play with this in sympy. So far there is only:
$ examples/qft.py
which calculates some simple scattering amplitude by direct
multiplications of the dirac matrices. Do you have some particular
example you'd like to play with? If so, I can try to implement that,
it should not be difficult.
Ondrej
I would like to implement the bosonic and fermionic operators (and
associated many particle fock states) for the many-particle
(non-relativistic) Schrodinger equation. More specifically, like
those in Fetter and Walecka many body quantum text. As a start I
would be happy with the operators for a single state. This would
involve three basic entities:
a - annihilation operator
adagger - creation operation
number state - |n>
along with either the commutation relations (bosons) or
anticommutation relations. I need to run right now, and I am just
starting to play with sympy, but there is a lot of potential for this
type of thing with sympy and it could have a significant impact.
Cheers,
Brian
> Ondrej
>
> >
>
Indeed. It should not be difficult, I'll try to give it a shot in the
comming days. If I get something, I'll let you know.
Ondrej
P.S. We are interested in any feedback you have when playing with SymPy.