Line data Source code
1 : MODULE m_clebsch
2 : CONTAINS
3 160734 : REAL FUNCTION clebsch(aj,bj,am,bm,cj,cm)
4 : ******************************************************************
5 : * Program calculates Clebsch-Gordan coefficients *
6 : * See: Landau and Lifshitz, Vol.3 *
7 : * cj,cm *
8 : * C *
9 : * aj,am,bj,bm *
10 : * Written by A.Soldatov (IAE) *
11 : ******************************************************************
12 : IMPLICIT NONE
13 :
14 : REAL, INTENT (IN) :: aj,bj,am,bm,cj,cm
15 :
16 : INTEGER n,k,i,i1,i2,i3,i4,i5,i6,i7,i8,i9,i10
17 : REAL x,s,c,e
18 : REAL f(100)
19 : INTRINSIC sqrt,exp,min0,max0
20 :
21 160734 : n = 100
22 160734 : k = 0
23 160734 : x = 2.0
24 160734 : f(1)=0.e0 ; f(2)=0.e0
25 :
26 : IF ( k <= 0 ) THEN
27 160734 : k=1
28 15912666 : DO i = 3,n
29 15751932 : f(i)=f(i-1)+log(x)
30 15912666 : x=x+1.e0
31 : ENDDO
32 : ENDIF
33 :
34 160734 : i = am + bm - cm + .1e0
35 160734 : IF (( i < 0 ).OR.( i > 0 )) THEN
36 160734 : clebsch = 0.e0
37 : RETURN
38 : ENDIF
39 :
40 160734 : i1 = aj + bj - cj + 1.1e0
41 160734 : IF ( i1 <= 0 ) THEN
42 160734 : clebsch = 0.e0
43 : RETURN
44 : ENDIF
45 :
46 160734 : i2 = aj - bj + cj + 1.1e0
47 160734 : IF ( i2 <= 0 ) THEN
48 160734 : clebsch = 0.e0
49 : RETURN
50 : ENDIF
51 :
52 160734 : i3 = bj + cj - aj + 1.1e0
53 160734 : IF ( i3 <= 0 ) THEN
54 160734 : clebsch = 0.e0
55 : RETURN
56 : ENDIF
57 :
58 160734 : x = aj + bj + cj + 2.1e0
59 160734 : i4 = x
60 160734 : i = x + .6e0
61 160734 : i = i4 - i
62 160734 : IF (( i < 0 ).OR.( i > 0 )) THEN
63 160734 : clebsch = 0.e0
64 : RETURN
65 : ENDIF
66 :
67 160734 : x = aj + am + 1.1e0
68 160734 : i5 = x
69 160734 : IF (i5 <= 0 ) THEN
70 160734 : clebsch = 0.e0
71 : RETURN
72 : ENDIF
73 :
74 160734 : i = x + .6e0
75 160734 : i = i - I5
76 160734 : IF (( i < 0 ).OR.( i > 0 )) THEN
77 160734 : clebsch = 0.e0
78 : RETURN
79 : ENDIF
80 :
81 160734 : i6 = aj - am + 1.1e0
82 160734 : IF (i6 <= 0 ) THEN
83 160734 : clebsch = 0.e0
84 : RETURN
85 : ENDIF
86 :
87 160734 : x = bj + bm + 1.1e0
88 160734 : i7=X
89 160734 : IF (i7 <= 0 ) THEN
90 160734 : clebsch = 0.e0
91 : RETURN
92 : ENDIF
93 :
94 160734 : i = x + .6e0
95 160734 : i = i - i7
96 160734 : IF (( i < 0 ).OR.( i > 0 )) THEN
97 160734 : clebsch = 0.e0
98 : RETURN
99 : ENDIF
100 :
101 160734 : i8 = bj - bm + 1.1e0
102 160734 : IF (i8 <= 0 ) THEN
103 160734 : clebsch = 0.e0
104 : RETURN
105 : ENDIF
106 :
107 160734 : x = cj + cm + 1.1e0
108 160734 : i9 = x
109 160734 : IF (i9 <= 0 ) THEN
110 160734 : clebsch = 0.e0
111 : RETURN
112 : ENDIF
113 :
114 160734 : i = x + .6e0
115 160734 : i = i - i9
116 160734 : IF (( i < 0 ).OR.( i > 0 )) THEN
117 160734 : clebsch = 0.e0
118 : RETURN
119 : ENDIF
120 :
121 160734 : i10 = cj - cm + 1.1e0
122 160734 : IF (i10 <= 0 ) THEN
123 160734 : clebsch = 0.e0
124 : RETURN
125 : ENDIF
126 :
127 160734 : x = f(i1) + f(i2) + f(i3) - f(i4)
128 160734 : i = i5 - i6
129 160734 : IF ( i == 0 ) THEN
130 98540 : i = i7 - i8
131 98540 : IF ( i == 0 ) THEN
132 84480 : i = i4/2
133 84480 : i5 = i4*0.5e0 + 0.6e0
134 84480 : i = i - i5
135 84480 : IF (( i < 0 ).OR.( i > 0 )) THEN
136 160734 : clebsch = 0.e0
137 : RETURN
138 : ENDIF
139 :
140 84480 : i6 = i5 - i6 + 1
141 84480 : i7 = i5 - i8 + 1
142 84480 : i8 = i5 - i10 + 1
143 84480 : s = x*0.5e0 + f(i5) - f(i6) - f(i7) - f(i8)
144 84480 : s = exp(s)
145 84480 : i5 = i8/2
146 84480 : i6 = i8*0.5e0 + 0.6e0
147 84480 : i5 = i5 - I6
148 84480 : IF ( i5 == 0 ) THEN
149 43840 : s = 1.e0 - s - 1.e0
150 : ENDIF
151 :
152 84480 : clebsch = s*SQRT( cj + cj + 1.e0 )
153 84480 : return
154 : ENDIF
155 : ENDIF
156 :
157 76254 : x = x + f(i5) + f(i6) + f(i7) + f(i8) + f(i9) + f(i10)
158 76254 : x = x*0.5e0
159 76254 : i10 = MIN0(i1,i6,i7)
160 76254 : i2 = i1 - i5
161 76254 : i3 = i1 - i8
162 76254 : i9 = MAX0(0,i2,i3) + 1
163 76254 : i1 = i1 + 1
164 76254 : i6 = i6 + 1
165 76254 : i7 = i7 + 1
166 76254 : i8 = i9/2
167 76254 : e = 1.e0
168 76254 : i5 = i9*0.5e0 + 0.6e0
169 76254 : i8 = i8 - i5
170 :
171 76254 : IF ( i8 == 0 ) THEN
172 22866 : e = -1.e0
173 : ENDIF
174 : S = 0.e0
175 164988 : DO i = i9, i10
176 88734 : c = x-f(i)-f(i1-i)-f(i6-i)-f(i7-i)-f(i-i2)-f(i-i3)
177 88734 : s = s + e*exp(c)
178 164988 : e = 1.e0 - e - 1.e0
179 : ENDDO
180 76254 : clebsch = SQRT( cj + cj + 1.e0 )*s
181 76254 : RETURN
182 :
183 :
184 :
185 : END FUNCTION clebsch
186 : END MODULE m_clebsch
|