User Functions
All of these functions are exported from ThreeBodyTB for your convenience.
Crystal / Energy
ThreeBodyTB.CrystalMod.makecrys
— Functionmakecrys(A,coords,types; units=missing)
Return a crystal object from 3×3 lattice, nat × 3 coords in crystal units, and nat element strings.
Will use units set by set_units
, with default to Ang. Can override with units
.
Note: also export-ed directly from ThreeBodyTB for convenience
julia> makecrys([10.0 0 0; 0 10.0 0; 0 0 10.0], [0.0 0.0 0.0], ["H"])
A1= 10.00000 0.00000 0.00000
A2= 0.00000 10.00000 0.00000
A3= 0.00000 0.00000 10.00000
H 0.00000 0.00000 0.00000
makecrys(filename::String)
Read filename, return crystal object. File can be POSCAR, or simple quantum espresso inputfile
The entire file can be in the string instead of the filename. The code decides which is which by looking for newlines"
makecrys(lines::Array{String,1})
Read string array, return crystal object. string array can be POSCAR, or simple quantum espresso inputfile
ThreeBodyTB.scf_energy
— Functionscf_energy(c::crystal)
Calculate energy, force, and stress for a crystal. Does self-consistent-field (SCF) calculation if using self-consistent electrostatics.
returns energy, tight-binding-crystal-object, error-flag
Arguments
c::crystal
: the structure to calculate on. Only required argument.database=missing
: Source of coeficients. Will be loaded from pre-fit coefficients if missing.smearing=0.01
: Gaussian smearing temperature, in Ryd. Usually can leave as default.grid=missing
: k-point grid, e.g. [10,10,10], default chosen automaticallyconv_thr = 1e-5
: SCF convergence threshold (Ryd).iter = 75
: number of iterations before switch to more conservative settings.mix = -1.0
: initial mixing. -1.0 means use default mixing. Will automagically adjust mixing if SCF is failing to converge.mixing_mode =:pulay
: default is Pulay mixing (DIIS). Other option is :simple, for simple linear mixing of old and new electron-density. Will automatically switch to simple if Pulay fails.
scf_energy(d::dftout)
SCF energy using crystal structure from DFT object.
scf_energy(tbc::tbc_crys)
SCF energy using crystal structure from TBC object.
ThreeBodyTB.scf_energy_force_stress
— Functionscf_energy_force_stress(c::crystal; database = missing, smearing = 0.01, grid = missing)
Calculate energy, force, and stress for a crystal.
return energy, force, stress, tight_binding_crystal_object
Arguments
c::crystal
: the structure to calculate on. Only required argument.database=missing
: Source of coeficients. Will be loaded from pre-fit coefficients if missing.smearing=0.01
: Gaussian smearing temperature, in Ryd. Usually can leave as default.grid=missing
: k-point grid, e.g. [10,10,10], default chosen automatically
scf_energy_force_stress(tbc::tb_crys; database = missing, smearing = 0.01, grid = missing)
Calculate energy, force, and stress for a tight binding crystal object. This allows the calculation to run without re-doing the SCF calculation. Assumes SCF already done!
returns energy, force, stress, tightbindingcrystal_object
ThreeBodyTB.relax_structure
— Functionrelax_structure(c::crystal; mode="vc-relax")
Find the lowest energy atomic configuration of crystal c
.
Arguments
c::crystal
: the structure to relax, only required argumentmode="vc-relax"
: Default (variable-cell relax) will relax structure and cell, anything else will relax structure only.database=missing
: coefficent database, default is to use the pre-fit pbesol databasesmearing=0.01
: smearing temperature (ryd), default = 0.01grid=missing
: k-point grid, e.g. [10,10,10], default chosen automaticallynsteps=100
: maximum iterationsupdate_grid=true
: update automatic k-point grid during relaxationconv_thr = 2e-3
: Convergence threshold for gradientenergy_conv_thr = 2e-4
: Convergence threshold for energy in Ryd
ThreeBodyTB.TB.read_tb_crys
— Functionfunction read_tb_crys(filename, tbc::tb_crys)
Reads and returns from filename
a tb_crys
object. See write_tb_crys
If cannot find "filename"
, will look for "filename.xml"
, "filename.gz"
, "filename.xml.gz"
Can read gzipped files directly.
ThreeBodyTB.TB.write_tb_crys
— Functionfunction write_tb_crys(filename, tbc::tb_crys)
Writes to filename
a tb_crys
object, using xml formatting. See read_tb_crys
Plotting
ThreeBodyTB.BandStruct.plot_bandstr
— Functionfunction plot_bandstr(h::tb_crys; kpath, names = missing, proj_types=missing, proj_orbs = missing, proj_nums=missing)
Plot the band structure of a tb_crys
object. Can also perform a projected band structure if you specify at least one of proj_types
, proj_orbs
, proj_nums
.
k-path specified by a kpath array and names.
Must do scf calculation before plotting.
Arguments
h::tb_crys
- The tight-biding object we want to plot bands from. Only required argument.kpath=[0.5 0 0 ; 0 0 0; 0.5 0.5 0.5; 0 0.5 0.5; 0 0 0 ;0 0 0.5]
-nk
× 3 array k-point path (high symmetry points).npts=30,
- number of points between high-symmetry k-points.names=missing
-nk
string array. Names of the high-symmetry k-pointsproj_types=missing
- types to project onto. Eitherproj_types="H"
orproj_types=["H", "O"]
are valid.proj_orbs=missing
- orbitals to project onto. eitherproj_orbs=:s
orproj_orbs=[:s, :p]
.proj_nums=missing
- atom numbers to project onto. Eitherproj_nums=1
orproj_nums=[1, 2]
efermi=missing
- allows you to specify fermi energy. Default is to take fromh
color="blue"
- specify line colorMarkerSize=missing"
- specify markersizeyrange=missing"
- specify y-range. e.g.yrange=[-0.7, 0.3]
plot_hk=false
- plot things besides the normal band structure. Can be one of:Seig, :Heig, :Hreal, :Himag, :Sreal, :Simag
to plot H or S eigvals or components. Primarily for debugging.align="vbm"
- default or"valence"
is to align valence band max to zero energy. Can also be"min"
, which aligns on the minimum eigenvalue, or"fermi"
or"ef"
, which align on the Fermi level,clear_pervious=true
- clears the plot before adding new stuff.do_display=true
- display the plot. Iffalse
, can be used with display-less nodes. You can still usesavefig
fromPlots
to produce saved images.
function plot_bandstr(h::tb)
Plots using tb
ThreeBodyTB.DOS.dos
— Functionfunction dos(tbc::tb_crys; grid=missing, npts=missing, proj_type=missing, do_display=true)
DOS, using tetrahedral integration
grid
is the k-point grid. Defaults to 1.6 times the default grid for energy integration.npts
is number of energiesproj_type
can be"none"
,"atomic"
, or"orbital"
. Defaults toatomic
if more than one atom type.do_display=false
will suppress plotting
return energies, dos, projected_dos, pdos_names
ThreeBodyTB.BandStruct.plot_compare_dft
— Functionfunction plot_compare_dft(tbc, bs; tbc2=missing)
Plots a band structure comparison between a tight-binding crystal object (tb_crys
) and a band structure directly from dft (either a dftout
or bs
object).
The k-points are fixed by the bs
object.
tbc2
is an optional second tbc_crys
.
ThreeBodyTB.BandStruct.plot_compare_tb
— Functionfunction plot_compare_tb(h1::tb_crys, h2::tb_crys; h3=missing)
Plot a comparison between different tight binding objects h1
, h2
, and optionally h3
. Options similar to plot_bandstr
but more limited.
function plot_compare_tb(h1::tb, h2::tb; h3=missing)
ThreeBodyTB.BandStruct.band_summary
— Functionfunction band_summary(tbc, kgrid, fermi=missing)
Produces summary of band structure. See below functions for more specific versions of function that automatically generate the k-points.
Note: gaps are not well-defined for non-magnetic systems with odd numbers of electrons, as they are required to be metals.
Returns direct_gap, indirect_gap, gaptype, bandwidth
-direct_gap
: minimum gap at one k-point between nominally filled and empty bands. Can be non-zero in metals. -indirect_gap
: LUMO - HOMO. Can be negative if material has a direct gap everywhere, but the conduction band at some k-point is below the valence band at a different k-point. Physically these are indirect gap semimetals. -gaptype
: is :metal
for all metals, :direct
or :indirect
for insulators. -bandwidth
: HOMO - minimumbandenergy. Included semicore states if they are in the TB calculation.
function band_summary(tbc::tb_crys; kgrid=missing, kpts=missing)
Will automatically generate standard k-grid by default.
function band_summary(tbc::tb_crys_kspace)
Will use internal k-points by default.
Utility
ThreeBodyTB.set_units
— Functionfunction set_units(;energy=missing, length=missing, both=missing)
Set global units for energy/length. Run with no arguments to check/return current units.
- Default units are
"eV"
and"Angstrom"
(or"Å"
or"Ang"
). - Choose atomic units by setting
energy="Ryd"
andlength="Bohr"
. - Set both at the same time with
both="atomic"
orboth="eVAng"
- Internally, all units are atomic. Only main public facing functions actually change units.
ThreeBodyTB.TB.Hk
— Functionfunction Hk(h::tbcryskspace, kpoint)
Calculate band structure at a k-point from a tb_crys_kspace
object. Note, can only return precalculated k-points. Need real-space version to get arbitrary k-points.
#Returns
vect
- Eigenvectors numwan × numwan complex matrix at kpointvals
- Eigenvalues (num_wan)hk
- Hamiltonian at kpointsk
- Overlap matrix at kpointvals0
- <vect | Hk0 | vect> where Hk0 is the non-scf part of the Hamiltonian.
#Arguments
h::tb_crys_kspace
- tbcryskspace objectkpoint
- e.g.[0.0,0.0,0.0]
scf=missing
- default is to take SCF from h.
function Hk(h::tb_k, kpoint)
Calculate band structure at a k-point from tb_k
. Must be pre-calculated k-point.
function Hk(hk,sk, h::tb, kpoint)
Hk function with pre-allocated memory hk, sk
function Hk(h::tb, kpoint)
Calculate band structure at a k-point from tb
function Hk(h::tb_crys, kpoint)
Calculate band structure at a k-point from tb_crys
ThreeBodyTB.TB.calc_bands
— Functionfunction calc_bands(tbc::tb_crys, kpoints::Array{Float64,2})
Calculate bandstructure for k-points from k-point array. Returns eigenvalues.
Arguments
tbc::tb_crys
- The tight binding objectkpoints::Array{Float64,2}
- k-point array. e.g.[0.0 0.0 0.0; 0.0 0.0 0.1]
function calc_bands(h, kpoints::Array{Float64,2})
Calculate bandstructure for k-points from k-point array. h is a tb
or tb_k
object.
ThreeBodyTB.CalcTB.calc_tb_fast
— Functioncalc_tb_fast(crys::crystal; database=missing, use_threebody=true, use_threebody_onsite=true)
Construct tb_crys
from crystal stucture, but does not solve. This is usually called internally by functions like scf_energy
, but you can use it directly if you want. Until you do a SCF energy calculation, the electron density and Fermi level will be wrong.
Arguments
crys::crystal
- Required crystal structuredatabase=missing
- Source of coefficients. Will load from default source if not specified.use_threebody=true
- Use three-body off-site interactions. Only turn off for testing purposes.use_threebody_onsite=true
- Use three-body on-site interactions. Only turn off for testing purposes.verbose=true
- set to false for less output.var_type=missing
- variable type oftb_crys
. Default isFloat64
.