フォーマット文にあたるものがよくわからなくて整形出力できない・・・IPython では TeX 表記もできる様だから、数学記号的に出力させられるのではないかとも思うのだが、よく分からんｗ
In [9]: chk3(4,4,4, 4, 0, -4)
0.104297703129124 sqrt(2002)/429 ({2: 1, 7: 1}, {11: 1, 3: 2, 13: 1})
In [10]: chk3(4,4,2, 1, -1, 0)
0.144400037573994 17*sqrt(385)/2310 ({17: 2}, {11: 1, 2: 2, 3: 2, 5: 1, 7: 1})
In [11]: chk3(3.5,3,1.5, 1.5, -1, -0.5)
0 0 ({0: 1}, {})
余りチェックしていない・・・
pylab 環境下で動作。
import sympy import numpy def fa(n): if n < 0: return 0 else: return math.factorial(n) def threeJ(j1, j2, j3, m1, m2, m3): j1 = sympy.Rational(j1) j2 = sympy.Rational(j2) j3 = sympy.Rational(j3) m1 = sympy.Rational(m1) m2 = sympy.Rational(m2) m3 = sympy.Rational(m3) if m1 + m2 + m3 <> 0: return 0 if j1 + j2 < j3: return 0 if j2 + j3 < j1: return 0 if j3 + j1 < j2: return 0 J = j1 + j2 + j3 + 1 jj3 = j1 + j2 - j3 jj2 = j1 - j2 + j3 jj1 = -j1 + j2 + j3 jm1 = j1 - m1 jp1 = j1 + m1 jm2 = j2 - m2 jp2 = j2 + m2 jm3 = j3 - m3 jp3 = j3 + m3 phase = (-1)**(j1 - j2 - m3) factor2 = sympy.Rational(fa(jj3) * fa(jj2) * fa(jj1) * fa(jp1) * fa(jm1) * fa(jp2) * fa(jm2) * fa(jp3) * fa(jm3), fa(J)) factor = sympy.sqrt(factor2) # jm231 = j2 - j3 - m1 jm132 = j1 - j3 + m2 kmin = max(0, jm231, jm132) kmax = min(jj3, jm1, jp2) s = 0 for k in arange(kmin, kmax + 1): s = s + sympy.Rational( (-1)**k, fa(k) * fa(jj3 - k) * fa(jm1 - k) * fa(jp2 - k) * fa(k - jm231) * fa(k - jm132) ) return phase * factor * s def threeJ0(j1, j2, j3): J2 = j1 + j2 + j3 if J2 % 2 == 1: return 0 J = J2 / 2 k1 = J - j1 k2 = J - j2 k3 = J - j3 s = sympy.Rational(fa(2 * k1) * fa(2 * k2) * fa(2 * k3), fa(J2 + 1)) t = sympy.Rational(fa(J), fact(k1) * fa(k2) * fa(k3)) return (-1)**J * sympy.sqrt(s) * t def dicsort(d): return sorted(d.items(), key=lambda x: x[0]) def pr(i): a = sympy.factorint((i * i).p) b = sympy.factorint((i * i).q) return a, b #dicsort(a), dicsort(b) def chk0(j1, j2, j3): x = threeJ0(j1, j2, j3) print x.evalf(), x, pr(x) def chk3(j1, j2, j3, m1, m2, m3): x = threeJ(j1, j2, j3, m1, m2, m3) print x.evalf(), x, pr(x) def chk03(j1, j2, j3): print chk0(j1, j2, j3) print chk3( j1, j2, j3, 0, 0, 0) return
