Wannierize

GU_to_GX_GY(G, X, Y, frozen)

Compute dΩ/dX and dΩ/dY from dΩ/dU.

Acutally they are the conjugate gradients, e.g., $\frac{d \Omega}{d U^*}$.

Arguments

• G: n_bands * n_wann * n_kpts array for gradient dΩ/dU
• X: n_wann * n_wann * n_kpts array for X
• Y: n_bands * n_wann * n_kpts array for Y
• frozen: n_bands * n_kpts BitMatrix for frozen bands

U_to_X_Y(U::Array{T,3}, frozen::BitMatrix) where {T<:Complex}

Convert the U layout to the (X, Y) layout.

See also X_Y_to_U.

Arguments

• U: n_bands * n_wann * n_kpts
• frozen: n_bands * n_kpts

XY_to_X_Y(XY::Matrix{T}, n_bands::Int, n_wann::Int)

Convert the XY layout to the (X, Y) layout.

See also X_Y_to_U.

Arguments

• XY: n_bands * n_wann * n_kpts contiguous array
• n_bands: number of bands, to be used to reshape XY
• n_wann: number of wannier functions, to be used to reshape XY

X_Y_to_U(X::Array{T,3}, Y::Array{T,3})

Convert the (X, Y) layout to the U layout.

There are three formats: U, (X, Y), and XY stored contiguously in memory. For each kpoint,

• U: size(U) = (n_bands, n_wann), the format used in the rest of the code
• (X, Y): size(X) = (n_wann, n_wann), size(Y) = (n_bands, n_wann), intermediate format
• XY: this is the format used in the optimizer

X_Y_to_XY(X::Array{T,3}, Y::Array{T,3}) where {T<:Complex}

Convert the (X, Y) layout to the XY layout.

See also X_Y_to_U.

check_frozen_bands(frozen_bands, n_wann)

Sanity check that the number of frozen bands at each kpoint <= n_wann.

Arguments

• frozen_bands: the BitMatrix of frozen bands
• n_wann: the number of wannier functions

disentangle(model; random_gauge=false, f_tol=1e-7, g_tol=1e-5, max_iter=200, history_size=3)

Run disentangle on the Model.

Arguments

• model: model

Keyword arguments

• random_gauge: use random U matrices as initial guess
• f_tol: tolerance for spread convergence
• g_tol: tolerance for gradient convergence
• max_iter: maximum number of iterations
• history_size: history size of LBFGS

get_fg!_disentangle(model::Model)

Return a tuple of two functions (f, g!) for spread and gradient, respectively.

get_frozen_bands(E, dis_froz_max, dis_froz_min)

Generate a BitMatrix of frozen bands by checking the two frozen windows.

Arguments

• E: the energy eigenvalues of the Hamiltonian
• dis_froz_max: the upper bound of the frozen window
• dis_froz_min: the lower bound of the frozen window

get_frozen_proj(E, U, dis_proj_max)

Get frozen bands according to band projectability.

Arguments

• E: the energy eigenvalues of the Hamiltonian
• U: the gauge rotation matrices
• dis_proj_max: the upper bound projectability. Bands with projectability >= dis_proj_max are frozen.

omega(bvectors, M, X, Y)

Compute WF spread in the (X, Y) layout.

Arguments

• bvectors: bvecoters
• M: n_bands * n_bands * * n_bvecs * n_kpts overlap array
• X: n_wann * n_wann * n_kpts array
• Y: n_bands * n_wann * n_kpts array

omega(model, X, Y)

Compute WF spread in the (X, Y) layout.

Arguments

• model: Model
• X: n_wann * n_wann * n_kpts array
• Y: n_bands * n_wann * n_kpts array

omega_grad(bvectors, M, X, Y, frozen)

Compute gradient of WF spread in the (X, Y) layout.

Arguments

