Book a Demo!
CoCalc Logo Icon
StoreFeaturesDocsShareSupportNewsAboutPoliciesSign UpSign In
Download
6042 views
Kernel: Python 3
import sys sys.path.append('../code') from init_mooc_nb import * init_notebook() import wraparound from matplotlib import cm from matplotlib.colors import hsv_to_rgb pi_ticks = [(-np.pi, r'$-\pi$'), (0, '$0$'), (np.pi, r'$\pi$')] def d_wave(w=None, direction=None): """Creates a d-wave system. Parameters: ----------- w : int Width of the system, if None the system is infinite direction : str Direction of translation symmetry, if None it's an infinite system in x and y. """ def hopx(site1, site2, p): return -p.t * pauli.sz - p.delta * pauli.sx def hopy(site1, site2, p): return -p.t * pauli.sz + p.delta * pauli.sx def onsite(site, p): return (4 * p.t - p.mu) * pauli.sz lat = kwant.lattice.square() if not w: def ribbon_shape(pos): (x, y) = pos return True sym = kwant.TranslationalSymmetry(*lat.prim_vecs) else: if direction == 'topo': def ribbon_shape(pos): (x, y) = pos return (0 <= y - x < w) sym = kwant.TranslationalSymmetry((1, 1)) elif direction == 'triv': def ribbon_shape(pos): (x, y) = pos return (0 <= y < w) sym = kwant.TranslationalSymmetry((1, 0)) sys = kwant.Builder(sym) sys[lat.shape(ribbon_shape, (0, 0))] = onsite sys[kwant.HoppingKind((1, 0), lat)] = hopx sys[kwant.HoppingKind((0, 1), lat)] = hopy return sys def graphene_infinite(): lat = kwant.lattice.honeycomb() a, b = lat.sublattices sys = kwant.Builder(kwant.TranslationalSymmetry(*lat.prim_vecs)) sys[lat.shape(lambda pos: True, (0, 0))] = 0 sys[kwant.builder.HoppingKind((0, 0), a, b)] = lambda site1, site2, p: p.t_1 sys[kwant.builder.HoppingKind((0, 1), a, b)] = lambda site1, site2, p: p.t_23 sys[kwant.builder.HoppingKind((-1, 1), a, b)] = lambda site1, site2, p: p.t_23 return sys def plot_dets(sys, p, ks, chiral=False): B = np.array(sys.symmetry.periods).T A = B.dot(np.linalg.inv(B.T.dot(B))) def momentum_to_lattice(k): k, residuals = np.linalg.lstsq(A, k)[:2] return list(k) sys = wraparound.wraparound(sys).finalized() kys, kxs = np.meshgrid(ks, ks) dets = np.zeros_like(kxs, dtype=complex) for i, kx in enumerate(ks): for j, ky in enumerate(ks): kx, ky = momentum_to_lattice([kx, ky]) ham = sys.hamiltonian_submatrix(args=(p, kx, ky)) if chiral: # Bring the chiral symmetric Hamiltonian in offdiagonal form U = (pauli.s0 + 1j * pauli.sx) / np.sqrt(2) ham = U.dot(ham).dot(U.T.conjugate()) dets[i, j] = ham[1, 0] H = np.angle(dets) / (2 * np.pi) V = np.abs(dets) H = np.mod(H, 1) V /= np.max(V) V = 1 - V**2 S = np.ones_like(H) HSV = np.dstack((H, S, V)) RGB = hsv_to_rgb(HSV) bounds = (ks.min(), ks.min(), ks.max(), ks.max()) pl = holoviews.RGB(RGB, bounds=bounds, label=r'$\det(h)$', kdims=['$k_x$', '$k_y$']) return pl(plot={'xticks': pi_ticks, 'yticks': pi_ticks})(style={'interpolation': None}) def Weyl_slab(w=5): def hopx(site1, site2, p): return 0.5j * p.t * pauli.sx - p.t * pauli.sz def hopy(site1, site2, p): return - p.t * pauli.sz def hopz(site1, site2, p): return 0.5j * p.t * pauli.sy - p.t * pauli.sz def onsite(site, p): return 6 * p.t * pauli.sz - p.mu * pauli.sz lat = kwant.lattice.general(np.eye(3)) sys = kwant.Builder(kwant.TranslationalSymmetry([1, 0, 0], [0, 1, 0])) def shape(pos): (x, y, z) = pos return (0 <= z < w) sys[lat.shape(shape, (0, 0, 0))] = onsite sys[kwant.HoppingKind((1, 0, 0), lat)] = hopx sys[kwant.HoppingKind((0, 1, 0), lat)] = hopy sys[kwant.HoppingKind((0, 0, 1), lat)] = hopz return sys
Populated the namespace with: np, matplotlib, kwant, holoviews, init_notebook, interact, display_html, plt, pf, SimpleNamespace, pprint_matrix, scientific_number, pretty_fmt_complex from code/edx_components: MoocVideo, PreprintReference, MoocDiscussion, MoocCheckboxesAssessment, MoocMultipleChoiceAssessment, MoocPeerAssessment, MoocSelfAssessment from code/functions: spectrum, hamiltonian_array, h_k, pauli

