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).sparse = :auto
: Default is to use dense matricies fornat < 100
. Can betrue
orfalse
to force choice.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. Starting default is smaller for larger unit cells.mixing_mode = :simple
: default is simple. Other options are:simple
and:DIIS
/:pulay
(direct inversion of iterative subspace). Will automatically switch to simple if Pulay fails.
scf_energy(d::dftout)
SCF energy using crystal structure from DFT object.
scf_energy(tbc::tb_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 automaticallysparse = :auto
: Default is to use dense matricies fornat < 100
. Can betrue
orfalse
to force choice.
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 Rydsparse = :auto
: Default is to use dense matricies fornat < 100
. Can betrue
orfalse
to force choice.
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
ThreeBodyTB.Symmetry.get_symmetry
— Functionfunction get_symmetry(c::crystal; verbose=true, sym_prec = 5e-4, magmoms=missing)
Return space group number, print other symmetry information if verbose=true
. Based on Spglib.jl
ThreeBodyTB.Symmetry.get_standard_crys
— Functionfunction get_standard_crys(c::crystal; sym_prec = 5e-4, magmoms=missing, to_primitive=true)
Return standardized crystal structure, per spglib convention. Check if space-group is correct.
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.BandStruct.plot_bandstr_sym
— Functionfunction plot_bandstr_sym(c::crystal; sym_prec = 5e-4, npts=30, efermi = missing, color="blue", MarkerSize=missing, yrange=missing, plot_hk=false, align = "vbm", proj_types = missing, proj_orbs = missing, proj_nums=missing, clear_previous=true, do_display=true, color_spin = ["green", "orange"], spin = :both, nspin = 1, database=missing)
Plots the band structure using a set of K-points determined by the space group. Will standardize the crystal structure and run the SCF calculation automatically. Conventions similar to Setyawan and Curtarolo CMS 2010 as well as the jarvis-tools package. Other options similar to plot_bandstr
function plot_bandstr_sym(tbc::tb_crys; sym_prec = 5e-4, npts=30, efermi = missing, color="blue", MarkerSize=missing, yrange=missing, plot_hk=false, align = "vbm", proj_types = missing, proj_orbs = missing, proj_nums=missing, clear_previous=true, do_display=true, color_spin = ["green", "orange"], spin = :both, nspin = 1, database=missing)
Plots the band structure using a set of K-points determined by the space group. This version takes in a tight-binding crystal object from a pervious SCF calculation. If the crystal is in the standardized unit cell, then it will not repeat the SCF calculation, but otherwise it will.
ThreeBodyTB.BandStruct.plot_bandstr_dos
— Functionfunction plot_bandstr_dos(c::crystal; sym_prec = 5e-4, npts=30, efermi = missing,color="blue", MarkerSize=missing,yrange=missing,plot_hk=false, align = "fermi", proj_types = missing, proj_orbs = missing, proj_nums=missing, clear_previous=true, do_display=true, color_spin = ["green", "orange"], spin = :both, nspin = 1, database=missing, smearing = 0.025)
Run SCF in standardized crystal structure and calculate band structure along symmetry-derived k-lines, as well as the DOS. Then plot them layed out side-by-side.
See plot_bandstr_sym
and dos
function plot_bandstr_dos(tbc::tb_crys)
Plot symmetry-derived k-lines and DOS layout out together. This version takes in tb_crys
object from a previous SCF calculation. Will need re-run SCF if tbc.crys
is not standardized.
See plot_bandstr_sym
and dos
ThreeBodyTB.DOS.dos
— Functionfunction dos(tbc::tb_crys; grid=missing, smearing=0.04, npts=missing, proj_type=missing, do_display=true)
DOS, mostly for testing.
npts
is number of energiesproj_type
can be"none"
,"atomic"
, or"orbital"
do_display=false
will suppress the actual plot
The combination of smearing
and grid
are important to get converged results.
See also dos
return energies, dos, projected_dos, pdos_names
ThreeBodyTB.DOS.dos_realspace
— Functionfunction dos_realspace(tbc; direction=3, grid=missing, smearing=0.005, npts=missing, do_display=true, use_sym=false, energy_grid=100, width=missing, energy_lims=missing, scissors_shift=0.0, scissors_shift_atoms=[])
Makes an atom and spatially resoloved DOS-style plot. Along direction=1,2,3
(lattice vector 1,2,3) the atoms are plotted, using transparancy to denote DOS amplitude. Try it out to see. Potentially useful for interfaces. Example usage: c = makecrys([8.0 0 0; 0 8.0 0; 0 0 10.0], [0 0 0; 0 0 0.5], [:Li, :Cl]) #quasi 1D LiCl chain
en, tbc, flag = scf_energy(c*[1,1,10]) #ten unit cells
dos_realspace(tbc, direction=3)
Adjust the energy_grid
and width
to adjust plotting parameters. Scissors shift in energy unit to certain atoms will add scissors shift to open band gaps is desired.
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; scf=missing, spin=1)
Calculate band structure at a k-point from tb_k
. Must be pre-calculated k-point.
function Hk(hk,sk, h::tb, kpoint; spin=1)
Hk function with pre-allocated memory hk, sk
function Hk(h::tb, kpoint; spin=1)
Calculate band structure at a k-point from tb
function Hk(h::tb_crys, kpoint; spin=1)
Calculate band structure at a k-point from tb_crys
ThreeBodyTB.TB.calc_bands
— Functionfunction calc_bands(tbc::tb_crys, kpoints::Array{Float64,2}; spin=1)
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]
spin::Int
- spin component (1 or 2) if magnetic, only 1 is non-sp.
function calc_bands(h, kpoints::Array{Float64,2}; spin=1)
Calculate bandstructure for k-points from k-point array. h is a tb
or tb_k
object.
ThreeBodyTB.set_bin_dirs
— Functionfunction set_bin_dirs(;qe=missing, mpi=missing, pseudodir=missing, templatedir=missing, wannier=missing)
Set directories where things like quantum espresso bins are located. If run with everything missing, instead print current dirs.
qe
- set bin directory of quantum espresso. Needs pw and pp installed. No useful default.mpi
- set mpi command. something like "mpirun -np "pseudo
- set directory for pseudopotentials. Default is to use ../pseudo/gbrv_pbesol/ as distributed with this codetemplate
- set directory for quantum espresso template files. Uses ../template_inputs/ by default.