-
Notifications
You must be signed in to change notification settings - Fork 42
/
wignerD.m
executable file
·80 lines (70 loc) · 3.72 KB
/
wignerD.m
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
function [D_l, d_l, R_l] = wignerD(l, alpha, beta, gamma)
%WIGNERD Returns the Wigner D-matrix and d-matrix for SH rotation
%
% WIGNERD computes the rotation matrices for rotation of functions in the
% spherical harmonic domain, as given directly by Wigner for complex SHs.
% Since it much faster to compute the rotation matrices ith recursive
% formulas, these are included here mostly for comparison.
%
% Main reference:
% D. A. Varshalovich, A. N. Moskalev, and V. K. Khersonskii.
% Quantum theory of angular momentum. World Scientific Pub., 1988.
% p.77 - 4.3(5)
%
% Inputs: l degree (or band) index
% alpha, beta, gamma (z, y', z'')-rotation Euler angles
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Archontis Politis, 6/5/2015
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% analytic expressions of small-d for order l = 1, 2, included for
% comparison (taken from Varshalovich etal., "Quantum Theory of Angular
% Momentum", 1988, pg.119
% note that the row and column order is inverted from -m:m compared to the
% book (m:-1:-m)
d1 = @(bta) [ (1+cos(bta))/2 -sin(bta)/sqrt(2) (1-cos(bta))/2;
sin(bta)/sqrt(2) cos(bta) -sin(bta)/sqrt(2);
(1-cos(bta))/2 sin(bta)/sqrt(2) (1+cos(bta))/2];
d2 = @(bta) [ (1+cos(bta))^2/4 -(1+cos(bta))*sin(bta)/2 sqrt(6)*(sin(bta))^2/4 -(1-cos(bta))*sin(bta)/2 (1-cos(bta))^2/4;
(1+cos(bta))*sin(bta)/2 (1+cos(bta))*(2*cos(bta)-1)/2 -sqrt(3/2)*cos(bta)*sin(bta) (1-cos(bta))*(2*cos(bta)+1)/2 -(1-cos(bta))*sin(bta)/2;
sqrt(6)*(sin(bta))^2/4 sqrt(3/2)*sin(bta)*cos(bta) 3*(cos(bta))^2/2-1/2 -sqrt(3/2)*cos(bta)*sin(bta) sqrt(6)*(sin(bta))^2/4;
(1-cos(bta))*sin(bta)/2 (1-cos(bta))*(2*cos(bta)+1)/2 sqrt(3/2)*cos(bta)*sin(bta) (1+cos(bta))*(2*cos(bta)-1)/2 -(1+cos(bta))*sin(bta)/2;
(1-cos(bta))^2/4 (1-cos(bta))*sin(bta)/2 sqrt(6)*(sin(bta))^2/4 (1+cos(bta))*sin(bta)/2 (1+cos(bta))^2/4];
% precompute factorials
fctrl = factorial(0:2*l).';
% if beta close to zero force small d to identity matrix
D_l = zeros(2*l+1);
if beta <= eps;
d_l = eye(2*l+1);
D_l = diag(exp(-1i*(alpha+gamma)*(-l:l)));
else
d_l = zeros(2*l+1);
for m = -l:l
for n = -l:l
% Varshalovich, Eq.4.3.1(5)
k = (max(0,n-m):min(l+n,l-m)).';
m1_t = (-1).^k;
fact_t = sqrt(fctrl(l+n+1)*fctrl(l-n+1)*fctrl(l+m+1)*fctrl(l-m+1)) ./ ...
(fctrl(l+n-k+1).*fctrl(l-m-k+1).*fctrl(k+m-n+1).*fctrl(k+1));
cos_beta = cos(beta/2).^(2*l+n-m-2*k);
sin_beta = sin(beta/2).^(2*k+m-n);
d_l_mn = (-1)^(m-n) * sum(m1_t.*fact_t.*cos_beta.*sin_beta);
d_l(n+l+1,m+l+1) = d_l_mn;
D_l(n+l+1,m+l+1) = exp(-1i*alpha*m)*d_l_mn*exp(-1i*gamma*n);
end
end
end
% Interestingly, if the Wigner-D are used directly they result in an active
% body rotation , instead of rotation of the coordinate system (!). Since
% most operations here are defined in terms of coordinate system rotations,
% the proper result is given by the conjugate of the Wigner-D. Compare for
% example with the rotation matrices using the recursive formula function.
% This is also the same as using e^(ima) * d^l_mn(b) * e^(ing) for the
% entries of the Wigner-D, something that can be found in the literature as
% well as the above formula.
% Compute rotation matrices up to l
R_l = conj(D_l);
end