Introduction

Ashvin Vishwanath from the University of California, Berkeley will introduce Weyl semimetals and other examples of gapless, yet topological, systems.

MoocVideo("MAWwa4r1qIc", src_location='10.1-intro')

Topological invariants of Fermi surfaces

The idea that leads us to the topology of gapless systems is extremely simple. It is:

If we consider momentum as an external conserved parameter, we can study topological closings of the gap in momentum space.

Let's consider the simplest type of topological invariant, one we learned about at the very beginning of this course. Remember the simplest topological invariant of a 0D Hamiltonian, the number of filled energy levels? What if we take two points in momentum space, k1\mathbf{k}_1 and k2\mathbf{k}_2, and consider a Hamiltonian such that the number of filled states changes by nn between these two points? We can conclude that there are at least nn Fermi surfaces that lie on every path between k1\mathbf{k}_1 and k2\mathbf{k}_2 in momentum space.

Now we just need to take this idea and apply it to more interesting systems and topological invariants!

What types of topological invariants are relevant? Aside from special circumstances, we cannot make use of time-reversal or particle-hole symmetries: in momentum space these only have an immediate effect in isolated k\mathbf{k}-points, where every momentum component is either 00 or π\pi. So there are no paths in momentum space for which either of the symmetries is effective in each point.

So we are left with only two symmetry classes: A and AIII (no symmetry at all or sublattice/chiral symmetry), and with only two invariants: if there is a sublattice symmetry, a winding number can be defined, and without it there's a Chern number.

Graphene and protected Dirac cones

We've already analysed the 0D Chern number that stabilizes the usual Fermi surfaces. Let's go one dimension higher, and study winding numbers in systems with sublattice symmetry around 1D loops.

For a winding number to be nonzero, we need to consider 1D loops in momentum space. As a reminder, with sublattice symmetry the Hamiltonian can always be brought to the form

H=(0h(k)h(k)0)H = \begin{pmatrix} 0 & h(\mathbf{k}) \\ h^\dagger(\mathbf{k}) & 0 \end{pmatrix}

The topological invariant is a nonzero winding of deth(k)\det h(\mathbf{k}) when k\mathbf{k} goes around some contour. Since h(k)h(\mathbf{k}) is continuous, this means that its determinant also has to vanish somewhere inside this contour.

To study a particular example where this appears, let's return to graphene, which we studied as a simple limit of Haldane model. For graphene we have the Hamiltonian

h(kx,ky)=t1eika1+t2eika2+t3eika3,h(k_x, k_y) = t_1 e^{i \mathbf{k} \cdot \mathbf{a_1}} + t_2 e^{i \mathbf{k} \cdot \mathbf{a_2}} + t_3 e^{i \mathbf{k} \cdot \mathbf{a_3}},

where t1,t2,t3t_1, t_2, t_3 are the three hoppings connecting a site in one of the two graphene sublattices, and a1,a2,a3a_1, a_2, a_3 are the lattice vectors connecting one unit cell to its neighbors.

To consider something specific, let's take t2=t3=tt_2 = t_3 = t and vary t1t_1. This is what the band structure and deth\det h look like:

%%output size=150 p = SimpleNamespace(t_1=1.0, t_23=1.0) sys = graphene_infinite() ks = np.sqrt(3) * np.linspace(-np.pi, np.pi, 80) kwargs = dict(title=lambda p: r'Graphene, $t_1 = {:.2} \times t$'.format(p.t_1), zticks=3) ts = np.linspace(1, 2.4, 8) (holoviews.HoloMap({p.t_1: spectrum(sys, p, k_x=ks, k_y=ks, **kwargs) for p.t_1 in ts}, kdims=['$t_1$']) + holoviews.HoloMap({p.t_1: plot_dets(sys, p, ks) for p.t_1 in ts}, kdims=['$t_1$']))

The left panel shows the band structure, and you see that it has gapless points. The right panel shows deth\det h by using hue as a phase and intensity as magnitude (so white is deth=0\det h = 0). There are two Dirac points (you see 6, but this is because we plot more than one Brillouin zone).

