PtH2Pt nanojunction

Experimental Motivation

• R. H. M. Smit, Y. Noat, C. Untiedt, N. D. Lang, M. C. van Hemert and J. M. van Ruitenbeek, Measurement of the conductance of a hydrogen molecule, Nature 419, 906 (2002). [PDF]
• M. Kiguchi, R. Stadler, I. S. Kristensen, D. Djukic, and J. M. van Ruitenbeek, Evidence for a single hydrogen molecule connected by an atomic chain, Phys. Rev. Lett. 98, 146802 (2007). [PDF]

NEGF+Tight-Binding modeling of electronic transport

• GPAW scripts

NEGF+DFT for zero-bias transmission function

• Calculate zero-bias transmission function using pt_h2_trans.py script:

#Transport calculations using transport object following https://wiki.fysik.dtu.dk/gpaw/documentation/transport/negftransport.html

from ase import Atoms
from gpaw.transport.calculator import Transport 
from gpaw.atom.basis import BasisMaker
from gpaw.occupations import FermiDirac
from gpaw.poisson import PoissonSolver
from gpaw.mixer import Mixer
from ase.visualize import view

a = 2.41  # Pt binding lenght
b = 0.90  # H2 binding lenght
c = 1.70  # Pt-H binding lenght
L = 7.00  # width of unit cell

#  L-Lead         scat region       R-Lead
#------------ ------------------- -----------                     
# Pt--Pt--Pt-|-Pt--Pt-H-H-Pt--Pt-|-Pt--Pt--Pt
#  0   1   2    3  4  5 6  7   8   9   10  11
#------------ ------------------- -----------                     

atoms = Atoms('Pt5H2Pt5', pbc=( 0, 0,1), cell=[ L, L, 9 * a + b + 2 * c])
atoms.positions[:5,  2] = [i * a for i in range(5)]
atoms.positions[-5:, 2] = [i * a + b + 2 * c for i in range(4, 9)]
atoms.positions[5:7, 2] = [4 * a + c, 4 * a + c + b]
atoms.positions[:, 0:1] = L / 2.
atoms.center()

# setup leads
pl_atoms1 = range(3)    #  3 atoms 0~2  is L-Lead
pl_atoms2 = range(9,12) #  3 atoms 9~11 is R-Lead
pl_cell1 = (L, L, 3 * a) # cell size of lead  3 Pt bonds 
pl_cell2 = pl_cell1
# visulize device with ag
view(atoms)

t = Transport(h=0.3,
              xc='PBE',
              basis='szp(dzp)',
              kpts=(1,1,1),  # 
              occupations=FermiDirac(0.1),
              mode='lcao',
              poissonsolver=PoissonSolver(nn=2, relax='GS'),
              txt='pt_h2_trans.txt',
              mixer=Mixer(0.1, 5, weight=100.0),
              pl_atoms=[pl_atoms1, pl_atoms2],
              pl_cells=[pl_cell1, pl_cell2],
              pl_kpts=(1,1,10), # lead is periodic along transport direction 
              plot_energy_range = [-4.,4.],    # min and max energy of transmission
              plot_energy_point_num = 201,     # Number of energy points
              #edge_atoms=[[0, 2], [0, 5]],    # edge and mol_atoms should be 
              #mol_atoms=range(1,5),           # specified to be able to restart.
              #analysis_mode=True,             # for restarting jobs 
              #scat_restart=True,              # need to specify mol_atom
              #lead_restart=True,
              #guess_steps=1,       
              non_sc=True,                    #True = Normal DFT (default) for zero bias transmission, False = NEGF-DFT
              )

atoms.set_calculator(t)
t.calculate_iv() # for zero-bias transmission, do not put any arguments

• Plot zero-bias transmission function using plot_trans.py script:

from gpaw.transport.analysor import Transport_Plotter
import numpy as np
import sys
from pylab import *

fd=0
plotter=Transport_Plotter(fd)
plotter.plot_setup()
# the following Emax , nE  should be the same as specified in the transport calculations
Emax=3.0
Emin=-Emax
nE=201   
bias_step=0   

