On Thu, Jun 16, 2011 at 11:00 PM, Ondrej Certik <ondrej.cer...@gmail.com> wrote:
> On Thu, Jun 16, 2011 at 8:34 PM, Ondrej Certik <ondrej.cer...@gmail.com> 
> wrote:
>> On Thu, Jun 16, 2011 at 6:38 PM, Ondrej Certik <ondrej.cer...@gmail.com> 
>> wrote:
>>> On Thu, Jun 16, 2011 at 12:04 AM, Sean Vig <sean.v....@gmail.com> 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
wouldn't it be easier to use cos(pi/2 - beta)?

Aaron Meurer

>
> Now when you want to calculate the wigner d function for parameters
> outside of this domain, you need to use the symmetries, see the
> section 4.4.
>
> Example: calculate d(1, 1, 0, -pi/2), as you needed above. Then we use
> this symmetry: d(j, m, mp, -beta) = d(j, mp, m, beta), and write:
>
> d(1, 1, 0, -pi/2) = d(1, 0, 1, pi/2) = sin(pi/2)/sqrt(2) = 1/sqrt(2)
>
> and that's it.
>
> All should be clear now.
>
> Ondrej
>
> --
> You received this message because you are subscribed to the Google Groups 
> "sympy" group.
> To post to this group, send email to sympy@googlegroups.com.
> To unsubscribe from this group, send email to 
> sympy+unsubscr...@googlegroups.com.
> For more options, visit this group at 
> http://groups.google.com/group/sympy?hl=en.
>
>

-- 
You received this message because you are subscribed to the Google Groups 
"sympy" group.
To post to this group, send email to sympy@googlegroups.com.
To unsubscribe from this group, send email to 
sympy+unsubscr...@googlegroups.com.
For more options, visit this group at 
http://groups.google.com/group/sympy?hl=en.

Reply via email to