Coverage for /builds/ase/ase/ase/transport/selfenergy.py: 82.14%

1# fmt: off 

2 

3import numpy as np 

4 

5 

6class LeadSelfEnergy: 

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

8 

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

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

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

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

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

14 self.eta = eta 

15 self.energy = None 

16 self.bias = 0 

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

18 

19 def retarded(self, energy): 

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

21 if energy != self.energy: 

22 self.energy = energy 

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

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

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

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

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

28 

29 return self.sigma_mm 

30 

31 def set_bias(self, bias): 

32 self.bias = bias 

33 

34 def get_lambda(self, energy): 

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

36 

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

38 conjugate. 

39 """ 

40 sigma_mm = self.retarded(energy) 

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

42 

43 def get_sgfinv(self, energy): 

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

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

46 

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

48 v_11 = v_00.copy() 

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

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

51 

52 delta = self.conv + 1 

53 while delta > self.conv: 

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

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

56 v_01_dot_b = np.dot(v_01, b) 

57 v_00 -= v_01_dot_b 

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

59 v_11 -= v_01_dot_b 

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

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

62 delta = abs(v_01).max() 

63 

64 return v_00 

65 

66 

67class BoxProbe: 

68 """Box shaped Buttinger probe. 

69 

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

71 """ 

72 

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

74 from Transport.Hilbert import hilbert 

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

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

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

78 se.real = hilbert(se.imag) 

79 se.imag -= 1 

80 self.selfenergy_e = eta * se 

81 self.energies = energies 

82 self.S = S 

83 

84 def retarded(self, energy): 

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