I'm not 100% sure if I should be posting this in the physics stack exchange but, I am pretty new to Mathematica, just installed it five weeks ago, and I am having a little trouble implementing the method in M P Lopez Sancho et al 1984 J. Phys. F: Met. Phys. 14 1205 to calculate the surface spectral function of my model. I need to calculate the surface (zeroth-layer) Greens function given by:

$$G_{00}(\omega) = (\omega - H_{00} - H_{01} T(\omega))^{-1}\,.$$

This quickly gives us the surface spectral function via:

$$A_{0}^0(\boldsymbol{k_{\parallel}}\epsilon) = -\frac{1}{\pi}Im\sum_{\alpha} G_{00}^{0 \alpha, 0 \alpha}(\boldsymbol{k_{\parallel}}\epsilon + i \eta)\,.$$

The transfer matrix $T(\omega)$ is built from the following: $$T(\omega) = t_0 + \tilde{t}_0 t_1 + \ldots + \tilde{t}_0 \tilde{t}_1 \ldots \tilde{t}_{n-1} t_{n}$$ where $t_i$ and $\tilde{t}_i$ are given by: $$ t_{i} = (I_2 - t_{i-1}\tilde{t}_{i-1} - \tilde{t}_{i-1}t_{i-1})^{-1} t_{i-1}^2$$ $$ \tilde{t}_{i} = (I_2 - t_{i-1}\tilde{t}_{i-1} - \tilde{t}_{i-1}t_{i-1})^{-1} \tilde{t}_{i-1}^2$$ $$ t_0 = (\omega - H_{00})^{-1}H_{01}^{\dagger}$$ $$ \tilde{t}_0 = (\omega - H_{00})^{-1}H_{01}$$

Where $H_{00}$ and $H_{01}$ are $2\times2$ matrices that describe the parts of the Hamiltonian that correspond to no hopping in the $kz$ direction and nearest-neighbor hopping in the $kz$ direction, respectively, so far, I have written this code up. Still, there is something flawed with it, and I am unable to pinpoint what precisely is wrong with it. I'd appreciate it greatly if anyone could help me out here.

My Attempt

First I put in my tight-binding model that I have already calculated:

q = {kx, ky, kz};
a1 = {1, 0, 0};
a2 = {-(1/2), Sqrt[3]/2, 0};
a3 = {-(1/2), -Sqrt[3]/2, 0};
(* This is the in plane momentum Subscript[k, x],Subscript[k, y] *)
k = {q.a1, q.a2, q.a3};
(* Tight binding model *)
d2[kx_, ky_, A1_, A2_, A3_, A4_] :=
   A1 Sin[kz] + A2 Sin[2 kz] + A3 Sum[Cos[k[[j]]], {j, 1, 3}] Sin[kz]
    + A4 (Sin[k[[1]] - k[[2]]] + Sin[k[[2]] - k[[3]]] + Sin[k[[3]] - k[[1]]]);
d3[kx_, ky_, B1_, B2_, B3_] := 
  B1 + B2 (3 - Sum[Cos[k[[j]]], {j, 1, 3}]) + B3 (1 - Cos[kz]);
f = Sum[Cos[k[[j]]], {j, 1, 3}];
g = Sin[k[[1]] - k[[2]]] + Sin[k[[2]] - k[[3]]] + Sin[k[[3]] - k[[1]]];

d[kx_, ky_, kz_, A1_, A2_, A3_, A4_, B1_, B2_, B3_] :=
  {0, d2[kx, ky, kz, A1, A2, A3, A4], d3[kx, ky, kz, B1, B2, B3]};
Energy[kx_, ky_, kz_] := Norm[d[kx, ky, kz, 1, 1, 1, 0.5, -0.8, 1, 1]];    

I then define the matrices $H_{00}$ and $H_{01}$:

Subscript[h, 00] = {{B3/2, -A1/2 - (A3 f)/2}, {A1/2 + (A3 f)/2, -B3/2}};
Subscript[h, 01] = {{B1 + B2 (3 - f), -I A4 g}, {I A4 g, -B1 - B2 (3 - f)}};

Then make I define the path in the Brillouin zone to travel through, in this case $\Gamma \rightarrow K \rightarrow M$ and omega values to sum over later:

kpath = {{0, {0, 0, 0, ω}}, {1, {4 π/3, 0, 0, ω}}, {2, {π, -π/Sqrt[3], 0, ω}}};
if1 = Interpolation[kpath, InterpolationOrder -> 1];
(*Small number for use with ω*)
δ = 0.01;
omega = Range[-6, 6, 0.1] + I δ;

Then I define the surface Greens function and the matrices $t_0$ and $\tilde{t}_0$:

(* Definition of Subscript[G, 00] *)
surfaceGreens[kx_, ky_, kz_, \[Omega]_, A1_, A2_, A3_, A4_, B1_, B2_, 
   B3_] = Inverse[\[Omega] *IdentityMatrix[2] - Subscript[h, 00] - 
    Subscript[h, 01]. transferMatrix];
sg = {};

Subscript[t, 0][\[Omega]_, A1_, A3_, A4_, B1_, B2_, B3_] := 
  Inverse[\[Omega] I - Subscript[h, 00]].ConjugateTranspose[Subscript[
    h, 01]];
