5  入力ファイルを物理の供述書として読む

入力ファイルの一行一行は、計算機への命令であると同時に、こちらが置いた物理の仮定の供述書でもある。
その供述を、読み落とさずに読めるか。

この章のゴール:Co の SCF 入力ファイルを冒頭から末尾まで一行ずつ「物理の供述書」として読み、&control&system&electrons の各行、PAW 擬ポテンシャル、半芯電子、ecutwfc/ecutrhoK_POINTSoccupations/smearing/degaussnspin/starting_magnetization が、それぞれどんな仮定を置いているのかを自分の言葉で述べられるようになる。深掘りは各章に預け、ここでは全体を「構造・電子・磁性・数値条件」の四色で見渡す。

前章で nspin を 2 に上げ、コバルトに磁石になりうるかを問えるようにした。だがその一行は、入力ファイルという供述書の一節にすぎない。この章では、smoke 束の強磁性(FM)入力 co.scf.in を冒頭から末尾まで読み下す。見取り図を先に渡しておく。入力ファイルのどの行も、必ず構造・電子・磁性・数値条件のどれか一つの役割を担っている。この四色の地図を握れば、入力ファイルは命令の羅列ではなく、「私はコバルトをこう仮定して計算した」という調書として読めるようになる。

まず全体を眺める。これが smoke 束で実際に走らせた FM 入力である。

&control
   calculation = 'scf'
   prefix      = 'co_hcp_fm'
   outdir      = './tmp'
   pseudo_dir  = './pseudo'
/
&system
   ibrav       = 4
   celldm(1)   = 4.738
   celldm(3)   = 1.6228
   nat         = 2
   ntyp        = 1
   ecutwfc     = 80
   ecutrho     = 640
   nspin       = 2
   starting_magnetization(1) = 0.7
   occupations = 'smearing'
   smearing    = 'mv'
   degauss     = 0.02
/
&electrons
   conv_thr    = 1.0d-8
   mixing_beta = 0.3
/
ATOMIC_SPECIES
   Co  58.933  Co.pbe-spn-kjpaw_psl.0.3.1.UPF
ATOMIC_POSITIONS (crystal)
   Co  0.000000  0.000000  0.000000
   Co  0.333333  0.666667  0.500000
K_POINTS automatic
   12 12 8  0 0 0

三つのブロック(&control, &system, &electrons)と、三つのカード(ATOMIC_SPECIES, ATOMIC_POSITIONS, K_POINTS)でできている。ブロックは & で開き / で閉じる入れ子で、上から順に「何を計算するか」「どんな系か」「どう解くか」を述べる。順に読んでいこう。

1 &control:何を計算するか

最初のブロックは、計算の種類と入出力の段取りを宣言する。

calculation = 'scf' が、このブロックでいちばん重い一行である。SCF(self-consistent field、自己無撞着場)とは、与えた原子配置に対して、電子状態と全エネルギーを「仮定と結果がつじつまの合う」ところまで反復して求める計算を指す。つまりこの一行で、私たちは「いまから、固定したこの原子配置での電子状態と全エネルギーを求める」と宣言している。原子を動かして構造を最適化する 'relax' でも、状態密度の前処理 'nscf' でもない。本書の計算の大半は、この 'scf' を土台に積み上がっている。

残る prefixoutdirpseudo_dir は、出力ファイルにつける名札と、作業ファイル・擬ポテンシャルの置き場所を指すだけで、物理の仮定はほとんど含まない。&control は段取りのブロックであって、供述の中身はまだ薄い。供述の中核は、次のブロックにある。

2 &system:構造・電子・磁性が同居する中核

&system は、この入力で最も多くの仮定が同居する場所である。結晶の形、電子の解像度、金属の扱い、磁性の自由度――性質の異なる四種類の前提が、一つのブロックに肩を並べる。だからこそ、ここを役割ごとに仕分けて読むことが、入力を読む技術の中心になる。実値で順に読む。

構造を述べる行ibrav=4 は六方(hexagonal)格子を選び、celldm(1)=4.738(bohr 単位=2.507 Å=実験値)が底面の辺長 \(a\) を、celldm(3)=1.6228 が軸比 \(c/a\) を与える。nat=2 は単位胞内の原子数、ntyp=1 は原子種の数である。これらは第2章で手渡した hcp 構造そのものを、数件の数として「置く」行にあたる。

電子の解像度を述べる行ecutwfc=80(Ry)と ecutrho=640(Ry)は、電子の波をどこまで細かく表すかの上限を決める。ecutwfc が波動関数の、ecutrho が電荷密度のカットオフである。ここで一つ目を引くのは、ecutrhoecutwfc のちょうど 8 倍になっていることだ。この倍率には理由があり、章末のコラムで扱う。値の意味と収束の見方は次章(第5章)が本丸である。

