分子軌道法 ー 計算代数幾何学からのアプローチ(4)

ただ一度の対角化で、反復計算を行わず、HF自己無同着計算を行う方法を紹介する。

環RにイデアルIが定義されている場合、その剰余環R/Iを定義することができる。前の三回で紹介した「代数的量子化学モデル」の場合にも、当然ながらこれは可能である。

例えば、HeH+のRHFモデルは多項式環$Z[x,y,e,R]$で初回に解説したイデアルIによって定義されている。R/Iの基底は単項式によって構成される。HeH+の場合、これは以下のようになる。

Q_{base}=\left(
\begin{array}{c}
ye^{3}, \\
ye^{2}, \\
ye, \\
y, \\
e^{3}, \\
e^{2}, \\
e, \\
1
\end{array}
\right)

この種の単項式基底の計算は、イデアルIのグレブナー基底を用いることで効率的に計算できる。

ここで R/Iにおいて$Q_{base}$に多項式環Rの任意の変数$p$を掛けて見よう。この結果は線形変換$M_p$で表すことができる。

p Q_{base}=M_p Q_{base}

Rに属する任意の変数p,qに対して上で定義される線形変換は可換である。

M_p M_q=M_q M_p

そして、$M_p$の固有値は、変数$p$が連立多項式I=0の根として取りうる値$W_p$である。

M_p V =W_p V

可換性から$V$は他の$M_q$の固有ベクトルでもある。連立多項式I=0の一組の根において、他の変数$q$の取る値$w_q$を、ベクトル$V$の固有値として取得できる。

w_q = (V|M_p|V)/(V|V)

HeH+のRHFモデルにおいて、$Q_{base}$と変数x,y,e,Rの積を表現する変換行列

$M_x,M_y,M_e,M_R$

は以下の通りである。これらは８×８行列である。Pythonのリスト形式で書いてある。

Mx=
[[0, 0, 0, 0, 
  11.080151325865552, -4.533316123788358,
  -0.35836894442630574, 3.2388337744224454]
 [0, 0, 0, 0,
  21.07390644331965, -5.5752716552661035,
  -5.8499651704764695, 11.541120304847897]
 [0, 0, 0, 0,
  13.19994598967836, -0.9827794471874843,
  -7.318902812741645, 9.908463596300761]
 [0, 0, 0, 0, 
  2.811521428790274, 0.22225715991443362,
  -2.008695249222711, 1.953024347868306]
 [8.931697287643605, -8.397156206184581,
  6.545956632156977, -3.5365217362474493,
  0, 0, 0, 0]
 [24.310030931798412, -21.919490812550833,
  15.652716678834729, -6.447239472143048,
  0, 0, 0, 0]
 [19.979101497085832, -16.546029346336013,
  9.929624789762936, -1.5541021159188364, 
  0, 0, 0, 0]
 [5.207839905693149, -4.059742766816225,
  2.1933186156303637, -0.19450451855497125,
  0, 0, 0, 0]]

My=
[[0, 0, 0, 0, 1, 0, 0, 0]
 [0, 0, 0, 0, 0, 1, 0, 0]
 [0, 0, 0, 0, 0, 0, 1, 0]
 [0, 0, 0, 0, 0, 0, 0, 1]
 [-5.535885815356438, 8.47906706551838, 
  -11.030456182591696, 13.225760703827165,
  0, 0, 0, 0]
 [-17.8798032915008, 25.616243135133775,
  -32.0468764996916, 37.56104565721137,
  0, 0, 0, 0]
 [-17.43234485713924, 23.374791314708975,
  -28.052043082927877, 32.302584588309806,
   0, 0, 0, 0]
 [-5.258640275660543, 6.840988600806554,
  -8.202496498256208, 9.80982277074971,
  0, 0, 0, 0]]

Me=
[[-3.6740043108251648, 1, 0, 0, 0, 0, 0, 0]
 [-4.865463887410308, 0, 1, 0, 0, 0, 0, 0]
 [-2.8627363447398095, 0, 0, 1, 0, 0, 0, 0]
 [-0.6201909048515171, 0, 0, 0, 0, 0, 0, 0]
 [0, 0, 0, 0, -3.6740043108251648, 1, 0, 0]
 [0, 0, 0, 0, -4.865463887410308, 0, 1, 0]
 [0, 0, 0, 0, -2.8627363447398095, 0, 0, 1]
 [0, 0, 0, 0, -0.6201909048515171, 0, 0, 0]]

MR=
[[1.46, 0, 0, 0, 0, 0, 0, 0]
 [0, 1.46, 0, 0, 0, 0, 0, 0]
 [0, 0, 1.46, 0, 0, 0, 0, 0]
 [0, 0, 0, 1.46, 0, 0, 0, 0]
 [0, 0, 0, 0, 1.46, 0, 0, 0]
 [0, 0, 0, 0, 0, 1.46, 0, 0]
 [0, 0, 0, 0, 0, 0, 1.46, 0]
 [0, 0, 0, 0, 0, 0, 0, 1.46]]

下の表に、$M_x$の異なる固有ベクトル$V_i, i=1,2,3,4$に対して変換行列
$M_p, p=x,y,e,R$の期待値

\frac{(V_i,M_p \cdot V_i)}{(V_i,V_i)}

を計算した結果を示す。ここでは無意味な複素数値解は省略する。前回までに解説した結果と比較して、同様に良い精度で解が求められていることがわかるだろう。

	$M_x$	$M_y$	$M_e$	$M_R$
$V_1$	-0.8015	-0.3371	-1.5997	1.46
$V_3$	0.8015	0.3371	-1.5997	1.46
$V_4$	0.6039	-1.1152	-0.5379	1.46

繰り返していうと、反復計算を使わず、ただ一回の行列対角化で自己無同着計算が行える、という点が、今回紹介したアルゴリズムの特徴である。