Tight Binding

WHERE CAN ONE CALCULATE AN EIGENVECTOR FOR USE IN THE PLOT AND EXPORT FUNCTIONS?

Features

The wannier90 module has the following features:

  1. Read output files from the VASP and wannier90 program
  2. Read Slater-Koster nearest-neighbour parameter lists (“standard” tight-binding, like 1st-nearest-neighbour approximation)
  3. Change or drop input parameters
  4. Create unit cells or supercells from input parameters
  5. Create finite structures (ribbons, dots) from input parameters
  6. Merge parameters from several input files (e.g. bulk & defect parameters)
  7. Calculate eigenvalue problems and bandstructures
  8. Apply magnetic field
  9. Plot eigenvectors (pseudo-real space representation - proper real space representation using Wannier orbitals to come).
  10. Export or use any data easily

Load Hamiltonian parameters

The Hamiltonian class is the central class of the module. At the moment, it contains too many features and should be split.

Anyway, after loading the parameters from a file, you have a Hamiltonian object which provides you with a lot of functions you can apply on the data.

Load the module:

from wannier90.w90hamiltonian import Hamiltonian

You can load data from a Slater-Koster nearest-neighbour parameter file:

ham=Hamiltonian.from_nth_nn_list("/path/to/nearestneighbourfile.dat")

Or from a wannier90 & VASP output:

ham=Hamiltonian.from_file(wannier90hr_graphene,poscarfile,wannier90woutfile)

You can find example data to play with in exampledata in the envTB directory.

Plot a bandstructure

Before you can calculate a bandstructure, you need a path in the Brillouin zone. The function standard_paths() gives you the standard paths for some common crystal structures:

path=ham.standard_paths('hexagonal',100)[2]

Alternatively, you can use point_path() to create a path between vertices you provide:

path = ham.point_path([[0,0,0],[0.5,0,0]],100)

Then, you can plot the bandstructure to a file:

ham.plot_bandstructure(path,'/tmp/myplot.png','d')

Alternatively, you can do this:

ham.plot_bandstructure('hexagonal','/tmp/myplot.png')

Or even shorter, if you only want to display the plot (you have to admit that this is really simple):

ham.plot_bandstructure('hexagonal')

Or store the data in a variable and save it to a file:

data=ham.bandstructure_data(path, 'd')
numpy.savetxt('/tmp/bs.dat', numpy.real(data), fmt="%12.6G")

You can calculate the Bloch eigenvalues of a specific point using bloch_eigenvalues().

Modify the Hamiltonian

After all, tight-binding is about using the parameters of the infinite crystal lattice for something different. The functions create_supercell_hamiltonian() and create_modified_hamiltonian() (only a wrapper for the first function, actually) give you that feature.

Drop orbitals:

If you can drop orbitals with a good conscience (e.g. the sigma system of graphene), use:

new_ham=ham.create_modified_hamiltonian(usedorbitals=(0,1))

A new Hamiltonian called new_ham is created, where only the first two orbitals are used. new_ham has all the functionalities of the original Hamiltonian ham (plot bandstructure, calculate eigenvectors, modify once more).

Drop hopping parameters:

If you create parameters with wannier90, you probably don’t want to use all hopping parameters. If you only want to keep hopping parameters to a chosen set of neighbour unit cells, use the parameter usedhoppingcells:

cells_to_keep=ham.unitcells_within_zone(2,'d',numpy.inf)
new_ham=ham.create_modified_hamiltonian(usedhoppingcells=cells_to_keep)

In this example, cells_to_keep contains all cells up to the second row of adjacent cells.

Create a supercell

A supercell is a cell that contains of several copies of the original cell. It is defined by

  • the coordinates of the unit cell copies
  • the new lattice vectors.

Both are given in the basis of the current lattice vectors.

For example, the following command creates a rectangular unit cell consisting of two hexagonal unit cells:

new_ham=ham.create_supercell_hamiltonian(cellcoordinates=[[0,0,0],[1,0,0]],latticevecs=[[1,-1,0],[1,1,0],[0,0,1]])

Plot cell diagram

This command creates a hexagonal unit cell consisting of four smaller hexagonal cells. The lattice vectors are, obviously, twice as long as the current ones:

new_ham=ham.create_supercell_hamiltonian(cellcoordinates=[[0,0,0],[1,0,0],[0,1,0],[1,1,0]],latticevecs=[[2,0,0],[0,2,0],[0,0,1]])

Add cell diagram and bandstructure plot

Create a ribbon