金属の扱いを述べる行occupations='smearing', smearing='mv', degauss=0.02(Ry)の三行は、金属を解くための仕掛けである。金属では電子がフェルミ準位まで切れ目なく詰まり、占有・非占有の境目が鋭く切り立つ。これをそのまま計算すると数値が暴れるので、境目をわずかにぼかして滑らかに扱う。smearing='mv' はそのぼかし方(Methfessel-Paxton 系の手法)を、degauss はぼかしの幅を指す。これは物性の値ではなく、金属を安定に解くための数値の道具である。幅の選び方が結果に与える影響は第7章で正面から扱う。

磁性を述べる行nspin=2 がスピン分極を許す宣言で、これがなければ(nspin=1 なら)電子は上向き・下向きが同数に固定され、磁石にはなりえない。starting_magnetization(1)=0.7 は、最初の反復に与える磁化の傾き、すなわち初期条件である。これは答えではなく出発点にすぎず、SCF が収束すれば終状態は初期値に依らない――この事実は第8章で確かめる。

3 &electrons:どう解くか、その作法

SCF は一度の計算では解けない。仮の電荷密度から出発し、そこから電子状態を求め、得た状態から新しい密度を作り、また状態を求める――この往復を、入口と出口がそろうまで繰り返す。&electrons は、その反復の作法を述べるブロックである。

conv_thr = 1.0d-8 は、収束の判定基準である。反復のあいだに全エネルギー(厳密には密度由来の量)の変化がこの値(Ry 単位)を下回ったら、「つじつまが合った」とみなして反復を止める。小さくするほど厳しく、止まるまでの反復は増える。

mixing_beta = 0.3 は、反復で新しい密度と古い密度をどの割合で混ぜるかを決める。新しい密度をそのまま全部採用すると、金属では往復が振動して収束しにくい。そこで新密度を 3 割だけ混ぜ、残り 7 割は古い密度を残して、ゆっくり近づける。金属では小さめに取って安定させるのが定石である。

この二行は、構造でも磁性でもなく、解き方の作法に属する。同じ系を解くなら、conv_thr を厳しくしても緩めても、収束しさえすれば同じ答えに行き着く。物理の前提を変える行ではない点が、&system の各行と性格が違う。SCF の中身そのものは第8章で扱う。

4 原子の中身、置き場所、波数の標本化

ブロックを抜けると、三枚のカードが続く。

ATOMIC_SPECIES は、原子種ごとに質量と擬ポテンシャルのファイルを結びつける。ここで指定した Co.pbe-spn-kjpaw_psl.0.3.1.UPF(pslibrary (Dal Corso 2014) 由来)が、本書の計算の土台を一手に背負う。これは PAW(projector augmented-wave)法、交換相関に PBE(Perdew–Burke–Ernzerhof (1996) の頭文字。電子どうしが避け合う効果の近似のしかたで、GGA と呼ばれる型の標準的なもの。詳しくは付録の用語集)、スカラー相対論(すなわちスピン軌道相互作用なし)の擬ポテンシャルで、価電子を 17 個――3d⁷4s² の 9 個に加え、ふつうは芯に押し込める 3s²3p⁶ の 8 個――まで含めている。この 3s・3p を半芯(semicore)電子と呼ぶ。半芯まで価電子に含めるのは、コバルトのような 3d 金属の磁性をきちんと表すために効いてくる選択である。擬ポテンシャルの違いが磁気モーメントの傾向に与える影響は、第8章で別の擬ポテンシャルと比べて確かめる。

ATOMIC_POSITIONS (crystal) は、celldm で決めた格子のなかに、原子をどこへ置くかを分率座標で書く。Co(0,0,0) と Co(1/3,2/3,1/2) の二原子で、これが hcp の最小単位を作る。(crystal) は座標を格子ベクトルの分数で読む指定である。

K_POINTS automatic 12 12 8 は、無限に続く周期構造を有限のデータで扱うために、ブリュアンゾーン(逆空間の単位胞)を \(12\times12\times8\) の格子で標本化する指定である。底面に対して \(c\) 軸方向が短いので、その方向の刻みも少なくしてある。末尾の 0 0 0 は標本格子の原点をずらさない指定である。周期の標本化をどこまで細かくすべきかは第6章で扱い、そこで磁気モーメント \(M\) が k 点で揺れることまで見る。

5 四つの色で全体を見渡す

ここまで一行ずつ読んできた供述書を、役割の色で塗り直すと、入力ファイルの骨格が一目で見える。図 1 のように、すべての行は構造・電子・磁性・数値条件の四つに振り分けられる。

flowchart TB
    A["co.scf.in"] --> B["構造<br>ibrav=4, celldm(1), celldm(3)<br>nat, ATOMIC_POSITIONS"]
    A --> C["電子<br>ecutwfc, ecutrho<br>Co PAW 擬ポテンシャル"]
    A --> D["磁性<br>nspin=2<br>starting_magnetization"]
    A --> E["数値条件<br>smearing, degauss, conv_thr<br>mixing_beta, K_POINTS"]
    style B fill:#d4edda,stroke:#28a745
    style C fill:#d1ecf1,stroke:#17a2b8
    style D fill:#fff3cd,stroke:#ffc107
    style E fill:#f8d7da,stroke:#dc3545
