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).

Aaron Meurer


-- 
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.

Reply via email to