@@ -66,9 +66,6 @@ def timer(label):
66
66
options = {
67
67
"nsteps": 15000,
68
68
"store_states": True,
69
- "rtol": 1e-12,
70
- "atol": 1e-12,
71
- "min_step": 1e-18,
72
69
"method": "vern9",
73
70
"progress_bar": "enhanced",
74
71
}
@@ -83,9 +80,9 @@ And let us set up the system Hamiltonian and bath parameters:
83
80
#
84
81
# We use the Hamiltonian employed in
85
82
# https://www.pnas.org/content/106/41/17255 and operate
86
- # in units of Hz :
83
+ # in units of GHz :
87
84
88
- Hsys = 3e10 * 2 * np.pi * Qobj([
85
+ Hsys = 30 * 2 * np.pi * Qobj([
89
86
[200, -87.7, 5.5, -5.9, 6.7, -13.7, -9.9],
90
87
[-87.7, 320, 30.8, 8.2, 0.7, 11.8, 4.3],
91
88
[5.5, 30.8, 0, -53.5, -2.2, -9.6, 6.0],
@@ -99,9 +96,9 @@ Hsys = 3e10 * 2 * np.pi * Qobj([
99
96
``` {code-cell} ipython3
100
97
# Bath parameters
101
98
102
- lam = 35 * 3e10 * 2 * np.pi
103
- gamma = 1 / 166e-15
104
- T = 300 * 0.6949 * 3e10 * 2 * np.pi
99
+ lam = 35 * 30 * 2 * np.pi
100
+ gamma = 1 / 166e-6
101
+ T = 300 * 0.6949 * 30 * 2 * np.pi
105
102
beta = 1 / T
106
103
```
107
104
@@ -114,18 +111,18 @@ env = DrudeLorentzEnvironment(T=T, lam=lam, gamma=gamma)
114
111
```
115
112
116
113
``` {code-cell} ipython3
117
- wlist = np.linspace(0, 200 * 3e10 * 2 * np.pi, 100)
118
- tlist = np.linspace(0, 1e-12 , 1000)
114
+ wlist = np.linspace(0, 200 * 30 * 2 * np.pi, 100)
115
+ tlist = np.linspace(0, 1e-3 , 1000)
119
116
120
- J = env.spectral_density(wlist) / (3e10 * 2 * np.pi)
117
+ J = env.spectral_density(wlist) / (30 * 2 * np.pi)
121
118
122
119
fig, axes = plt.subplots(1, 2, sharex=False, figsize=(10, 3))
123
120
124
121
fig.subplots_adjust(hspace=0.1) # reduce space between plots
125
122
126
123
# Spectral density plot:
127
124
128
- axes[0].plot(wlist / (3e10 * 2 * np.pi), J, color='r', ls='--', label="J(w)")
125
+ axes[0].plot(wlist / (30 * 2 * np.pi), J, color='r', ls='--', label="J(w)")
129
126
axes[0].set_xlabel(r'$\omega$ (cm$^{-1}$)', fontsize=20)
130
127
axes[0].set_ylabel(r"$J(\omega)$ (cm$^{-1}$)", fontsize=16)
131
128
axes[0].legend()
@@ -193,7 +190,7 @@ linestyles = [
193
190
for m in range(7):
194
191
Q = basis(7, m) * basis(7, m).dag()
195
192
axes.plot(
196
- np.array(tlist) * 1e15 ,
193
+ np.array(tlist) * 1e6 ,
197
194
np.real(expect(outputFMO_HEOM.states, Q)),
198
195
label=m + 1,
199
196
color=colors[m % len(colors)],
@@ -232,7 +229,7 @@ And now let's plot the Bloch-Redfield solver results:
232
229
fig, axes = plt.subplots(1, 1, figsize=(12, 8))
233
230
234
231
for m, Q in enumerate(Q_list):
235
- axes.plot(tlist * 1e15 , expect(outputFMO_BR.states, Q), label=m + 1)
232
+ axes.plot(tlist * 1e6 , expect(outputFMO_BR.states, Q), label=m + 1)
236
233
237
234
axes.set_xlabel(r'$t$ (fs)', fontsize=30)
238
235
axes.set_ylabel(r"Population", fontsize=30)
@@ -334,7 +331,7 @@ with timer("ME ODE solver"):
334
331
fig, axes = plt.subplots(1, 1, figsize=(12, 8))
335
332
336
333
for m, Q in enumerate(Q_list):
337
- axes.plot(tlist * 1e15 , expect(outputFMO_ME.states, Q), label=m + 1)
334
+ axes.plot(tlist * 1e6 , expect(outputFMO_ME.states, Q), label=m + 1)
338
335
339
336
axes.set_xlabel(r'$t$', fontsize=20)
340
337
axes.set_ylabel(r"Population", fontsize=16)
@@ -363,7 +360,7 @@ with timer("ME ODE solver"):
363
360
fig, axes = plt.subplots(1, 1, figsize=(12, 8))
364
361
for m, Q in enumerate(Q_list):
365
362
axes.plot(
366
- tlist * 1e15 ,
363
+ tlist * 1e6 ,
367
364
expect(outputFMO_ME_nodephase.states, Q),
368
365
label=m + 1,
369
366
)
0 commit comments