Topological Hall effect in four terminal devices: Difference between revisions

From phys824
Jump to navigationJump to search
(Created page with "zgnr_adatoms_4t.py <pre> from math import sqrt import random import itertools as it import tinyarray as ta import numpy as np import matplotlib.pyplot as plt import kwant cla...")
 
No edit summary
 
(5 intermediate revisions by the same user not shown)
Line 1: Line 1:
zgnr_adatoms_4t.py
KWANT script below describes a skyrmion and computes the Hall Resistance and Spin Hall angle of the four-terminal device.
 
<pre>
<pre>
from math import sqrt
from types import SimpleNamespace
import random
from matplotlib import pyplot
import itertools as it
from math import cos, sin, pi
import tinyarray as ta
import numpy as np
import numpy as np
import matplotlib.pyplot as plt
import scipy.stats as reg
import kwant
import kwant
from sympy import S, Eq, solve, Symbol


class Honeycomb(kwant.lattice.Polyatomic):
lat = kwant.lattice.square()
    """Honeycomb lattice with methods for dealing with hexagons"""


    def __init__(self, name=''):
s_0 = np.identity(2)
        prim_vecs = [[0.5, sqrt(3)/2], [1, 0]] # bravais lattice vectors
s_z = np.array([[1, 0], [0, -1]])
        # offset the lattice so that it is symmetric around x and y axes
s_x = np.array([[0, 1], [1, 0]])
        basis_vecs = [[-0.5, -1/sqrt(12)], [-0.5, 1/sqrt(12)]]
s_y = np.array([[0, -1j], [1j, 0]])
        super(Honeycomb, self).__init__(prim_vecs, basis_vecs, name)
        self.a, self.b = self.sublattices


    def hexagon(self, tag):
def HedgeHog(site,ps):
         """ Get sites belonging to hexagon with the given tag.
         x,y = site.pos
            Returns sites in counter-clockwise order starting
         r = ( x**2 + y**2 )**0.5
            from the lower-left site.
         theta = (np.pi/2)*(np.tanh((ps.r0 - r)/ps.delta) + 1)
         """
         if r != 0:
         tag = ta.array(tag)
            Ex = (x/r)*np.sin(theta)*s_x + (y/r)*np.sin(theta)*s_y + np.cos(theta)*s_z
         #        a-sites b-sites
         else:
        deltas = [(0, 0), (0, 0),
            Ex = s_z
                  (1, 0), (0, 1),
         return 4*s_0 + ps.Ex * Ex
                  (0, 1), (-1, 1)]
     
         lats = it.cycle(self.sublattices)
         return (lat(*(tag + delta)) for lat, delta in zip(lats, deltas))


    def hexagon_neighbors(self, tag, n=1):
def Lead_Pot(site,ps):
        """ Get n'th nearest neighbor hoppings within the hexagon with
    return 4*s_0 + ps.Ex * s_z
            the given tag.
        """
        hex_sites = list(self.hexagon(tag))
        return ((hex_sites[(i+n)%6], hex_sites[i%6]) for i in range(6))


def random_placement(builder, lattice, density):
def MakeSystem(ps, show = False):
     """ Randomly selects hexagon tags where adatoms can be placed.
     H  = kwant.Builder()
        This avoids the edge case where adatoms would otherwise
        be placed on incomplete hexagons at the system boundaries.
    """
    assert 0 <= density <= 1
    system_sites = builder.sites()


     def hexagon_in_system(tag):
     def shape_2DEG(pos):
         return all(site in system_sites for site in lattice.hexagon(tag))
        x,y = pos
         return ( (abs(x) < ps.L) and (abs(y) < ps.W) ) or (
            (abs(x) < ps.W) and (abs(y) < ps.L))
   
    H[lat.shape(shape_2DEG,(0,0))] = HedgeHog
    H[lat.neighbors()] = -s_0
   
    # Leads
    sym_x = kwant.TranslationalSymmetry((-1,0))
    H_lead_x = kwant.Builder(sym_x)
    shape_x = lambda pos: abs(pos[1])<ps.W and pos[0]==0
    H_lead_x[lat.shape(shape_x,(0,0))] = Lead_Pot
    H_lead_x[lat.neighbors()] = -s_0
 
    sym_y = kwant.TranslationalSymmetry((0,-1))
    H_lead_y = kwant.Builder(sym_y)
    shape_y = lambda pos: abs(pos[0])<ps.W and pos[1]==0
    H_lead_y[lat.shape(shape_y,(0,0))] = Lead_Pot
    H_lead_y[lat.neighbors()] = -s_0
   
    H.attach_lead(H_lead_x)
    H.attach_lead(H_lead_y)
    H.attach_lead(H_lead_y.reversed())
    H.attach_lead(H_lead_x.reversed())
   
    if show:
        kwant.plot(H)


     # get set of tags for sites in system (i.e. tags from only one sublattice)
     return H
    system_tags = (s.tag for s in system_sites if s.family == lattice.a)
    # only allow tags which have complete hexagons in the system
    valid_tags = list(filter(hexagon_in_system, system_tags))
    N = int(density * len(valid_tags))
    total_hexagons=len(valid_tags)
    valuef=random.sample(valid_tags, N)
    return valuef


