master-thesis/python/energy_flow_proper/02_finite_temperature/notebook.org

• Configuration and Setup
• Using the Data
  • Jupyter Setup
  • Load the Data
  • Calculate System Energy
  • Calculate the Heat Flow
  • Calculate the Interaction Energy
• Close the Data File

Configuration and Setup

This will be tangled into the config file that can be used with the HOPS cli.

  from hops.core.hierarchy_parameters import HIParams, HiP, IntP, SysP, ResultType
  from hops.core.hierarchyLib import HI
  from hops.util.bcf_fits import get_ohm_g_w
  from hops.util.truncation_schemes import TruncationScheme_Power_multi
  import hops.util.bcf
  import numpy as np
  import hops.util.matrixLib as ml
  from stocproc import StocProc_FFT

  wc = 5
  s = 1

  # The BCF fit
  bcf_terms = 5
  g, w = get_ohm_g_w(bcf_terms, s, wc)

  integration = IntP(t_max=30, t_steps=int(30 // 0.001))
  system = SysP(
      H_sys=0.5 * np.array([[-1, 0], [0, 1]]),
      L=0.5 * np.array([[0, 1], [1, 0]]),
      psi0=np.array([0, 1]),
      g=g,
      w=w,
      bcf_scale=0.5,
      T=0.5,
  )

  params = HIParams(
      SysP=system,
      IntP=integration,
      HiP=HiP(
          nonlinear=True,
          normalized_by_hand=True,
          result_type=ResultType.ZEROTH_AND_FIRST_ORDER,
          truncation_scheme=TruncationScheme_Power_multi.from_g_w(
              g=.5 * g, w=w, p=1, q=0.5, kfac=1.2
          ),
          save_therm_rng_seed=True,
      ),
      Eta=StocProc_FFT(
          spectral_density=hops.util.bcf.OhmicSD_zeroTemp(
              s,
              1,
              wc,
          ),
          alpha=hops.util.bcf.OhmicBCF_zeroTemp(
              s,
              1,
              wc,
          ),
          t_max=integration.t_max,
          intgr_tol=1e-4,
          intpl_tol=1e-4,
          negative_frequencies=False,
      ),
      EtaTherm=StocProc_FFT(
        spectral_density=hops.util.bcf.Ohmic_StochasticPotentialDensity(
            s, 1, wc, beta=1 / system.__non_key__["T"]
        ),
        alpha=hops.util.bcf.Ohmic_StochasticPotentialCorrelations(
            s, 1, wc, beta=1 / system.__non_key__["T"]
        ),
        t_max=integration.t_max,
        intgr_tol=1e-4,
        intpl_tol=1e-4,
        negative_frequencies=False,
    ),
  )

Using the Data

Jupyter Setup

  import numpy as np
  import matplotlib.pyplot as plt
  import utilities as ut
  import figsaver as fs

Load the Data

  from hopsflow import hopsflow, util
  from hops.core.hierarchyData import HIMetaData

Now we read the trajectory data.

  class result:
      hd = HIMetaData("data", ".").get_HIData(params, read_only=True)
      N = hd.samples
      τ = hd.get_time()
      ψ_1 = hd.aux_states
      ψ = hd.stoc_traj
      seeds = hd.rng_seed

  result.N

10000

Calculate System Energy

Simple sanity check.

  _, e_sys, σ_e_sys = util.operator_expectation_ensemble(
      iter(result.ψ),
      system.H_sys,
      result.N,
      params.HiP.nonlinear,
      save="./results/new_energy.npy",
  )

  with fs.hiro_style():
      plt.gcf().set_size_inches(4, 3)
      plt.errorbar(result.τ, e_sys.real, yerr=σ_e_sys.real, ecolor="yellow")
      plt.ylabel(r"$\langle H_S\rangle$")
      plt.xlabel(r"$τ$")
      fs.export_fig("system_energy")

100% 9999/9999 [00:18<00:00, 537.35it/s]

Calculate the Heat Flow

Now let's calculate the heatflow. In this simple case it is engouh to know the first hierarchy states.

First we set up some parameter objects for the alogrithm.

  hf_system = hopsflow.SystemParams(
      system.L, system.g, system.w, system.bcf_scale, params.HiP.nonlinear
  )

  η = params.Eta
  ξ = params.EtaTherm
  ξ.calc_deriv = True

  hf_therm = hopsflow.ThermalParams(ξ=ξ, τ=result.τ, num_deriv=False)

Now we can apply our tooling to one trajectory for testing.

  hf_sample_run = hopsflow.HOPSRun(result.ψ[0], result.ψ_1[0], hf_system)
  hf_sample_run_therm = hopsflow.ThermalRunParams(hf_therm, result.seeds[0])
  first_flow = hopsflow.flow_trajectory_coupling(hf_sample_run, hf_system)
  first_flow_therm = hopsflow.flow_trajectory_therm(hf_sample_run, hf_sample_run_therm)
  plt.plot(result.τ, first_flow)
  plt.plot(result.τ, first_flow_therm)

<matplotlib.lines.Line2D

at

0x7eff709eb220>

And now for all trajectories.

  full_flow = hopsflow.heat_flow_ensemble(
      iter(result.ψ),
      iter(result.ψ_1),
      hf_system,
      result.N,
      (iter(result.seeds), hf_therm),
      every=result.N // 10,
      save="results/flow.npy",
  )

  with ut.hiro_style():
      _, ax = ut.plot_convergence(result.τ, full_flow, transform=lambda y: -y)
      ax.legend()

We can integrate the energy change in the bath:

  for i in range(len(full_flow)):
      e_bath = util.integrate_array(-full_flow[i][1], result.τ)
      plt.plot(result.τ, e_bath)

Calculate the Interaction Energy

First we calculate it from energy conservation.

  e_int = (1/2 - e_sys - e_bath).real
  plt.plot(result.τ, e_int)

<matplotlib.lines.Line2D

at

0x7eff7033e2e0>

And then from first principles:

  _, e_int_ex, _ = hopsflow.interaction_energy_ensemble(
      result.ψ,
      result.ψ_1,
      hf_system,
      result.N,
      (result.seeds, hf_therm),
      save="results/interaction_energy.npy",
  )

And both together:

  plt.plot(result.τ, e_int, label="integrated")
  plt.plot(result.τ, e_int_ex, label="exact")
  plt.legend()

<matplotlib.legend.Legend at 0x7eff70317df0>

Seems to work :P.

Close the Data File

We need to release the hold on the file.

 result.hd.close()