Coverage for /builds/ase/ase/ase/vibrations/placzek.py: 96.88%

1# fmt: off 

2 

3import numpy as np 

4 

5import ase.units as u 

6from ase.calculators.excitation_list import polarizability 

7from ase.vibrations.raman import Raman, RamanPhonons 

8from ase.vibrations.resonant_raman import ResonantRaman 

9 

10 

11class Placzek(ResonantRaman): 

12 """Raman spectra within the Placzek approximation.""" 

13 

14 def __init__(self, *args, **kwargs): 

15 self._approx = 'PlaczekAlpha' 

16 ResonantRaman.__init__(self, *args, **kwargs) 

17 

18 def set_approximation(self, value): 

19 raise ValueError('Approximation can not be set.') 

20 

21 def _signed_disps(self, sign): 

22 for a, i in zip(self.myindices, self.myxyz): 

23 yield self._disp(a, i, sign) 

24 

25 def _read_exobjs(self, sign): 

26 return [disp.read_exobj() for disp in self._signed_disps(sign)] 

27 

28 def read_excitations(self): 

29 """Read excitations from files written""" 

30 self.ex0E_p = None # mark as read 

31 self.exm_r = self._read_exobjs(sign=-1) 

32 self.exp_r = self._read_exobjs(sign=1) 

33 

34 def electronic_me_Qcc(self, omega, gamma=0): 

35 self.calculate_energies_and_modes() 

36 

37 V_rcc = np.zeros((self.ndof, 3, 3), dtype=complex) 

38 pre = 1. / (2 * self.delta) 

39 pre *= u.Hartree * u.Bohr # e^2Angstrom^2 / eV -> Angstrom^3 

40 

41 om = omega 

42 if gamma: 

43 om += 1j * gamma 

44 

45 for i, r in enumerate(self.myr): 

46 V_rcc[r] = pre * ( 

47 polarizability(self.exp_r[i], om, 

48 form=self.dipole_form, tensor=True) - 

49 polarizability(self.exm_r[i], om, 

50 form=self.dipole_form, tensor=True)) 

51 self.comm.sum(V_rcc) 

52 

53 return self.map_to_modes(V_rcc) 

54 

55 

56class PlaczekStatic(Raman): 

57 def read_excitations(self): 

58 """Read excitations from files written""" 

59 self.al0_rr = None # mark as read 

60 self.alm_rr = [] 

61 self.alp_rr = [] 

62 for a, i in zip(self.myindices, self.myxyz): 

63 for sign, al_rr in zip([-1, 1], [self.alm_rr, self.alp_rr]): 

64 disp = self._disp(a, i, sign) 

65 al_rr.append(disp.load_static_polarizability()) 

66 

67 def electronic_me_Qcc(self): 

68 self.calculate_energies_and_modes() 

69 

70 V_rcc = np.zeros((self.ndof, 3, 3), dtype=complex) 

71 pre = 1. / (2 * self.delta) 

72 pre *= u.Hartree * u.Bohr # e^2Angstrom^2 / eV -> Angstrom^3 

73 

74 for i, r in enumerate(self.myr): 

75 V_rcc[r] = pre * (self.alp_rr[i] - self.alm_rr[i]) 

76 self.comm.sum(V_rcc) 

77 

78 return self.map_to_modes(V_rcc) 

79 

80 

81class PlaczekStaticPhonons(RamanPhonons, PlaczekStatic): 

82 pass 

83 

84 

85class Profeta(ResonantRaman): 

86 """Profeta type approximations. 

87 

88 Reference 

89 --------- 

90 Mickael Profeta and Francesco Mauri 

91 Phys. Rev. B 63 (2000) 245415 

92 """ 

93 

94 def __init__(self, *args, **kwargs): 

95 self.set_approximation(kwargs.pop('approximation', 'Profeta')) 

96 self.nonresonant = kwargs.pop('nonresonant', True) 

97 ResonantRaman.__init__(self, *args, **kwargs) 

98 

99 def set_approximation(self, value): 

100 approx = value.lower() 

101 if approx in ['profeta', 'placzek', 'p-p']: 

102 self._approx = value 

103 else: 

104 raise ValueError('Please use "Profeta", "Placzek" or "P-P".') 

105 

106 def electronic_me_profeta_rcc(self, omega, gamma=0.1, 

107 energy_derivative=False): 

108 """Raman spectra in Profeta and Mauri approximation 

109 

110 Returns 

111 ------- 

112 Electronic matrix element, unit Angstrom^2 

113 """ 

114 self.calculate_energies_and_modes() 

115 

116 V_rcc = np.zeros((self.ndof, 3, 3), dtype=complex) 

117 pre = 1. / (2 * self.delta) 

118 pre *= u.Hartree * u.Bohr # e^2Angstrom^2 / eV -> Angstrom^3 

119 

120 def kappa_cc(me_pc, e_p, omega, gamma, form='v'): 

121 """Kappa tensor after Profeta and Mauri 

122 PRB 63 (2001) 245415""" 

123 k_cc = np.zeros((3, 3), dtype=complex) 

124 for p, me_c in enumerate(me_pc): 

125 me_cc = np.outer(me_c, me_c.conj()) 

126 k_cc += me_cc / (e_p[p] - omega - 1j * gamma) 

127 if self.nonresonant: 

128 k_cc += me_cc.conj() / (e_p[p] + omega + 1j * gamma) 

129 return k_cc 

130 

131 mr = 0 

132 for a, i, r in zip(self.myindices, self.myxyz, self.myr): 

133 if energy_derivative >= 0: 

134 V_rcc[r] += pre * ( 

135 kappa_cc(self.expm_rpc[mr], self.ex0E_p, 

136 omega, gamma, self.dipole_form) - 

137 kappa_cc(self.exmm_rpc[mr], self.ex0E_p, 

138 omega, gamma, self.dipole_form)) 

139 if energy_derivative: 

140 V_rcc[r] += pre * ( 

141 kappa_cc(self.ex0m_pc, self.expE_rp[mr], 

142 omega, gamma, self.dipole_form) - 

143 kappa_cc(self.ex0m_pc, self.exmE_rp[mr], 

144 omega, gamma, self.dipole_form)) 

145 mr += 1 

146 self.comm.sum(V_rcc) 

147 

148 return V_rcc 

149 

150 def electronic_me_Qcc(self, omega, gamma): 

151 self.read() 

152 Vel_rcc = np.zeros((self.ndof, 3, 3), dtype=complex) 

153 approximation = self.approximation.lower() 

154 if approximation == 'profeta': 

155 Vel_rcc += self.electronic_me_profeta_rcc(omega, gamma) 

156 elif approximation == 'placzek': 

157 Vel_rcc += self.electronic_me_profeta_rcc(omega, gamma, True) 

158 elif approximation == 'p-p': 

159 Vel_rcc += self.electronic_me_profeta_rcc(omega, gamma, -1) 

160 else: 

161 raise RuntimeError( 

162 'Bug: call with {} should not happen!'.format( 

163 self.approximation)) 

164 

165 return self.map_to_modes(Vel_rcc)