def cross(W, L):
def Transport(Hf,EE,ps):
     def shape(pos):
     smatrix = kwant.smatrix(Hf, energy=EE, args=[ps])
         return ((-W <= pos[0] <= W and -L <= pos[1] <= L) or  # vertical
    G=np.zeros((4,4))
                (-W <= pos[1] <= W and -L <= pos[0] <= L))   # horizontal
    for i in range(4):
    return shape
         a=0
        for j in range(4):
            G[i,j] = smatrix.transmission(i, j)
            if i != j:
                a += G[i,j]
        G[i,i] = -a
       
    V = np.linalg.solve(G[:3,:3], [1.,0,0])
    Hall = V[2] - V[1]
   
    return G, Hall
ps = SimpleNamespace(L=45, W=40, delta=10, r0=20, Ex=1.)


## Pauli matrices ##
H = MakeSystem(ps, show=True)
sigma_0 = ta.array([[1, 0], [0, 1]])  # identity
Hf = H.finalized()
sigma_x = ta.array([[0, 1], [1, 0]])
sigma_y = ta.array([[0, -1j], [1j, 0]])
sigma_z = ta.array([[1, 0], [0, -1]])


## Hamiltonian value functions ##
def Vz(site):
    Hd = HedgeHog(site,ps)
    return (Hd[0,0] - Hd[1,1]).real


def onsite_potential(site, params):
def Vy(site):
     return params['ep'] * sigma_0
    Hd = HedgeHog(site, ps)
     return Hd[0,1].imag


def potential_shift(site, params):
kwant.plotter.map(H, Vz);
    return params['mu'] * sigma_0
kwant.plotter.map(H, Vy);


def kinetic(site_i, site_j, params):
# Hall Resistance
    return -params['gamma'] * sigma_0
ps = SimpleNamespace(L=20, W=15, delta=3, r0=6, Ex=1.)


def rashba(site_i, site_j, params):
H = MakeSystem(ps, show=False)
    d_ij = site_j.pos - site_i.pos
Es = np.linspace(0.1,5.,50)
    rashba = 1j * params['V_R'] * (sigma_x * d_ij[1] - sigma_y * d_ij[0])
Hf = H.finalized()
    return rashba + kinetic(site_i, site_j, params)
dataG , dataHall = [],[]


def spin_orbit(site_i, site_j, params):
for EE in Es:
     so = (1j/3.) * params['V_I'] * sigma_z
    ps.delta = EE
    return so
    energy = 2.
    G,Hall = Transport(Hf, energy, ps)
     dataHall.append(Hall)


lat = Honeycomb()
pyplot.plot(Es, dataHall, 'o-', label="Skyrmion")
pv1, pv2 = lat.prim_vecs
pyplot.xlabel('Domain width $\delta$')
ysym = kwant.TranslationalSymmetry(pv2 - 2*pv1) # lattice symmetry in -y direction
pyplot.ylabel('Hall Resistance')
xsym = kwant.TranslationalSymmetry(-pv2) # lattice symmetry in -x direction
pyplot.title('Topologial Hall  Resistance')
pyplot.legend();
pyplot.show()


