In [1]:
import sys

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

init_notebook()

import scipy
import scipy.linalg as sla
from matplotlib import cm


import warnings

warnings.simplefilter("ignore", UserWarning)


def bhz(L, W, H, system_type):
    """A cuboid region of BHZ material with two leads attached.

    parameters for leads and scattering region can be defined separately
    """
    # Onsite and hoppings matrices used for building BHZ model
    def onsite(site, p):
        return (p.C + 2 * p.D1 + 4 * p.D2) * pauli.s0s0 + (
            p.M + 2 * p.B1 + 4 * p.B2
        ) * pauli.s0sz

    def hopx(site1, site2, p):
        return -p.D2 * pauli.s0s0 - p.B2 * pauli.s0sz + p.A2 * 0.5j * pauli.sxsx

    def hopy(site1, site2, p):
        return -p.D2 * pauli.s0s0 - p.B2 * pauli.s0sz + p.A2 * 0.5j * pauli.sysx

    def hopz(site1, site2, p):
        return -p.D1 * pauli.s0s0 - p.B1 * pauli.s0sz + p.A1 * 0.5j * pauli.szsx

    def shape_lead(pos):
        (x, y, z) = pos
        return (0 <= z < H) and (0 <= y < W)

    def shape_syst(pos):
        (x, y, z) = pos
        return (0 <= z < H) and (0 <= y < W) and (0 <= x < L)

    def hopx_phase(site1, site2, p):
        x1, y1, z1 = site1.pos
        x2, y2, z2 = site2.pos
        return hopx(site1, site2, p) * np.exp(-0.5j * p.Bz * (x1 - x2) * (y1 + y2))

    lat = kwant.lattice.general(np.identity(3))

    if system_type == "sys":
        syst = kwant.Builder()
        syst[lat.shape(shape_syst, (0, 0, 0))] = lambda site, p: onsite(
            site, p
        ) - p.mu_scat * np.eye(4)
    elif system_type == "lead":
        sym = kwant.TranslationalSymmetry((1, 0, 0))
        syst = kwant.Builder(sym)
        syst[lat.shape(shape_lead, (0, 0, 0))] = lambda site, p: onsite(
            site, p
        ) - p.mu_lead * np.eye(4)
    elif system_type == "infinite":
        syst = kwant.Builder(kwant.TranslationalSymmetry(*lat.prim_vecs))
        syst[lat.shape(lambda pos: True, (0, 0))] = lambda site, p: onsite(
            site, p
        ) - p.mu_lead * np.eye(4)

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


def bhz_scatter(L, W, H):
    syst = bhz(L, W, H, "sys")
    lead = bhz(L, W, H, "lead")
    syst.attach_lead(lead)
    syst.attach_lead(lead.reversed())
    return syst


def cond_mu(p, L, W, H):
    p.mu_lead = 0.7
    syst = bhz_scatter(L, W, H)
    sys_leads_fixed = syst.finalized().precalculate(energy=0, params=dict(p=p))
    mus = np.linspace(-0.4, 0.4, 40)
    cond = [
        kwant.smatrix(
            sys_leads_fixed, energy=0, params=dict(p=p.update(mu_scat=mu))
        ).transmission(1, 0)
        for mu in mus
    ]
    return np.array(cond), mus


def plot_cond_mu(cond, mus):
    xdim, ydim = [r"$\mu$", r"$G\,[e^2/h]$"]
    kwargs = {"kdims": [xdim, ydim]}
    plot = holoviews.Path((mus, cond), **kwargs).opts(
        plot={"xticks": 3, "yticks": [0, 2, 4, 6, 8]}, style={"color": "r"}
    )
    return plot.redim.range(**{xdim: (-0.4, 0.4), ydim: (0, 8)}).relabel("Conductance")


def plot_bands(p, L, W, H):
    lead = bhz(L, W, H, "lead")
    kwargs = {
        "k_x": np.linspace(-np.pi / 3, np.pi / 3, 101),
        "ylims": [-1, 1],
        "yticks": 5,
        "xticks": [(-np.pi / 3, r"$-\pi/3$"), (0, r"$0$"), (np.pi / 3, r"$\pi/3$")],
    }
    p.mu_lead = 0
    return spectrum(lead, p, **kwargs)


def plot_cond_spect(mu, cond_plot, bands_plot):
    return cond_plot * holoviews.VLine(mu).opts(
        style={"color": "b"}
    ) + bands_plot.relabel("Spectrum") * holoviews.HLine(mu).opts(style={"color": "b"})


def plot_warping(A=1.2, B=1.8, C=1.5, Kmax=1.0):
    def evaluate_on_grid(X, Y, func):
        """ X, Y should be in np.meshgrid form. It's enough for func to work on floats. """
        data = []
        for xx, yy in zip(X, Y):
            data.append([func(i, j) for i, j in zip(xx, yy)])
        data = np.array(data)
        return data

    def get_energy_function(A, B, C):
        """ Used for plotting of hexagonal warping. """

        def func(kx, ky):
            matrix = (
                A * (kx ** 2 + ky ** 2) * pauli.s0
                + B * (kx * pauli.sy - ky * pauli.sx)
                + C * 0.5 * ((kx + 1j * ky) ** 3 + (kx - 1j * ky) ** 3) * pauli.sz
            )
            return sla.eigh(matrix)[0]

        return func

    zmin, zmax = -1.0, 3.5
    xylims = (-1.2, 1.2)
    zlims = (-1.0, 3.5)
    kdims = [r"$k_x$", r"$k_y$"]
    vdims = [r"E"]
    # Generate a circular mesh
    N = 100
    r = np.linspace(0, Kmax, N)
    p = np.linspace(0, 2 * np.pi, N)
    r, p = np.meshgrid(r, p)
    x, y = r * np.cos(p), r * np.sin(p)
    energies = evaluate_on_grid(x, y, func=get_energy_function(A, B, C))

    xy_ticks = [-1.2, 0, 1.2]
    zticks = [-1.0, 0.0, 1.0, 2.0, 3.0]
    style = {"xticks": xy_ticks, "yticks": xy_ticks, "zticks": zticks}
    kwargs = {
        "extents": (xylims[0], xylims[0], zlims[0], xylims[1], xylims[1], zlims[1]),
        "kdims": kdims,
        "vdims": vdims,
    }

    # hex_cmap colormap is defined below.
    plot = holoviews.Overlay(
        [
            holoviews.TriSurface(
                (x.flat, y.flat, energies[:, :, i].flat), **kwargs
            ).opts(style=dict(cmap=hex_cmap, linewidth=0), plot=style)
            for i in range(energies.shape[-1])
        ]
    )
    return plot.opts(plot={"Overlay": {"fig_size": 350}})


# Custom colormap for the hexagonal warping plot
cmap_list = [
    ((value + 1) / 4.0, colour)
    for value, colour in zip([-1.0, 0.0, 3.0], ["Blue", "White", "Red"])
]
hex_cmap = matplotlib.colors.LinearSegmentedColormap.from_list("custom", cmap_list)
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:50.463904.