/tex/edu/kaku/numangm.tex

Clebsch-Gordan 係数その他の数値計算
ここでは、Clebsch-Gordan 係数と、Racah 係数の数値計算の手法を記録しておく。
一般的な参考書として、D.A. Varshalovich 達の ‘Quantum Theory of Angular momentum’ を挙げておこう。
１) Clebsch-Gordan 係数
j1 + j2 = j3 という角運動量の関係を想定する。方位量子数は、mi (i = 1, 2, 3) を出来るだけ使用する。
１-１) 最も単純な方法は、Wigner や Racah に依る一般式を直接計算する事である。この場合、階乗計算による数値の桁溢れを如何にして防ぐかが問題となる。
桁溢れが起こるとそこで計算をあきらめるとか、もっと桁数の大きな計算機 (ソフト) を利用するという以外に、次の３個の対策があり得る。

1-1-1 全ての数値を素数の積として表現する。やってみると、かなり大きな素数が登場するので、素数表の準備と維持が結構面倒臭い。大体は小さな素数の積で表現できるのに、例外的に大きな素数が登場してピンチに陥る事があるという意味である。

M. Rodberg 著 (だったかな？) の 3-j and 6-j symbols という本では、この例外的な素数を欄外かどこかへ移して記述をしていたと思う。

(j1 + j2 + j3 + 1)! が登場するから、この数値以下の全ての素数が登場するから当然の事である。
厳密な計算が出来ると言う特徴を有する。

1-1-2 階乗部分を、次の式で置き換える。
nX
n! = exp (log (i))
i=1

この手法は、かなり大きな自然数に対しても有効であるという長所がある。一方、精度が悪く、計算時間がかかるという短所がある。

1-1-3 階乗を以下の数値で置き換える。
n! ==> n!/B n
この置き換えで階乗の値は変化するが、Wigner や Racah の式を見ると、結果には定数を除いて影響が無い事が分かる。単精度では、B = 16、倍精度では B = 256 を採用するのが妥当だと思う。B をこのように選ぶと、階乗計算に際し、２進数を基調とする計算機では、有効数字が劣化しないという特徴がある。
階乗の計算を回避出来たとして、もう一度一般式を見ると、一般式の和の項は項毎に符号が反転すると言う欠点が目に付く。この事は、和に登場する最大値付近での相殺により有効数字の桁数が制限される事を意味する。従って、階乗に関する問題を乗り切ったからといって、この式を使うプログラムで精度が確保出来ていると言う保証は無い。
1

相隣る項の和 (実際は差) を明示的に計算した結果をプログラムに組み込むのが一つの対策と考えるかもしれないが、相隣る項の和の符号が一定であると言う保証は無い。
計算したいのは、ユニタリー変換の係数である事を思い出すと、その絶対値は１を越える事は無い。そのようなある意味で小さな数値を計算するのに、大きな数値が登場するのはどこかおかしいと思うのが自然であろう。
独立変数の広い範囲にわたり計算可能とする為には、漸化式を利用するのが賢明であるように思える。計算したい対象の性質を考えると、二つの手法があり得る。
１-２) m2 に関する漸化式を利用する
ここでの議論では、j1 , j2 , j3 , m3 は固定して考えるので、< j1 m2 j2 m2 |j3 m3 >= (m1 , m2 )という略号を使用する。

C 係数の定義式に
j23 = (j1 + j2 )2 = j21 + j2 + j1+ j2− + j1− j2+ + 2 j1z j2zを作用させる事により、以下の漸化式を得る。

{j3 (j3 +1)−j1 (j1 +1)−j2 (j2 +1)−2m1 m2 }(m1 , m2 ) = j1− j2+ (m1 +1, m2 −1)+j1+ j2− (m1 −1, m2 +1)pここで、j ± ≡ (j ∓ m)(j ± m + 1) を定義した。
補助的知識として、j2 の磁気量子数が最大である場合、(m1 , m2 = j2 ) > 0 である。C 係数の符号に関する Condon-Shortley の規約は、この１行に集約出来る様である。
更に、当然の知識として次の和則がある。

X(m1 , m2 )2 = 1

m1 +m2 =m3

もう一つの当然の知識として、(m1 , m2 > j2 ) = 0 がある。
基本的なアルゴリズムは、以下の様にまとめられる。

1-2-0 (m3 − (j2 + 1), j2 + 1) = 0, (m3 − j2 , j2 ) = 1 を仮定する。(m3 − j2 , j2 ) は正ではあるが、絶対値は当面不定である。

1-2-1 上で与えた漸化式を用い、m2 の大きな方から小さい方へ向けて、物理的に０とは異なる値を持ち得る全ての m2 に対する (m1 , m2 ) の値を計算する。
この値は、絶対値が正である不定定数を除いて、相対的には正しい。

1-2-2 ユニタリー性を利用して、不定定数の計算し、これまで計算した数値に掛ける。
これで原理的には任意の C 係数が計算できる訳だが、実用上の注意がある。この注意は、
次の 1-3) の事項に対しても該当する。
この計算を実行してみると、C 係数の絶対値は、計算の中間あたりで大きくなり、両端では小さな値をとる。即ち、計算の前半では精度が確保されるが、後半では精度は急速に落ち2て来る。前半を計算した後で、今度は (m3 + j2 , −j2 ) から、大きな m2 に向かって計算を開始し、中間点で両計算を接続する。
更に、注意を追加しよう。接続点での C 係数がたまたま０になる事がある。パリティーC 係数の存在を考えると、ありそうなトラブルである。接続は、必ず２点での C 係数の値を用いて行うべきである。
１-３) j3 に関する漸化式を利用する
ここでは、j1 , j2 , m1 , m2 , m3 は固定して考えるので、C 係数は C(j3 ) と略記しよう。

j1z の行列要素と Wigner-Eckart の定理を利用すると、以下の j3 に関する漸化式が導ける。

A(j3 ) C(j3 ) = B(j3 ) C(j3 − 1) + B(j3 + 1) C(j3 + 1)
A(j3 ) = m1 − m3

j1 (j1 + 1) − j2 (j2 + 1) + j3 (j3 + 1)
2 j3 (j3 + 1)

sB(j3 ) =(j32 − m23 ){j32 − (j1 − j2 )2 }{(j1 + j2 + 1)2 − j32 }(2j3 )2 {(2j3 )2 − 1}計算を開始するには、sign( C(j3 = j1 + j2 ) ) = (−1)j1 +j2 −j3 という知識が必要である。
後は、1-2) の後半に述べた実用上の注意を守れば良い。

1-2) と 1-3) のどちらを利用すべきか？ という事は、自分が利用したい C 係数の取り得る範囲を考えれば判断が出来るだろう。
もう一つの実用上のコメントを述べよう。磁気量子数が全て０であるような C 係数の計算は、< `f ||Yλ ||`i > という換算行列要素の計算上必要である。この場合、パリティーの保存則があるから、ここで与えた漸化式は比の形になる。即ち、桁落ちによる精度劣化という不安要素を無視して良い。次の計算手順が考えられる。

1-3-0 < j1 0, j1 0|0 0 >, < j1 0 (j1 − 1) 0|1 0 > の知識を仮定する。但し、j1 , j2 の内、大きい方を j1 と置き直しておく。
1-3-1 j3 を０又は１から、実際に計算したい j3 又はその最小値まで、漸化式を用いて増やす。

1-3-2 漸化式を用いて、j2 を実際に計算したい値まで小さくする。
1-3-3 漸化式を用いて、j3 を必要な値まで持ち上げる。
２) Racah 係数の数値計算

Igal Talmi に依ると、W 係数の産みの親である G. Racah さえも、W 係数ではなく、6-j 係数を使用しているようである。

Racah による、一般式を利用する場合の手法は、C 係数の計算と同じである。漸化式を利用する手法を記録しておく。
3

出発点となるのは、Biedenharn-Elliot の恒等式である。