• bvectors: bvecoters
• M: n_bands * n_bands * * n_bvecs * n_kpts overlap array
• X: n_wann * n_wann * n_kpts array
• Y: n_bands * n_wann * n_kpts array
• frozen: n_bands * n_kpts array for frozen bands

orthonorm_freeze(U, frozen)

Normalize and freeze a block of a matrix.

Conditions:

• Block form: U = vcat(Uf, Ur)
• Semiunitary
  • U' * U = I
  • Uf' * Uf + Ur' * Ur = I
• Frozen: Uf * Uf' = I
• Also: Uf * Ur' = 0

Strategy:

1. orthogonalize Uf
2. project Uf out of Ur
3. orthogonalize the range of Ur

Arguments

• U: the matrix to be orthonormalized and frozen
• frozen: the BitVector specifying which bands are frozen

set_frozen_degen!(frozen_bands, E, atol=1e-4)

Freeze bands which are degenerate with the frozen_bands.

In some cases, we might want to freeze the whole set of degenerated eigen vectors.

Arguments

• frozen_bands: the BitMatrix of frozen bands
• E: the energy eigenvalues of the Hamiltonian
• atol: the tolerance of degeneracy

set_frozen_proj!(model, dis_proj_max; degen=false, degen_atol=1e-4)

Set frozen bands of the Model according to projectability.

Arguments

• model: the Model to be set
• dis_proj_max: the upper bound projectability. Bands with projectability >= dis_proj_max are frozen.

Keyword Arguments

• degen: whether to freeze the whole set of degenerated eigen vectors
• degen_atol: the tolerance of degeneracy

set_frozen_win!(model, dis_froz_max, dis_froz_min=-Inf; degen=false, degen_atol=1e-4)

Set frozen bands of the Model according to two energy windows.

Arguments

• model: the Model to be set
• dis_froz_max: the upper bound of the frozen window
• dis_froz_min: the lower bound of the frozen window

Keyword Arguments

• degen: whether to freeze the whole set of degenerated eigen vectors
• degen_atol: the tolerance of degeneracy

zero_froz_grad!(G, frozen)

Set gradient of frozen bands to 0.

This is used in test.

Arguments

• G: gradient of the spread, in XY layout
• frozen: BitMatrix for frozen bands, n_bands * n_kpts

get_fg!_maxloc(model::Model)

Return a tuple of two functions (f, g!) for spread and gradient, respectively.

max_localize(model; f_tol=1e-7, g_tol=1e-5, max_iter=200, history_size=3)

Maximally localize spread functional w.r.t. all kpoints on a unitary matrix manifold.

Arguments

• model: model

Keyword arguments

• f_tol: tolerance for spread convergence
• g_tol: tolerance for gradient convergence
• max_iter: maximum number of iterations
• history_size: history size of LBFGS

compute_error(model, U::Array{Complex{T},3})

Compute the smoothness error of the gauge.

parallel_transport(model::Model{T}; use_U=false, log_interp=false)

Parallel transport the gauge from the first kpoint to all other kpoints.

Assumptions:

• the kpoints are contained in a N1 * N2 * N3 cartesian grid
• the neighbor list must contain the six cartesian neighbors along x, y, z directions

Arguments

• model: model

Keyword arguments

• use_U: use the gauge U instead of random matrix
• log_interp: use logarithmic interpolation method

Choose the column and the target point to contract the vector path of the columns of matrix_path

Create an interpolation path between x and y, unit vectors by normalizing the linear interpolation of parameter t[i].

Compute the parallel transport of the matrix path (matrixpath) along the homotopy of the frame (framepath) that starts with the columns (columns[]) of matrix_path.

size(framepath) = (ncol, ncol, nk, nt) size(matrixpath) = (ncol, ncol, n_k)

Contract the matrix_path to a single matrix point in the space of unitaries using parallel transport.

i.e., contract Obs(k) matrices to constant vectors.

size(matrixpath) = nwann x nwann x nk t: vector of kpoint indexes along a different k direction