\!\(\*OverscriptBox[\(t\), \(~\)]\), 0] [\[Omega]_, A1_, A3_, A4_, 
   B1_, B2_, B3_] := 
  Inverse[\[Omega] I - Subscript[h, 00]].Subscript[h, 01];

Here we loop over the path in momentum space and omega to build up the transfer matrix $T(\omega)$ and the surface greens function, this seems to be where my issue lies. I don't seem to be appending to my lists $t$ and $\tilde{t}$ in my while loop below:

       \[Omega] = omega[[j]];
       t = {};    
       \!\(\*OverscriptBox[\(t\), \(~\)]\) = {};
       (* This is the first element of the t and Overscript[t, ~] lists that will be used to calculate T(\[Omega]). *)
       Subscript[t, 0][\[Omega], 1, 1, 0.5, -0.8, 1, 1];
       Subscript[\!\(\*OverscriptBox[\(t\), \(~\)]\), 0][\[Omega], 1, 1, 0.5, -0.8, 1, 1]; 
       t = {Subscript[t, 0]};
       \!\(\*OverscriptBox[\(t\), \(~\)]\) = {Subscript[\!\(\*OverscriptBox[\(t\), \(~\)]\), 0]};

(* Build up the list of matrices t and Overscript[t, ~] to calculate T, this incredibly slow... *)
       n = 2;                                               
       While[n < 5 , AppendTo[t, Inverse[IdentityMatrix[2] - t[[n - 1]].
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[n - 1]] - 
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[n - 1]].t[[n - 1]]].t[[
      n - 1]]^2] ; AppendTo[
\!\(\*OverscriptBox[\(t\), \(~\)]\), 
   Inverse[IdentityMatrix[2] - t[[n - 1]].
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[n - 1]] - 
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[n - 1]].t[[n - 1]]].
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[n - 1]]^2]; n++];

     (* Calculation of the transfer matrix T(\[Omega]), placeholder for right now. 
Remember to put in more efficient expression to calculate T(\[Omega]) once above errors are fixed.*)
     transferMatrix = t[[1]] + 
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[1]].t[[2]] + 
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[1]].
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[2]].t[[3]] + 
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[1]].
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[2]].
\!\(\*OverscriptBox[\(t\), \(~\)]\)[[3]].t[[4]];

   (*Calculate the surface Green's function*)
   AppendTo[sg, surfaceGreens[kx, ky, kz, ω, 1, 1, 1, 0.5, -0.8, 1, 1]];,
  {j, 1, Length[omega]}],
 {i, 1, Length[kpath]}]

Finally, we plot the surface spectral function using DensityPlot:

greens = Total[sg];
spectralFunction[kx_, ky_, kz_, ω_] := -1/π Im[Tr[greens]];
DensityPlot[{spectralFunction @@ if1[i]}, {i, 0, 2}]
  • $\begingroup$ Welcome to Mathematica.SE! I suggest the following: 1) As you receive help, try to give it too, by answering questions in your area of expertise. 2) Take the tour and check the faqs! 3) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, [by clicking the checkmark sign](tinyurl.com/4srwe26 $\endgroup$
    – Dunlop
    Feb 28, 2020 at 4:36
  • 1
    $\begingroup$ Perhaps you can highlight what exactly doesn't work with your code. This may assist others in helping you out. $\endgroup$
    – Dunlop
    Feb 28, 2020 at 4:36
  • 3
    $\begingroup$ As @Dunlop mentioned, you should edit your post to include what isn't working, how you know it's not working, and what the expected output is. I'll point out some things you might want to double-check. You have && in the While part of your Do loop. This should probably be a semi-colon as && is a logical operator for boolean values. In your Do loop, you define the Greek letter omega as a function, but then use the omega symbol without calling the actual function, so omega has no value when you call t0 and l0. A4 seems to go undefined in t0 and l0. $\endgroup$
    – MassDefect
    Feb 28, 2020 at 4:44
  • 3
    $\begingroup$ You define Subscript[h, 00] and Subscript[h, 01] but don't actually use them anywhere that I can see. You do use the variable h in surfaceGreens but it doesn't have any definition. After running your code, try evaluating t, l, sg, and especially spectralFunction@@if1[0]. You'll find that the outputs are full of undefined variables like t0, l0, Greek omega, and h. If you're looking to plot, then these need to take on numerical values at some point. If you're using DensityPlot there can be two variables (x and y) that only take on values during the plotting, but no more. $\endgroup$
    – MassDefect
    Feb 28, 2020 at 4:50
  • $\begingroup$ Thanks for the help @MassDefect and @Dunlop, I really appreciate it. I've made some edits to my code above (replacing && with a semicolon and fixing the definitions of my functions to take in all the inputs needed). The issue seems to be with my While part of the Do loop, I don't seem to be appending to my lists $t$ and $\tilde{t}$ like I need to be doing to calculate the transfer matrix so when I try to calculate the third entry I see errors that say parts of the list do not exist. The error specifically is Part::partw: Part 2 of {Subscript[{}, 0]} does not exist. $\endgroup$ Feb 28, 2020 at 19:45


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Browse other questions tagged or ask your own question.