In the following example of a zigzag Graphene nanoribbon, one has to accomplish the following steps:

  1. Create a rectangular unit cell out of the hexagonal unit cell (see above):

    ham2=ham.create_supercell_hamiltonian([[0,0,0],[1,0,0]],[[1,-1,0],[1,1,0],[0,0,1]])
    
  2. Create a ribbon unit cell which has the width of the ribbon:

    ham3=ham2.create_supercell_hamiltonian([[0,i,0] for i in range(unitcells)],[[1,0,0],[0,unitcells,0],[0,0,1]])
    
  3. Remove all hoppings to neighbouring cells in y direction and drop the first and last orbital of the cell to make it a zigzag ribbon:

    ham4=ham3.create_modified_hamiltonian(ham3.drop_dimension_from_cell_list(1),usedorbitals=range(1,ham3.nrorbitals()-1))
    

Plot BS, unit cell geometry and neighbour cell geometry

Shift energy

Since it’s convenient to have the Fermi energy at 0 eV, but DFT software doesn’t automatically do that, you can shift the energy range:

new_ham=ham.create_modified_hamiltonian(energyshift=3.34)

Add magnetic field

new_ham=ham.create_modified_hamiltonian(magnetic_B=1,gauge_B='landau_x')

Depending on the symmetry of the system, you have to choose the gauge (see the function documentation).

Export modified Hamiltonian

blar

Plot Wannier orbitals

Es gibt doch irgendwo schon die lineare Interpolation und das Einlesen, oder?

Plot electron density

Add me, I’m pretty much the coolest feature, right?

Mixins

blarblar tbc blar

More examples

See Example section.

More features

You can find more documentation about the methods when you click on their names.

  • latticevectors(): lattice vectors of the system.
  • reciprocal_latticevectors(): reciprocal lattice vectors of the system.
  • orbitalspreads(): orbital spreads (=sizes) of the basis orbitals.
  • orbitalpositions(): orbital positions of the basis orbitals.
  • maincell_eigenvalues(): calculate eigenvalues of the main unit cell as if it were not in a periodic crystal. Use this to calculate the spectrum of finite structures.
  • bloch_eigenvalues(): calculate the eigenvalues for a single vector k.
  • create_orbital_vector_list(): Concatenate the amplitudes of a solution vector with information about their basis elements (coordinates, spread). Use this if you want to export more information about the eigenvectors than just the vector itself.
  • plot_vector(): plot an eigenvector by putting circles with a size proportional to the eigenvector amplitudes on the orbital positions. The function doesn’t use the real space probability density of the orbitals.
  • unitcells_within_zone(): Returns a list of unit cells within a certain area. Use it e.g. if you want to drop hopping parameters.
  • drop_dimension_from_cell_list(): Takes a list of unit cell coordinates and drops the x,y or/and z dimensionf rom the list - this way you can create a 2D material from a 3D material or a 1D material from a 2D material.
  • standard_paths(): Create standard paths within the Brillouin zones of the possible crystal lattices.
  • unitcellcoordinates(): Cartesian coordinates of the given unit cells.
  • drawunitcells(): Plot the main cell and the cells where there exist hopping matrix elements to.
  • point_path(): Create a path between given points. Use it to create a k-point path for the bandstructure.
  • apply_electrostatic_potential(): Apply an electrostatic potential.
  • shift_fermi_energy_to_zero(): Applies an energy shift so that the new Fermi energy is at 0.

Code reference

Alle Funktionen durchtesten!!!! Sparse matrix Umstellung hat vl noch Spuren hinterlassen

class envtb.wannier90.w90hamiltonian.BandstructurePlot[source]

Combine several bandstructure plots.

Call plot(kpoints,data) for every bandstructure plot. Then, call save(filename) to save to a file.

plot(kpoints, data)[source]

Add a bandstructure plot to the figure.

kpoints: list of kpoints data: list of eigenvalues for each kpoint.

class envtb.wannier90.w90hamiltonian.Hamiltonian[source]

TODO: durchschleifen der argumente bei bloch_eigenvalues etc. ist bloed. vl. argumente bei allen anderen mit *args und auf dokumentation von bloch_eigenvalues verweisen? TODO: unitcellcoordinates and numbers is used ambiguously TODO: performance: use in-place operations for numpy -=, +=, *= and consider numpy.fromfunction TODO: to ensure a variable is a numpy array: a = array(a, copy=False) TODO: REFACTOR REFACTOR REFACTOR TODO: sparse matrices: http://docs.scipy.org/doc/scipy/reference/sparse.html,