Propagate U, defined at the first kpt, to the given list of kpts. Those must be neighbors, and only the first kpoint is assumed to have been rotated.

Arguments

• U: n_wann * n_wann * n_kpts gauge matrices
• kpts: list of kpoints along which to propagate the gauge matrices
• dk: the distance between two consecutive kpoints, in fractional coordinates
• M: overlap matrices
• kpoints: BVectors.kpoints, kpoint coordinates of the grid
• kpb_k: BVectors.kpb_k
• kpb_b: BVectors.kpb_b

get_fg!_rotate(model::Model)

Return a tuple of two functions (f, g!) for spread and gradient, respectively.

opt_rotate(model; f_tol=1e-7, g_tol=1e-5, max_iter=200, history_size=3)

Maximally localize spread functional w.r.t. single unitary matrix W.

Arguments

• model: model

Keyword arguments

• f_tol: tolerance for spread convergence
• g_tol: tolerance for gradient convergence
• max_iter: maximum number of iterations
• history_size: history size of LBFGS

rotate_U(U::Array{T,3}, W::Matrix{T}) where {T<:Complex}

Rotate the U matrices at each kpoint by the same W matrix.

$\forall \bm{k}$, $U_{\bm{k}} W$

Useful once we have the optimal rotation matrix W, then update the initial U matrices by rotating them by W.

split_eig(E, U, eig_groups)

Split eigenvalues into several groups.

The separation is done by

1. construct Wannier gauge Hamiltonian, $H_{\bm{k}} = U_{\bm{k}}^\dagger [\epsilon_{n \bm{k}}] U_{\bm{k}}$
2. diagonalize the Hamiltonian, the eigenvalues are sorted in ascending order, and they are split into several groups according to the indices in eig_groups.

Arguments

• E: eigenvalues
• U: (semi-)Unitary matrices gauge transformation
• eig_groups: a Vector, in which each element is a Vector of indices of eigenvalues that belong to the same group.

split_eig(E, U, n_val)

Split eigenvalues into two groups.

The separation is done by

1. construct Wannier gauge Hamiltonian, $H_{\bm{k}} = U_{\bm{k}}^\dagger [\epsilon_{n \bm{k}}] U_{\bm{k}}$
2. diagonalize the Hamiltonian, the eigenvalues are sorted in ascending order, so that the first n_val eigenvalues are the occupied states, and the rest are the unoccupied states.

Arguments

• E: eigenvalues
• U: (semi-)Unitary matrices gauge transformation
• n_val: number of valence states

split_model(model, n_val)

Split the Model into two Models.

Arguments

• model: the Model to be split
• n_val: number of valence WFs

split_model(model, eig_groups)

Split the Model into several Models.

Arguments

• model: the Model to be split
• eig_groups: a Vector, in which each element is a Vector of indices of eigenvalues that belong to the same group.

split_unk(dir, Us, outdirs)

Rotate UNK files.

These are large matrices, so we read/write to disk for each kpoint sequentially, inside the function.

Arguments

• dir: directory where UNK files are stored
• Us: a Vector of Wannier (semi-)unitary gauge matrices for rotating band groups
• outdirs: a Vector of output directories for each band group

Keyword arguments

• binary: whether to write in Fortran binary format

split_unk(dir, Uv, Uc, outdir_val="val", outdir_cond="cond")

Rotate UNK files.

These are large matrices, so we read/write to disk for each kpoint sequentially, inside the function.

Arguments

• dir: directory where UNK files are stored
• Uv: the Wannier (semi-)unitary matrix for rotating valence bands
• Uc: the Wannier (semi-)unitary matrix for rotating conduction bands
• outdir_val: output directory for valence bands
• outdir_cond: output directory for conduction bands

Keyword arguments

• binary: whether to write in Fortran binary format

split_wannierize(model::Model, n_val::Int)

Split the model and run parallel transport to smoothen the gauge.

