In [1]:
import sys

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

init_notebook()
import scipy.sparse.linalg as sl

my = 0.5 * (pauli.sys0 + pauli.sysz)
s0s0sz = np.kron(pauli.s0s0, pauli.sz)
s0szsz = np.kron(pauli.s0sz, pauli.sz)
mys0 = np.kron(my, pauli.s0)
s0s0sx = np.kron(pauli.s0s0, pauli.sx)
s0s0sy = np.kron(pauli.s0s0, pauli.sy)
szsxsz = np.kron(pauli.szsx, pauli.sz)
s0sysz = np.kron(pauli.s0sy, pauli.sz)
sxsxsz = np.kron(pauli.sxsx, pauli.sz)
sysxsz = np.kron(pauli.sysx, pauli.sz)


def make_qshe_sc(l=40, w=10, lead=False):
    def shape(pos):
        (x, y) = pos
        return (1.0 * y ** 2 / l ** 2 + 1.0 * x ** 2 / w ** 2) <= 0.25

    def onsite(site, p):
        (x, y) = site.pos
        return (
            (p.M - 4 * p.B) * s0szsz
            - 4 * p.D * s0s0sz
            + p.gaps(x, y)[1] * mys0
            + p.gaps(x, y)[0] * s0s0sx
        )

    def hopx(site1, site2, p):
        return p.B * s0szsz + p.D * s0s0sz + 0.5j * p.A * szsxsz

    def hopy(site1, site2, p):
        return p.B * s0szsz + p.D * s0s0sz - 0.5j * p.A * s0sysz

    lat = kwant.lattice.square()
    syst = kwant.Builder()
    syst[lat.shape(shape, (0, 0))] = onsite
    syst[kwant.HoppingKind((1, 0), lat)] = hopx
    syst[kwant.HoppingKind((0, 1), lat)] = hopy

    if lead:
        sym = kwant.TranslationalSymmetry((0, 1))
        lead = kwant.Builder(sym)
        lead[lat(0, 0)] = 1.5 * p.B * s0szsz
        lead[kwant.HoppingKind((0, 1), lat)] = -p.B * s0szsz
        syst.attach_lead(lead)
        syst.attach_lead(lead.reversed())
    return syst.finalized()


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

    def onsite(site, p):
        (x, y) = site.pos
        return (
            (p.M - 4 * p.B) * s0szsz
            - 4 * p.D * s0s0sz
            + p.gaps(x, y)[1] * mys0
            + p.gaps(x, y)[0] * s0s0sx
        )

    def hopx(site1, site2, p):
        return p.B * s0szsz + p.D * s0s0sz + 0.5j * p.A * szsxsz

    def hopy(site1, site2, p):
        return p.B * s0szsz + p.D * s0s0sz - 0.5j * p.A * s0sysz

    lat = kwant.lattice.square()
    sym = kwant.TranslationalSymmetry((0, 1))
    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.finalized()


def make_2d_pwave(w, l):
    def shape(pos):
        (x, y) = pos
        return (1.0 * y ** 2 / l ** 2 + 1.0 * x ** 2 / w ** 2) <= 0.25

    def hopx(site1, site2, p):
        (x1, y1) = site1.pos
        (x2, y2) = site2.pos
        phi = p.phase(0.5 * (x1 + x2), 0.5 * (y1 + y2))
        return -p.t * pauli.sz + 1j * p.delta * (
            np.cos(phi) * pauli.sx + np.sin(phi) * pauli.sy
        )

    def hopy(site1, site2, p):
        (x1, y1) = site1.pos
        (x2, y2) = site2.pos
        phi = p.phase(0.5 * (x1 + x2), 0.5 * (y1 + y2))
        return -p.t * pauli.sz - 1j * p.delta * (
            np.cos(np.pi / 2 + phi) * pauli.sx + np.sin(np.pi / 2 + phi) * pauli.sy
        )

    def onsite(site1, p):
        return (4 * p.t - p.mu) * pauli.sz

    lat = kwant.lattice.square()
    syst = kwant.Builder()
    syst[lat.shape(shape, (w / 2 - 1, 0))] = onsite
    syst[kwant.HoppingKind((1, 0), lat)] = hopx
    syst[kwant.HoppingKind((0, 1), lat)] = hopy

    return syst.finalized()


def bhz_slab(l, w, h):
    lat = kwant.lattice.general(np.identity(3))
    syst = kwant.Builder()

    def shape(pos):
        (x, y, z) = pos
        return (0 <= z < h) and (1.0 * y ** 2 / l ** 2 + 1.0 * x ** 2 / w ** 2) <= 0.25

    def onsite(site, p):
        (x, y, z) = site.pos
        phi = p.phase(x, y)
        return (
            (p.C + 2 * p.D1 + 4 * p.D2) * s0s0sz
            + (p.M + 2 * p.B1 + 4 * p.B2) * s0szsz
            + p.delta * (np.cos(phi) * s0s0sx + np.sin(phi) * s0s0sy)
        )

    def hopx(site1, site2, p):
        return -p.D2 * s0s0sz - p.B2 * s0szsz + p.A2 * 0.5j * sxsxsz

    def hopy(site1, site2, p):
        return -p.D2 * s0s0sz - p.B2 * s0szsz + p.A2 * 0.5j * sysxsz

    def hopz(site1, site2, p):
        return -p.D1 * s0s0sz - p.B1 * s0szsz + p.A1 * 0.5j * szsxsz

    syst[lat.shape(shape, (0, 0, 0))] = onsite
    syst[kwant.HoppingKind((1, 0, 0), lat)] = hopx
    syst[kwant.HoppingKind((0, 1, 0), lat)] = hopy
    syst[kwant.HoppingKind((0, 0, 1), lat)] = hopz

    return syst.finalized()


def calc_energies(syst, p, num_orbitals, num_states):
    ham = syst.hamiltonian_submatrix(params=dict(p=p), sparse=True).tocsc()
    energies, states = sl.eigsh(ham, sigma=0, k=num_states)
    densities = (
        np.linalg.norm(states.reshape(-1, num_orbitals, num_states), axis=1) ** 2
    )
    return energies, states, densities
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 2020-06-07 at 08:23:44.345304.