図 1: 入力ファイルの各行を四つの役割に塗り分ける。構造(緑)・電子(青)・磁性(黄)・数値条件(赤)。どの行も必ずどれか一色に属する。

色とは別に、ブロックの入れ子も読めるようにしておきたい。図 2 のとおり、入力は &control(何を計算するか)→ &system(どんな系か)→ &electrons(どう解くか)の順で、外枠から中身へ降りていく。供述の重みは中央の &system に集まる。

flowchart LR
    A["&control<br>何を計算するか<br>calculation='scf'"] --> B["&system<br>どんな系か<br>構造・電子・磁性・金属の扱い"]
    B --> C["&electrons<br>どう解くか<br>conv_thr, mixing_beta"]
    C --> D["ATOMIC_SPECIES /<br>ATOMIC_POSITIONS /<br>K_POINTS"]
    style B fill:#fff3cd,stroke:#ffc107
図 2: 三つのブロックの入れ子と役割。外枠の段取りから系の中身へ、最後に解き方へと降りる。供述の中核は &system に集まる。

小課題 4.1

smoke 束の入力に現れる次の八行を、「構造」「電子」「磁性」「数値条件」の四つの役割に塗り分けよ。塗り分けたうえで、どれか一行について「迷いやすいが、なぜその色になるか」を一言で説明せよ。

ibrav=4 / ecutwfc=80 / nspin=2 / degauss=0.02 / celldm(3)=1.6228 / starting_magnetization(1)=0.7 / K_POINTS automatic 12 12 8 / ecutrho=640

解答例

  • 構造ibrav=4, celldm(3)=1.6228(六方格子の種類と軸比=結晶の形を決める)。
  • 電子ecutwfc=80, ecutrho=640(波の解像度=電子状態の表現力を決める)。
  • 磁性nspin=2, starting_magnetization(1)=0.7(スピン分極の自由度と、その初期条件)。
  • 数値条件degauss=0.02, K_POINTS automatic 12 12 8(フェルミ面のぼかし幅と、波数の標本化)。

迷いやすいのは degauss である。これは \(\mu_\text{B}\) や eV のような物性の値に見えて、実はそうではない。金属を安定に解くために占有の境目をぼかす幅、すなわち数値的な道具に属する。この区別を取り違えると、第7章で見る「degauss の選び方で \(M\) が揺れる」という現象を、物性そのものの変化と読み違えてしまう。

よくある誤解

starting_magnetization(1)=0.7 を見て、「初期に 0.7 と入れたから、出てくる磁気モーメントも 0.7 あたりに引っぱられるのではないか」と疑いたくなる。だが、これは初期条件であって答えではない。SCF はエネルギーの低い終状態を探す手続きで、出発点を変えても同じ谷に転がり込むなら、終状態は初期値に依らない。初期値を 0.1 から 1.0 まで変えても同じ終状態へ収束することは、第8章で実データを並べて確かめる。初期条件はゴールを決めず、せいぜいゴールまでの反復の歩数を変えるだけである。「入れた値が出てくる」ように見える行と、「出発点を指すだけ」の行を、混同しないことが大切である。

コラム:なぜ Co の PAW は ecutrhoecutwfc の 8 倍に取るのか

擬ポテンシャルによっては、電荷密度のカットオフ ecutrho を波動関数の ecutwfc の 4 倍に取れば足りる。ところがこの Co PAW では 8 倍にしている(smoke 束で 80→640 Ry)。理由は二つ重なっている。一つは、コバルトが 3d 電子をもつ金属で、電荷密度が芯の近くで鋭く立ち上がること。もう一つは、半芯の 3s・3p を価電子に含めているため、密度を表す波がさらに細かくなることだ。電荷密度はもともと波動関数の二乗で、二乗した分だけ細かい波を含むから、ecutrhoecutwfc より高く取る必要がある。コバルトではその細かさが一段増すぶん、倍率も 4 倍では足りず 8 倍に取る。ecutrhoecutwfc と一緒に決まる相棒であって、単独で吊り上げるつまみではない、と覚えておくとよい。

6 この章のまとめ

  • 入力ファイルは命令の羅列ではなく、置いた物理の仮定を述べる供述書である。どの行も構造・電子・磁性・数値条件のどれか一つの役割に属する。
  • ブロックは &control(何を計算するか=calculation='scf')→ &system(どんな系か)→ &electrons(どう解くか)の順で降りる。供述の中核は &system で、構造(ibrav/celldm)・解像度(ecutwfc/ecutrho)・金属の扱い(smearing/degauss)・磁性(nspin/starting_magnetization)が同居する。
  • &electronsconv_thrmixing_betadegauss は物性の値ではなく、解き方・解きやすさの作法に属する。物性を述べる行と取り違えない。
  • Co の PAW は半芯と 3d のために ecutrhoecutwfc の 8 倍に取る。starting_magnetization は初期条件であって答えではない。
  • 次章からは、この供述書で「値の意味」を先送りにした二つの解像度――ecutwfcK_POINTS――に順に向き合う。まず波の解像度、ecutwfc から。