Arguments

• model: the Model to be split
• n_val: number of valence WFs

split_wannierize(model::Model, eig_groups)

Split the model and run parallel transport to smoothen the gauge.

Arguments

• model: the Model to be split
• eig_groups: a Vector, in which each element is a Vector of indices of eigenvalues that belong to the same group.

Model for spin polarized system with constraint.

Traditionally, we run two independent Wannierizations for spin up and spin down. Here we add a constraint to maximally overlap the spin-up and spin-down WFs, so that they map one-by-one to each other.

disentangle(model; f_tol=1e-7, g_tol=1e-5, max_iter=200, history_size=3)

Run disentangle on a MagModel.

Arguments

• model: MagModel
• λ: Lagrange multiplier of the ↑↓ overlap term

Keyword arguments

• f_tol: tolerance for spread convergence
• g_tol: tolerance for gradient convergence
• max_iter: maximum number of iterations
• history_size: history size of LBFGS

get_fg!_disentangle(model::MagModel)

Return a tuple of two functions (f, g!) for spread and gradient, respectively.

omega_updn(M)

Compute QPPM Eq. 8.

Arguments

• M: the overlap matrix between up and down WFs, size: (nwann, nwann), should be the matrix returned from overlap_updn.

overlap_updn(up::Model{T}, dn::Model{T}, M::Array{Complex{T},3}) where {T<:Real}

Compute the overlap between up and down WFs.

Actually N - Ω↑↓, according to QPPM Eq. 8, where N = n_wann.

Arguments

• M: the MagModel.M matrices, size (nbands, nbands, n_kpts)
• Uup: the up gauge matrices, size: (nbands, nwann, n_kpts)
• Udn: the down gauge matrices, size: (nbands, nwann, n_kpts)

overlap_updn_grad(model::MagModel, Uup, Udn)

Compute gradients of overlap_updn.

$\frac{d \Omega}{d U^{\uparrow}}$ and $\frac{d \Omega}{d U^{\downarrow}}$.

TODO: this is actually the gradient of Tr[overlap_updn]

Arguments

• M: the MagModel.M matrices, size (nbands, nbands, n_kpts)
• Uup: the up gauge matrices, size: (nbands, nwann, n_kpts)
• Udn: the down gauge matrices, size: (nbands, nwann, n_kpts)

overlap_updn_grad(model::MagModel, Xup, Yup, Xdn, Ydn)

Compute gradient of overlap_updn.

$\frac{d \Omega}{d X^{\uparrow}}$, $\frac{d \Omega}{d Y^{\uparrow}}$, $\frac{d \Omega}{d X^{\downarrow}}$, $\frac{d \Omega}{d Y^{\downarrow}}$.

disentangle(model, r₀, λ=1.0; f_tol=1e-7, g_tol=1e-5, max_iter=200, history_size=3)

Run disentangle on the Model with center penalty.

Arguments

• model: model
• r₀: 3 * n_wann, WF centers in cartesian coordinates
• λ: penalty strength

Keyword arguments

• f_tol: tolerance for spread convergence
• g_tol: tolerance for gradient convergence
• max_iter: maximum number of iterations
• history_size: history size of LBFGS

get_fg!_center_disentangle(model, r₀, λ=1.0)

Return a tuple of two functions (f, g!) for spread and gradient, respectively.

omega_center(bvectors, M, X, Y, r₀, λ)

Compute WF spread with center penalty, in the (X, Y) layout.

Arguments

• bvectors: bvecoters
• M: n_bands * n_bands * * n_bvecs * n_kpts overlap array
• X: n_wann * n_wann * n_kpts array
• Y: n_bands * n_wann * n_kpts array
• r₀: 3 * n_wann, WF centers in cartesian coordinates
• λ: penalty strength

omega_center(mode, X, Y, r₀, λ)

Compute WF spread with center penalty, in the (X, Y) layout.

Arguments