TODO: plots konsistenter zur aussenwelt machen, dh steuerung von aussen erlauben (mehrere
plots uebereinander, nebeneinander etc.
TODO: alles sparse + abgeleitete hamiltonians auf die alten verweisen,
nicht explizit schreiben
TODO: moeglichkeit zu einem output-logfile mit versionsnummer, zB
mit globaler variable LOG
TODO; die reihenfolge in create_supercell ist komisch, zB beruecksichtigt energy_shift schon die vergroesserung
der zelle

TODO: apply_electrostatic_potential fuer alle fkten mit 1,2,3 argumenten

There are several ways to initialize the Wannier90 Hamiltonian: 1) Hamiltonian.from_file(wannier90filename,poscarfilename,wannier90woutfilename) 2) Hamiltonian.from_raw_data(unitcellmatrixblocks,unitcellnumbers,latticevecs,orbitalspreads,orbitalpositions) 3) Hamiltonian.from_nth_nn_list(nnfile,customhopping):

See the documentation of those methods.

apply_electrostatic_potential(potential)[source]

Apply an electrostatic potential to the system and return the new Hamiltonian.

potential: a LinearInterpolationNOGrid object. If the object is a 1D/2D interpolation, the y and z/z coordinate are not used.

Don’t forget to match the unit of length and potential!

Best function yet in OOP.

Very Pythonic: potential can be anything with ().

bandstructure_data(kpoints, basis='c', usedhoppingcells='all')[source]

Calculates the bandstructure for a given kpoint list. For direct plotting, use plot_bandstructure(kpoints,filename).

If kpoints is a string, this string will be interpreted as the name of the crystal structure (see standard_paths), and the crystal structure’s default kpoint path will be used.

usedhoppingcells: If you don’t want to use all hopping parameters, you can set them here (get the list of available cells with unitcellnumbers() and strip the list from unwanted cells). basis: ‘c’ or ‘d’. Determines if the kpoints are given in cartesian reciprocal coordinates or direct reciprocal coordinates.

Return: A list of eigenvalues for each kpoint is returned. To sort by band, use data.transpose().

If MPI is used, ONLY THE ROOT PROCESS returns the data, the others return None.

bloch_eigenvalues(k, basis='c', usedhoppingcells='all', return_evecs=False)[source]

Calculates the eigenvalues of the eigenvalue problem with Bloch boundary conditions for a given vector k.

The function uses a dense matrix eigenvalue solver because it returns all eigenvalues, so don’t let the matrices get too big.

usedhoppingcells: If you don’t want to use all hopping parameters, you can set them here (get the list of available cells with unitcellnumbers() and strip the list from unwanted cells). basis: ‘c’ or ‘d’. Determines if the kpoints are given in cartesian reciprocal coordinates or direct reciprocal coordinates. return_evecs: If True, evecs are also returned as the second return value.

create_modified_hamiltonian(usedhoppingcells='all', usedorbitals='all', energyshift=None, magnetic_B=None, gauge_B='landau_x', mixin_ham=None, mixin_hoppings=None, mixin_cells=None, mixin_assoc=None, onsite_potential=None)[source]

Creates a Hamiltonian with dropped orbitals or hopping cells. This is just a wrapper for create_supercell_hamiltonian().

create_orbital_vector_list(vector, include_third_dimension=False, include_spread=False)[source]

Create a list of orbital positions with given eigenvector amplitudes. Only the real part from the eigenvector is kept.

vector: Vector to connect to the orbital positions. include_third_dimension: Include the z position of the points. Default is False. include_spread: Include the spread of the orbital. Default is False.

Return: A matrix containing the following columns:
x y (z) (spread) value
create_supercell_hamiltonian(cellcoordinates, latticevecs, usedhoppingcells='all', usedorbitals='all', energyshift=None, magnetic_B=None, gauge_B='landau_x', mixin_ham=None, mixin_hoppings=None, mixin_cells=None, mixin_assoc=None, onsite_potential=None, output_maincell_only=False)[source]

Creates the matrix elements for a supercell containing several unit cells, e.g. a volume with one unit cell missing in the middle or one slice of a nanoribbon.

