xiaoqiugood的个人博客分享 http://blog.sciencenet.cn/u/xiaoqiugood

博文

pwscf计算费米面

已有 12171 次阅读 2013-9-16 07:27 |个人分类:Fermi surface et al|系统分类:科研笔记

例子:计算fcc Cu的费米面  【转自Valenhou+解读】
步骤:

1)先进行一次自洽计算,如前面的例子;

2)kvecs_FS.x产生要计算费米面时的密集网格k点,这时的k不仅包含了第一布里渊区的,还有其他区域的;

3)进行一次非自洽计算,计算这些k点的本征值;

4)然后采用bands_FS.x将计算的k点以及相应的本征值,转换成xcrysden软件的bxsf格式,以便采用xcrysden来画图。


一、自洽计算

cat > cu.scf.in << EOF


&control

calculation='scf'

restart_mode='from_scratch',

pseudo_dir = './',

outdir='./'

prefix='cu'

tstress = .true.

tprnfor = .true.

/

&system

ibrav = 2, celldm(1) =6.73, nat= 1, ntyp= 1,

ecutwfc = 25.0, ecutrho = 300.0

occupations='smearing', smearing='gaussian', degauss=0.02

/

&electrons

diagonalization='david'

conv_thr = 1.0e-8

mixing_beta = 0.7

/

ATOMIC_SPECIES

Cu 63.55 Cu.pz-d-rrkjus.UPF

ATOMIC_POSITIONS

Cu 0.0 0.0 0.0

K_POINTS (automatic)

8 8 8 0 0 0


pw.x < cu.scf.in>cu.scf.out


二、准备计算费米面要用到的密集网格k点
采用shell脚步来准备


Sysname='cu'

Calc_Type='FS'     #【这个参数似乎没用到,新版例子08中没有】
nabc=' 12 12 12 '
n_start=3
n_last=6
#
E_Fermi=`grep Fermi cu.scf0.out cut -c 26-36`   【侯博回答,此处为笔误,cu.scf0.out实际应改为cu.scf.out
a1=`grep 'b(1)' cu.scf.out cut -c 24-54`
a2=`grep 'b(2)' cu.scf.out cut -c 24-54`
a3=`grep 'b(3)' cu.scf.out cut -c 24-54`

cat > kvecs_FS.in <<EOF
$a1
$a2
$a3
$nabc
$Sysname
EOF

cat > kvecs_FS.in < kvecs_FS.out


kvecs_FS.x< kvecs_FS.in > kvecs_FS.out
 

注释:这里$Sysname用来标记所计算的体系,这个例子中给它赋值为'cu',$Calc_Type用来标记所进行的是Fermi surface的计算,可以赋值为'FS',$nabc用来设置在计算费米面时所要计算的k点网格,这里设置为'12x12x12'的k点网格。$n_start用来标记从哪个能带开始计算起,这里设置为3,$n_last用来设置要计算到哪根能带为止,这里设置为6,一般它们两个要设置在费米能级所在的带的上下,也就是它们的范围之内要包括了费米能级【如何有效地设置这两个参数啊?设置前怎么知道】。$a1, $a2,$a3用来标记所得到k点网格产生的k点坐标【12x12x12产生的各网格点的坐标?】$E_Fermi用来标记自洽计算中得到的费米能级值【即,从前一步自洽计算得到的结果(或自洽)文件中提取费米能 Ferimi energy】

  【推测,此步应该产生一个名为kves_cu包含了k点坐标的文件,见第二步:cat  kvecs_$Sysname >> cu.fs.in

【其中的a1,a2,a3其实就是你*.scf.out里的倒空间基矢b1,b2,b3复制粘贴就行了。

