HOPSFlow-Paper/subprojects/cycle_shift/finding_relaxation_time.org
2023-11-28 13:08:17 -05:00

6.9 KiB
Raw Permalink Blame History

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.

  def make_cycle(θ):
      (p_H, p_L) = timings(3. / θ, 3. / θ)

      return OttoEngine(
          δ=[.8, .8],
          ω_c=[2, 2],
          ψ_0=qt.basis([2], [1]),
          description=f"Classic Cycle",
          k_max=4,
          bcf_terms=[4] * 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 >

/hiro/HOPSFlow-Paper/media/branch/main/subprojects/cycle_shift/.ob-jupyter/e0c333c986674c4d0bcb9fafbd5cd47a6fa77210.svg

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.

We tried 100, now we try 60. 60 is a bit too short. So let's try 70

  long_cycle = make_cycle(45)
  ot.integrate_online(long_cycle, 100_000)

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,
  )

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 45.0

/hiro/HOPSFlow-Paper/media/branch/main/subprojects/cycle_shift/.ob-jupyter/3690e1e23aa7103ad49ab3612aa2b0baaa1e5dad.svg

/hiro/HOPSFlow-Paper/media/branch/main/subprojects/cycle_shift/.ob-jupyter/2c1a4d916249a5998d36181e93f93a3a46712b94.svg

  ot.plot_energy(long_cycle)
<Figure size 1200x400 with 1 Axes> <AxesSubplot: xlabel= $\tau$ ylabel= Energy >

/hiro/HOPSFlow-Paper/media/branch/main/subprojects/cycle_shift/.ob-jupyter/7fa0ef03b52d099133087281c56b939c001c3563.svg

pu.plot_with_σ(long_cycle.t, long_cycle.bath_energy_flow().sum_baths())
<Figure size 1200x400 with 1 Axes> <AxesSubplot: > ((<matplotlib.lines.Line2D at 0x7fb4f3b63bb0>) <matplotlib.collections.PolyCollection at 0x7fb4f3b73d00>)

/hiro/HOPSFlow-Paper/media/branch/main/subprojects/cycle_shift/.ob-jupyter/6b4e457a1dd6b1339a8fe78be12d02718f63d2d4.svg

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]

/hiro/HOPSFlow-Paper/media/branch/main/subprojects/cycle_shift/.ob-jupyter/ed0b6c9c04abab45ec849730cf1d7ed025ece4e8.svg