Functions#
Support functions for spinWannier.
Utility functions for the spinWannier package.
- spinWannier.wannier_utils.S_xyz_expectation_values(eigvecs)#
For a k-space Hamiltonian Hk (function of kx and ky), calculate the eigen-energies and expectation values of spin Sx, Sy, Sz, at all the k-points in ‘kpoints’.
- Parameters:
Hk (function) – A function returning a 2x2 FEG Hamiltonian.
kpoints (array of tuples) – Array of k-points.
- Returns:
array of expectation values for energy and spin expectation values
- Return type:
four arrays
- spinWannier.wannier_utils.W_gauge_to_H_gauge(kpoints, O_mn_k_W_matrices, U_mn_k_matrices=[], hamiltonian=True)#
- Transform from Wannier gauge to Hamiltonian gauge, i.e., for every k-point either diagonalize Hamiltonian
or use the matrix from previous Hamiltonian diagonalization (‘U_mn_k’) to transform some other operator to H gauge.
- Parameters:
O_mn_k_W (numpy array) – The array (at all k-points) of k-space operators in Wannier gauge.
U_mn_k (numpy array) – The unitary matrices which diagonalize the Hamiltonian for each k-point.
hamiltonian (bool) – If True, the operator is Hamiltonian and the U_mn_k will be determined; for any other operator they have to be provided.
- Returns:
The k-space operator dictionary in Hamiltonian gauge.
- Return type:
dict
- spinWannier.wannier_utils.band_with_spin_projection_under_threshold_for_kpoint(kpoint=(0.0, 0.0, 0.0), orbital_characters_considered={0: [4, 5, 6, 7, 8], 1: [1, 2, 3], 2: [1, 2, 3]}, threshold=0.25, PROCAR_file='PROCAR', n_ions=3, skip_lowest_bands=10)#
Return the lowest band with the spin projection below the threshold for ‘kpoint’.
- Parameters:
kpoint (tuple) – The k-point.
orbital_characters_considered (dict) – The orbital characters considered.
threshold (float) – The threshold.
PROCAR_file (str) – The PROCAR file.
n_ions (int) – The number of ions.
skip_lowest_bands (int) – The number of lowest bands to skip.
- Returns:
The lowest band with the spin projection below the threshold for ‘kpoint’.
- Return type:
int
- spinWannier.wannier_utils.check_file_exists(file_name)#
- Check if file exists. If yes, add him a number in parentheses that does not exist.
Return the new name.
- Parameters:
file_name (str) – The file name.
- Returns:
The new file name.
- Return type:
str
- spinWannier.wannier_utils.coerce_R_vectors_to_basic_supercell(R_tr=(-3, 2, 0), R_mesh_ijk=(5, 5, 1))#
wannier90_hr.dat has sometimes the R vectors shifted (probably to facilitate the ‘minimal replice’ interpolation. Therefore, given some uniform grid centered at zero (e.g. 5x5x1) and a translated R vectors from hr_dat file (e.g. (-3,2,0) ), we would like to get the R vector’s image in the uniform grid, to be able to compare with our results. For (-3,2,0) this would be (2,2,0) in case of (5,5,1) grid (adding 5 to the ‘x’ coordinate). -> Algorithm: translate the ‘R’ by a supercell vector (in negative and positive direction and for x, y, z separately -> 8 possible translations: so searching the immediate vicinity of the R vector. Return a dictionary from the ‘hr_dat grid’ to the uniform grid.
- Parameters:
R_tr (tuple of 3 ints) – The translated R vector.
R_mesh_ijk (tuple of 3 ints) – The uniform grid.
- Returns:
The translated R vector in the uniform grid.
- Return type:
tuple of 3 ints
- spinWannier.wannier_utils.coerce_to_positive_angles(angles)#
Coerce angles to be positive (i.e., in the range [0, 2pi]).
- Parameters:
angles (np.array) – The angles.
- Returns:
The coerced angles.
- Return type:
np.array
- spinWannier.wannier_utils.convert_wann_TB_model_for_QuT(folder_in='./', folder_out='./', seedname_out='system', winfile='wannier90.win', wann_centres_file='wannier90_centres.xyz', hr_R_file='hr_R_dict.pickle', spn_R_file='spn_R_dict.pickle')#
- Convert the outputs of my wannier to TB to Joaquin’s code.
- Generates files
{system_name}_Mx{Mx:d}My{My:d}Mz{Mz:d}_E{E:.2f}.uc …. unit cell {system_name}_Mx{Mx:d}My{My:d}Mz{Mz:d}_E{E:.2f}.xyz …. wannier centers in cartesian coordinates {system_name}_Mx{Mx:d}My{My:d}Mz{Mz:d}_E{E:.2f}_hr.dat …. real-space Hamiltonian {system_name}_Mx{Mx:d}My{My:d}Mz{Mz:d}_E{E:.2f}_sxhr.dat …. real-space spin S_x operator {system_name}_Mx{Mx:d}My{My:d}Mz{Mz:d}_E{E:.2f}_syhr.dat …. real-space spin S_y operator {system_name}_Mx{Mx:d}My{My:d}Mz{Mz:d}_E{E:.2f}_szhr.dat …. real-space spin S_z operator
- Parameters:
winfile (str, optional) – _description_. Defaults to ‘wannier90.win’.
hr_R_file (str, optional) – _description_. Defaults to ‘hr_R_dict.pickle’.
spn_R_file (str, optional) – _description_. Defaults to ‘spn_R_dict.pickle’.
- spinWannier.wannier_utils.dict_to_matrix(data_dict, num_wann, spin_index=False)#
Convert hr or spn dictionary to a matrix accepted by wannierBerri: first index is m, second n, third is the index of an R_vector from iRvec list.
- Parameters:
data_dict (dict) – The dictionary with R vectors as keys and hopping as a complex number matrix as the values.
num_wann (int) – The number of Wannier functions.
spin_index (bool) – If True, the dictionary has spin index.
- Returns:
The array of R vectors. numpy array: The data matrix.
- Return type:
numpy array
- spinWannier.wannier_utils.eigenval_dict(eigenval_file='wannier90.eig', win_file='wannier90.win')#
Return the eigenvalues as a dictionary with the keys being the k-point tuples.
- Parameters:
eigenval_file (str) – The name of the eigenvalues file.
win_file (str) – The name of the wannier90.win file.
- Returns:
The eigenvalues dictionary.
- Return type:
dict
- spinWannier.wannier_utils.eigenval_for_kpoint(kpoint=(0.0, 0.0, 0.0), band=-4, eigenval_file='wannier90.eig', win_file='wannier90.win')#
- Return the eigenvalue at ‘kpoint’ and ‘band’.
!! kpoint must be one from the kpoints in the wannier90.win file !!
- Parameters:
kpoint (tuple) – The k-point.
band (int) – The band.
eigenval_file (str) – The eigenvalue file.
win_file (str) – The win file.
- Returns:
The eigenvalue at ‘kpoint’ and ‘band’.
- Return type:
float
- spinWannier.wannier_utils.fermi_surface_spin_texture(kpoints2D, bands2D, Sx2D, Sy2D, Sz2D, ax=None, E=0, E_F=0, E_thr=0.01, fig_name=None, quiver_scale=1, scatter_for_quiver=True, scatter_size_quiver=1, scatter_size=0.8, reduce_by_factor=1, kmesh_limits=None, colorbar_Sx_lim=[-1, 1], colorbar_Sy_lim=[-1, 1], colorbar_Sz_lim=[-1, 1], n_points_for_one_angstrom_radius=120, ylim_margin=0.1, contour_for_quiver=True, contour_line_width=2.5, arrow_linewidth=0.005, arrow_head_width=3, quiver_angles='xy', quiver_scale_units='xy', inset_with_units_of_arrows=True, color_middle=(0.85, 0.85, 0.85, 1.0), savefig=True, showfig=True)#
Plot scatter points with a spin texture on a constant energy xy surface (probably at E=EF) if the energy difference of each given point is lower than some threshold. If there is more such points, grab the one with minimum difference from the energy surface.
E and E_thr in eV
kpoints2D: 2D array of shape (Nkpoints, 3) bands2D: 2D array of shape (Nkpoints, Nbands)
circumf_distance_of_arrows: distance between the arrows on the circumference of the circle = k (1/Angstrom) * phi (rad)
- Parameters:
kpoints2D (np.array) – The 2D k-points.
bands2D (np.array) – The 2D bands.
Sx2D (np.array) – The 2D Sx.
Sy2D (np.array) – The 2D Sy.
Sz2D (np.array) – The 2D Sz.
ax (plt.axis, optional) – The axis.
E (float, optional) – The energy.
E_F (float, optional) – The Fermi energy.
E_thr (float, optional) – The energy threshold.
fig_name (str, optional) – The figure name.
quiver_scale (float, optional) – The quiver scale.
scatter_for_quiver (bool, optional) – If True, scatter for quiver.
scatter_size_quiver (int, optional) – The scatter size for quiver.
scatter_size (float, optional) – The scatter size.
reduce_by_factor (int, optional) – The reduce by factor.
kmesh_limits (list, optional) – The k-mesh limits.
colorbar_Sx_lim (list, optional) – The colorbar limits for Sx.
colorbar_Sy_lim (list, optional) – The colorbar limits for Sy.
colorbar_Sz_lim (list, optional) – The colorbar limits for Sz.
n_points_for_one_angstrom_radius (int, optional) – The number of points for one Angstrom radius.
ylim_margin (float, optional) – The y-limit margin.
contour_for_quiver (bool, optional) – If True, contour for quiver.
contour_line_width (float, optional) – The contour line width.
arrow_linewidth (float, optional) – The arrow line width.
arrow_head_width (int, optional) – The arrow head width.
quiver_angles (str, optional) – The quiver angles.
quiver_scale_units (str, optional) – The quiver scale units.
inset_with_units_of_arrows (bool, optional) – If True, the inset with units of arrows.
color_middle (tuple, optional) – The middle color.
savefig (bool, optional) – If True, save the figure.
showfig (bool, optional) – If True, show the figure.
- spinWannier.wannier_utils.files_wann90_to_dict_pickle(model_dir='./', disentanglement=False)#
Convert wannier90 files to pickled dictionaries files.
- Parameters:
model_dir (str) – The model directory.
disentanglement (bool) – If True, disentanglement is considered.
- Returns:
The number of bands and the number of wannier functions.
- Return type:
tuple
- spinWannier.wannier_utils.get_2D_kpoint_mesh(G, limits=[-0.5, 0.5], Nk=100)#
Get a 2D mesh of k-points in the reciprocal space.
- Parameters:
G (numpy array) – The reciprocal lattice vectors.
limits (list of 2 floats) – The limits of the mesh.
Nk (int) – The number of k-points in each direction.
- Returns:
The k-points in the reciprocal space. list of tuples: The k-points in the Cartesian space.
- Return type:
list of tuples
- spinWannier.wannier_utils.get_DFT_kgrid(fin='wannier90.win')#
Get the ‘mp_grid’ from wannier90.win file which tells the k-grid used in the DFT calculation.
- Parameters:
fin (str) – The name of the wannier90.win file.
- Returns:
The k-grid.
- Return type:
tuple of 3 ints
- spinWannier.wannier_utils.get_kpoint_names(fwin='wannier90.win')#
Parse wannier90.win file to get k-point names as a list of tuples. The k-point names are in the ‘begin kpoints’ and ‘end kpoints’ section.
- Parameters:
fwin (str) – The name of the wannier90.win file.
- Returns:
The k-point names.
- Return type:
list of tuples
- spinWannier.wannier_utils.get_kpoint_path(kpoint_matrix, G, Nk)#
Interpolate between kpoints in the kpoint_matrix to obtain a k-point path with Nk points in each segment.
- Parameters:
kpoint_matrix (list of tuples) – The k-point matrix.
G (numpy array) – The reciprocal lattice vectors.
Nk (int) – The number of k-points in each segment.
- Returns:
The k-points in the reciprocal space. list of tuples: The k-points in the Cartesian space. list of floats: The k-point distance position.
- Return type:
list of tuples
- spinWannier.wannier_utils.get_skiprows_hr_dat(fin='wannier90_hr.dat')#
Get the number of skiprows in wannier90_hr.dat file.
- Parameters:
fin (str) – The name of the hr.dat file.
- Returns:
The number of skiprows.
- Return type:
int
- spinWannier.wannier_utils.hr_wann90_to_dict(fin_wannier90='wannier90_hr.dat')#
Convert wannier90 hr.dat to dictionary in the form that we are using: R vectors as keys and hopping as a complex number matrix as the values.
- Parameters:
fin_wannier90 (str) – The name of the hr.dat file.
- Returns:
The hopping dictionary.
- Return type:
dict
- spinWannier.wannier_utils.interpolate_operator(operator_dict, u_dis_dict, u_dict, latt_params, reciprocal_latt_params, R_grid, U_mn_k=None, hamiltonian=True, kpoints=[(0, 0, 0)], save_real_space=False, real_space_fname='hr_R_dict.dat', save_folder='./tb_model_wann90/', save_as_pickle=True, verbose=False)#
- Takes operator O (Hamiltonian, spin-operator …) evaluated on a coarse DFT k-mesh and Hamiltonian eigenstates onto ‘kpoints’ using the
(semi)unitary transformation defined by U_dis*U from the Hamiltonian gauge (eigenstate gauge) to the Wannier gauge,
performing Fourier transform to the real-space using a set of lattice vectors
inverse Fourier transforming to the k-space for all the k-points from ‘kpoints’.
- Parameters:
operator_dict (dict) – The operator dictionary.
u_dis_dict (dict) – The disentangled unitary dictionary.
u_dict (dict) – The unitary dictionary.
latt_params (np.array) – The lattice parameters.
reciprocal_latt_params (np.array) – The reciprocal lattice parameters.
R_grid (list) – The R grid.
U_mn_k (dict, optional) – The Hamiltonian eigenstates.
hamiltonian (bool, optional) – If True, the Hamiltonian is calculated.
kpoints (list, optional) – The k-points.
save_real_space (bool, optional) – If True, save the real space.
real_space_fname (str, optional) – The real space file name.
save_folder (str, optional) – The save folder.
save_as_pickle (bool, optional) – If True, save as pickle.
- Returns:
interpolated operator matrices (values) at the given k-points (keys) in the Hamiltonian gauge; unitary matrices which diagonalize the Hamiltonian (only if hamiltonian==True).
- Return type:
dict
- spinWannier.wannier_utils.load_dict(fin='spn_dict.pickle', text_file=False)#
If not text_file, than it’s binary.
- Parameters:
fin (str) – The name of the file to load: spn_dict.pickle / hr_R_dict.pickle / u_dict.pickle / u_dis_dict.pickle spn_dict.dat / hr_R_dict.dat / u_dict.dat / u_dis_dict.dat
text_file (bool) – If True, load as a text file.
- Returns:
The loaded dictionary.
- Return type:
dict
- spinWannier.wannier_utils.load_eigenvals(eigenval_file='wannier90.eig')#
Parse the wannier90 .eig file to get the k-point and band-resolved eigenvalues.
- Parameters:
eigenval_file (str) – The name of the eigenvalues file.
- Returns:
2D list (n_kpoints X n_bands) containing the eigenvalues.
- Return type:
list of lists
- spinWannier.wannier_utils.load_lattice_vectors(win_file='wannier90.win')#
Return a 3x3 matrix where 1st row is the 1st lattice vector etc.
- Parameters:
win_file (str) – The name of the wannier90.win file.
- Returns:
The lattice vectors matrix.
- Return type:
numpy array
- spinWannier.wannier_utils.magmom_OOP_or_IP(INCAR_path)#
Return True if MAGMOM in INCAR_path lies purely OOP, else return False.
- Parameters:
INCAR_path (str) – The path to the INCAR file.
- Returns:
True if MAGMOM in INCAR_path lies purely OOP, else False.
- Return type:
bool
- spinWannier.wannier_utils.magmom_direction(INCAR_path)#
Return 0 / 1 / 2 if spin is along x / y / z direction, respectively.
- Parameters:
INCAR_path (str) – The path to the INCAR file.
- Returns:
0 / 1 / 2 if spin is along x / y / z direction, respectively.
- Return type:
int
- spinWannier.wannier_utils.material_name_to_latex(material)#
Convert the material name to LaTeX format.
- Parameters:
material (str) – The material name.
- Returns:
The LaTeX formatted material name.
- Return type:
str
- spinWannier.wannier_utils.matrix_to_dict(data_matrix, Rvecs, spin_index=False, spin_names=['x', 'y', 'z'])#
- Convert data matrix (e.g. the Ham_R or SS_R file from symmetrization) into a dictionary with
lattice vectors as keys and if spin_index is True then also (‘x’, ‘y’, ‘z’) for spin.
- Parameters:
data_matrix (numpy array) – The data matrix.
Rvecs (numpy array) – The lattice vectors.
spin_index (bool) – If True, the dictionary has spin index.
spin_names (list of 3 str) – The spin names.
- Returns:
The dictionary with R vectors as keys and hopping as a complex number matrix as the values.
- Return type:
dict
- spinWannier.wannier_utils.operator_exp_values(eigvecs, Operator)#
Return the expectation values of Operator (an N x N matrix, where N is the size of the eigenvectors) for all eigenvectors.
- Parameters:
eigvecs (array of matrices) – Array of M square matrices. Each N x N square matrix contains N eigenvectors (columns of the matrix).
Operator (matrix) – N x N Hermition matrix acting on the eigenvectors.
- Returns:
M x N array of expectation values.
- spinWannier.wannier_utils.outer(s, M)#
Outer product between a 2x2 spin matrix and a general orbital nxn M matrix, so that the spin blocks are the big blocks.
- Parameters:
s (np.array) – Pauli (2 x 2) matrix sx, sy or sz
M (np.array) – general orbital (n x n) matrix
- Returns:
The resulting matrix.
- Return type:
np.array
- spinWannier.wannier_utils.outer2(s, M)#
Works for s matrix of any dimension, but slower than ‘outer(s, M)’ Outer product between a nxn spin matrix and a general orbital nxn M matrix, so that the spin blocks are the big blocks.
- Parameters:
s (np.array) – Pauli (2 x 2) matrix sx, sy or sz
M (np.array) – general orbital (n x n) matrix
- Returns:
The resulting matrix.
- Return type:
np.array
- spinWannier.wannier_utils.parse_KPOINTS_file(KPOINTS_file_path)#
Parse the KPOINTS file.
- Parameters:
KPOINTS_file_path (str) – The KPOINTS file path.
- Returns:
The k-point matrix, the number of k-points, and the k-path ticks.
- Return type:
tuple
- spinWannier.wannier_utils.plot_bands_spin_texture(kpoints, kpath, kpath_ticks, Eigs_k, S_mn_k_H_x, S_mn_k_H_y, S_mn_k_H_z, NW, E_F=0, fout='spin_texture_1D_home_made.jpg', fig_caption='Wannier interpolation', yaxis_lim=[-5, 5], savefig=True, showfig=True)#
Output a figure with Sx, Sy, and Sz-projected band structure.
- Parameters:
kpoints (list of tuples) – The k-points in the reciprocal space.
kpath (list of floats) – The k-point distance position.
kpath_ticks (list of floats) – The k-point ticks.
Eigs_k (dict) – The eigenvalues dictionary.
S_mn_k_H_x (dict) – The Sx operator dictionary in Hamiltonian gauge.
S_mn_k_H_y (dict) – The Sy operator dictionary in Hamiltonian gauge.
S_mn_k_H_z (dict) – The Sz operator dictionary in Hamiltonian gauge.
E_F (float) – The Fermi level in eV.
fout (str) – The name of the output file.
fig_caption (str) – The caption of the figure.
yaxis_lim (list of 2 floats) – The y-axis limits.
savefig (bool) – If True, save the figure.
showfig (bool) – If True, show the figure.
- spinWannier.wannier_utils.read_PROCAR_lines_for_kpoint(kpoint=(0.0, 0.0, 0.0), PROCAR_file='PROCAR')#
Get the section of PROCAR belonging to the first occurence of ‘kpoint’ in PROCAR_file.
- Parameters:
kpoint (tuple) – The k-point.
PROCAR_file (str) – The PROCAR file.
- Returns:
The lines of PROCAR belonging to the first occurence of ‘kpoint’.
- Return type:
list
- spinWannier.wannier_utils.real_space_grid_from_hr_dat(fname='wannier90_hr.dat')#
Get the real-space grid from the seedname_hr.dat file. See Pizzi 2020 Sec. 4.2.
- Parameters:
fname (str) – The name of the hr.dat file.
- Returns:
The real-space grid.
- Return type:
list of tuples
- spinWannier.wannier_utils.real_to_W_gauge(kpoints, O_mn_R_W)#
Perform inverse Fourier transformation from real-space operator to k-space in Wannier gauge.
- Parameters:
kpoints (list of tuples) – The k-points.
O_mn_R_W (dict) – The real-space operator dictionary.
- Returns:
The k-space operator dictionary.
- Return type:
dict
- spinWannier.wannier_utils.real_to_W_gauge_accelerated(kpoints, O_mn_R_W)#
Perform inverse Fourier transformation from real-space operator to k-space in Wannier gauge.
- Parameters:
kpoints (list of tuples) – The k-points.
O_mn_R_W (dict) – The real-space operator dictionary.
- Returns:
Array of O_mn matrices.
- Return type:
numpy array
- spinWannier.wannier_utils.reciprocal_lattice_vectors(real_space_lattice_vector_matrix)#
Return the reciprocal lattice vectors from the real-space lattice vectors.
- Parameters:
real_space_lattice_vector_matrix (numpy array) – The real-space lattice vectors matrix.
- Returns:
The reciprocal lattice vectors.
- Return type:
list of numpy arrays
- spinWannier.wannier_utils.replace_middle_of_cmap_with_custom_color(color_middle=(0.85, 0.85, 0.85, 1.0), middle_range=0.1)#
Replace the middle of the colormap with a custom color.
- Parameters:
color_middle (tuple) – The middle color.
middle_range (float) – The middle range.
- Returns:
The new colormap.
- Return type:
cmap
- spinWannier.wannier_utils.save_bands_and_spin_texture(kpoints_rec, kpoints_cart, kpath, Eigs_k, S_mn_k_H_x, S_mn_k_H_y, S_mn_k_H_z, kmesh_2D=False, fout='bands_spin.pickle', save_folder='./tb_model_wann90/')#
Save the bands and spin texture information for given kpoints.
- Parameters:
kpoints_rec (list of tuples) – The k-points in the reciprocal space.
kpoints_cart (list of tuples) – The k-points in the Cartesian space.
kpath (list of floats) – The k-point distance position.
Eigs_k (dict) – The eigenvalues dictionary.
S_mn_k_H_x (dict) – The Sx operator dictionary in Hamiltonian gauge.
S_mn_k_H_y (dict) – The Sy operator dictionary in Hamiltonian gauge.
S_mn_k_H_z (dict) – The Sz operator dictionary in Hamiltonian gauge.
kmesh_2D (bool) – If True, the k-points are in a 2D mesh.
fout (str) – The name of the output file.
save_folder (str) – The folder where to save the file
- spinWannier.wannier_utils.save_bands_and_spin_texture_old(kpoints_rec, kpoints_cart, kpath, Eigs_k, S_mn_k_H_x, S_mn_k_H_y, S_mn_k_H_z, kmesh_2D=False, fout='bands_spin.pickle', save_folder='./tb_model_wann90/')#
Save the bands and spin texture information for given kpoints.
- Parameters:
kpoints_rec (list of tuples) – The k-points in the reciprocal space.
kpoints_cart (list of tuples) – The k-points in the Cartesian space.
kpath (list of floats) – The k-point distance position.
Eigs_k (dict) – The eigenvalues dictionary.
S_mn_k_H_x (dict) – The Sx operator dictionary in Hamiltonian gauge.
S_mn_k_H_y (dict) – The Sy operator dictionary in Hamiltonian gauge.
S_mn_k_H_z (dict) – The Sz operator dictionary in Hamiltonian gauge.
kmesh_2D (bool) – If True, the k-points are in a 2D mesh.
fout (str) – The name of the output file.
save_folder (str) – The folder where to save the file
- spinWannier.wannier_utils.selected_band_plot(band=0)#
- spinWannier.wannier_utils.split_spn_dict(spn_dict, spin_names=['x', 'y', 'z'])#
- Split dictionary with keys as ((Rx, Ry, Rz), spin_component_name) to three dictionaries
with keys as (Rx, Ry, Rz).
- Parameters:
spn_dict (dict) – The dictionary with R vectors as keys and hopping as a complex number matrix as the values.
spin_names (list of 3 str) – The spin names.
- Returns:
The dictionary with R vectors as keys and hopping as a complex number matrix as the values. dict: The dictionary with R vectors as keys and hopping as a complex number matrix as the values. dict: The dictionary with R vectors as keys and hopping as a complex number matrix as the values.
- Return type:
dict
- spinWannier.wannier_utils.spn_to_dict(model_dir='./', fwin='wannier90.win', fin='wannier90.spn', formatted=False, fout='spn_dict.pickle', save_as_text=False)#
Convert wannier90.spn file to a dictionary object and save as pickle. If ‘text_file’ == True, then save as a human-readable text file.
- Parameters:
model_dir (str) – The directory of the model.
fwin (str) – The win file.
fin (str) – The spn file.
formatted (bool) – If True, the spn file is formatted.
fout (str) – The output file.
save_as_text (bool) – If True, save as a text file.
- Returns:
The spin-projection matrices.
- Return type:
np.array
- spinWannier.wannier_utils.u_to_dict(fin='wannier90_u.mat', fout='u_dict.pickle', text_file=False, write_sparse=False)#
Convert _u.mat from wannier90 to pickled python dictionary.
- Parameters:
fin (str) – The input file.
fout (str) – The output file.
text_file (bool) – If True, save as a text file.
write_sparse (bool) – If True, write the sparse matrix.
- Returns:
The number of bands and the number of wannier functions.
- Return type:
tuple
- spinWannier.wannier_utils.uniform_real_space_grid(R_mesh_ijk=(5, 5, 1))#
- From the maxima in each direction, make a regular mesh
(-R_max_i, +R_max_i) x (-R_max_j, +R_max_j) x (-R_max_k, +R_max_k).
- Parameters:
R_mesh_ijk (tuple of 3 ints) – The number of mesh points in each direction.
- Returns:
The real-space grid.
- Return type:
list of tuples
- spinWannier.wannier_utils.unite_spn_dict(spn_x_dict, spn_y_dict, spn_z_dict, spin_names=['x', 'y', 'z'])#
- Unite dictionary with keys as (Rx, Ry, Rz) to three dictionaries
with keys as ((Rx, Ry, Rz), spin_component_name).
- Parameters:
spn_x_dict (dict) – The dictionary with R vectors as keys and hopping as a complex number matrix as the values.
spn_y_dict (dict) – The dictionary with R vectors as keys and hopping as a complex number matrix as the values.
spn_z_dict (dict) – The dictionary with R vectors as keys and hopping as a complex number matrix as the values.
spin_names (list of 3 str) – The spin names.
- Returns:
The dictionary with R vectors as keys and hopping as a complex number matrix as the values.
- Return type:
dict
- spinWannier.wannier_utils.wannier_energy_windows(wann_bands_lims, eigenval_file='wannier90.eig')#
Get energy window limits for wannierization.
- Parameters:
wann_bands_lims (tuple of two int) – Zero-indexed indices of the minimum and maximum bands included in wannierization. E.g. (0, 21) if the bottom 22 bands should be wannierized.
E_F (float, optional) – The Fermi level in eV. Defaults to 0.0.
eigenval_file (str, optional) – Wannier90 .eig file name. Defaults to “wannier90.eig”.
- Returns:
the wannierization and frozen energy minima and maxima.
- Return type:
4 floats
Utility functions for the wannier_quality method of the spinWannier.WannierTBmodel.WannierTBmodel class.
- spinWannier.wannier_quality_utils.compare_eigs_bandstructure_at_exact_kpts(dft_bands, wann_bands, num_kpoints, num_wann, f_name_out='WannierBerri_quality_error_Fermi_corrected.dat')#
Compare the DFT and Wannierized band structures at the exact k-points.
- Parameters:
dft_bands (np.array) – DFT bands.
wann_bands (np.array) – Wannierized bands.
num_kpoints (int) – Number of k-points.
num_wann (int) – Number of Wannierized bands.
f_name_out (str, optional) – Output file name. Defaults to ‘WannierBerri_quality_error_Fermi_corrected.dat’.
- Returns:
Array with the DFT energy in the first column and the Wannierization error in the second column.
- Return type:
np.array
- spinWannier.wannier_quality_utils.duplicate_kpoints_for_home_made(data, NK)#
Duplicate also the last k-point (in dictionary the keys are unique, so actually the data in the dictionaries where keys are k-points contain only one of each k-point, so if k-path starts and ends with the same k-point, only the first one is recorded.
- Parameters:
data (np.array) – Data to duplicate.
NK (int) – Number of k-points.
- Returns:
Duplicated data.
- Return type:
np.array
- spinWannier.wannier_quality_utils.get_NKpoints(OUTCAR='OUTCAR')#
Return number of kpoints from bands calculation from OUTCAR.
- Parameters:
OUTCAR (str, optional) – OUTCAR path. Defaults to ‘OUTCAR’.
- Returns:
number of k-points stated in the OUTCAR
- Return type:
int
- spinWannier.wannier_quality_utils.get_band_at_kpoint_from_EIGENVAL(EIGENVAL_path='./EIGENVAL', target_band=1, target_kpoint_string='0.0000000E+00 0.0000000E+00 0.0000000E+00')#
Get the energy of the target_band at the target_kpoint from the EIGENVAL file.
- Parameters:
EIGENVAL_path (str, optional) – Path to the EIGENVAL file. Defaults to ‘./EIGENVAL’.
target_band (int, optional) – Target band. Defaults to 1.
target_kpoint_string (str, optional) – Target k-point string. Defaults to ‘0.0000000E+00 0.0000000E+00 0.0000000E+00’.
- Returns:
Energy of the target_band at the target_kpoint.
- Return type:
float
- spinWannier.wannier_quality_utils.get_fermi(path='.')#
Extract Fermi energy from the DOSCAR file.
- Parameters:
path (str, optional) – Path to the directory with the DOSCAR file. Defaults to “.”.
- Returns:
Fermi energy in eV.
- Return type:
float
- spinWannier.wannier_quality_utils.get_fermi_corrected_by_matching_bands(nsc_calculation_path='../0_nsc_for_wann_25x25_frozmaxmargin_0.2eV', corrected_at_kpoint='0.0000000E+00 0.0000000E+00 0.0000000E+00', corrected_at_band=11, sc_calculation_path='../sc', fout_name='FERMI_ENERGY_corrected.in')#
Get the Fermi energy from the self-consistent calculation and correct it so that the band at the target_kpoint and target_band is at the same energy in the non-self-consistent calculation.
- Parameters:
path (str, optional) – Path to the directory with the DOSCAR file. Defaults to “.”.
nsc_calculation_path (str, optional) – Path to the non-self-consistent calculation directory. Defaults to ‘../0_nsc_for_wann_25x25_frozmaxmargin_0.2eV’.
corrected_at_kpoint (str, optional) – Target k-point string. Defaults to ‘0.0000000E+00 0.0000000E+00 0.0000000E+00’.
corrected_at_band (int, optional) – Target band. Defaults to 11.
sc_calculation_path (str, optional) – Path to the self-consistent calculation directory. Defaults to “../sc”.
fout_name (str, optional) – Output file name. Defaults to “FERMI_ENERGY_corrected.in”.
- Returns:
Corrected Fermi energy in eV.
- Return type:
float
- spinWannier.wannier_quality_utils.get_frozen_window_min_max(wannier90winfile='wannier90.win')#
Get the frozen window min and max from the wannier90.win file.
- Parameters:
wannier90winfile (str, optional) – Path to the wannier90.win file. Defaults to ‘wannier90.win’.
- Returns:
Frozen window min. float: Frozen window max.
- Return type:
float
- spinWannier.wannier_quality_utils.integrate_error(error_by_energy, E_min=-1000.0, E_max=1000.0)#
Integrate the error in ‘f_name_in’ in the energy range [E_min, E_max] included.
- Parameters:
error_by_energy (np.array) – Error vs. energy.
E_min (float, optional) – Minimum energy. Defaults to -1e3.
E_max (float, optional) – Maximum energy. Defaults to 1e3.
- Returns:
Array with the integrated error.
- Return type:
np.array
- spinWannier.wannier_quality_utils.parse_eigenval_file(fin, spin=0)#
Parse the EIGENVAL file and return the kpoints, bands, number of kpoints, and number of bands.
- Parameters:
fin (str) – Path to the EIGENVAL file.
spin (int, optional) – Spin index. Defaults to 0.
- Returns:
kpoints. np.array: bands. int: number of kpoints. int: number of bands.
- Return type:
np.array
- spinWannier.wannier_quality_utils.plot_err_vs_bands(kpoints, kpath, kpath_ticks, Eigs_k, E_diff, S_diff, NW, fout='ERRORS_ALL_band_structure.jpg', yaxis_lim=None, savefig=True, showfig=True)#
Output a figure with RMSE_E, RMSE_Sx, RMSE_Sy, and RMSE_Sz-projected band structure.
- Parameters:
kpoints (np.array) – kpoints.
kpath (np.array) – kpath.
kpath_ticks (list) – kpath ticks.
Eigs_k (dict) – Eigs_k.
E_diff (np.array) – E_diff.
S_diff (np.array) – S_diff.
fout (str, optional) – Output file name. Defaults to ‘ERRORS_ALL_band_structure.jpg’.
yaxis_lim (list, optional) – y-axis limits. Defaults to None.
savefig (bool, optional) – Save the figure. Defaults to True.
showfig (bool, optional) – Show the figure. Defaults to True.
- spinWannier.wannier_quality_utils.plot_err_vs_energy(error_by_energy, Ef, title='Wannierization RMS error vs. energy', fig_name_out='wannier_quality_error_by_energy.png', savefig=True, showfig=True)#
Plot the error vs. energy.
- Parameters:
error_by_energy (np.array) – Error vs. energy.
Ef (float) – Fermi energy.
title (str, optional) – Title of the plot. Defaults to “Wannierization RMS error vs. energy”.
fig_name_out (str, optional) – Output file name. Defaults to “wannier_quality_error_by_energy.png”.
savefig (bool, optional) – Save the figure. Defaults to True.
showfig (bool, optional) – Show the figure. Defaults to True.
- spinWannier.wannier_quality_utils.vasp_calc_collinear(EIGENVAL_path='./EIGENVAL')#
Get the N_eig from the EIGENVAL file.
- Parameters:
EIGENVAL_path (str, optional) – Path to the EIGENVAL file. Defaults to ‘./EIGENVAL’.
- Returns:
N_eig (1 - non-spin-polarized, 2 - spin-polarized).
- Return type:
int