【附nband 参数的英文释义: nbnd INTEGER Default: for an insulator, nbnd = number of valence bands (nbnd = # of electrons /2); for a metal, 20% more (minimum 4 more) number of electronic states (bands) to be calculated. Note that in spin-polarized calculations the number of k-point, not the number of bands per k-point, is doubled 】

 
二、进行非自洽的计算,得到这些k点的本征值:  【计算费米面之前有必要进行能带结构计算吗】
cat > cu.fs.in << EOF
&control  
calculation='nscf'  
prefix='cu',  
pseudo_dir = './',  
outdir='./'
/
&system  
ibrav=2, celldm(1) =6.73, nat=1, ntyp=1,  
ecutwfc = 25.0, ecutrho = 300, nbnd=8
/
&electrons  
conv_thr = 1.0e-8    mixing_beta = 0.7  
/
ATOMIC_SPECIES
Cu 63.55 Cu.pz-d-rrkjus.UPF
ATOMIC_POSITIONS
Cu 0.0 0.0 0.0
K_POINTS

EOF


cat  kvecs_$Sysname >> cu.fs.in   【侯博回答:kvecs_$Sysname ,即是kvecs_cu;表示将所产生的K点坐标追加到文件cu.fs.in中吗?】

 

pw.x <cu.fs.in> cu.fs.out


三、准备将计算的本征值和对应的k点写成xcrysden的bxsf格式
mv cu.fs.out  Bands.out
cat > input_FS <<EOF
$n_start
$n_last
$E_Fermi
$Sysname
$Calc_Type
$nabc
$a1
$a2
$a3
EOF
 

bands_FS.x <bands.out>bands_FS.out


mv Bands.bxsf  cu.fs.bxsf




仿照上述脚本写得计算金属Sc费米面的完整脚步如下:

   

#!/bin/sh

####################################################################

PW_ROOT=/home/pwscf/pwscf/espresso-4.0.5/bin/

PSEUDO_DIR=/home/workplace/pwscf-phononpy/pseudo

TMP_DIR=/home/workplace/pwscf-phononpy/tmp

#export PARA_PREFIX='mpirun -np 2'

# or   export PARA_PREFIX='mpirun' ,export  PARA_POSTFIX= -np 3

export PW_ROOT PSEUDO_DIR TMP_DIR

export name='Sc'    #【定义环境变量等】


#######################【自洽计算,产生$name.scf.out_$a文件,即Sc.scf.out_0文件】################################ self-consistent calculation -Non-Spin-Polarised case  #【为什么要进行自洽计算,目的是.....】


for a in 0

do

cat > $name.scf.in_$a << EOF

&control

   calculation = 'scf'

   restart_mode='from_scratch',

   prefix='$name',

   pseudo_dir = '$PSEUDO_DIR/',

   outdir='$TMP_DIR/'

   tstress=.t.,

   tprnfor=.t.

/

&system

 ibrav=0,

 nat=2,

 ntyp=1,

 ecutwfc=90.0, ecutrho = 360.0

 occupations='smearing',

 smearing='methfessel-paxton',

 degauss=0.02

/

&electrons

conv_thr = 1.0d-8  

mixing_beta = 0.7

/

#&CELL

#cell_dynamics= 'damp-pr'

# press=$a

#/


ATOMIC_SPECIES

Sc    44.9559    Sc.pbe-nsp-van.UPF


CELL_PARAMETERS {bohr}

6.254820382  0.00000000  0.00000000

-3.127410191  5.41683334  0.00000000

0.000000000  0.00000000  9.958412336


ATOMIC_POSITIONS {crystal}

Sc 0.3333333333333286 0.6666666666666714 0.2500000000000000

Sc 0.6666666666666714 0.3333333333333286 0.7500000000000000


K_POINTS {automatic}

8 8 4 0 0 0

EOF

$PW_ROOT/pw.x  < $name.scf.in_$a >$name.scf.out_$a

done


#########【准备并用kvecs_FS.x程序产生计算费米面要用到的密集网格k点,各k点坐标存储于文件kves_$Sysname】#############

#kves_$Sysname的定义见keves_FS.f源程序



# prepare input file $name.fs.in

#

Sysname='$name'

Calc_Type='FS'

nabc=' 16 16  8'     #【k点分割数与一般的计算相同吗】

n_start=3

n_last=6

#

E_Fermi=`grep Fermi $name.scf.out_$a| cut-c 26-36`    #$name.scf.out_$a,由前一步自洽计算产生

a1=`grep 'b(1)' $name.scf.out_$a| cut -c 24-54`

a2=`grep 'b(2)' $name.scf.out_$a| cut -c 24-54`

a3=`grep 'b(3)' $name.scf.out_$a| cut -c 24-54`

 

cat > kvecs_FS.in <<EOF

$a1                            #$a1具有xx,xx,xx的数字格式,相当于一维向量,与nabc格式类似,见上

$a2

$a3

$nabc

$Sysname

EOF

 

kvecs_Fs.x < kvecs_FS.in >kvecs_FS.out



#########################【非自洽计算,产生文件$name.nscf.out】########################

# 进行一次非自洽计算,计算上一步产生的各k点的本征值;


cat > $name.nscf.in << EOF

&control

  calculation='nscf',

  prefix='$name'

  pseudo_dir = '$PSEUDO_DIR/',

  outdir='$TMP_DIR/'

/

&system

   ibrav = 0,

  nat=  2,

  ntyp= 1,

  ecutwfc = 90,

  ecutrho = 400.00,

#  occupations='smearing',   #【非自洽计算这些参数均无必要】

 # smearing='methfessel-paxton',

 # degauss=0.02

  nbnd=8  #【一定要定义nbnd吗?默认值是.....】

/

&electrons

 diagonalization='cg'  

  mixing_beta = 0.7

  conv_thr =  1.0d-10          #【1.0d-10与1.0e-10含义一样吗】

/

#【diagonalization='cg' 英文注释:conjugate-gradient-like band-by-band diagonalization. Typically slower than 'david' but it uses less memory     and is more robust (it seldom fails)】


#&CELL

#cell_dynamics= 'damp-pr'

# press=$a

#/


ATOMIC_SPECIES

Sc    44.9559    Sc.pbe-nsp-van.UPF


CELL_PARAMETERS {bohr}

6.254820382  0.00000000  0.00000000

-3.127410191  5.41683334  0.00000000

0.000000000  0.00000000  9.958412336


ATOMIC_POSITIONS {crystal}

Sc 0.3333333333333286 0.6666666666666714 0.2500000000000000

Sc 0.6666666666666714 0.3333333333333286 0.7500000000000000

K_POINTS    #【K_POINTS已在上一步产生】

EOF

 

cat kvecs_$Sysname >>$name.nscf.in 【将上一步产生的各k点坐标,附加到$name.nscf.in 文件后面

 

$PW_ROOT/pw.x< $name.nscf.in >$name.nscf.out


##########【采用bands_FS.x程序将上一步得到的k点及相应本征值转化成xcrysden能读取的bxsf格式】##########################


#【下面是即将见证奇迹的时刻-采用bands_FS.x将计算的k点以及相应的本征值,转换成xcrysden软件的bxsf格式,以便采用xcrysden来画图。 】


# prepare input data (input_FS, Bands.out)for bands_FS

 

mv $name.nscf.out  Bands_NSP.out

   #【注意这利用了非自洽计算产生的$name.nscf.out文件,该步的意义可能是bands_FS.x陈旭默认读取文件名为 Bands_NSP.out?】

 

cat > input_FS <<EOF

$n_start $n_last

$E_Fermi

$Sysname

$nabc

$a1

$a2

$a3

EOF

 

$PW_ROOT/bands_FS.x < Bands_NSP.out >bands_fs.out    #【这一操作将产生Bands_FS.bxsf文件】

mv Bands_FS.bxsf  $name.fs_NSP.bxsf    #【这一操作意义不大,是为了记忆方便】

 

$ECHO "  Fermi surface plot: use 'xcrysden --bxsf $name.fs_NSP.bxsf' to plot ...\c"

$ECHO" done"


#######################################################





ESPRESSO 5.0.2 pp/examples/example02解读

   This example shows how to use pw.x to calculate the DOS of Ni     and how to plot the Fermi Surface using XCrysDen


The calculation proceeds as follows (for the meaning of the cited input variables see the appropriate INPUT_* file)


1) make a self-consistent calculation for Ni (like in example 1).    (input=ni.scf.in, output=ni.scf.out)