cellcoordinates: unit cells contained in the supercell [[0,0,0],[1,0,0],...]. Use integergrid3d() to easily create a unit cell grid (which you can modify). latticevectors: lattice vectors of the new unit cell in the basis of the old unit cells (!!). E.g. a supercell of four graphene unit cells could have latticevecs=[[2,0,0],[0,2,0],[0,0,1]]. If you want to create a ribbon or a molecule, use a high value in one of the coordinates (e.g. a long y lattice vector). usedhoppingcells: If you don’t want to use all hopping parameters, you can set the cells to “hop” to here (list of cell coordinates). usedorbitals: a list of orbitals to use. Default is ‘all’. Note: this only makes sense if the selected orbitals don’t interact with other orbitals. energyshift: Shift energy scale by energyshift, e.g. to shift the Fermi energy to 0. Also shifts the Fermi energy variable of the Hamiltonian. magnetic_B: Magnetic field in perpendicular direction (in T) gauge_B: ‘landau_x’: Landau gauge for systems with x periodicity (A=(-By,0,0))

‘landau_y:’: Landau gauge for systems with y periodicity (A=(0,Bx,0)) ‘symmetric’: Symmetric gauge for systems with x and y periodicity (A=1/2(-By,Bx,0))
mixin_ham: A “mix-in” Hamiltonian from which some matrix elements are used. Default is None.
The mixin is done before the other modifications (magnetic field, energyshift)
mixin_hoppings: A list of matrix elements from the main Hamiltonian that should be substituted with matrix elements from mixin_ham.
The conjugate hopping is generated automatically (i.e. (0,1) will be automatically expanded to (0,1),(1,0) ). Example: Substitute all matrix elements from orbitals 0,1 to orbitals 0,1,2,3: mixin_hoppings=[(0,0),(0,1),(0,2),(0,3),(1,0),(1,1),(1,2),(1,3)]
mixin_cells: List of unit cells. Substitution can be restricted to specific unit cells.
E.g. the main cell: mixin_cells=[[0,0,0]] E.g. the main cell and the ones left and right: mixin_cells=[[-1,0,0],[0,0,0],[1,0,0]] Default is None, which means that all cells that exist in both Hamiltonians will be substituted. Mind that the cell coordinates have to be lists, not tuples!
mixin_assoc: Association list for orbitals in the current Hamiltonian and the mixin Hamiltonian. Type is dictionary.

If the mixin Hamiltonian describes a different system, the orbital numbers may not be the same. Only the relevant orbitals (those mentioned in mixin_hoppings) have to be here. Default is None, which means that the orbital numbers are assumed to be identical. Example: Main Hamiltonian Mixin Hamiltonian

0 0 1 1 2 2 3 3 4 4 100 15 101 16 102 17 103 18 104 19

mixin_assoc={0:0,1:1,2:2,3:3,4:4,100:15,101:16,102:17,103:18,104:19}

onsite_potential: List of numbers. The values will be added to the diagonal (=onsite) matrix elements of the main cell. This approximates an electrostatic potential in the system. Default is None.

output_maincell_only: If True, only the main cell matrix block will be calculated. This makes sense for a big system of which you want to calculate the eigenvalues, not the bandstructure. Default is false

Return: New Hamiltonian with the new properties.

drawunitcells(unitcellnumbers='all')[source]

Create a plot of a list of unit cells.

unitcellnumbers: Numbers of unit cells to plot. Default value is ‘all’, then unitcellnumbers() is used.

drop_dimension_from_cell_list(dimension, unitcellnumbers='all')[source]

Takes a list of unit cell coordinates and drops the x,y or/and z dimension from the list - this way you can create a 2D material from a 3D material or a 1D material from a 2D material. The function deletes all unit cell numbers that have a nonzero entry in that dimension.

dimension: the dimension to drop: 0, 1 or 2. Can also be a list, e.g. (0,1) drops the first and second dimension. unitcellnumbers: List of unit cell numbers. Default value is ‘all’.

fermi_energy()[source]

Return the system’s Fermi energy.

classmethod from_file(wannier90filename, poscarfilename, wannier90woutfilename, outcarfilename)[source]

A constructor to create an object based on data from files. wannier90filename: Path to the wannier90_hr.dat file poscarfilename: Path to the VASP POSCAR file wannier90woutfilename: Path to the wannier90.wout file outcarfilename: Path to the VASP OUTCAR file

classmethod from_nth_nn_list(nnfile, customhopping=None)[source]

A constructor to create a nth-nearest-neighbour Hamiltonian.

nnfile: File containing the system information (see example data) customhopping: Dictionary, containing hopping parameters overriding those in nnfile.

Example: {0:ONSITE,1:1STNN,2:2NDNN}
classmethod from_raw_data(unitcellmatrixblocks, unitcellnumbers, latticevecs, orbitalspreads, orbitalpositions, fermi_energy)[source]