We also see that the winding numbers around these two Dirac points have opposite signs (because by going around them clockwise you encounter red, blue and green colors in opposite orders). This must always be the case since the winding number around the edges of the complete Brillouin zone must vanish - as you walk down every edge of the Brillouin zone twice in opposite directions, their contributions always cancel.

As t1t_1 increases, the two poles move towards each other, eventually annihilating and leaving a completely gapped dispersion relation. Let's now try to obtain an effective model for the dispersion at each pole and at the phase transition point.

We know that deth\det h has to vanish next to some point k0\mathbf{k}_0. We can expand the Hamiltonian to a linear order next to this point, which immediately leaves us with a Hamiltonian

H(k)=(0eiα(vxδkx+ivyδky)eiα(vxδkxivyδky)0),H(\mathbf{k}) = \begin{pmatrix} 0 & e^{i\alpha} (v_x \delta k_x + i v_y \delta k_y) \\ e^{-i\alpha} (v_x \delta k_x - i v_y \delta k_y) & 0 \end{pmatrix},

where δk\mathbf{\delta k} is of course the difference between k\mathbf{k} and the Dirac point momentum. Of course this is the 2D Dirac equation, which should be very familiar now.

At the phase transition where the two Dirac points annihilate, we can also quickly guess that the correct dispersion should be a quadratic function along the axis connecting the two Dirac points, and linear along the other axis (this is also what we see in the plot). We thus have

H(k)=(0eiα(βδk12+m+iv2δk2)eiα(βδk12+miv2δk2)0),H(\mathbf{k}) = \begin{pmatrix} 0 & e^{i\alpha} (\beta \delta k_1^2 + m + i v_2 \delta k_2) \\ e^{-i\alpha} (\beta \delta k_1^2 + m - i v_2 \delta k_2) & 0 \end{pmatrix},

such that for m>0m>0 we have a fully gapped Hamiltonian, and for m<0m<0 there are two Dirac points.

dd-wave superconductors and edge states

Gapless points with Dirac dispersion were known for quite some time before graphene. They exist in the cuprate family of high temperature superconductors, known to have a dd-wave order parameter. These materials are layered, with weak couplings between the layers, so in the study of these complicated systems, one often starts with a simplified two-dimensional Hamiltonian.

This Hamiltonian just has the usual kinetic energy term of a single particle band and a superconducting pairing proportional to kx2ky2k_x^2 - k_y^2, namely

H=(k2/2mμΔ(kx2ky2)Δ(kx2ky2)μk2/2m).H = \begin{pmatrix} k^2/2m -\mu & \Delta (k_x^2 - k_y^2) \\ \Delta (k_x^2 - k_y^2) & \mu-k^2/2m \end{pmatrix}.

There is no spin-orbit coupling here, so the Hamiltonian has a spinless time-reversal symmetry H=HH = H^*. It also has a particle-hole symmetry H=τyHτyH= - \tau_y H^* \tau_y. Their product, the chiral symmetry H=τyHτyH = -\tau_y H \tau_y allows the Hamiltonian to have gapless points where both the single-particle dispersion and the pairing vanish.

Difference between sublattice symmetries

Time-reversal symmetry ensures that the winding points come in pairs at opposite momenta, just like in graphene. In graphene however, the chiral symmetry operator σz\sigma_z commutes with the time-reversal symmetry operator. This means that applying time-reversal symmetry changes the direction of a loop in momentum space, but leaves the winding number invariant. In dd-wave superconductors on the other hand, the chiral symmetry operator τy\tau_y is odd under time-reversal (i.e. the operators anticommute), and the winding is invariant under it.

This means that a Dirac point at momentum kk and positive winding must come together with a Dirac point at k-k and also positive winding. Since the total winding over the Brillouin zone must be 0, this means that in superconducting systems the Dirac points come in quadruplets: two with positive winding and two with negative winding.

The dd-wave superconductor Hamiltonian gives just that: there are 4 Dirac points at kx=ky=kF/2|k_x| = |k_y| = k_F / \sqrt{2}.

question = r"What happens if you make the 2D $d$-wave Hamiltonian 3D, by adding coupling between 2D layers?" answers = ["The Dirac points couple and gap out.", "In 3D you cannot have a $d$-wave pairing.", "There will remain isolated gapless points in the larger 3D Brillouin zone.", "You get a closed 1D Dirac line of gap closings in the 3D Brillouin zone."] explanation = (r"The real and imaginary parts of the solutions of $\det h(\mathbf{k})=0$ form two surfaces " r"in the Brillouin zone. The intersection of these two surfaces is a line.") MoocMultipleChoiceAssessment(question=question, answers=answers, correct_answer=3, explanation=explanation)