# adatom lattice, for visualization only
# Spin Hall angle
adatom_lat = kwant.lattice.Monatomic(lat.prim_vecs, name='adatom')
def solvematrix4(G,I0,I2,V1,V3):
    # I=GxV
    I1,I3,V0,V2=S('I1 I3 V0 V2'.split())
    equations = [Eq(I0, G[0][0]*V0+G[0][1]*V1+G[0][2]*V2+G[0][3]*V3), \
                Eq(I1, G[1][0]*V0+G[1][1]*V1+G[1][2]*V2+G[1][3]*V3), \
                Eq(I2, G[2][0]*V0+G[2][1]*V1+G[2][2]*V2+G[2][3]*V3), \
                Eq(I3, G[3][0]*V0+G[3][1]*V1+G[3][2]*V2+G[3][3]*V3)]
    return solve(equations, [I1,I3,V0,V2], set=True, dict=False, rational=None, manual=True, minimal=True, particular=True, implicit=False, quick=True)


def site_size(s):                                        
def plot_SH_angle(fsys, energies, ps):
     return 0.25 if s.family in lat.sublattices else 0.35
     data = []


def site_color(s):
     I0=0
     return '#000000' if s.family in lat.sublattices else '#808080'
     I2=0
 
     V1=1.0
def create_lead_h(W, symmetry, axis=(0, 0)):
     V3=0
    lead = kwant.Builder(symmetry)
     lead[lat.wire(axis, W)] = 0. * sigma_0
    lead[lat.neighbors(1)] = kinetic
    return lead
 
def create_lead_v(L,W, symmetry):
    axis=(0, 0)
    lead = kwant.Builder(symmetry)
    lead[lat.wire(axis, W)] = 0. * sigma_0
    lead[lat.neighbors(1)] = kinetic
    return lead
 
def create_hall_cross(W, L, adatom_density=0.2, show_adatoms=False):
    ## scattering region ##
    sys = kwant.Builder()
    sys[lat.shape(cross(W, L), (0, 0))] = onsite_potential
    sys[lat.neighbors(1)] = kinetic
 
    adatoms = random_placement(sys, lat, adatom_density)
 
    sys[(lat.hexagon(a) for a in adatoms)] = potential_shift
    sys[(lat.hexagon_neighbors(a, 1) for a in adatoms)] = rashba
    sys[(lat.hexagon_neighbors(a, 2) for a in adatoms)] = spin_orbit
    if show_adatoms:
        # no hoppings are added so these won't affect the dynamics; purely cosmetic
        sys[(adatom_lat(*a) for a in adatoms)] = 0
     ## leads ##
    leads = [create_lead_h(W, xsym)]
    leads += [create_lead_v(L,W, ysym)]
    leads += [lead.reversed() for lead in leads]  # right and bottom leads
    for lead in leads:
        sys.attach_lead(lead)
    return sys
 
def plot_bands(sys, params):
    fsys = sys.finalized()
    fig, (ax1, ax2) = plt.subplots(1, 2)
    ax1.set_title('left lead band structure')
    ax2.set_title('top lead band structure')
    for ax in (ax1, ax2):
        ax.set_xlabel('$k$', size=20)
        ax.set_ylabel('$E$', size=20)
    kwant.plotter.bands(fsys.leads[0], ax=ax1, args=(params,))
    kwant.plotter.bands(fsys.leads[1], ax=ax2, args=(params,))
    fig.tight_layout()
    plt.show()
 
def plot_ldos(sys, Es, params):
    fsys = sys.finalized()
    ldos = kwant.ldos(fsys, energy=Es, args=(params,))
     ldos = ldos[::2] + ldos[1::2] # sum spins
    kwant.plotter.map(fsys, ldos, vmax=0.1)
 
def plot_conductance(sys, energies,lead_i,lead_j, params):
    fsys = sys.finalized()
    data = []


     for energy in energies:
     for energy in energies:
         smatrix = kwant.smatrix(fsys, energy, args=(params,))
         smatrix = kwant.greens_function(fsys, energy, args=[ps])
         data.append(smatrix.transmission(lead_i,lead_j))
         G=np.zeros((4,4))
 
        for i in range(4):
    plt.figure()
            a=0
    plt.plot(energies, data)
            for j in range(4):
    plt.xlabel("energy (t)")
                G[i,j]=smatrix.transmission(i, j)
    plt.ylabel("conductance (e^2/h)")
                if i != j: a+=G[i,j]
    plt.show()
            G[i,i]=-a
 
        a=solvematrix4(G,I0,I2,V1,V3)[1]
if __name__ == '__main__':
        result=list(a)
#    params = dict(gamma=1., ep=0, mu=0., V_I=0.017, V_R=0)      # Tl adatoms
    params = dict(gamma=1., ep=0, mu=0, V_I=0.0032, V_R=0)       # In adatoms


    W=16
        I1=result[0][0]
    L=32
        I3=result[0][1]
    adatom_density=0.3
        V0=result[0][2]
        V2=result[0][3]
        print(I1,I3,V0,V2)


    sys = create_hall_cross(W, L, adatom_density, show_adatoms=True)
        Gspin=np.zeros((4,4))
    kwant.plot(sys, site_color=site_color, site_size=site_size)
        for i in range(4):
            a=0
            for j in range(4):
                ttdagger = smatrix._a_ttdagger_a_inv(i, j)
                sigma_z_matrix = np.kron(np.eye(len(ttdagger)/2),s_z)   
                Gspin[i,j]=np.trace(sigma_z_matrix.dot(ttdagger)).real
                if i != j: a+=Gspin[i,j]
            Gspin[i,i]=-a


    plot_bands(sys, params)
        IS0=Gspin[1,0]*(V1-V0)+Gspin[2,0]*(V2-V0)+Gspin[3,0]*(V3-V0)
        data.append(IS0/I1)


     sys = create_hall_cross(W, L, adatom_density, show_adatoms=False)
     pyplot.figure()
     plot_ldos(sys, 1e-5, params)
     pyplot.plot(energies, data)
      
     pyplot.xlabel("energy (t)")
     Es = np.linspace(-0.2, 0.2, 100)
     pyplot.ylabel("spin hall angle")
     plot_conductance(sys, Es, 1, 0, params)   # T01
     pyplot.show()
    plot_conductance(sys, Es, 2, 0, params)   # T02
plot_SH_angle(Hf, Es, ps)
</pre>
</pre>

Latest revision as of 17:30, 7 December 2016

KWANT script below describes a skyrmion and computes the Hall Resistance and Spin Hall angle of the four-terminal device. 

from types import SimpleNamespace
from matplotlib import pyplot
from math import cos, sin, pi
import numpy as np
import scipy.stats as reg
import kwant
from sympy import S, Eq, solve, Symbol

lat = kwant.lattice.square()

s_0 = np.identity(2)
s_z = np.array([[1, 0], [0, -1]])
s_x = np.array([[0, 1], [1, 0]])
s_y = np.array([[0, -1j], [1j, 0]])

def HedgeHog(site,ps):
        x,y = site.pos
        r = ( x**2 + y**2 )**0.5
        theta = (np.pi/2)*(np.tanh((ps.r0 - r)/ps.delta) + 1)
        if r != 0:
            Ex = (x/r)*np.sin(theta)*s_x + (y/r)*np.sin(theta)*s_y + np.cos(theta)*s_z
        else:
            Ex = s_z
        return 4*s_0 + ps.Ex * Ex
       

def Lead_Pot(site,ps):
    return  4*s_0 +  ps.Ex * s_z

def MakeSystem(ps, show = False):
    H  = kwant.Builder()

    def shape_2DEG(pos):
        x,y = pos
        return  ( (abs(x) < ps.L) and (abs(y) < ps.W) ) or ( 
            (abs(x) < ps.W) and (abs(y) < ps.L))
    
    H[lat.shape(shape_2DEG,(0,0))] = HedgeHog
    H[lat.neighbors()] = -s_0
    
    # Leads
    sym_x = kwant.TranslationalSymmetry((-1,0))
    H_lead_x = kwant.Builder(sym_x)
    shape_x = lambda pos: abs(pos[1])<ps.W and pos[0]==0 
    H_lead_x[lat.shape(shape_x,(0,0))] = Lead_Pot
    H_lead_x[lat.neighbors()] = -s_0
   
    sym_y = kwant.TranslationalSymmetry((0,-1))
    H_lead_y = kwant.Builder(sym_y)
    shape_y = lambda pos: abs(pos[0])<ps.W and pos[1]==0 
    H_lead_y[lat.shape(shape_y,(0,0))] = Lead_Pot 
    H_lead_y[lat.neighbors()] = -s_0
    
    H.attach_lead(H_lead_x)
    H.attach_lead(H_lead_y)
    H.attach_lead(H_lead_y.reversed())
    H.attach_lead(H_lead_x.reversed())
    
    if show:
        kwant.plot(H)

    return H
 

def Transport(Hf,EE,ps):
    smatrix = kwant.smatrix(Hf, energy=EE, args=[ps])
    G=np.zeros((4,4))
    for i in range(4):
        a=0
        for j in range(4): 
            G[i,j] = smatrix.transmission(i, j)
            if i != j:
                a += G[i,j]
        G[i,i] = -a 
        
    V = np.linalg.solve(G[:3,:3], [1.,0,0])
    Hall = V[2] - V[1]
    
    return G, Hall
	
ps = SimpleNamespace(L=45, W=40, delta=10, r0=20, Ex=1.)

H = MakeSystem(ps, show=True)
Hf = H.finalized()

def Vz(site):
    Hd = HedgeHog(site,ps)
    return (Hd[0,0] - Hd[1,1]).real 

def Vy(site):
    Hd = HedgeHog(site, ps)
    return Hd[0,1].imag 

kwant.plotter.map(H, Vz);
kwant.plotter.map(H, Vy);

# Hall Resistance
ps = SimpleNamespace(L=20, W=15, delta=3, r0=6, Ex=1.)

H = MakeSystem(ps, show=False)
Es = np.linspace(0.1,5.,50)
Hf = H.finalized()
dataG , dataHall = [],[]

for EE in Es: 
    ps.delta = EE
    energy = 2.
    G,Hall = Transport(Hf, energy, ps)
    dataHall.append(Hall)

pyplot.plot(Es, dataHall, 'o-', label="Skyrmion")
pyplot.xlabel('Domain width $\delta$')
pyplot.ylabel('Hall Resistance')
pyplot.title('Topologial Hall  Resistance')
pyplot.legend();
pyplot.show()

# Spin Hall angle
def solvematrix4(G,I0,I2,V1,V3):
    # I=GxV
    I1,I3,V0,V2=S('I1 I3 V0 V2'.split())
    equations = [Eq(I0, G[0][0]*V0+G[0][1]*V1+G[0][2]*V2+G[0][3]*V3), \
                 Eq(I1, G[1][0]*V0+G[1][1]*V1+G[1][2]*V2+G[1][3]*V3), \
                 Eq(I2, G[2][0]*V0+G[2][1]*V1+G[2][2]*V2+G[2][3]*V3), \
                 Eq(I3, G[3][0]*V0+G[3][1]*V1+G[3][2]*V2+G[3][3]*V3)]
    return solve(equations, [I1,I3,V0,V2], set=True, dict=False, rational=None, manual=True, minimal=True, particular=True, implicit=False, quick=True)

def plot_SH_angle(fsys, energies, ps):
    data = []

    I0=0
    I2=0
    V1=1.0
    V3=0

    for energy in energies:
        smatrix = kwant.greens_function(fsys, energy, args=[ps])
        G=np.zeros((4,4))
        for i in range(4):
            a=0
            for j in range(4):
                G[i,j]=smatrix.transmission(i, j)
                if i != j: a+=G[i,j]
            G[i,i]=-a
        a=solvematrix4(G,I0,I2,V1,V3)[1]
        result=list(a)

        I1=result[0][0]
        I3=result[0][1]
        V0=result[0][2]
        V2=result[0][3]
        print(I1,I3,V0,V2)

        Gspin=np.zeros((4,4))
        for i in range(4):
            a=0
            for j in range(4):
                ttdagger = smatrix._a_ttdagger_a_inv(i, j)
                sigma_z_matrix = np.kron(np.eye(len(ttdagger)/2),s_z)    
                Gspin[i,j]=np.trace(sigma_z_matrix.dot(ttdagger)).real
                if i != j: a+=Gspin[i,j]
            Gspin[i,i]=-a

        IS0=Gspin[1,0]*(V1-V0)+Gspin[2,0]*(V2-V0)+Gspin[3,0]*(V3-V0)
        data.append(IS0/I1)

    pyplot.figure()
    pyplot.plot(energies, data)
    pyplot.xlabel("energy (t)")
    pyplot.ylabel("spin hall angle")
    pyplot.show()
plot_SH_angle(Hf, Es, ps)
Retrieved from "https://wiki.physics.udel.edu/wiki_phys824/index.php?title=Topological_Hall_effect_in_four_terminal_devices&oldid=2445"