A constructor used to create a custom Hamiltonian. unitcellmatrixblocks: Hopping elements, arranged by unit cells. unitcellnumbers: coordinates of the unit cells “hopped” to. latticevecs: lattice vectors. orbitalspreads: spreads of the orbitals orbitalpositions: positions of the orbitals

hermitian_hoppinglist(unitcellnumbers)[source]

The function removes unit cells from a list of unit cells whose “parity partners” are missing to ensure a Hermitian Bloch matrix.

Return: kept, removed

kept: Kept unit cell numbers removed: removed unit cell numbers (just for control purposes)

If hopping to a specific unit cell is not used, one has to make sure that the parity inversed unit cell (=the cell with the “negative” coordinates”) is also dropped. That’s because the matrix elements of the bloch matrix look like this:

... + gamma_i e^ikR + gamma_i e^-ikR + ...

The sum of the two terms is cos(ikR) and real.

–> The function drops the terms which miss their partner and thus won’t become real.

Note: It makes sense to remove not only the “parity partner”, but all unit cells which are identical due to symmetry.

integergrid3d(i, j, k)[source]

Creates an integer grid with the dimensions i,j,k

latticevectors()[source]

Return the system’s lattice vectors.

maincell_eigenvalues(solver='dense', return_evecs=False, **kwargs)[source]

Calculates the eigenvalues of the main cell (no hopping to adjacent unit cells).

solver: eigenvalue solver. There are:
‘dense’: Assuming a dense matrix; returns all eigenvalues. Uses scipy.linalg.eig. E.g. >>> evals=ham.maincell_eigenvalues() ‘scipy_arpack’: find a given number of eigenvalues and eigenvectors of a BIG, SPARSE matrix (including shift-invert). It can never give you all eigenvalues. Uses ARPACK through scipy.sparse.linalg.eigsh. You have to supply additional parameters for eigsh using **kwargs, e.g. >>> ham.maincell_eigenvalues(‘arpack’,k=10,sigma=0.0,ncv=100) See http://docs.scipy.org/doc/scipy/reference/generated/scipy.sparse.linalg.eigsh.html for the available parameters. You will probably need k,sigma, and maybe nvc, which. Consider using which=’SM’ if E_F=0.

return_evecs: Also return eigenvectors.

maincell_hamiltonian_matrix()[source]

Returns the Hamiltonian matrix for the main cell, without hopping parameters to other cells. This is the matrix whose eigenvalues you can calculate using maincell_eigenvalues().

nrorbitals()[source]

Returns the number of orbitals/bands.

orbitalpositions()[source]

Return the wannier90 orbital positions.

orbitalspreads()[source]

Return the wannier90 orbital spreads.

plot_bandstructure(kpoints, filename=None, basis='c', usedhoppingcells='all', mark_reclattice_points=False, mark_fermi_energy=False)[source]

Calculate the bandstructure at the points kpoints (given in cartesian reciprocal coordinates - use direct_to_cartesian_reciprocal(k) if you want to use direct coordinates) and save the plot to filename. The ending of filename determines the file format. If filename=None (default), the plot will not be saved (you can display it using pyplot.show() ).

If kpoints is a string, this string will be interpreted as the name of the crystal structure (see standard_paths), and the crystal structure’s default kpoint path will be used.

usedhoppingcells: If you don’t want to use all hopping parameters, you can set them here (get the list of available cells with unitcellnumbers() and strip the list from unwanted cells). basis: ‘c’ or ‘d’. Determines if the kpoints are given in cartesian reciprocal coordinates or direct reciprocal coordinates. mark_reclattice_points: You can mark important reciprocal lattice points, like Gamma or K. This variable can be (i) True if you use a string for kpoints (ii) a list which contains the names of the points and the points: mark_reclattice_points=[names,points]. The points have to be in cartesian coordinates. Default is False. mark_fermi_energy: If you supply the Fermi energy here, a line will be drawn. If True, the Fermi energy will be taken from fermi_energy(). Default is False.

If MPI is used, ONLY THE ROOT PROCESS plots. This coincides with bandstructure_data, where also only the root process returns all the bandstructure data.

Return: lines: List of matplotlib.lines.Line2D objects that were drawn. You can change the style, color etc., like:

for line in lines:
line.set_color(‘red’)

fermi_energy_line: The fermi energy mark Line2D object. lattice_point_lines: The lattice point marks Line2D object.

plot_orbital_positions()[source]

