Coverage for local_installation/dynasor/tools/damped_harmonic_oscillator.py: 100%

30 statements  

« prev     ^ index     » next       coverage.py v7.3.2, created at 2024-12-21 12:02 +0000

1import numpy as np 

2from numpy.typing import NDArray 

3 

4 

5def acf_position_dho(t: NDArray[float], w0: float, gamma: float, A: float = 1.0): 

6 r""" 

7 The damped damped harmonic oscillator (DHO) autocorrelation function for the position. 

8 The definition of this function can be found in the `dynasor documentation 

9 <dynasor.materialsmodeling.org/theory.html#damped-harmonic-oscillator-model>_`. 

10 

11 Parameters 

12 ---------- 

13 t 

14 Time, usually as an array. 

15 w0 

16 Natural angular frequency of the DHO. 

17 gamma 

18 Damping of DHO. 

19 A 

20 Amplitude of the DHO. 

21 """ 

22 

23 t = np.abs(t) 

24 

25 if 2 * w0 > gamma: # underdamped 

26 we = np.sqrt(w0**2 - gamma**2 / 4.0) 

27 return A * np.exp(-gamma * t / 2.0) * ( 

28 np.cos(we * t) + 0.5 * gamma / we * np.sin(we * t)) 

29 elif 2 * w0 < gamma: # overdamped 

30 tau = 2 / gamma 

31 tau_S = tau / (1 + np.sqrt(1 - (w0 * tau)**2)) 

32 tau_L = tau / (1 - np.sqrt(1 - (w0 * tau)**2)) 

33 return A / (tau_L - tau_S) * (tau_L * np.exp(-t/tau_L) - tau_S * np.exp(-t/tau_S)) 

34 else: 

35 tau = 2 / gamma 

36 return A * np.exp(-t/tau) * (1 + t / tau) 

37 

38 

39def acf_velocity_dho(t: NDArray[float], w0: float, gamma: float, A: float = 1.0): 

40 r""" 

41 The damped damped harmonic oscillator (DHO) autocorrelation function for the velocity. 

42 The definition of this function can be found in the `dynasor documentation 

43 <dynasor.materialsmodeling.org/theory.html#damped-harmonic-oscillator-model>_`. 

44 

45 Parameters 

46 ---------- 

47 t 

48 Time, usually as an array. 

49 w0 

50 Natural angular frequency of the DHO. 

51 gamma 

52 Damping of DHO. 

53 A 

54 Amplitude of the DHO. 

55 """ 

56 

57 t = np.abs(t) 

58 

59 if 2 * w0 > gamma: # underdamped 

60 we = np.sqrt(w0**2 - gamma**2 / 4.0) 

61 return A * w0**2 * np.exp(-gamma * t / 2.0) * ( 

62 np.cos(we * t) - 0.5 * gamma / we * np.sin(we * t)) 

63 elif 2 * w0 < gamma: # overdamped 

64 tau = 2 / gamma 

65 tau_S = tau / (1 + np.sqrt(1 - (w0 * tau)**2)) 

66 tau_L = tau / (1 - np.sqrt(1 - (w0 * tau)**2)) 

67 return A / (tau_L - tau_S) * (np.exp(-t/tau_S)/tau_S - np.exp(-t/tau_L)/tau_L) 

68 else: 

69 tau = 2 / gamma 

70 return A * np.exp(-t/tau) * (1 - t / tau) 

71 

72 

73def psd_position_dho(w: NDArray[float], w0: float, gamma: float, A: float = 1.0): 

74 r""" 

75 The power spectral density (PSD) function for the damped harmonic oscillator (DHO) 

76 (i.e., the Fourier transform of the autocorrelation function) for the position. 

77 

78 The definition of this function can be found in the `dynasor documentation 

79 <dynasor.materialsmodeling.org/theory.html#damped-harmonic-oscillator-model>_`. 

80 

81 Parameters 

82 ---------- 

83 w 

84 Angular frequency, usually as an array. 

85 w0 

86 Natural angular frequency of the DHO. 

87 gamma 

88 Damping of DHO. 

89 A 

90 Amplitude of the DHO. 

91 """ 

92 return 2 * w0**2 * A * gamma / ((w**2 - w0**2)**2 + (w * gamma)**2) 

93 

94 

95def psd_velocity_dho(w: NDArray[float], w0: float, gamma: float, A: float = 1.0): 

96 r""" 

97 The power spectral density (PSD) function for the damped harmonic oscillator (DHO) 

98 (i.e., the Fourier transform of the autocorrelation function) for the position. 

99 

100 The definition of this function can be found in the `dynasor documentation 

101 <dynasor.materialsmodeling.org/theory.html#damped-harmonic-oscillator-model>_`. 

102 

103 Parameters 

104 ---------- 

105 w 

106 Angular frequency, usually as an array. 

107 w0 

108 Natural angular frequency of the DHO. 

109 gamma 

110 Damping of DHO. 

111 A 

112 Amplitude of the DHO. 

113 """ 

114 return w**2 * psd_position_dho(w, w0, gamma, A)