()()
j1 j2 j3
j1 j2 j3
`1 `2 `3
`01 `02 `03()()()
000000Xj```j```j123231312=(−)φ(λ)
0λ```λ```32121 λ3λここで位相は、φ(λ) = j1 + j2 + j3 + `1 + `2 + `3 + `01 + `02 + `03 + λBiedenharn と Rose はこの式が W 係数の漸化式の一般形であると述べ、`01 = 1/2, `02 =`2 + 1/2, `03 = `3 + 1/2 の場合に W 係数が満足すべき漸化式を与えている。
一つの変数が１/２である場合の次の 6-j 係数の知識が前提となる。
s()(j1 + j2 + j3 + 2)(j1 + j2 − j3 + 1)
j1j2
j3= (−1)j1 +j2 +j3 +1(2j1 + 1)(2j1 + 2)(2j2 + 1)(2j2 + 2)
j2 + 1/2 j1 + 1/2 1/2(j1j2
j3
j2 + 1/2 j1 − 1/2 1/2s)= (−1)j1 +j2 +j3(j3 + j1 − j2 )(j2 + j3 − j1 + 1)
2j1 (2j1 + 1)(2j2 + 1)(2j2 + 2)

Biedenharn-Elliot の恒等式で、`01 = 1/2, `02 = j3 ± 1/2, `03 = j2 ± 1/2 を代入すると、
λ = `1 ± 1/2 という２項だけの和となり、４個の漸化式が登場する。パラメータの値を調節し、j2 , j3 , `1 , `2 , `3 を固定し、j1 だけの漸化式を導く事が出来るようである。
j1 E(j1 + 1){j1 + 1} + F (j1 ){j1 } + (j1 + 1)E(j1 ){j1 − 1} = 0
ここで、({j1 } ≡j1 j2 j3
`1 `2 `3)©
ªE(j1 ) ≡ [j12 − (j2 − j3 )2 ][(j2 + j3 + 1)2 − j12 ][j12 − (`2 − `3 )2 ][(`2 + `3 + 1)2 − j12 ]F (j1 ) ≡ (2j1 + 1) {P (j1 )[−P (j1 ) + P (j2 ) + P (j3 )] + P (`1 )[P (j1 ) + P (j2 ) − P (j3 )]+P (`3 )[P (j1 ) − P (j2 ) + P (j3 )] − 2P (j1 )P (`1 )}ここで、P (j) ≡ j(j + 1)
符号の規則は、以下の通りの様である。

sign({j1 }) = (−1)j2 +j3 +`2 +`3
実用上の注意を２つ。

4

j1 の最大値と最小値の双方から、漸化式を用いて j1 の各種の値に対して j1 を計算し、中央付近で規格化するのは、C 係数の場合と同じである。E(j), F (j) には大きな数の２乗の差が登場し得る。ここで情報が相殺する可能性がある。

j(j + 1) − `(` + 1) = (j − `)(j + ` + 1)と変形し、右辺の和や差の部分を整数演算し、その結果を浮動小数点数に変換しておくのが良いだろう。

３) 9-j 係数の数値計算

9-j 係数の漸化式を用いた計算の経験は無いので、書く事は出来ない。次の情報だけを残しておこう。参考書を見ると、次の式が書いてある。

jjJ1212j3 j4 J34 JJ 
13 J24=X(2J 0(−1)

0(2J + 1)

J0=

X(0(−1)2J (2J 0 + 1)

J0=

X

(2J 0(−1)(2J 0 + 1)

J0

j1 j3 J13J24 J J 0
J12 J34 J
J13 J24 J 0
j2 j4 J24
J J13 J 0)()()(

j2 j4 J24
j3 J 0 J34
j1 j3 J13J34 J 0 j4
J12 J34 J
j4 J 0 j3)(

J12 J34 J
J 0 j1 j2)()(

j2 j4 J24J 0 J12 j1
j1 j3 J13J 0 j2 J12)))２) に登場した W 係数の漸化式を利用すると、一連の W 係数が同時に計算出来るから、9-j 係数の計算には有利に思える。
実際に計算する場合には、どの式を利用するのが有利なのだろう？
和の範囲が最も狭い式を採用すると言う指針でプログラムを組んでみると、一つの式を利用する場合よりもそれなりに計算時間を節約できる。この場合、計算精度での優劣は経験していない。

5