2) make a band structure calculation for Ni (input=ni.dos.in,   output=ni.dos.out) on a uniform k-point grid (automatically

  generated).

   In this example the Fermi level is calculated with the   tetrahedra method (not in the actual band structure calculation but in

  the subsequent DOS calculation). If preferred, a gaussian broadening  may be specified in this or in the subsequent step.


3) the program dos.x reads file filpun (ni.pun) and calculates the DOS on a  uniform grid of energies from Emin to Emax, with grid step Delta E.

  The output DOS is in file ni.dos, ready for plotting. 【如何设定uniform grid?为什么要求uniform grid】


4) the program projfwc.x projects the crystal wavefunctions on an    orthogonalized basis set of atomic orbitals, calculates the Loewdin

  charges, spilling parameter, and the projected DOS (total DOS in file    'ni.pdos_tot', s and d componentin files 'ni.pdos_atm#1(Ni)_wfc#1(s)'

  and 'ni.pdos_atm#1(Ni)_wfc#2(d)' respectively). (input=ni.pdos.in,   output=ni.pdos.in)


5) Fermi Surface plot, courtesy of    Eyvaz Isaev

     Theoretical Physics Department;   Moscow State Institute of Steel and Alloys;   (eyvaz_isaev@yahoo.com, e.isaev@misis.ru)

     【斯人已去】