Edge states

Now let's see how bulk-edge correspondence can be made to work for gapless systems. The idea here is to consider the projection of the wave vector parallel to a continuous sample boundary ParseError: KaTeX parse error: Got function '\mathbin' with no arguments as subscript at position 3: k_\̲m̲a̲t̲h̲b̲i̲n̲{||} as a parameter, and to apply the bulk-edge correspondence to the remaining lower-dimensional Hamiltonian.

Whenever the line corresponding to a constant ParseError: KaTeX parse error: Got function '\mathbin' with no arguments as subscript at position 3: k_\̲m̲a̲t̲h̲b̲i̲n̲{||} crosses a Dirac point, the winding number of the Hamiltonian ParseError: KaTeX parse error: Got function '\mathbin' with no arguments as subscript at position 5: H(k_\̲m̲a̲t̲h̲b̲i̲n̲{||}) changes by the winding of the Dirac point. This means that for certain values of momentum parallel to the boundary, a zero energy edge state will appear.

For a dd-wave superconductor this will only happen for some crystalline orientations, as you can see for yourself:

%%opts VLine (color='k') Curve (ls='--') p = SimpleNamespace(mu=2.0, t=1.0, delta=1.0) k = np.arccos(1-p.mu/p.t/4) ks = np.linspace(-np.pi, np.pi, 50) sys0 = d_wave() sys1 = d_wave(50, direction='topo') sys2 = d_wave(50, direction='triv') det_plot = plot_dets(sys0, p, ks, chiral=True) det_plot1= (det_plot * holoviews.Curve(([-np.pi, np.pi], [np.pi, -np.pi])) * holoviews.Curve(([-np.pi, np.pi-2*k], [np.pi-2*k, -np.pi])) * holoviews.Curve(([-np.pi+2*k, np.pi], [np.pi, -np.pi+2*k]))).relabel('$\det(h)$') det_plot2 = det_plot * holoviews.VLine(k) * holoviews.VLine(-k) kwargs = dict(k_x=ks, ylims=[-2, 2], xticks=pi_ticks, yticks=3) (spectrum(sys1, p, title='Ribbon with edge states', **kwargs) * holoviews.VLine(-2*k) * holoviews.VLine(2*k) + det_plot1 + spectrum(sys2, p, title='Ribbon without edge states', **kwargs) * holoviews.VLine(-k) * holoviews.VLine(k) + det_plot2).cols(2)

On the right panels you once again see deth\det h, with added lines denoting different the values of kk_{\mathbin{||}} crossing the Dirac points. If the sample boundary is along the (1, 0) axis, the Dirac points have coinciding kk_{\mathbin{||}}, and their windings cancel, so that no single value of kk_{\mathbin{||}} carries an edge state.

On the other hand, the crystal boundary (1, 1), which lies at an angle π/4\pi/4 with respect to the crystallographic axes, has a total winding of +2 at k=0k_{\mathbin{||}}=0 and a winding of −1 for k=±kFk_{\mathbin{||}}=\pm k_F. In this case, each k<kF|k_{\mathbin{||}}|<k_F carries a single edge state.

These edge states were known to exist long before the discovery of topological insulators, and it is fascinating to see how they perfectly fit to the theory of topological phenomena!

Weyl points

So far we've seen two examples of Dirac points in two dimensions, the surface of a 3D topological insulator and graphene. You might wonder, why don't we have such cones in three dimensions? These do indeed exist and are called Weyl points instead of Dirac points. The reason is historical - Dirac's equation for the electron (which is in 3D) involves states with four components, two for the electron and two for the hole. The direct generalization of graphene to 3D3D that we will discuss involves states with two electron component. Such electron states with linear dispersion were studied first by Weyl, and have strange properties as we will illustrate below.

Let us start by writing the low-energy Hamiltonian for the three dimensional generalization of graphene:

H(k)=(σxkx+σyky+σzkz).H({\bf k})=(\sigma_x k_x+\sigma_y k_y+\sigma_z k_z).

Here you might think of σx,y,z\sigma_{x,y,z} as the spin of the electron (just as on the surface of a topological insulator).

Next we try the usual thing we would do with a two-dimensional Dirac cone - namely see what happens when we gap it out by applying a magnetic field σB\bf\sigma\cdot B. Adding such a term, we find that the Hamiltonian transforms as follows:

H(k)H(k)+σB=σ(k+B).H({\bf k})\rightarrow H({\bf k})+{\bf\sigma\cdot B}={\bf\sigma\cdot (k+B)}.

The key observation here is that the addition of a magnetic field effectively shifts the wave-vector as

kk~=k+B.{\bf k}\rightarrow \tilde{\bf k}={\bf k+ B}.

So applying the most general perturbation we can think of does not gap out the Weyl point where the energy vanishes. Instead, the perturbation only shifts the Weyl point around in momentum space. This feels like some kind of topological protection.

%%opts Surface [azimuth=45] sys = Weyl_slab(w=10) p = SimpleNamespace(t=1.0, mu=None) mus = np.linspace(-0.4, 2, 13) kwargs = dict(k_x=np.linspace(-np.pi, 0), k_y=np.linspace(-np.pi, np.pi), title=lambda p: 'Weyl semimetal, $\mu = {:.2}$'.format(p.mu), num_bands=4) holoviews.HoloMap({p.mu: spectrum(sys, p, **kwargs) for p.mu in mus}, kdims=[r'$\mu$'])

Is there a sense in which Weyl points are "topological"? They are clearly protected, but is there some topological reason for the protection? As in the rest of this section, the topology of gapless system becomes apparent by looking at the Hamiltonian in lower dimensional subspaces of momentum space. For the case of Weyl, the momentum space is three dimensional, so let us look at two dimensional subspaces of momentum space.

A natural subspace to choose is to fix kz=mk_z=m. The Weyl Hamiltonian then becomes that of a massive 2D Dirac cone

H2D,Dirac(kx,ky;m)H(kx,ky,kz=m)=(σxkx+σyky+mσz).H_{2D,Dirac}(k_x,k_y;m)\equiv H(k_x,k_y,k_z=m)=(\sigma_x k_x+\sigma_y k_y+m\sigma_z).

As we talked about in week 4 with Chern insulators, the massive Dirac model has a Chern number, which changes by 11 if mm changes sign.

So we can think of the Weyl Hamiltonian in the momentum planes at fixed kzk_z as Chern insulators with Chern numbers nCh=0n_{Ch}=0 (i.e. trivial) if kz<0k_z < 0 and nCh=1n_{Ch}=1 (topological) if kz>0k_z > 0. The Hamiltonian at kz=0k_z=0 is at the phase transition point of the Chern insulator, which supports a gapless Dirac point.

Systems with Weyl points are known as Weyl semimetals. Just like other topological phases, Weyl semimetals have an interesting surface spectrum. We can understand this easily by viewing the Weyl point as a stack of Chern insulators in momentum space. For any surface in a plane that contains the zz-axis, we can treat kzk_z as a conserved quantity. At this kz=mk_z=m, the Hamiltonian is just that of a Chern insulator with an appropriate Chern number. For the range of kzk_z where the Chern number nCh(kz)=1n_{Ch}(k_z)=1, the surface spectrum supports chiral edge states with an energy approximated at low energy by

E(kx,kz)v(kz)kx.E(k_x,k_z)\approx v(k_z)k_x.

We can consider the edge states over a range of kzk_z together to visualize the "surface states".

The unique property of the surface states is that if we set kx=0k_x=0 then the energy vanishes on a line in the surface spectrum. This line actually terminates at kz=0k_z=0, where the Chern number changes. Such lines, which are referred to as "Fermi arcs," are the unique bounday properties (hence the bulk-boundary correspondence) of Weyl semimetals.

At large enough kzk_z, the two dimensional Hamiltonian H2D,Dirac(kx,ky;kz)H_{2D,Dirac}(k_x,k_y;k_z) becomes trivial i.e. nCh(kz)=0n_{Ch}(|k_z|\rightarrow \infty)=0. This means that if the Chern number is nCh=1n_{Ch}=1 in a range of kzk_z, then nCh(kz)n_{Ch}(k_z) must change twice resulting in two Weyl points. So Weyl points come in pairs. These points map onto the ends of the Fermi arcs on the surface.

question = r"What protects the surface state of Weyl semi-metals from scattering inside the bulk Weyl point?" answers = ["Chiral symmetry.", "The energy gap in the bulk.", "Absence of scattering.", "The non-zero Chern number of the bulk."] explanation = (r"The bulk has gapless states due to the Weyl point. " "Therefore, only momentum conservation protects surface states from going into the bulk.") MoocMultipleChoiceAssessment(question=question, answers=answers, correct_answer=2, explanation=explanation)

Questions about what you just learned? Ask them below!

MoocDiscussion("Questions", "Topology in gapless systems")