Plot the positions of the orbitals in the unit cell.

plot_vector(vector, scale=1, figsize=None)[source]

Plot a vector with geometry by putting circles on the positions of the orbitals. The size of the circles corresponds to the absolute square, the color to the sign.

vector: vector to plot scale: scale factor for the circles. figsize: w,h tuple in inches

point_path(corner_points, nrpointspersegment)[source]

Generates a path connecting the corner_points with nrpointspersegment points per segment (excluding the next point), resulting in sum(nrpointspersegment)+1 points. The points in corner_points can have any dimension. nrpointspersegment is a list with one element less than corner_points. If nrpointspersegment is an integer, it is assumed to apply to each segment.

Example: my_hamiltonian.point_path([[0,0],[1,1],[2,2]],[2,2]) gives [[0.0, 0.0], [0.5, 0.5], [1.0, 1.0], [1.5, 1.5], [2, 2]] (note: those are 5 points, which is sum([2,2])+1)

or equivalently: my_hamiltonian.point_path([[0,0],[1,1],[2,2]],2)

reciprocal_latticevectors()[source]

Return the system’s reciprocal lattice vectors.

shift_fermi_energy_to_zero()[source]

Applies an energy shift so that the new Fermi energy is at 0.

The new Hamiltonian is returned.

standard_paths(name, nrpointspersegment=100)[source]

Gives the standard path for a Bravais lattice in direct reciprocal coordinates.

At the moment, there are ‘hexagonal’, ‘fcc’, ‘1D’ and ‘1D-symmetric’.

name: Name of the lattice nrpointspersegment: optional; if > 1, a list of intermediate points connecting the main points is also returned and can be used for a bandstructure path (nrpointspersegment points per segment). Default value: 100

Return: points,names(,path)

points: points in the path names: names of the points (path: path with intermediate points. Only returned if nrpointspersegment is > 1)

unitcellcoordinates(unitcellnumbers='all')[source]

Cartesian coordinates of the given unit cells.

unitcellnumbers: a list of the unit cell numbers. Default value is ‘all’, then unitcellnumbers() is used.

unitcellnumbers()[source]

Returns the numbers of the unit cells supplied in the wannier90_hr.dat file.

unitcells_within_zone(zone, basis='c', norm_order=2)[source]

Returns a list of unit cells within a certain area. The function is comparing the same point in each cell (e.g. always the bottom left end).

zone: can be a number or a tuple:
number: radius to include cells within. tuple: area to include cells within, in the sense of distance from the origin along a direction.

basis: determines if zone is given in cartesian (‘c’) or direct (‘d’) coordinates. IMPORTANT: If direct coordinates are used, use integers for zone, not float!

norm_order: if zone is a number (=radius), norm_order is the norm to use (mathematical definition, see http://docs.scipy.org/doc/numpy/reference/generated/numpy.linalg.norm.html). Default is 2 (=Euclidean norm) Short version: 2 gives you a “circle”, numpy.inf a “square”.

Examples: Cells within 30 Angstrom: unitcells_within_zone(30) Cells within a 6x8x1 Angstrom cuboid: unitcells_within_zone((3.0,4.0,0.5)) Cells within a 4x4x4 block in direct coordinates: unitcells_within_zone((2,2,2),’d’)

write_matrix_elements(outputfile, usedhoppingcells='all', usedorbitals='all')[source]

Write the wannier90 matrix elements to a file (*.wetb) readable by Florian’s code. Information contained in the file: - lattice vectors - orbital positions - orbital spreads - matrix elements of the chosen unit cells and orbitals

Description of the file format: The file contains three sections, divided by an arbitrary number of blank lines (one at least, obviously). Lines starting with # are comments and are ignored. Anything in a line after the data is also ignored (i.e. you can write anything in the same line after the data, with or without #). First block: lattice vectors in rows Second block: orbital spread (first column) and position (other columns) of every orbital Third block: Matrix elements.

Column 1-3: Unit cell number Column 4: Orbital number in main unit cell Column 5: Orbital number in other unit cell (the one the electron “hops” to) Column 6&7: Real & imaginary part of the matrix element

outputfile: Name of the output file (*.wetb - Wannier90-Environmental-dependent-Tight-Binding) usedhoppingcells: If you don’t want to use all hopping parameters, you can set them here (get the list of available cells with unitcellnumbers() and strip the list from unwanted cells). usedorbitals: a list of used orbitals to use. Default is ‘all’. Note: this only makes sense if the selected orbitals don’t interact with other orbitals.