mirror of
https://github.com/vale981/HOPSFlow-Paper
synced 2025-03-04 17:21:42 -05:00
12 KiB
12 KiB
In here, we'll try to find out how long the cycle has to be, to ensure complete thermalization.
Boilerplate
import figsaver as fs
import plot_utils as pu
from hiro_models.one_qubit_model import StocProcTolerances
from hiro_models.otto_cycle import OttoEngine
import hiro_models.model_auxiliary as aux
import numpy as np
import qutip as qt
import utilities as ut
import stocproc
import matplotlib.pyplot as plt
import otto_utilities as ot
import ray
ray.shutdown()
#ray.init(address='auto')
ray.init()
from hops.util.logging_setup import logging_setup
import logging
logging_setup(logging.INFO)
plt.rcParams['figure.figsize'] = (12,4)Cycle
This model reflects the currently "best-working" configuration. The only parameter we'll vary is the cycle time θ. We'll be interested in whether the system thermalizes.
Therefore we need a way to keep the coupling switching speed constant.
def timings(τ_c, τ_i):
τ_th = (1 - 2 * τ_c) / 2
τ_i_on = τ_th - 2 * τ_i
timings_H = (0, τ_c, τ_c + τ_th, 2 * τ_c + τ_th)
timings_L_hot = (τ_c, τ_c + τ_i, τ_c + τ_i + τ_i_on, τ_c + 2 * τ_i + τ_i_on)
timings_L_cold = tuple(time + timings_H[2] for time in timings_L_hot)
return timings_H, (timings_L_cold, timings_L_hot)Now we define the prototype. The numeric accuracy is jacked down, as we don't need precision. We only use two cycles because of the long cycle time. With δ=.4 θ=70 seemed ok. So let's try with .6 and 50 ok. 1 and 40 is ok too. Let's try 1 and 45. Seems OK, but let's ease off on the coupling strength and try .8. That doesn't seem to be long enough. Lets try 55. .7 has better convergence, so let's try that.
def make_cycle(θ):
(p_H, p_L) = timings(3. / θ, 3. / θ)
return OttoEngine(
δ=[.7, .7],
ω_c=[1, 1],
ψ_0=qt.basis([2], [1]),
description=f"Classic Cycle",
k_max=4,
bcf_terms=[5] * 2,
truncation_scheme="simplex",
driving_process_tolerances=[StocProcTolerances(1e-3, 1e-3)] * 2,
thermal_process_tolerances=[StocProcTolerances(1e-3, 1e-3)] * 2,
T=[0.5, 4],
therm_methods=["tanhsinh", "tanhsinh"],
Δ=1,
num_cycles=2,
Θ=θ,
dt=0.001,
timings_H=p_H,
timings_L=p_L,
streaming_mode=True,
shift_to_resonance=(False, False),
L_shift=(0, 0),
) ot.plot_cycle(make_cycle(100))| <Figure | size | 1200x400 | with | 1 | Axes> | <AxesSubplot: | xlabel= | $\tau$ | ylabel= | Operator Norm | > |
Cursory Scanning
We can now test the model to find wether it allows enough time for complete thermalization. We'll start with a really long cycle.
long_cycle = make_cycle(60) ot.integrate_online(long_cycle, 100000)We can also import from taurus.
from hiro_models.model_auxiliary import import_results
import_results(
other_data_path="./taurus/.data",
other_results_path="./taurus/results",
interactive=False,
models_to_import=[long_cycle],
force=True,
)[INFO root 62054] Skipping cb958c7faea73706443c8601e965ef4e3c3dbd120dc415df21e77d8c98b80abd. [WARNING root 62054] Importing taurus/.data/fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a/_7/fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a_78e537f7d2c8ad853e73ab4421b7e329_1.h5 to .data/fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a/_7/fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a_78e537f7d2c8ad853e73ab4421b7e329_1.h5. [WARNING root 62054] The model description is 'Classic Cycle'. [WARNING root 62054] Importing taurus/results/flow_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz to results/flow_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz. [WARNING root 62054] Importing taurus/results/interaction_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz to results/interaction_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz. [WARNING root 62054] Importing taurus/results/interaction_power_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz to results/interaction_power_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz. [WARNING root 62054] Importing taurus/results/system_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz to results/system_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz. [WARNING root 62054] Importing taurus/results/system_power_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz to results/system_power_fe0aece8736438482a61340ae25a9b715d7c29d9347391915858eae7f449ee0a.npz. [INFO root 62054] Skipping 4f0fc3dd9c5abe8846ace8fe3d4aae731be36651bb8a350f565e151067be0e1f.
Now we look at the system energy.
f, a, *_ = pu.plot_with_σ(long_cycle.t, long_cycle.system_energy())
a.set_xlim(0, long_cycle.Θ)| 0.0 | 60.0 |
ot.plot_energy(long_cycle)| <Figure | size | 1200x400 | with | 1 | Axes> | <AxesSubplot: | xlabel= | $\tau$ | ylabel= | Energy | > |
| <Figure | size | 1200x400 | with | 1 | Axes> | <AxesSubplot: | xlabel= | $\tau$ | ylabel= | Energy | > |
pu.plot_with_σ(long_cycle.t, long_cycle.bath_energy_flow().sum_baths())| <Figure | size | 1200x400 | with | 1 | Axes> | <AxesSubplot: | > | ((<matplotlib.lines.Line2D at 0x7f65f7866d90>) <matplotlib.collections.PolyCollection at 0x7f65f7866fd0>) |
We would like to know how far away from the thermal state the system is.
def thermal_state(Ω, T):
ρ = np.array([[np.exp(-Ω/T), 0], [0, 1]])
ρ /= np.sum(np.diag(ρ))
return ρ import hops.util.utilities
from hopsflow.util import EnsembleValue
with aux.get_data(long_cycle) as data:
trace_dist_c = hops.util.utilities.trace_distance(data, relative_to=thermal_state(long_cycle.T[0], long_cycle.energy_gaps[0]))
trace_dist_h = hops.util.utilities.trace_distance(data, relative_to=thermal_state(long_cycle.T[1], long_cycle.energy_gaps[1]))
f, a = plt.subplots()
pu.plot_with_σ(long_cycle.t, EnsembleValue(trace_dist_c), ax=a, label=r"$||\rho(\tau)-\rho_c||$")
pu.plot_with_σ(long_cycle.t, EnsembleValue(trace_dist_h), ax=a, label=r"$||\rho(\tau)-\rho_h||$")
a.plot(long_cycle.t, (long_cycle.H(long_cycle.t)[:, 0, 0] - 1)/2, label="H Modulation")
a.set_xlabel(r"$\tau$")
#a.set_xlim(155)
a.legend()
fs.export_fig("thermalization")/nix/store/vkzza81mzwyk5br1c6cm67g48xycvmvl-python3-3.9.15-env/lib/python3.9/site-packages/matplotlib/cbook/__init__.py:1369: ComplexWarning: Casting complex values to real discards the imaginary part return np.asarray(x, float) /nix/store/vkzza81mzwyk5br1c6cm67g48xycvmvl-python3-3.9.15-env/lib/python3.9/site-packages/matplotlib/axes/_axes.py:5340: ComplexWarning: Casting complex values to real discards the imaginary part pts[0] = start /nix/store/vkzza81mzwyk5br1c6cm67g48xycvmvl-python3-3.9.15-env/lib/python3.9/site-packages/matplotlib/axes/_axes.py:5341: ComplexWarning: Casting complex values to real discards the imaginary part pts[N + 1] = end /nix/store/vkzza81mzwyk5br1c6cm67g48xycvmvl-python3-3.9.15-env/lib/python3.9/site-packages/matplotlib/axes/_axes.py:5344: ComplexWarning: Casting complex values to real discards the imaginary part pts[1:N+1, 1] = dep1slice /nix/store/vkzza81mzwyk5br1c6cm67g48xycvmvl-python3-3.9.15-env/lib/python3.9/site-packages/matplotlib/axes/_axes.py:5346: ComplexWarning: Casting complex values to real discards the imaginary part pts[N+2:, 1] = dep2slice[::-1]
Convergence
Let's see what hierarchy depth is really required.
plt.cla()
long_cycle = make_cycle(55)
ref = None
ks = (sorted([5,6,7]))[::-1]
for k_max in ks:
cyc = long_cycle.copy()
cyc.num_cycles = 1
cyc.ω_c = [.5, .5]
cyc.k_max = k_max
cyc.streaming_mode = False
cyc.δ = [.7,.7]
ot.integrate_online(cyc, 16*2)
with aux.get_data(cyc) as data:
acc = data.stoc_traj[:]
if k_max == max(ks):
ref = acc
else:
plt.plot(cyc.t, np.linalg.norm((np.abs(ref-acc)).sum(axis=0),axis=1) / np.linalg.norm((np.abs(ref)).sum(axis=0),axis=1) * 100, label=k_max)
plt.legend() [INFO hops.core.integration 78407] Choosing the nonlinear integrator.
[INFO root 78407] Starting analysis process.
[INFO root 78407] Started analysis process with pid 79308.
[INFO hops.core.hierarchy_data 78407] Creating the streaming fifo at: /home/hiro/Documents/Projects/UNI/master/eflow_paper/python/otto_motor/subprojects/relaxation/results_cf28536b03d6a994d7d9761ccb3171527d6b2888be6ec8db22e2504b479d93f8.fifo
[INFO hops.core.integration 78407] Using 16 integrators.
[INFO hops.core.integration 78407] Some 0 trajectories have to be integrated.
[INFO hops.core.integration 78407] Using 19448 hierarchy states.
0it [00:00, ?it/s]
[INFO hops.core.integration 78407] Choosing the nonlinear integrator.
[INFO root 78407] Starting analysis process.
[INFO root 78407] Started analysis process with pid 79318.
[INFO hops.core.hierarchy_data 78407] Creating the streaming fifo at: /home/hiro/Documents/Projects/UNI/master/eflow_paper/python/otto_motor/subprojects/relaxation/results_7ce1482692af554bfd6e3fa80425b8aad003eb240166c0d1e0595ebc60feaf04.fifo
[INFO hops.core.integration 78407] Using 16 integrators.
[INFO hops.core.integration 78407] Some 0 trajectories have to be integrated.
[INFO hops.core.integration 78407] Using 8008 hierarchy states.
0it [00:00, ?it/s]
[INFO hops.core.integration 78407] Choosing the nonlinear integrator.
[INFO root 78407] Starting analysis process.
[INFO root 78407] Started analysis process with pid 79323.
[INFO hops.core.hierarchy_data 78407] Creating the streaming fifo at: /home/hiro/Documents/Projects/UNI/master/eflow_paper/python/otto_motor/subprojects/relaxation/results_9f7e27b79a81cf1c8cb244929fff927fc21dd5733dd0adbc968a20df514ce80b.fifo
[INFO hops.core.integration 78407] Using 16 integrators.
[INFO hops.core.integration 78407] Some 0 trajectories have to be integrated.
[INFO hops.core.integration 78407] Using 3003 hierarchy states.
0it [00:00, ?it/s]
<matplotlib.legend.Legend at 0x7fa4ffac1e80>
-
For ω_c=2:
- For δ=.8 i found great variations from k=5 to k=6. This is BAD. For .4 it is ok. For .7 too at least for 16 trajectories. For .7 k=4 seems OK.
- same for ω_c=1,.5