Coverage for /builds/kinetik161/ase/ase/transport/greenfunction.py: 82.50%

40 statements  

« prev     ^ index     » next       coverage.py v7.2.7, created at 2023-12-10 11:04 +0000

1import numpy as np 

2 

3 

4class GreenFunction: 

5 """Equilibrium retarded Green function.""" 

6 

7 def __init__(self, H, S=None, selfenergies=[], eta=1e-4): 

8 self.H = H 

9 self.S = S 

10 self.selfenergies = selfenergies 

11 self.eta = eta 

12 self.energy = None 

13 self.Ginv = np.empty(H.shape, complex) 

14 

15 def retarded(self, energy, inverse=False): 

16 """Get retarded Green function at specified energy. 

17 

18 If 'inverse' is True, the inverse Green function is returned (faster). 

19 """ 

20 if energy != self.energy: 

21 self.energy = energy 

22 z = energy + self.eta * 1.j 

23 

24 if self.S is None: 

25 self.Ginv[:] = 0.0 

26 self.Ginv.flat[:: len(self.S) + 1] = z 

27 else: 

28 self.Ginv[:] = z 

29 self.Ginv *= self.S 

30 self.Ginv -= self.H 

31 

32 for selfenergy in self.selfenergies: 

33 self.Ginv -= selfenergy.retarded(energy) 

34 

35 if inverse: 

36 return self.Ginv 

37 else: 

38 return np.linalg.inv(self.Ginv) 

39 

40 def calculate(self, energy, sigma): 

41 """XXX is this really needed""" 

42 ginv = energy * self.S - self.H - sigma 

43 return np.linalg.inv(ginv) 

44 

45 def apply_retarded(self, energy, X): 

46 """Apply retarded Green function to X. 

47 

48 Returns the matrix product G^r(e) . X 

49 """ 

50 return np.linalg.solve(self.retarded(energy, inverse=True), X) 

51 

52 def dos(self, energy): 

53 """Total density of states -1/pi Im(Tr(GS))""" 

54 if self.S is None: 

55 return -self.retarded(energy).imag.trace() / np.pi 

56 else: 

57 GS = self.apply_retarded(energy, self.S) 

58 return -GS.imag.trace() / np.pi 

59 

60 def pdos(self, energy): 

61 """Projected density of states -1/pi Im(SGS/S)""" 

62 if self.S is None: 

63 return -self.retarded(energy).imag.diagonal() / np.pi 

64 else: 

65 S = self.S 

66 SGS = np.dot(S, self.apply_retarded(energy, S)) 

67 return -(SGS.diagonal() / S.diagonal()).imag / np.pi