在凝聚态物理与材料计算领域，传统的密度泛函理论（DFT）在处理包含部分填充 d 或 f 电子的强关联体系（如过渡金属氧化物、镧系/锕系元素化合物）时，往往面临严重失效。动态平均场理论（DMFT）通过引入局域多体相互作用的动态效应，补足了这一短板。
本文将以SrVO3为例介绍在ABINIT中开展DFT+DMFT计算的核心逻辑与关键控制参数。
教程参考：
https://docs.abinit.org/tutorial/dmft/
环境配置
ABINIT安装，参考：Abinit-10.4.7安装教程
环境变量填写：
export ABI_HOME=Replace_with_absolute_path_to_abinit_top_level_dir # Change this lineexport PATH=$ABI_HOME/src/98_main/:$PATH # Do not change this line: path to executableexport ABI_TESTS=$ABI_HOME/tests/ # Do not change this line: path to tests direxport ABI_PSPDIR=$ABI_TESTS/Pspdir/ # Do not change this line: path to pseudos dir
本篇推送所涉及计算文件路径：ABINIT安装包/tests/tutoparal/Input/tdmft*.abi
ABINIT默认赝势位置：
ABINIT安装包/tests/Pspdir/
可前往https://www.pseudo-dojo.org/ 下载
OmegaMaxent安装
下载:
https://github.com/dbergeron1/OmegaMaxEnt
https://www.physique.usherbrooke.ca/MaxEnt/index.php/Main_Page
具体编译设置根据实际环境修改
conda install -c conda-forge fftw openblas liblapackcd /public/software/abinit-10.4.7/OmegaMaxEnt/OmegaMaxEnt/buildexport LIBRARY_PATH=$CONDA_PREFIX/lib:$LIBRARY_PATHexport LD_LIBRARY_PATH=$CONDA_PREFIX/lib:$LD_LIBRARY_PATHcmake .. \-DCMAKE_PREFIX_PATH=$CONDA_PREFIX \-DCMAKE_EXE_LINKER_FLAGS="-L$CONDA_PREFIX/lib -Wl,-rpath,$CONDA_PREFIX/lib" \-DCMAKE_SHARED_LINKER_FLAGS="-L$CONDA_PREFIX/lib"make -j
LDA 基态预计算与能带特征分析
教程输入文件tdmft_1.abi包含两个数据集
数据集 1：LDA 基态自洽计算，核心设置为nnsclo=3、nline=3，保证空带波函数的对角化精度，同时采用 PAW 赝势、合适的 k 点网格与平面波截断能，完成体系的基态收敛。
数据集 2：能带结构与 fatbands 计算，通过pawfatbnd=1开启轨道分辨的能带权重输出，可分别获得 V 的 d 轨道、O 的 p 轨道在各能带中的权重占比，明确关联轨道对应的能带范围。
# SrVO3, a simple archetypal correlated material## In this test, we compute the DFT electronic structure of SrVO3 in the first# dataset and then the orbital character of the bands for the t2g and eg# orbitals of the vanadium atom and the p orbitals of the oygen atoms.#Multi-dataset parametersndtset 2 # Two datasetsjdtset 12getden -1 # Second one takes the density of the first#Definition of the unit cellacell 3*7.2605 # Cubic cell withrprim 1.00.00.0 # real space primitive translation vectors0.01.00.00.00.01.0#Definition of the atom types and pseudopotentialsnatom 5 # Five atomsntypat 3 # Three typesznucl 23388 # First atom type should be the correlated one V# then, we have Sr and O.typat 12333 # V Sr O3xred # This keyword indicates that the location of the atoms# will follow, one triplet of number for each atom.# We use relative distance, along the translational# lattice vectors.0.000.000.000.500.500.500.500.000.000.000.500.000.000.000.50pp_dirpath "$ABI_PSPDIR" # This is the path to the directory where# pseudopotentials for tests are storedpseudos "Psdj_paw_pw_std/V.xml, Psdj_paw_pw_std/Sr.xml, Psdj_paw_pw_std/O.xml" # Name and location of the pseudopotentials#Planewave basis set, number of bands and occupationsecut 12.0 # Maximal plane-wave kinetic energy cut-off, in Hartreepawecutdg 20.0 # PAW: Energy Cutoff for the Double Gridnband 30 # Number of bandsoccopt 3 # Occupation option for metaltsmear 1200 K # Temperature of smearingpawprtvol 3 # Printing additional information#First dataset specific parametersnstep 40 # Number of iterations for the DFT convergencenline 3 # Number of line minimisationsnnsclo 3 # Number of non-self consistent loopstolwfr 1.0d-15 # Tolerance of the DFT convergence#K point gridngkpt 333 # Reciprocal space vectors are built from# the rprim parameters. This is the size of the# reciprocal k-points.nshiftk 4 # Convergence of the density with regular shifts.shiftk 0.50.50.50.50.00.00.00.50.00.00.00.5istwfk *1#Second dataset specific parametersnbandkss2 2kssform2 3pawfatbnd2 1 # Fatbands, only resolved for total angular momentaiscf2 -2 #kptopt2 -4 # Bandstructure: 4 segmentsndivk2 18202028 # Lengths of each segmentskptbounds2 # The four segments of the band structure are# connecting 5 points of the Brillouin zone,# given here.1/41/41/4 # R'0.00.00.0 # Gamma1/20.00.0 # X1/21/20.0 # M0.00.00.0 # Gamma
计算完成后，可通过 xmgrace 工具绘制 fatbands 图：
-
V 的 d 轨道（l=2）fatbands：
xmgrace tdmft_1o_DS2_FATBANDS_at0001_V_is1_l0002 -par ../Input/tdmft_fatband.par须注意tdmft_fatband.par文件路径
-
O 的 p 轨道（l=1）fatbands：
xmgrace tdmft_1o_DS2_FATBANDS_at0003_O_is1_l0001 tdmft_1o_DS2_FATBANDS_at0004_O_is1_l0001 tdmft_1o_DS2_FATBANDS_at0005_O_is1_l0001 -par ../Input/tdmft_fatband.par
分析结果可明确：SrVO3中 V 的 t2g 轨道主要对应 21-23 号能带，eg 轨道对应 24-25 号能带，O 的 p 轨道对应 12-20 号能带，费米能级位于 t2g 能带中间，且 t2g 轨道与 O 的 p 轨道存在显著杂化。
DFT+DMFT 自洽计算
基于收敛的 LDA 基态，可正式开展 DFT+DMFT 自洽计算，案例中存在演示窄窗口（仅 t2g 能带）与宽窗口（t2g+O-p 能带）两种关联轨道构建方案，对应教程中的tdmft_2.abi输入文件。
输入文件核心架构输入文件包含两个数据集：数据集 1 为高精度的 LDA 基态计算，与预计算流程保持一致，确保基态波函数的收敛性；数据集 2 为 DFT+DMFT 计算
核心参数设置如下：
窄窗口方案：dmftbandi=21、dmftbandf=23，仅采用 t2g 对应的 3 条能带构建 Wannier 轨道，配合dmft_t2g=1限定关联作用于 t2g 子空间；
宽窗口方案：dmftbandi=12、dmftbandf=23，包含 O-p 能带与 t2g 能带共 12 条能带，构建的 Wannier 轨道具备更强的局域性，可更完整地捕捉轨道杂化效应；
求解器与相互作用设置：
dmft_solv=5启用 CTQMC 求解器，
upawu=3.13333 eV、
jpawu=0.75833 eV设置库仑相互作用，
dmft_dc=1采用 FLL 双计数方案；
数值网格设置：匹配 CTQMC 求解器的要求，设置合理的dmft_nwli、dmft_nwlo、dmftqmc_l、dmftqmc_n参数，平衡计算精度与成本。
tdmft_2&3.abi
# SrVO3, a simple archetypal correlated material## In this test, we compute the DFT electronic structure of SrVO3 in the first# dataset and then run a DFT+DMFT calculation from the DFT one.## DATASET 1: DFT# DATASET 2: DFT+DMFT# Multi-dataset parametersndtset 2jdtset 12getwfk -1#Definition of the unit cellacell 3*7.2605 # Cubic cell withrprim 1.00.00.0# real space primitive translation vectors0.01.00.00.00.01.0#Definition of the atom types and pseudopotentialsnatom 5 # Five atomsntypat 3 # Three typesznucl 23388 # First atom type should be the correlated on V# then, we have Sr and O.typat 12333 # V Sr O3xred 0.000.000.00 # This keyword indicates that the location of the atoms# will follow, one triplet of number for each atom.# We use relative distance, along the translational# lattice vectors.0.500.500.500.500.000.000.000.500.000.000.000.50pp_dirpath "$ABI_PSPDIR" # This is the path to the directory where# pseudopotentials for tests are storedpseudos "Psdj_paw_pw_std/V.xml, Psdj_paw_pw_std/Sr.xml, Psdj_paw_pw_std/O.xml" # Name and location of the pseudopotentials#Planewave basis set, number of bands and occupationsecut 12.0 # Maximal plane-wave kinetic energy cut-off, in Hartreepawecutdg 20.0 # PAW: Energy Cutoff for the Double Gridnband 30 # Number of bandsoccopt 3 # Occupation option for metaltsmear 1200 K # Temperature of smearingpawprtvol 3 # Printing additional informationprtvol 4#First dataset specific parametersnstep1 30 # Number of iterations for the DFT convergencenline1 5 # Number of line minimisationsnnsclo1 5 # Number of non-self consistent loopstolvrs 1.0d-7#K point gridngkpt 333 # Reciprocal space vectors are built from# the rprim parameters. This is the size of the# reciprocal k-points.nshiftk 4 # Convergence of the density with regular shifts.shiftk 1/21/21/21/20.00.00.01/20.00.00.01/2istwfk *1#DFT aloneusedmft1 0#Second dataset specific parametersnstep2 10 # Number of iterations for the DFT+DMFT convergencenline2 10 # Number of line minimisationsnnsclo2 10 # Number of non-self consistent loops#DFT+DMFTusedmft2 1 # Active DMFTdmftbandi 21 # First band included in the projection. Initialdmftbandf 23# and final bands.#根据实际计算需求修改选定能带dmft_nwlo 100 # Logarythmic frequency meshdmft_nwli 100000 # Linear freqeuncy meshdmft_iter 1 # Number of iterations of the DMFT part.# We often use single-shot, since anyway the charge density# changes through the DFT+DMFT anyway.dmftcheck 0dmft_rslf 1 # Read self-energy, ifnothing(like here) initialize.dmft_mxsf 0.7 # Mixing of the old andnew self-energy at every iterations.dmft_dc 1 # Double counting type. 1 is Fully Localized Limit(FLL)dmft_t2g 1 # Special value for t2g only calculation.#CTQMCdmft_solv 5 # Choice of solver: Internal CT-SEG.dmftqmc_l 50 # Number of time slices forG(tau).dmftqmc_n 3.d6 # Number of QMC sweepsdmftqmc_therm 10000 # Thermalizationdmftctqmc_gmove 0 # Global move occurence in QMCdmftctqmc_order 50 # Perturbation orderdmftctqmc_chains 1 # force 1 chain per MPI task for testing purposes# auto-setting match OpenMP threads and alter convergence#DFT+Uusepawu1 1usepawu 10dmatpuopt 1 # The density matrix: the simplest expression.lpawu 2 -1 -1 # Angular momentum for the projected Hamiltonianf4of2_sla3 0.0 0.0 0.0upawu1 0.00 0.0 0.0 eVupawu2 3.1333333333333333 0.0 0.0 eV # Values of U for each angular momentumjpawu2 0.7583333333333333 0.0 0.0 eV # Values of J
自能与准粒子重整化权重 Z 计算
准粒子重整化权重 Z 是衡量体系关联强度的核心物理量，Z 值越小，关联效应越强，其取值可通过零频附近的自能虚部计算得到。
自能文件位于tdmft_2o_DS2Self-omega_iatom0001_isppol1，包含松原频率、自能实部、自能虚部等核心数据；
head -n 8 tdmft_2o_DS2Self-omega_iatom0001_isppol1 > self.datxmgrace -block self.dat -bxy 1:3
虚时间格林函数 G (τ) 与统计噪声分析
局域虚时间格林函数输出至Gtau.dat文件，是解析延拓获得实频率谱函数的核心依据，同时可用于验证计算的统计收敛性。
xmgrace -nxy Gtau.dat
局域谱函数
实频率谱函数需通过最大熵方法完成从虚时间格林函数到实频率谱函数的解析延拓，推荐采用 OmegaMaxEnt 代码完成该计算。
生成 OmegaMaxEnt 默认参数模板
OmegaMaxEnt# 在Gtau.dat所在的工作目录执行，生成默认输入参数文件
修改 OmegaMaxEnt 输入参数文件
# 核心数据文件：必须是你的DMFT输出的Gtau.datdata file: Gtau.dat# ---------------------- 预处理器参数 ----------------------bosonic data(yes/[no]): noimaginary time data(yes/[no]): yestemperature(in energy units, k_B=1): 0.001finite value at infinite frequency(yes/[no]): n
执行解析延拓计算
OmegaMaxEnt# 直接在当前目录执行，会自动读取修改好的参数文件
并绘图
xmgrace OmegaMaxEnt_final_result/optimal_spectral_function_*.dat
K空间谱函数
输出和实验 ARPES 完全对应的动量 - 能量分辨谱函数 A (k,ω)，
程核心逻辑
用高精度 DMFT 计算输出收敛、低噪声的虚频自能 Σ(iω)
用 OmegaMaxEnt 把虚频自能解析延拓到实频，得到 Σ(ω)
把实频自能传回 ABINIT，计算每个 k 点、每个能量点的谱函数 A (k,ω)
输出 k 分辨谱函数文件，画图得到关联能带图。
用 tdmft_5.abi 做更高精度的 DMFT 计算，降低 QMC 噪声，保证解析延拓的准确性（只运行数据集1 和数据集2）
# Testing CTQMC options## == Convergency and starting# DATASET 1: LDA# DATASET 2: DMFTndtset 2jdtset 12prtvol 4pawprtvol 3getwfk2 1nline2 10nnsclo2 10getden3 1##### CONVERGENCE PARAMETERSnstep1 30nstep2 1nstep3 30ecut 20pawecutdg 60tolvrs 1.0d-10nband 30occopt 3 tsmear 1200 K##### PHYSICAL PARAMETERSnatom 5 ntypat 3 typat 12333znucl 23.038.08.0# V Sr O*3xred 0.000.000.00 #vectors (X) of atom positions in REDuced coordinates0.500.500.500.500.000.000.000.500.000.000.000.50acell 3*7.2605rprim 1.00.00.0 #Real space PRIMitive translations0.01.00.00.00.01.0# == Points k and symetrieskptopt 1ngkpt 666nshiftk 4shiftk 1/21/21/21/20.00.00.01/20.00.00.01/2istwfk *1# == LDA+DMFTusedmft1 0usedmft2 1usedmft3 1dmftbandi 21dmftbandf 23dmft_nwlo 1600dmft_nwli 100000dmft_iter 10dmftcheck 0dmft_rslf 0dmft_mxsf 0.8dmft_dc 1dmft_t2g 1# == CTQMCdmft_solv 5 # CTQMCdmftqmc_l 800dmftqmc_n 1.d8dmftqmc_therm 10000# In general the correct value for dmftctqmc_basis is 1 (the default)dmftctqmc_basis 2 # to preserve the test:dmftctqmc_check 0 # check calculationsdmftctqmc_correl 1 # correlationsdmftctqmc_grnns 0 # green noisedmftctqmc_meas 1 # modulo de mesure Edmftctqmc_mrka 0 # markov analysisdmftctqmc_mov 0 # moviedmftctqmc_order 050 # perturbationdmftctqmc_chains 1 # force 1 chain per MPI task for testing purposes# auto-setting match OpenMP threads and alter convergence# == DFT+Uusepawu1 1usepawu 10dmatpuopt 1 # The density matrix: the simplest expression.lpawu 2-1-1f4of2_sla 0.0 0.0 0.0upawu1 0.0 0.0 0.0 eVjpawu1 0.0 0.0 0.0 eVupawu2 3.1333333333333333 0.0 0.0 eVjpawu2 0.7583333333333333 0.0 0.0 eVupawu3 0.0000000000000000 0.0 0.0 eVjpawu3 0.0000000000000000 0.0 0.0 eV#upawu3 3.1333333333333333 0.0 0.0 eV#jpawu3 0.7583333333333333 0.0 0.0 eV##jpawu 1.1666666666666663 0.0 0.0 eV#################################BANDSTRUCTURE################################dmft_rslf3 1tolwfr3 1.0d-12nbandkss3 20kssform3 3pawfatbnd3 1dmft_kspectralfunc3 1#Parameters (to uncomment) for bands structureiscf3 -2kptopt3 -4#ndivk3 4557ndivk3 90505070#ndivk3 9557kptbounds3 1/21/21/2 #R'0.00.00.0 #Gamma1/20.00.0 #X1/21/20.0 #M0.00.00.0 #Gamma
确认自能文件
tdmft_5o_DS2Selfrotformaxent0001_isppol1_iflavor0001存在
虚频自能提取与 OmegaMaxEnt 解析延拓（得到实频自能 Σ(ω)）
# 1. 新建专门的解析延拓目录，避免文件混乱mkdir Spectralcd Spectral# 2. 提取低噪声的低频虚频自能head -n 26 ../tdmft_5o_DS2Selfrotformaxent0001_isppol1_iflavor0001 > self.dat# 3. 生成OmegaMaxEnt默认参数文件OmegaMaxEnt
修改OmegaMaxEnt_input_params.dat
data file: self.datOPTIONAL PREPROCESSING TIME PARAMETERSDATA PARAMETERSbosonic data(yes/[no]):imaginary time data(yes/[no]):temperature(in energy units, k_B=1):finite value at infinite frequency(yes/[no]): yes
执行自能的解析延拓计算
OmegaMaxEnt用实频率绘制自能量的虚部
xmgrace OmegaMaxEnt_final_result/optimal_spectral_function.dat
ABINIT k 分辨谱函数计算所需文件准备
cd ..# 2. 把延拓好的实频自能复制到主目录，命名为self_ra.datcp Spectral/OmegaMaxEnt_final_result/optimal_spectral_function.dat self_ra.dat# 3. 生成ABINIT需要的实频网格文件# 第一步：统计self_ra.dat的行数，写入网格文件wc -l self_ra.dat > tdmft_5o_DS3_spectralfunction_realgrid# 第二步：把自能数据追加到网格文件末尾cat self_ra.dat >> tdmft_5o_DS3_spectralfunction_realgrid# 4. 生成多轨道实频自能文件（SrVO3的3个t2g轨道简并，复制3次）cp self_ra.dat tdmft_5i_DS3Self_ra-omega_iatom0001_isppol1# 追加2次，让文件里有3份完全相同的自能数据，对应3个t2g轨道cat self_ra.dat >> tdmft_5i_DS3Self_ra-omega_iatom0001_isppol1cat self_ra.dat >> tdmft_5i_DS3Self_ra-omega_iatom0001_isppol1# 5. 复制自能的酉矩阵文件（用于非立方体系的轨道旋转，SrVO3是立方晶系，直接复制即可）cp tdmft_5o_DS2.UnitaryMatrix_iatom0001 tdmft_5i_DS3.UnitaryMatrix_iatom0001# 6. 备份虚频自能文件，用于ABINIT重启校验cp tdmft_5o_DS2Self-omega_iatom0001_isppol1 tdmft_5o_DS3Self-omega_iatom0001_isppol
修改 tdmft_5.abi，添加 DS3 数据集，执行 k 分辨谱函数计算
# 原来的ndtset 2 改成下面两行ndtset 3jdtset 3
执行 k 分辨谱函数ABINIT程序的计算（注意资源分配）
k 分辨关联能带画图
1. 提取画图所需数据
cp tdmft_5o_DS3_DFTDMFT_SpectralFunction_kres bands_dmft.datgrep " BAND" -A 261 tdmft_5o_DS3_FATBANDS_at0001_V_is1_l0001 | grep -v BAND > bands_dft.dat
gnuplot 画图脚本
#Graphset nokeyset tics font ", 30"set xzeroaxis lw 1 lt 3 lc rgb "white"set tics textcolor rgb "black"setxtics("R'"1.000000,"{/Symbol G}"90.00000000,"X"140.000000000,"M"190.00000000,"{/Symbol G}"260.0000000)set ytics 1.0set ylabel "Energy(eV)"set yrange[-3.8:3.8]set zrange[0:100]#Colorsset palette defined(0.0"black", 0.2"dark-blue",2.0"cyan", 5"green", 7"yellow")set mouseset pm3d map# PLOTsplot "bands_dmft.dat" u 3:1:2 with pm3d,\"bands_dft.dat" using($1+1):2:(0.0) with points pt 3 ps 0.3 linecolor rgb "white"# OUPUTset term postscript colour eps enhancedset tics font ", 10"set output "bands.eps"replot
