1import numpy as np 

2 

3 

4class LeadSelfEnergy: 

5 conv = 1e-8 # Convergence criteria for surface Green function 

6 

7 def __init__(self, hs_dii, hs_dij, hs_dim, eta=1e-4): 

8 self.h_ii, self.s_ii = hs_dii # onsite principal layer 

9 self.h_ij, self.s_ij = hs_dij # coupling between principal layers 

10 self.h_im, self.s_im = hs_dim # coupling to the central region 

11 self.nbf = self.h_im.shape[1] # nbf for the scattering region 

12 self.eta = eta 

13 self.energy = None 

14 self.bias = 0 

15 self.sigma_mm = np.empty((self.nbf, self.nbf), complex) 

16 

17 def retarded(self, energy): 

18 """Return self-energy (sigma) evaluated at specified energy.""" 

19 if energy != self.energy: 

20 self.energy = energy 

21 z = energy - self.bias + self.eta * 1.j 

22 tau_im = z * self.s_im - self.h_im 

23 a_im = np.linalg.solve(self.get_sgfinv(energy), tau_im) 

24 tau_mi = z * self.s_im.T.conj() - self.h_im.T.conj() 

25 self.sigma_mm[:] = np.dot(tau_mi, a_im) 

26 

27 return self.sigma_mm 

28 

29 def set_bias(self, bias): 

30 self.bias = bias 

31 

32 def get_lambda(self, energy): 

33 """Return the lambda (aka Gamma) defined by i(S-S^d). 

34 

35 Here S is the retarded selfenergy, and d denotes the hermitian 

36 conjugate. 

37 """ 

38 sigma_mm = self.retarded(energy) 

39 return 1.j * (sigma_mm - sigma_mm.T.conj()) 

40 

41 def get_sgfinv(self, energy): 

42 """The inverse of the retarded surface Green function""" 

43 z = energy - self.bias + self.eta * 1.j 

44 

45 v_00 = z * self.s_ii.T.conj() - self.h_ii.T.conj() 

46 v_11 = v_00.copy() 

47 v_10 = z * self.s_ij - self.h_ij 

48 v_01 = z * self.s_ij.T.conj() - self.h_ij.T.conj() 

49 

50 delta = self.conv + 1 

51 while delta > self.conv: 

52 a = np.linalg.solve(v_11, v_01) 

53 b = np.linalg.solve(v_11, v_10) 

54 v_01_dot_b = np.dot(v_01, b) 

55 v_00 -= v_01_dot_b 

56 v_11 -= np.dot(v_10, a) 

57 v_11 -= v_01_dot_b 

58 v_01 = -np.dot(v_01, a) 

59 v_10 = -np.dot(v_10, b) 

60 delta = abs(v_01).max() 

61 

62 return v_00 

63 

64 

65class BoxProbe: 

66 """Box shaped Buttinger probe. 

67  

68 Kramers-kroning: real = H(imag); imag = -H(real) 

69 """ 

70 def __init__(self, eta, a, b, energies, S, T=0.3): 

71 from Transport.Hilbert import hilbert 

72 se = np.empty(len(energies), complex) 

73 se.imag = .5 * (np.tanh(.5 * (energies - a) / T) - 

74 np.tanh(.5 * (energies - b) / T)) 

75 se.real = hilbert(se.imag) 

76 se.imag -= 1 

77 self.selfenergy_e = eta * se 

78 self.energies = energies 

79 self.S = S 

80 

81 def retarded(self, energy): 

82 return self.selfenergy_e[self.energies.searchsorted(energy)] * self.S