ecaljdoc

Appearance

Dielectric function: epsPP0 ​

epsPP0(2025-5-8 version)

During QSGW, we calculate for in the GW calculation.

However, we need to calculate with larger number of k points to compare experiments.

Based on the no local-field correction equation , we can calculate interband and intraband contributions separately. The intraband contribution is the Drude weight at the limit . We use the command epsPP0 now.

We can include local-field corrections (in fact, we do include the local-field correction is QSGW iteration), but computationally expensive and not open to users currently.

Computational steps ​

Examples of epsPP0 ​

It is instructive to learn things from samples as

  • ecalj/Samples/EPS/GaAsEps
  • ecalj/Samples/EPS/CuEpsPP0

Follow job files in these directories(bash job should work). We explain the steps in job in the following steps.

You can start from only ctrl, GWinput. Then

gnuplot -p epsinter.glt 
gnuplot -p epsintra.glt 
gnuplot -p epsall.glt

Note that we usually need many k points if you like to have 0.1 level of error for dielectric constants (at the limit of ) (e.g, reasonable results for GaAs may require n1n2n3 20 20 20 in GWinput.)

step1: lmf収束計算行う。 ​

sigmファイルがある場合は, QSGW計算による固有値固有関数からの誘電関数を計算することができる。

step2: GWinputのパラメータを設定 ​

  • GWinput内でQforEPSで囲まれた箇所を探してください。 // defaultではこのように書かれている

    QforEPSau on
<QforEPS>
0 0 0.00050
0 0 0.00100
0 0 0.00200
</QforEPS>

    • この部分は誘電関数を計算するときのqベクトルを設定しています。 この値を変えることで誘電関数のq方向の依存性を調べることができます。

    • QforEPSau onがあるので、単位はa.u.になっています。すなわち、0 0 0.001であれば = (0,0,0.001) bohr^{-1}です

    • もしQforEPSau onがない場合、qの単位は です。は=alatでctrlで定義したものです。あるいはSiteInfo.chkに表示されています。

    光学応答の場合はq=0としたいですが、そうすると今のコードでは数値的に不安定です。 なので適宜小さくとってください。この不安定性は以下の図、ecalj/Samples/EPS/GaAsEpsでgnuplot -p epsinter.gltで100eVs弱のところに現れています。なのでたとえばGaAsEpsだと0,0,0.00050(EPS0001に対応）だとすこしにおおきくなっておりよろしくない、ということになります。gnuplotファイルepsinter.gltを見てください。

  • エネルギーメッシュの取り方
    誘電関数を計算する際のメッシュは対数メッシュでとられており、GWinput内で
    HistBin_dw    2d-3 
    HistBin_ratio 1.08

の値を変えることで、誘電関数のエネルギーメッシュを増やす(減らす)ことができます。メッシュを細かくとれば、誘電関数の構造がより現れるようになります。ただし、細かすぎると計算コストがすこし増えます。計算法はまずは虚部のウエイトをヒストグラム的に蓄積したあと、実部はヒルベルト変換の方法で求めています。なので数値的にはまあまあ安定です。

  • 計算メモリ量をへらしたい場合はGWinputに
KeepEigen False
KeepPbp False

等のコマンドを適宜書き込んで下さい。（we have to explain details...）

  • QforEPSau on これはQforEPSのqをa.u.ではかるという意味になります。 たとえば<QforEPS>に0d0 0d0 0.00005とあれば、これは q=(0d0 0d0 0.00005)/bohrと読んでください。 （過去のQforEPSunita on と同じ意味）。 確認するには,2pi/alatをEPSファイル最初の3つの数字に乗じて <QforEPS>にかかれているqになるのをみます。

  • 以前の--zmel0のオプションは2025-5-8で廃止しました。 これはGramSchmidt2をsugw.f90に挿入して すこし正確な直行化が可能になったためです。

  • epsPP0では最後にreadeps.pyを読んでて答えをまとめ上げてます。

  • Cuの場合だと、qがあまりにゼロに近いと計算が不安定です。たとえば

QforEPSau on
<QforEPS>
 0d0 0d0 0.00005
 0d0 0d0 0.0001
 0d0 0d0 0.0002
 0d0 0d0 0.0004
 0d0 0d0 0.0008
 0d0 0d0 0.0016
 0d0 0d0 0.0032
</QforEPS>

などとすると、 0d0 0d0 0.00005のグラフがおかしいです。 それ以外のグラフはほぼ重なります。だた、0.0032ぐらいになるとomega=0での実部がいくらかずれてきます。

  • CuepsPP0ではフェルミ面の寄与もあり、
gnuplot -p epsintra.glt

を実行してフェルミ面の寄与が見れます。This is Drude weight、q分解でみれています。 omega=0近傍でのふるまいはFetter-Waleckaにあるようにふるまっています。

gnuplot -p epsintra.glt
gnuplot -p epsall.glt
  • epsPP0の計算が終わるとEPS000*.nlfc.datというファイルがGWinputで設定したq点の数分作られます。この中は
    // example EPS0001.nlfc.dat
    q(1:3)   w(Ry)   eps    epsi  --- NO LFC
  0.01000000  0.00000000  0.00000000    0.0000E+00  0.276888086703565E+02  0.211273167660642E-17    0.361156744555281E-01 -0.275572453667495E-20
  0.01000000  0.00000000  0.00000000    0.1015E-04  0.276873286605807E+02  0.424191214559347E-01    0.361175202203266E-01 -0.553348246663633E-04
  0.01000000  0.00000000  0.00000000    0.3076E-04  0.276716446748960E+02  0.128515093142673E+00    0.361372966009293E-01 -0.167832020581202E-03
  0.01000000  0.00000000  0.00000000    0.5200E-04  0.276450859807482E+02  0.197440420455938E+00    0.361709489903162E-01 -0.258331892398948E-03
  0.01000000  0.00000000  0.00000000    0.7389E-04  0.276274571460234E+02  0.249411594405801E+00    0.361929258459201E-01 -0.326737828014011E-03

1~3行目にq点座標、4行目にエネルギー(Ry単位)、5・6行目が誘電関数の実部・虚部、7・8行目が誘電関数の逆数の実部・虚部 が書かれています。

** （改良計画）誘電関数計算にはバンド吸収端で(E-E0)**.5でキレイに立ち上がらない(GaAsなどの場合）、 という数値計算上の問題があります。メッシュを細かく取るときれいにできることは示したことがあります。 変な内挿法よりも必要なことだけ細かく取るという技法が有効だと思います。汎用化しないといけないです。

**（改良計画）q=0でのinterbandの寄与は<u|u>行列を使えば先に1/q^2の割り算ができるので もっと直接的にできるはずです。

(Here we keep old documents, but commented out...)