a.  First, one generates a  grid of k-points (all of them, not only those  in the Irreducible Brilloin Zone) using auxiliary code kvecs_FS.x

 

b. Then, the non-scf calculation is performed

c.   Then, auxiliary code bands_FS.x collects the data and produces a   file ni.fs.bxsf that can be read by XCrySDen (www.xcrysden.org) as:

     xcrysden --bxsf ni.fs.bxsf



  Additional info for customization of the script: 【脚本中的用户定义部分】

#

# A user has to edit so-called "user part" in order to define some required parameters.

#

#(1)  Sysname   - a nickname for your system

#(2)  Calc_Type - The Fermi Surface calculations (FS) or band-structure  calculations (Band) which will be included later.

#             Presently band-structure calculations could be carried out    by means ofplotband.xfrom PP (postprocessing) directory     or a little package #             distributed by E.Isaev (posted to pw_forum).

#(3)  nabc      - a number for dividing of each edge of a parallelepiped【平行六面体】.      Be careful, the total number of generated k-points is

#                     (na+1)*(nb+1)*(nc+1), i.e. including Gamma-point. 及如果定义为16*16*8=2048个,则实际为 (16+1)*(16+1)*(8+1)=2601个

# (4) n_start   - starting band's number for the Fermi Surface calculations.    It is obvious, we have to deal with the bands crossing the  Fermi level.

#                      【怎样确定该值】

# (5) n_last    - last band's number for FS calculations  【怎样确定该值】

#

# That's all!!! Present values in the script(may be edited):

Sysname='ni'

Calc_Type='FS'

nabc=' 16 16 16 '

n_start=2

n_last=5

#

#

# Nota Bene : You cantake morebands andthen choosefrom a XCrySDen menu   only those bands which cross the Fermi level  

# Nota Bene : If you have mistaken choosing bands to be considered for the FS construction, you do not need to restart all calculations.

#            Just edit "bands_FS" file and restart "bands_FS.x" manually. 【设置错了也没有关系,只需修改bands_FS 文件(是源文件bands_#     FS.f90?)即可】    It will read     Bands.out and result Bands.bxsf which you can rename as you like.

#  




run_exampe 脚本

侯博补充解释

这个错误来自bands_FS.x 在读入前一步的本征值计算输出文件Bands_NSP.out时所导致的。
   你检查该文件里总的k点数是不是超过了max_kpoints=100000(这个由于bands_FS.f90的开头的定义所限制的)。


     另外可能最关键的原因是你所采用的Sc.pbe-nsp-van.UPF里的Z valence=11,也就是它的价电子为11(包含了半芯态的3s2 3p6),因此你在做pw.x的计算中,Sc的原胞里
2个原子,那总的价电子数为22因此nbnd必须设置的要比22大,比如大10相应地,bands_FS.x的输入文件input_FS中的值n_startn_last也需做调整,将它们设置Fermi energy所在的能带(你的例子为22)的值附近




https://wap.sciencenet.cn/blog-567091-725266.html

上一篇:PWSCF计算电声耦合常数
下一篇:PWSCF计算电子态密度
收藏 IP: 128.84.125.*| 热度|

0

该博文允许注册用户评论 请点击登录 评论 (0 个评论)

数据加载中...
扫一扫,分享此博文

Archiver|手机版|科学网 ( 京ICP备07017567号-12 )

GMT+8, 2024-5-13 06:48

Powered by ScienceNet.cn

Copyright © 2007- 中国科学报社

返回顶部