On Thu, Jun 16, 2011 at 10:20 PM, Aaron S. Meurer <[email protected]> wrote: > > On Jun 16, 2011, at 11:17 PM, Ondrej Certik wrote: > >> 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? >> >> Ondřej >> > > I'm just noting that cos(pi/2 - beta) is simpler than sqrt(1 - cos(beta)**2).
Ah, I see. I am now putting all the equations into my notes, so that I can reference them. The thing is, that the original equation is a symbolic derivative equation with cos(beta) in it to various powers. After everything is simplified, you get sqrt(1 - cos(beta)**2) for one special case of the parameters (j=1,m=1,m'=0), and the canonical way to write it is using sin(beta). The original formula is actually one of many, how to calculate it (there are also other approaches, which contain different functions). The end result should however be the same always. 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.
