On Thu, Jun 16, 2011 at 10:06 PM, Aaron Meurer <[email protected]> wrote: > On Thu, Jun 16, 2011 at 11:00 PM, Ondrej Certik <[email protected]> > wrote: >> On Thu, Jun 16, 2011 at 8:34 PM, Ondrej Certik <[email protected]> >> wrote: >>> On Thu, Jun 16, 2011 at 6:38 PM, Ondrej Certik <[email protected]> >>> wrote: >>>> On Thu, Jun 16, 2011 at 12:04 AM, Sean Vig <[email protected]> wrote: >>>>> Hi all, >>>>> In working on some stuff with spin states, I ran into some problems with >>>>> the >>>>> current implementation of the Wigner small-d matrix, Rotation.d in >>>>> sympy.physics.quantum.spin. I had written methods to change bases using >>>>> the >>>>> Wigner D-function [0] and in testing decided to try >>>>>>>>qapply(JzBra(1,1)*JzKet(1,1).rewrite('Jx')) >>>>> which should be 1, as the spin ket is rewritten in the Jx basis and then >>>>> back to the Jz basis to apply the innerproduct, but I have that it gives >>>>> 0. >>>>> I traced this back to a bug in the Rotation.d function, which currently >>>>> has >>>>> an open issue [1]. For the Wigner D-function and the small d-matrix, the >>>>> conventions laid out in Varshalovich "Quantum Theory of Angular Momentum". >>>>> It seems the d-matrix is fine for positive values of angle argument, but >>>>> does not obey the symmetry d(j,m,mp,-beta)=(-1)**(m-mp)*d(j,m,mp,beta) and >>>>> does not agree with the tables in Varshalovich. Those terms that fail for >>>>> the j=1 case are in an XFAIL test in my branch. What is odd is that when I >>>> >>>> >>>> I think there is a bug in the code. See my comment here: >>>> >>>> http://code.google.com/p/sympy/issues/detail?id=2423#c3 >>>> >>>> The correct results for general beta are: >>>> >>>> d(1, 0, 1, beta) = sin(beta)/sqrt(2) >>>> d(1, 1, 0, beta) = -sin(beta)/sqrt(2) >>>> >>>> However, sympy gives: >>>> >>>>>>> Rotation.d(1,0,1,beta) >>>> ⎽⎽⎽ >>>> ╲╱ 2 ⋅(2⋅cos(β) + 2) >>>> ──────────────────── >>>> 4 >>>>>>> Rotation.d(1,1,0,beta) >>>> ⎽⎽⎽ >>>> -╲╱ 2 ⋅(cos(β) + 1) >>>> ─────────────────── >>>> 2 >>>> >>>> Which is wrong (it looks quite ok, that it changes sign, as it should, >>>> but something is wrong with the cos(beta) thing). The sympy code >>>> implements the "Eq. 7 in Section 4.3.2 of Varshalovich." >>>> >>>> >>>>> ran the equation used to define the d-matrix through Mathematica, I got >>>>> results that agreed with the sympy output, so the problem may be in the >>>>> equation and not a bug in the code. If anyone could take a look at that, >>>>> I'd >>>>> appreciate it. >>>> >>>> Which *exact* equation did you run through Mathematica? Eq. 7 in >>>> section 4.3.2? What *exactly* did you get? Did you get an expression >>>> involving cos(beta), just like sympy above? Can you paste here the >>>> Mathematica code? I'll run it with my Mathematica, to verify, that we >>>> didn't make a mistake. >>>> >>>> >>>> Now we just need to systematically look at it, and nail it down. We >>>> need to get the correct expressions: >>>> >>>> d(1, 0, 1, beta) = sin(beta)/sqrt(2) >>>> d(1, 1, 0, beta) = -sin(beta)/sqrt(2) >>>> >>>> one way or the other, for general beta. Then things will start to work. >>>> >>>> Sean, let me know if you have any questions to the above. >>>> >>>> Once we fix this, we'll move on to the other problems you raised. >>>> >>>> I am CCing Brian, who implemented that code in SymPy. However, we >>>> should be able to fix this ourselves Sean ---- all we need to do is to >>>> take the eq. (7), and see what expression we get for J=1, M=1, M'=0, >>>> just put it there by hand (don't use mathematica) and see what you >>>> get. >>>> >>>> Post here your results, and I'll verify them with my independent >>>> calculation and we'll nail it down. >>> >>> Ok, I think that I have nailed it down. There are actually several problems: >>> >>> 1) First of all, this is the range for beta: >>> >>> 0 <= beta <= pi >>> >>> see Varshalovich, page 74. >>> >>> 2) See the attached screenshot of my calculation, which calculates >>> d(j, 0, 1, beta) from the equation (7) on page 75, and shows, that it >>> is equal to >>> sin(beta)/sqrt(2), consistent with Varshalovich result in the table 4.4. >>> >>> >>> As such, from 2) it follows, that the sympy result (see my previous >>> email) is *wrong*, as can be checked by substituting beta = pi/3 and >>> checking against sin(beta)/sqrt(2). For beta=pi/2, it happens to be >>> equal, but that is pure accident. >>> >>> From 1) it follows, that you can't call it for beta=-pi/2, you need to >>> adjust alpha and gamma instead. >> >> Actually, I wasn't very clear here. You use the formula (7) to >> calculate the wigner d function for *any* J, M, M', for 0<=beta<=pi. >> In particular, you get the formulas: >> >> d(1, 1, 0, beta) = -sin(beta)/sqrt(2) >> d(1, 0, 1, beta) = +sin(beta)/sqrt(2) >> >> And I am stressing here that this is *only* valid for 0 <= beta <= pi. >> You can equivalently write the sin(beta) using cos as sin(beta) = >> sqrt(1-cos**2(beta)), which is valid in this whole domain. You >> actually only get expressions involving cos(beta) from (7), but since >> we are on this restricted domain, it doesn't matter. > > I admit I don't have a clue what you guys are talking about, but
I wasn't super confident in this either, but I am becoming now. I will put my understanding of everything + derivations into my http://theoretical-physics.net/ notes. So far there is some info here: http://en.wikipedia.org/wiki/Wigner_D-matrix but it's not very informative, at least I didn't understand from it what is going on when I read it couple weeks ago. > wouldn't it be easier to use cos(pi/2 - beta)? cos(pi/2 - beta) = sin(beta) but how would it make things easier? 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.
