Code examples

The easiest way to start out is to look at examples. You can find them in the simple module. Run them by including it:

from envtb import simple
simple.NameOfExample()

Example input files

The folder exampledata contains example data from third-party applications that you can use.

01_graphene_vasp_wannier90

A DFT-LDA calculation of Graphene, as well as the Wannier orbital matrix elements resulting from it. The data was calculated using VASP and wannier90.

VASP files:

  • EIGENVAL - kpoints and energies of the bandstructure
  • POSCAR - unit cell configuration (lattice vectors, atom positions)
  • PROCAR - electronic density of the ground state, projected on atomic orbitals
  • OUTCAR - output protocol, containing the Fermi energy

The easiest way to use the files is probably to start python in the example data directory. Then, you don’t have to supply an absolute path to the files:

cd /path/to/envtb/exampledata/01_graphene_vasp_wannier90
python
>>> from wannier90.w90hamiltonian import Hamiltonian
>>> ham=Hamiltonian.from_file("wannier90_hr.dat","POSCAR","wannier90.wout","OUTCAR")

If you run the IPython web interface, you can start it anyway you like and then execute the above cd command in IPython.

wannier90 files:

  • wannier90_hr.dat - matrix elements in the Wannier-orbital basis
  • wannier90.wout - output protocol, containing the orbital positions and spreads

02_graphene_3rdnn

Nth-nearest-neighbor parameterization files for the pi system of Graphene. It contains valid parameters for 1st-NN and 3rd-NN. If you have own parameters, you can use up to 10th-NN with the supplied neighbor lists.

Electrostatics

envtb.simple.SimpleElectrostaticProblem()[source]

Calculates the electrostatic potential of a 1D stripe of 1nm width and 30 nm height with a gridsize of 1 nm. At 11 nm, the potential is fixed to 8V.

Boundary conditions:

  • left, right: Periodic boundary conditions (using PeriodicContainer).
  • top: Neumann boundary condition, slope=1e18
  • bottom: Dirichlet boundary condition, potential=0

(default boundary condition)

envtb.simple.PotentialOfSimpleConductor2D()[source]

Rectangle with 5V conducting plate at the bottom.

Boundary conditions:

  • Left, right, top: default boundary condition 0V
  • bottom: Dirichlet BC, 5V
envtb.simple.PotentialOfGluedRectangles2D()[source]

Two rectangles with capacitor plates, glued together with Container. There are four capacitor plates. The boundary conditions are the default BC.

The relative position of the two rectangles to each other is given by the parameters position and offset of the connect() function (see documentation).

Quantum capacitance

envtb.simple.GrapheneQuantumCapacitance(equation='potential')[source]

Calculate the quantum capacitance of Graphene on SiO2, including temperature.

equation: ‘charge’ or ‘potential’, see documentation of QuantumCapacitanceSolver.

Boundary conditions:

  • left, right: periodic BC
  • top: Neumann BC, slope=0
  • bottom: backgate capacitor plate

Please read the code comments for further documentation.

Return:

voltages: The voltages where the capacitance is calculated (x values) capacitance: the quantum capacitance (y(x) values)

envtb.calculations.graphenequantumcapacitance.LoopQuantumCapacitanceWithSidegatesFixedSystem(sidegatevoltages, graphenesidegate_behavior='metal')[source]

The function calculates the quantum capacitance of a GNR, embedded in the following system:

  • height: 600nm
  • width: 400nm
  • gridsize: 1nm
  • GNR + sidegate vertical position: 300nm
  • GNR horizontal position: centered
  • GNR width: 80nm
  • distance sidegate-GNR: 30nm
  • backgate position: bottom
  • material between backgate + GNR/sidegates: SiO2, epsilon=3.9
  • left/right boundary condition: periodic
  • top BC: neumann, slope=0

The backgate voltage is varied from -60V to 60V in 0.5V.

sidegatevoltages: A list of side gate voltages which will be calculated
symmetrically and asymetrically.
graphenesidegate_behaviour: can have the values ‘metal’ and ‘graphene’,
the sidegates will behave accordingly.

Return:

voltages: the voltages where
the charge was calculated.
parameters: list of graphene
widths and sidegate voltages
capacitancequantumlist: list of
capacitances.

Tight-Binding

