#+PROPERTY: header-args :session 01_zero_temp :kernel python :pandoc t :async yes * Configuration and Setup This will be tangled into the [[file:config.py][config file]] that can be used with the HOPS cli. #+begin_src jupyter-python :results none :tangle config.py 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 = 6 g, w = get_ohm_g_w(bcf_terms, s, wc) integration = IntP(t_max=30, t_steps=int(30 // 0.01)) 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.8, T=0, ) 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=g, w=w, p=1, q=0.5, kfac=1.4 ), ), 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-3, intpl_tol=1e-3, negative_frequencies=False, ), EtaTherm=None, ) #+end_src * Hops Integration We can use multiple avenues. ** Local Integration #+begin_src vterm :term-name integration hi 500 integrate #+end_src #+RESULTS: :RESULTS: Linux ArLeenUX 5.15.11-zen1 x86_64 14:45:17 up 2 days 18:10, 1 user, load average: 1.55, 1.48, 1.44 impure  ~/D/P/U/m/m/p/e/01_zero_temperature  hi 500 integrate Loading the configuration from config.py. This might take a while... / JobManager started on ArLeenUX:36931 (bytearray(b'HOPS36931')) [TET 6.88ms [0.0c/s] TTG -- 0.0% ETA -- ORT --] /nix/store/h9pzzn81vrj2vzhavbjgv101yc52qpx0-python3-3.9.9-env/lib/python3.9/site-packages/jobmanager/jobmanager.py:130: UserWa rning: num_threads could not be set, MKL / openblas not found warnings.warn("num_threads could not be set, MKL / openblas not found") /nix/store/h9pzzn81vrj2vzhavbjgv101yc52qpx0-python3-3.9.9-env/lib/python3.9/site-packages/jobmanager/jobmanager.py:130: UserWa rning: num_threads could not be set, MKL / openblas not found warnings.warn("num_threads could not be set, MKL / openblas not found") /nix/store/h9pzzn81vrj2vzhavbjgv101yc52qpx0-python3-3.9.9-env/lib/python3.9/site-packages/jobmanager/jobmanager.py:130: UserWa rning: num_threads could not be set, MKL / openblas not found warnings.warn("num_threads could not be set, MKL / openblas not found") /nix/store/h9pzzn81vrj2vzhavbjgv101yc52qpx0-python3-3.9.9-env/lib/python3.9/site-packages/jobmanager/jobmanager.py:130: UserWa [TET 1.01s [0.0c/s] TTG -- 0.0% ETA -- ORT --] res_q #0 0/s 0kB|rem.:492, done:0, failed:0, prog.:0 w1:00:00:14 [5.3c/min] #1 - 2.72s [82.5c/s] [==============> ] TTG 00:00:10 w2:00:00:14 [5.0c/min] #1 - 2.09s [82.4c/s] [===========> ] TTG 00:00:11 w3:00:00:14 [5.1c/min] #1 - 2.20s [82.8c/s] [============> ] TTG 00:00:10 w4:00:00:14 [5.0c/min] #1 - 1.97s [82.9c/s] [==========> ] TTG 00:00:11 local res_q 0 285.1GB/s [TET 00:00:14 [18.4c/min] TTG 00:26:29 0.8% ETA 20220117_15:12:02 ORT 00:26:43] res_q #0 2.122GB/s 28.72GB|rem.:484, done:4, failed:0, prog.:4 ^Ccapi_return is NULL Call-back cb_f_in_zvode__user__routines failed. capi_return is NULL Call-back cb_f_in_zvode__user__routines failed. capi_return is NULL Call-back cb_f_in_zvode__user__routines failed. SystemExit, quit processing, reinsert current argument, please wait SystemExit, quit processing, reinsert current argument, please wait capi_return is NULL Call-back cb_f_in_zvode__user__routines failed. [TET 00:00:15 [18.4c/min] TTG 00:26:29 0.8% ETA 20220117_15:12:03 ORT 00:26:44] res_q #0 1.976GB/s 28.72GB|rem.:484, done:4, failed:0, prog.:4 ############## in JM SERVER EXIT HIServer start at 2022-01-17 14:45:19.308956 | runtime 1.500e+01s HIServer total number of jobs : 492 w1:00:00:15 [5.3c/min] #1 - [TET-3.72s----[74.9c/s]-TTG>00:00:10 27.9% ETA 20220117_14:45:45 ORT 00:00:13] w2:00:00:15 [5.0c/min] #1 - [TET-3.10s----[72.0c/s] TTG 00:00:11 22.3% ETA 20220117_14:45:46 ORT 00:00:14] w3:00:00:15 [5.1c/min] #1 - [TET-3.21s----[73.9c/s]>TTG 00:00:11 23.7% ETA 20220117_14:45:46 ORT 00:00:14] w4:00:00:15 [5.0c/min] #1 - [TET-2.97s----[73.6c/s] TTG 00:00:11 21.9% ETA 20220117_14:45:46 ORT 00:00:13] local res_q 0 285.1GB/s impure  ~/D/P/U/m/m/p/e/01_zero_temperature  exit 18.9s :END: And there we go. It is better to run the above command in a vterm-session. ** Remote/Distributed Integration We start the server locally. #+begin_src vterm :term-name local-server hi 1000 start-server #+end_src #+RESULTS: :RESULTS: Linux ArLeenUX 5.15.11-zen1 x86_64 16:34:02 up 2 days 19:58, 1 user, load average: 1.08, 1.67, 2.52 impure  ~/D/P/U/m/m/p/e/01_zero_temperature  hi 1000 start-server Loading the configuration from config.py. This might take a while... / JobManager started on ArLeenUX:42524 (bytearray(b'hierarchy')) [TET-00:12:05--[43.6c/min]-TTG-0.00ms-------------------------100%-------------------------ETA-20220117_16:46:11-ORT-00:12:05] res_q #0 14.84GB/s 10.51TB|rem.:0, done:500, failed:0, prog.:0 ############## in JM SERVER EXIT HIServer start at 2022-01-17 16:34:05.075876 | runtime 7.280e+02s HIServer total number of jobs : 500 | processed : 500 | succeeded : 500 | failed : 0 | timing in sec: min 1.386e+01 | max 4.145e+01 | avr 2.478e+01 | not processed : 0 | queried : 0 | not queried yet : 0 :END: And jack in with a remote client. In this case my box at home. ** Client Starting a client is trivial. #+begin_src vterm :term-name local-client client localhost #+end_src #+RESULTS: :RESULTS: Linux ArLeenUX 5.15.11-zen1 x86_64 16:34:08 up 2 days 19:58, 1 user, load average: 1.07, 1.66, 2.51 impure  ~/D/P/U/m/m/p/e/01_zero_temperature  client localhost /nix/store/zwwf4gkrcx2ly273ivhn7a0bwwl0r9ki-python3-3.9.9-env/lib/python3.9/site-packages/jobmanager/jobmanager.py:130: UserWa rning: num_threads could not be set, MKL / openblas not found warnings.warn("num_threads could not be set, MKL / openblas not found") /nix/store/zwwf4gkrcx2ly273ivhn7a0bwwl0r9ki-python3-3.9.9-env/lib/python3.9/site-packages/jobmanager/jobmanager.py:130: UserWa rning: num_threads could not be set, MKL / openblas not found warnings.warn("num_threads could not be set, MKL / openblas not found") /nix/store/zwwf4gkrcx2ly273ivhn7a0bwwl0r9ki-python3-3.9.9-env/lib/python3.9/site-packages/jobmanager/jobmanager.py:130: UserWa rning: num_threads could not be set, MKL / openblas not found warnings.warn("num_threads could not be set, MKL / openblas not found") /nix/store/zwwf4gkrcx2ly273ivhn7a0bwwl0r9ki-python3-3.9.9-env/lib/python3.9/site-packages/jobmanager/jobmanager.py:130: UserWa rning: num_threads could not be set, MKL / openblas not found warnings.warn("num_threads could not be set, MKL / openblas not found") w1:00:11:55 [2.0c/min] #23 - [TET 7.38s [0.0c/s] TTG -- 0.0% ETA -- ORT --] w2:00:11:55 [2.0c/min] #23 - [TET 4.46s [0.0c/s] TTG -- 0.0% ETA -- ORT --] w3:00:11:55 [2.1c/min] #23 - [TET 5.82s [0.0c/s] TTG -- 0.0% ETA -- ORT --] w4:00:11:55 [2.1c/min] #23 - [TET 9.32s [0.0c/s] TTG -- 0.0% ETA -- ORT --] local res_q 0 416.7GB/s :END: * Using the Data ** Jupyter Setup #+begin_src jupyter-python :results none import matplotlib.pyplot as plt import numpy as np import utilities as ut #+end_src ** Load the Data #+begin_src jupyter-python :results none from hopsflow import hopsflow, util from hops.core.hierarchyLib import HI #+end_src Now we read the trajectory data. #+begin_src jupyter-python class result: with HI(params, 500).get_data(read_only=True) as hd: N = hd.samples τ = hd.get_time() ρ = hd.get_rho_t() ψ_1 = np.array(hd.aux_states)[0:N] ψ = np.array(hd.stoc_traj)[0:N] result.N #+end_src #+RESULTS: : 1000 ** Calculate System Energy Simple sanity check. #+begin_src jupyter-python _, e_sys, σ_e_sys = util.operator_expectation_ensemble( iter(result.ψ), system.H_sys, result.N, params.HiP.nonlinear, save="./results/energy.npy" ) plt.errorbar(result.τ, e_sys.real, yerr=σ_e_sys.real, ecolor="yellow") #+end_src #+RESULTS: :RESULTS: : 100% 999/999 [00:00<00:00, 1645.43it/s] : [[file:./.ob-jupyter/bd7fe1b1f16f9711d71689389818fa54480dba99.svg]] :END: The energy bleeds out of the system. We don't reach the steady state yet. Also we don't loose all the energy. The energy eigenvalues of the system are. #+begin_src jupyter-python np.linalg.eig(system.H_sys)[0] #+end_src #+RESULTS: : array([-0.5, 0.5]) The begin and and values of the system energy expectation are. #+begin_src jupyter-python e_sys[0].real, e_sys[-1].real #+end_src #+RESULTS: | 0.5 | -0.44770384926040235 | And the total energy lost is: #+begin_src jupyter-python e_sys[0].real - e_sys[-1].real #+end_src #+RESULTS: : 0.9477038492604024 We do start in the state. #+begin_src jupyter-python system.psi0 #+end_src #+RESULTS: : array([0, 1]) ** 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. #+begin_src jupyter-python :results none hf_system = hopsflow.SystemParams( system.L, system.g, system.w, system.bcf_scale, params.HiP.nonlinear ) #+end_src Now we can apply our tooling to one trajectory for testing. #+begin_src jupyter-python hf_sample_run = hopsflow.HOPSRun(result.ψ[0], result.ψ_1[0], hf_system) first_flow = hopsflow.flow_trajectory_coupling(hf_sample_run, hf_system) with ut.hiro_style(): plt.plot(result.τ, first_flow) #+end_src #+RESULTS: [[file:./.ob-jupyter/c7984b39bab2dfd43f7a9f95bdf0badd202394a0.svg]] And now for all trajectories. #+begin_src jupyter-python full_flow = hopsflow.heat_flow_ensemble( result.ψ, result.ψ_1, hf_system, result.N, 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() #+end_src #+RESULTS: :RESULTS: : 100% 999/999 [00:01<00:00, 503.21it/s] [[file:./.ob-jupyter/c998467651c2660a92b64ec40f17672cbd6f16aa.svg]] :END: We can integrate the energy change in the bath: #+begin_src jupyter-python e_bath = util.integrate_array(-full_flow[-1][1], result.τ) plt.plot(result.τ, e_bath) #+end_src #+RESULTS: :RESULTS: | | [[file:./.ob-jupyter/1ef90eca2ed89295be26af97745a08d32a9aa037.svg]] :END: ** Calculate the Interaction Energy First we calculate it from energy conservation. #+begin_src jupyter-python e_int = (1/2 - e_sys - e_bath).real with ut.hiro_style(): plt.plot(result.τ, e_int) #+end_src #+RESULTS: [[file:./.ob-jupyter/10f25f6b2f73d65ae8fb7775fef8dce0379fd08e.svg]] And then from first principles: #+begin_src jupyter-python _, e_int_ex, σ_e_int_ex = hopsflow.interaction_energy_ensemble(result.ψ, result.ψ_1, hf_system, result.N) with ut.hiro_style(): plt.errorbar(result.τ, e_int_ex, yerr=σ_e_int_ex, ecolor="yellow") #+end_src #+RESULTS: :RESULTS: : 100% 999/999 [00:02<00:00, 406.62it/s] [[file:./.ob-jupyter/d14fe1b04f0527eec44c69990195936cdcbbd965.svg]] :END: And both together: #+begin_src jupyter-python with ut.hiro_style(): plt.errorbar(result.τ, e_int_ex, yerr=σ_e_int_ex, ecolor="yellow") plt.plot(result.τ, e_int) #+end_src #+RESULTS: [[file:./.ob-jupyter/7cebf22e2e864368752a1919b342b685fcdd1303.svg]]