Conductance of single-molecule: Difference between revisions

From phys824
Jump to navigationJump to search
No edit summary
 
(One intermediate revision by the same user not shown)
Line 1: Line 1:
===Using Tight Binding Method===
===Using Tight Binding Method===
Lead Hamiltonians:
 
<pre>
<pre>
import numpy as np
from ase.transport.calculators import TransportCalculator
import pylab
# onsite energies 0.0, nearest neighbor hopping -1.0, and
# second nearest neighbor hopping 0.2
H_lead = np.array([[ 0. , -1. ,  0.2,  0. ],
H_lead = np.array([[ 0. , -1. ,  0.2,  0. ],
                  [-1. ,  0. , -1. ,  0.2],
                    [-1. ,  0. , -1. ,  0.2],
                  [ 0.2, -1. ,  0. , -1. ],
                    [ 0.2, -1. ,  0. , -1. ],
                  [ 0. ,  0.2, -1. ,  0. ]])
                    [ 0. ,  0.2, -1. ,  0. ]])
</pre>
 
Scattering Hamiltonian:
H_scat = np.zeros((6, 6))
<pre>
 
H_scat = np.array([[ 0., -1 ,  0,   00, 0. ],
#Principal layers on either side of S
                  [-1 ,  0 , 0.2,  0, 0, 0. ],
H_scat[:2, :2] = H_scat[-2:, -2:] = H_lead[:2, :2]
                  [ 0 ,0.2 ,   0,-0.8, 0, 0. ],
 
                  [ 0 ,  0 ,-0.8,   0,0.2, 0. ],
#Scatering region (hydrogen molecule) - onsite 0.0 and hopping -0.8
                  [ 0 , 0 ,   0, 0.2, 0,-1 ],
H_scat[2:4, 2:4] = [[0.0, -0.8], [-0.8, 0.0]]
                  [ 0 , 0 ,   0,  0, -1, 0 ]])
 
#coupling to the leads - nearest neighbor only
H_scat[1, 2] = H_scat[2, 1] = H_scat[3, 4] = H_scat[4, 3] = 0.2
 
tcalc = TransportCalculator(h=H_scat,  # Scattering Hamiltonian
                            h1=H_lead, # Lead 1 (left)
                            h2=H_lead, # Lead 2  (right)
                            energies=np.arange(-3, 3, 0.02))
 
T_e = tcalc.get_transmission()
pylab.plot(tcalc.energies, T_e)
pylab.title('Transmission function')
pylab.show()
 
tcalc.set(pdos=[2, 3])
pdos_ne = tcalc.get_pdos()
pylab.plot(tcalc.energies, pdos_ne[0], ':')
pylab.plot(tcalc.energies, pdos_ne[1], '--')
pylab.title('Projected density of states')
pylab.show()
 
h_rot, s_rot, eps_n, vec_nn = tcalc.subdiagonalize_bfs([2, 3])
tcalc.set(h=h_rot, s=s_rot) # Set the rotated matrices
for n in range(2):
    print "eigenvalue, eigenvector:", eps_n[n],',', vec_nn[:, n]
 
pdos_rot_ne = tcalc.get_pdos()
pylab.plot(tcalc.energies, pdos_rot_ne[0], ':')
pylab.plot(tcalc.energies, pdos_rot_ne[1], '--')
pylab.title('Projected density of states (rotated)')
pylab.show()
 
h_cut, s_cut = tcalc.cutcoupling_bfs([2])
tcalc.set(h=h_cut, s=s_cut)
T_cut_bonding_e = tcalc.get_transmission()
pylab.plot(tcalc.energies, T_cut_bonding_e)
pylab.title('Transmission (bonding orbital cut)')
pylab.show()
</pre>
</pre>

Latest revision as of 22:22, 29 November 2012

Using Tight Binding Method

import numpy as np
from ase.transport.calculators import TransportCalculator
import pylab

# onsite energies 0.0, nearest neighbor hopping -1.0, and
# second nearest neighbor hopping 0.2
H_lead = np.array([[ 0. , -1. ,  0.2,  0. ],
                    [-1. ,  0. , -1. ,  0.2],
                    [ 0.2, -1. ,  0. , -1. ],
                    [ 0. ,  0.2, -1. ,  0. ]])

H_scat = np.zeros((6, 6))

#Principal layers on either side of S
H_scat[:2, :2] = H_scat[-2:, -2:] = H_lead[:2, :2]

#Scatering region (hydrogen molecule) - onsite 0.0 and hopping -0.8
H_scat[2:4, 2:4] = [[0.0, -0.8], [-0.8, 0.0]]

#coupling to the leads - nearest neighbor only
H_scat[1, 2] = H_scat[2, 1] = H_scat[3, 4] = H_scat[4, 3] = 0.2

tcalc = TransportCalculator(h=H_scat,  # Scattering Hamiltonian
                            h1=H_lead, # Lead 1 (left)
                            h2=H_lead, # Lead 2  (right)
                            energies=np.arange(-3, 3, 0.02))

T_e = tcalc.get_transmission()
pylab.plot(tcalc.energies, T_e)
pylab.title('Transmission function')
pylab.show()

tcalc.set(pdos=[2, 3])
pdos_ne = tcalc.get_pdos()
pylab.plot(tcalc.energies, pdos_ne[0], ':')
pylab.plot(tcalc.energies, pdos_ne[1], '--')
pylab.title('Projected density of states')
pylab.show()

h_rot, s_rot, eps_n, vec_nn = tcalc.subdiagonalize_bfs([2, 3])
tcalc.set(h=h_rot, s=s_rot) # Set the rotated matrices
for n in range(2):
    print "eigenvalue, eigenvector:", eps_n[n],',', vec_nn[:, n]

pdos_rot_ne = tcalc.get_pdos()
pylab.plot(tcalc.energies, pdos_rot_ne[0], ':')
pylab.plot(tcalc.energies, pdos_rot_ne[1], '--')
pylab.title('Projected density of states (rotated)')
pylab.show()

h_cut, s_cut = tcalc.cutcoupling_bfs([2])
tcalc.set(h=h_cut, s=s_cut)
T_cut_bonding_e = tcalc.get_transmission()
pylab.plot(tcalc.energies, T_cut_bonding_e)
pylab.title('Transmission (bonding orbital cut)')
pylab.show()