In [ ]:
from __future__ import print_function
import sisl
import numpy as np
import matplotlib.pyplot as plt
%matplotlib inline
In this example you will familiarize your-self to the concept of buffer atoms. A buffer atom is an atom that is completely neglected in the self-consistent calculation, but is used as an initialization for the bulk electrode regions.
Here a pristine graphene flake will be constructed and subsequently a Carbon chain will act as an STM-like tip to simulate STM experiments.
As the Carbon chain is terminated to vacuum, the dangling bonds will create spurious effects very different from a pristine, bulk chain. To understand why it is necessary to add buffer atoms it is useful to understand the TranSiesta method. Any TranSiesta calculation starts with calculating an initial guess for the Hamiltonian as input for the Green function method: \begin{equation} \mathbf G^{-1}(E) = \mathbf S (E+i\eta) - \mathbf H - \sum_i\boldsymbol\Sigma_i \end{equation} If the initial $\mathbf H'$ represents a Hamiltonian close to the open-boundary problem $\mathbf H$; it will converge with a higher probability, and in much less time. Improving the initial guess Hamiltonian is time well-worth spent as TranSiesta is, typically, more difficult to converge.
As an example consider the Hamiltonian for the chain:
<vacuum> C -- C -- C -- C -- C -- C ...
It is clear that the atom closest to the vacuum region resides in a very different chemical and potential landscape than an atom in the middle of the chain. If TranSiesta uses the initial Hamiltonian for the chain electrode as the atom closest to the vacuum region it will be very far from the potential landscape of a bulk electrode. So to mitigate this one can specify:
%block TBT.Atoms.Buffer
atom [ 1 -- 4 ]
%endblock
to remove the first 4 atoms from the calculation (note that negative indices counts from the end). Then the electrode will begin from the 5th atom which is relatively far from the dangling bond. This will be a much better initial guess for the Hamiltonian.
In [ ]:
graphene = sisl.geom.graphene(1.44)
elec = graphene.tile(2, axis=0)
elec.write('ELEC_GRAPHENE.fdf')
elec.write('ELEC_GRAPHENE.xyz')
In [ ]:
C1d = sisl.Geometry([[0,0,0]], graphene.atom[0], [10, 10, 1.4])
elec_chain = C1d.tile(4, axis=2)
elec_chain.write('ELEC_CHAIN.fdf')
elec_chain.write('ELEC_CHAIN.xyz')
chain = elec_chain.tile(4, axis=2)
In [ ]:
device = elec.repeat(5, axis=1).tile(4, axis=0)
# Attach the chain on-top of an atom
# First find an atom in the middle of the device
idx = device.close(device.center(what='xyz'), R=1.45)[1]
# Attach the chain at a distance of 2.25 along the third lattice vector
device = device.attach(idx, chain, 0, dist=2.25, axis=2)
# Add vacuum along chain, we really no not care how much vacuum, but it
# is costly on memory, not so much on performance.
device = device.add_vacuum(15, axis=2)
device.write('DEVICE.fdf')
device.write('DEVICE.xyz')
Create a new directory for a different range of buffer atoms, from 0 to 4, start by using 4 buffer atoms.
How does convergence behave for different number of buffer atoms?
REMARK there are 2 places in the fdf file you should change when changing the number of atoms (the electrode atom specification and the buffer atoms).
In [ ]:
tbt = sisl.get_sile('siesta.TBT.nc')
In [ ]: