In [1]:
import sys

sys.path.append("../../code")
from init_mooc_nb import *

init_notebook()
%output size = 150

import warnings

warnings.simplefilter("ignore", UserWarning)

bhz_parameters = {
    "topo": {"A": 0.5, "B": 1.00, "D": 0.0, "M": 0.2},
    "topo2": {"A": 0.5, "B": 1.00, "D": 0.0, "M": 1.0},
    "triv": {"A": 0.5, "B": 1.00, "D": 0.0, "M": -0.2},
    "lead": {"A_lead": 1.5, "B_lead": 1.00, "D_lead": 0.0, "M_lead": 0.0},
}

# Onsite and hopping functions for the BHZ model.
# Sometimes, we use different BHZ parameters in the
# scattering region and leads, so we treat them
# separately.
def _onsite(site, p, is_lead=False):
    if is_lead:
        B, D, M = p.B_lead, p.D_lead, p.M_lead
    else:
        B, D, M = p.B, p.D, p.M
    return (
        (M - 4 * B) * pauli.s0sz
        - 4 * D * pauli.s0s0
        + p.ez_y * np.kron(pauli.sy, (pauli.s0 + pauli.sz) / 2)
    )


def onsite(site, p):
    return _onsite(site, p, is_lead=False)


def onsite_lead(site, p):
    return _onsite(site, p, is_lead=True)


def _hopx(site1, site2, p, is_lead=False):
    if is_lead:
        A, B, D = p.A_lead, p.B_lead, p.D_lead
    else:
        A, B, D = p.A, p.B, p.D
    return B * pauli.s0sz + D * pauli.s0s0 + 1j * A * pauli.szsx


def hopx(site1, site2, p):
    return _hopx(site1, site2, p, is_lead=False)


def hopx_lead(site1, site2, p):
    return _hopx(site1, site2, p, is_lead=True)


def _hopy(site1, site2, p, is_lead=False):
    if is_lead:
        A, B, D = p.A_lead, p.B_lead, p.D_lead
    else:
        A, B, D = p.A, p.B, p.D
    return B * pauli.s0sz + D * pauli.s0s0 - 1j * A * pauli.s0sy


def hopy(site1, site2, p):
    return _hopy(site1, site2, p, is_lead=False)


def hopy_lead(site1, site2, p):
    return _hopy(site1, site2, p, is_lead=True)


def two_terminal(L, w):
    """ Make a two terminal system with the BHZ model. """

    def shape(pos):
        (x, y) = pos
        return 0 <= y < w and 0 <= x < L

    def lead_shape(pos):
        (x, y) = pos
        return 0 <= y < w

    lat = kwant.lattice.square()
    syst = kwant.Builder()

    def onsite_two_term(site, p):
        return onsite(site, p) - p.mu * np.eye(4)

    def onsite_two_term_lead(site, p):
        return onsite_lead(site, p) - p.mu_lead * np.eye(4)

    # Scattering region
    syst[lat.shape(shape, (0, 0))] = onsite_two_term
    syst[kwant.HoppingKind((1, 0), lat)] = hopx
    syst[kwant.HoppingKind((0, 1), lat)] = hopy

    # Leads
    lead = kwant.Builder(kwant.TranslationalSymmetry((-1, 0)))
    lead[lat.shape(lead_shape, (0, 0))] = onsite_two_term_lead
    lead[kwant.HoppingKind((1, 0), lat)] = hopx_lead
    lead[kwant.HoppingKind((0, 1), lat)] = hopy_lead

    # Attach leads
    syst.attach_lead(lead)
    syst.attach_lead(lead.reversed())

    return syst


def bhz(w=None):
    """Translationally invariant BHZ system with a infinite or fixed width w."""
    lat = kwant.lattice.square()
    if w is None:
        syst = kwant.Builder(kwant.TranslationalSymmetry(*lat.prim_vecs))
        syst[lat.shape(lambda pos: True, (0, 0))] = onsite
    else:

        def ribbon_shape(pos):
            (x, y) = pos
            return 0 <= y < w

        sym = kwant.TranslationalSymmetry((1, 0))
        syst = kwant.Builder(sym)
        syst[lat.shape(ribbon_shape, (0, 0))] = onsite

    syst[kwant.HoppingKind((1, 0), lat)] = hopx
    syst[kwant.HoppingKind((0, 1), lat)] = hopy
    return syst


def G_mu_plot(p, mus, color):
    syst = two_terminal(40, 40).finalized()
    G = [
        kwant.smatrix(syst, energy=0.0, params=dict(p=p)).transmission(1, 0)
        for p.mu in mus
    ]
    ydim = r"G $[e^2/h]$"
    kdims = [r"$\mu$", ydim]
    plot = holoviews.Path((mus, np.array(G)), kdims=kdims, label="Conductance")
    ticks = {"xticks": [-0.8, -0.4, 0, 0.4, 0.8], "yticks": [0, 2, 4, 6, 8, 10]}
    return plot.redim.range(**{ydim: (0, 10)}).opts(plot=ticks, style={"color": color})


def G_Ez_plot(p, E_zs):
    syst = two_terminal(40, 20).finalized()
    G = [
        kwant.smatrix(syst, energy=0.0, params=dict(p=p)).transmission(1, 0)
        for p.ez_y in E_zs
    ]
    ydim = r"G $[e^2/h]$"
    kdims = [r"$E_z$", ydim]
    plot = holoviews.Path((E_zs, np.array(G)), kdims=kdims, label="Conductance")
    ticks = {"xticks": [0, 0.05, 0.10, 0.15], "yticks": [0, 0.5, 1.0, 1.5, 2.0]}
    return plot.redim.range(**{ydim: (0, 2)}).opts(plot=ticks)
Populated the namespace with:
np, matplotlib, kwant, holoviews, init_notebook, SimpleNamespace, pprint_matrix, scientific_number, pretty_fmt_complex, plt, pf, display_html
from code/edx_components:
MoocVideo, MoocDiscussion, MoocCheckboxesAssessment, MoocMultipleChoiceAssessment, MoocSelfAssessment
from code/functions:
spectrum, hamiltonian_array, h_k, pauli
Using kwant 1.4.2 and holoviews 1.13.2
Executed on 2021-03-06 at 14:55:56.487230.