• model: Model
• X: n_wann * n_wann * n_kpts array
• Y: n_bands * n_wann * n_kpts array
• r₀: 3 * n_wann, WF centers in cartesian coordinates
• λ: penalty strength

omega_center_grad(bvectors, M, X, Y, frozen, r₀, λ)

Compute gradient of WF spread with center penalty, in the (X, Y) layout.

Arguments

• bvectors: bvecoters
• M: n_bands * n_bands * * n_bvecs * n_kpts overlap array
• X: n_wann * n_wann * n_kpts array
• Y: n_bands * n_wann * n_kpts array
• frozen: n_bands * n_kpts array for frozen bands
• r₀: 3 * n_wann, WF centers in cartesian coordinates
• λ: penalty strength

struct SpreadCenter

A struct containing both Spread and WF center penalty.

get_fg!_center_maxloc(model, r₀, λ=1.0)

Return a tuple of two functions (f, g!) for spread and gradient, respectively.

max_localize_center(model, r₀, λ=1.0; f_tol=1e-7, g_tol=1e-5, max_iter=200, history_size=3)

Maximally localize spread functional with center constraint on a unitary matrix manifold.

Arguments

• model: model
• r₀: 3 * n_wann, WF centers in cartesian coordinates
• λ: penalty strength

Keyword arguments

• f_tol: tolerance for spread convergence
• g_tol: tolerance for gradient convergence
• max_iter: maximum number of iterations
• history_size: history size of LBFGS

omega_center(bvectors, M, U, r₀, λ)

Compute WF spread with center penalty, for maximal localization.

Arguments

• bvectors: bvecoters
• M: n_bands * n_bands * * n_bvecs * n_kpts overlap array
• U: n_wann * n_wann * n_kpts array
• r₀: 3 * n_wann, WF centers in cartesian coordinates
• λ: penalty strength

omega_center(model, U, r₀, λ)

Compute WF spread with center penalty, for maximal localization.

Arguments

• model: Model
• U: n_wann * n_wann * n_kpts array
• r₀: 3 * n_wann, WF centers in cartesian coordinates
• λ: penalty strength

omega_center(model, r₀, λ)

Compute WF spread with center penalty, for maximal localization.

Arguments

• model: Model
• r₀: 3 * n_wann, WF centers in cartesian coordinates
• λ: penalty strength

omega_center_grad(bvectors, M, U, r₀, λ)

Compute gradient of WF spread with center penalty, for maximal localization.

Arguments

• bvectors: bvecoters
• M: n_bands * n_bands * * n_bvecs * n_kpts overlap array
• U: n_wann * n_wann * n_kpts array
• r₀: 3 * n_wann, the target WF centers in cartesian coordinates
• λ: penalty strength

omega_center_grad(bvectors, M, U, r, r₀, λ)

Compute gradient of WF spread with center penalty, for maximal localization.

Arguments

• bvectors: bvecoters
• M: n_bands * n_bands * * n_bvecs * n_kpts overlap array
• U: n_wann * n_wann * n_kpts array
• r: 3 * n_wann, the current WF centers in cartesian coordinates
• r₀: 3 * n_wann, the target WF centers in cartesian coordinates
• λ: penalty strength

disentangle(model; f_tol=1e-7, g_tol=1e-5, max_iter=200, history_size=3)

Run disentangle on a MagModel, with center constraints.

Arguments

• model: MagModel
• r₀: 3 * n_wann, WF centers in cartesian coordinates
• λc: Lagrange multiplier of the center term
• λs: Lagrange multiplier of the spin-up and spin-down overlap term

Keyword arguments

• f_tol: tolerance for spread convergence
• g_tol: tolerance for gradient convergence
• max_iter: maximum number of iterations
• history_size: history size of LBFGS

get_fg!_center_disentangle(model::MagModel, r₀, λc, λs)

Return a tuple of two functions (f, g!) for spread and gradient, respectively.