tc = plotter.tc(bias_step) # transmission 
ee=np.linspace(Emin,Emax,nE) # Energy range specified in transport calculation
plot(ee, tc, 'b-o')

xlabel('Energy (eV)')
ylabel('Transmission')
show()

NEGF+DFT for I-V characteristics

• Calculate transmission function at finite bias voltage using pt_h2_iv.py script:

#Transport calculations using transport object following https://wiki.fysik.dtu.dk/gpaw/documentation/transport/negftransport.html

from ase import Atoms
from gpaw.transport.calculator import Transport 
from gpaw.atom.basis import BasisMaker
from gpaw.occupations import FermiDirac
from gpaw.poisson import PoissonSolver
from gpaw.mixer import Mixer
from ase.visualize import view

a = 2.41  # Pt binding lenght
b = 0.90  # H2 binding lenght
c = 1.70  # Pt-H binding lenght
L = 7.00  # width of unit cell

#  L-Lead         scat region       R-Lead
#------------ ------------------- -----------                     
# Pt--Pt--Pt-|-Pt--Pt-H-H-Pt--Pt-|-Pt--Pt--Pt
#  0   1   2    3  4  5 6  7   8   9   10  11
#  |                                       |  
# C_i                                     C_f
#------------ ------------------- -----------                     
#  0   1   2                       0   1   2
#  |                                       |
# L_i                                     R_f
atoms = Atoms('Pt5H2Pt5', pbc=( 0, 0,1), cell=[ L, L, 9 * a + b + 2 * c])
atoms.positions[:5,  2] = [i * a for i in range(5)]
atoms.positions[-5:, 2] = [i * a + b + 2 * c for i in range(4, 9)]
atoms.positions[5:7, 2] = [4 * a + c, 4 * a + c + b]
atoms.positions[:, 0:1] = L / 2.
atoms.center()

# setup leads
pl_atoms1 = range(3)    #  3 atoms 0~2  is L-Lead
pl_atoms2 = range(9,12) #  3 atoms 9~11 is R-Lead
pl_cell1 = (L, L, 3 * a) # cell size of lead  3 Pt bonds 
pl_cell2 = pl_cell1
# visulize device with ag
view(atoms)

t = Transport(h=0.3,
              xc='PBE',
              basis='szp(dzp)',
              kpts=(1,1,1),  # 
              occupations=FermiDirac(0.1),
              mode='lcao',
              poissonsolver=PoissonSolver(nn=2, relax='GS'),
              txt='pt_h2_iv.txt',
              mixer=Mixer(0.1, 5, weight=100.0),
              pl_atoms=[pl_atoms1, pl_atoms2],
              pl_cells=[pl_cell1, pl_cell2],
              pl_kpts=(1,1,10), # lead is periodic along transport direction 
              plot_energy_range = [-4.,4.],    # min and max energy of transmission
              plot_energy_point_num = 201,     # Number of energy points
              edge_atoms=[[0, 2], [0, 11]],    # Defining edge atoms in the form [[L_i, R_f][C_i, C_f]] 
              mol_atoms=range(3,9),            # specified to be able to restart.
              #analysis_mode=True,             # for restarting jobs 
              #scat_restart=True,              # need to specify mol_atom
              #lead_restart=True,
              guess_steps=1,       
              #non_sc=False,                   #True = Normal DFT (default)for zero bias transmission, False = NEGF-DFT
              fixed_boundary=False
              )

atoms.set_calculator(t)
t.calculate_iv(0.5, 3)  # bias voltage set in the interval 0 to 0.5 eV using 3 steps: 0, 0.25 and 0.5 V
                        # could have 3rd argument specified the nth point of bias voltage to begin with

• Plot I-V characteristics using plot_iv.py script:

from gpaw.transport.analysor import Transport_Plotter
import numpy as np
from pylab import *
import sys
plotter = Transport_Plotter()
plotter.plot_setup()

nbias=2   #number of bias points to plot  

bias, current = plotter.iv(nbias)
bias=np.abs(bias)
plot(bias, current, 'r-o')
xlabel('Bias Voltage (V)')
ylabel('Current (microA)')
show()