envtb.simple.PlotGraphenenthNNBandstructure(nnfile, outfile)[source]

Plot bandstructure of Graphene along the default hexagonal path.

nnfile: path to the nearest neighbour parameter file, e.g. /path/to/envtb/exampledata/02_graphene_3rdnn/graphene3rdnnlist.dat outfile: plot image output path, e.g. /tmp/bandstructure.png

envtb.simple.PlotGrapheneWannierBandstructure(wannier90hr_graphene, poscarfile, wannier90woutfile, outfile)[source]

Plot bandstructure of Graphene along the default hexagonal path from Wannier90 data. You can use the files in /path/to/envtb/exampledata/01_graphene_vasp_wannier90

wannier90hr_graphene: path to the wannier90_hr.dat file containing the graphene bulk calculation poscarfile: path to the VASP POSCAR file of the graphene bulk calculation wannier90woutfile: path to the wannier90.wout file output: path to the output image file.

envtb.simple.plot_zigzag_graphene_nanoribbon_pz_bandstructure(wannier90hr_graphene, poscarfile, wannier90woutfile, outcarfile, width, output=None, usedhoppingcells_rings='all')[source]

Plot the pi bandstructure of a zigzag graphene nanoribbon based on a wannier90 calculation of bulk graphene. The pz orbitals have to be the first two orbitals in the wannier90 file. A .dat file with the bandstructure data is also saved as output.dat (each line representing one k-point).

wannier90hr_graphene: path to the wannier90_hr.dat file containing the graphene bulk calculation poscarfile: path to the VASP POSCAR file of the graphene bulk calculation outcarfile: path to the VASP OUTCAR file width: width of the ribbon (number of rings). output: path to the output image file. If None, nothing will be saved. Default is None. usedhoppingcells_rings: If you don’t want to use all hopping parameters, you can set the number of ‘rings’ surrounding the main cell here. If it is a list (e.g. range(5)), several plots are created. The default value is ‘all’.

Return: Hamiltonian: the generated GNR Hamiltonian. data: bandstructure data. path: path through reciprocal space.

envtb.simple.plot_armchair_graphene_nanoribbon_pz_bandstructure(wannier90hr_graphene, poscarfile, wannier90woutfile, width, output, usedhoppingcells_rings='all')[source]

Plot the pi bandstructure of an armchair graphene nanoribbon based on a wannier90 calculation of bulk graphene. The pz orbitals have to be the first two orbitals in the wannier90 file. In this example function, width is the number of rings in transversal direction. A .dat file with the bandstructure data is also saved as output.dat (each line representing one k-point).

wannier90hr_graphene: path to the wannier90_hr.dat file containing the graphene bulk calculation poscarfile: path to the VASP POSCAR file of the graphene bulk calculation width: width of the ribbon (number of rings). Must be an even number. output: path to the output image file. usedhoppingcells_rings: If you don’t want to use all hopping parameters, you can set the number of ‘rings’ surrounding the main cell here. If it is a list (e.g. range(5)), several plots are created. The default value is ‘all’.

envtb.simple.plot_zigzag_graphene_nanoribbon_pz_bandstructure_nn(nnfile, width, output, length, magnetic_B=None)[source]

Plot the pi bandstructure of a zigzag graphene nanoribbon based on a n-th nearest neighbour parameterization of bulk graphene. A .dat file with the bandstructure data is also saved as output.dat (each line representing one k-point).

nnfile: path to the nearest-neighbour input file (see example files) width: width of the ribbon (number of rings). output: path to the output image file. length: length in x-direction of the returned (not the bandstructure plot!) unit cell (in rings)

Return: Hamiltonian of the unit cell with given length

envtb.simple.plot_armchair_graphene_nanoribbon_pz_bandstructure_nn(nnfile, width, output, magnetic_B=None)[source]

Plot the pi bandstructure of an armchair graphene nanoribbon based on a n-th nearest neighbour parameterization of bulk graphene. In this example function, width is the number of rings in transversal direction. A .dat file with the bandstructure data is also saved as output.dat (each line representing one k-point).

nnfile: path to the nearest-neighbour input file (see example files) width: width of the ribbon (number of rings). Must be an even number. output: path to the output image file.

Table Of Contents

Previous topic

Quantum capacitance

Next topic

Tight Binding

This Page