add finite T to 11

This commit is contained in:
Valentin Boettcher 2022-07-22 19:25:26 +02:00
parent 9e6ef7cfad
commit 8a222a7ac9
2 changed files with 541 additions and 0 deletions

View file

@ -0,0 +1,217 @@
import figsaver as fs
import plot_utils as pu
import hops.util.bcf
import hops.util.bcf_fits
from hops.core.hierarchy_parameters import *
from stocproc import StocProc_FFT, StocProc_TanhSinh
from hops.util.abstract_truncation_scheme import TruncationScheme_Simplex
from hops.core.integration import HOPSSupervisor
import qutip
def ho_zero(max_HO_level=15, wc=2, s=1, bcf_terms=7, k_max=5, bcf_scale=1, T=1):
hops_bcf = hops.util.bcf.OhmicBCF_zeroTemp(
s,
1,
wc,
normed=True
)
g, w = hops_bcf.exponential_coefficients(bcf_terms)
integration = IntP(t_steps=(40, 1000))
q = qutip.operators.create(max_HO_level) + qutip.operators.destroy(max_HO_level)
p = (
qutip.operators.destroy(max_HO_level) - qutip.operators.create(max_HO_level)
) / 1j
system = SysP(
H_sys=0.25 * (p ** 2 + q ** 2).full(),
L=[0.5 * q.full()],
psi0=qutip.states.fock(max_HO_level, n=0).full().flatten(),
g=[g],
w=[w],
bcf_scale=[bcf_scale],
T=[.5],
description="A single HO at zero Temperature",
)
params = HIParams(
SysP=system,
IntP=integration,
HiP=HiP(
seed=0,
nonlinear=True,
result_type=ResultType.ZEROTH_AND_FIRST_ORDER,
truncation_scheme=TruncationScheme_Simplex(k_max),
save_therm_rng_seed=True,
auto_normalize=True,
accum_only=False,
rand_skip=None,
),
Eta=[
StocProc_FFT(
spectral_density=hops.util.bcf.OhmicSD_zeroTemp(
s,
1,
wc,
normed=True
),
alpha=hops_bcf,
t_max=integration.t_max,
intgr_tol=1e-5,
intpl_tol=1e-5,
negative_frequencies=False,
)
],
EtaTherm=[StocProc_TanhSinh(
spectral_density=hops.util.bcf.Ohmic_StochasticPotentialDensity(
s,
1,
wc,
beta=1/system.__non_key__["T"][0],
normed=True
),
alpha=hops.util.bcf.Ohmic_StochasticPotentialCorrelations(
s,
1,
wc,
beta=1/system.__non_key__["T"][0],
normed=True
),
t_max=integration.t_max,
intgr_tol=1e-5,
intpl_tol=1e-5,
negative_frequencies=False,
)],
)
return params
model_keys = [dict(bcf_scale=1, wc=ω_c) for ω_c in [1, 2]] + [dict(bcf_scale=s, wc=2) for s in [.5, 1, 2]]
multi_params = [ho_zero(**keys) for keys in model_keys]
multi_params = multi_params
import ray
ray.shutdown()
ray.init()
supervisors = []
for params in multi_params:
supervisor = HOPSSupervisor(
params,
10_000,
data_path="ho_data",
data_name="finite_t",
)
supervisor.integrate()
supervisors.append(supervisor)
from hopsflow import hopsflow
flow_hops = []
for supervisor in supervisors:
hf_sys = hopsflow.SystemParams.from_hi_params(supervisor.params)
hf_therm = hopsflow.ThermalParams.from_hi_params(supervisor.params)
flow_hops.append(
hopsflow.heat_flow_from_data(
supervisor.get_data(True),
hf_sys,
thermal_params=hf_therm,
every=1000,
save=f"flow_finite",
overwrite_cache=True
)
)
fig, ax = plt.subplots()
for params, flow, keys in zip(multi_params, flow_hops, model_keys):
pu.plot_with_σ(
params.IntP.t,
-1 * flow.for_bath(0),
ax=ax,
label=f"{len(params.SysP.g[0])}, {params.HiP.truncation_scheme.kmax}",
)
ax.legend()
from hopsflow import gaussflow as gf
exact_flows = []
for params, keys in zip(multi_params, model_keys):
α = gf.BCF(
params.IntP.t_max,
hops.util.bcf.OhmicBCF_nonZeroTemp(
s=1,
eta=params.SysP.bcf_scale[0],
w_c=keys["wc"],
beta=1 / params.SysP.__non_key__["T"][0],
),
num_terms=6,
resolution=0.01,
)
α_0 = gf.BCF(
params.IntP.t_max,
# factors=params.SysP.g[0],
# exponents=params.SysP.w[0],
hops.util.bcf.OhmicBCF_zeroTemp(
s=1,
eta=params.SysP.bcf_scale[0],
w_c=keys["wc"],
normed=True,
),
# hops_bcf,
num_terms=6,
resolution=0.001,
)
α_0_dot = gf.BCF(
params.IntP.t_max,
factors=-α_0.factors * α_0.exponents,
exponents=α_0.exponents,
# lambda t: 2 / (1j * np.pi) * (wc / (1 + 1j * wc * t)) ** 3,
# num_terms=6,
# resolution=0.001,
)
sys = gf.SystemParams(Ω=1, η=1, α_0=α_0)
flow = gf.Flow(sys, α, α_0_dot, n=0)
flow_τ = flow(params.IntP.t)
exact_flows.append(flow_τ)
_tt = np.linspace(0, 4, 1000)
fig, ax = pu.plot_complex(_tt, α(_tt)- α.approx(_tt), label="finite")
pu.plot_complex(_tt, α_0(_tt)- α_0.approx(_tt), label="finite", ax=ax)
# ut.plot_complex(result.τ, α(result.τ)- α.approx(result.τ), label="finite", ax=ax)
# ut.plot_complex(result.τ, α_0_dot(result.τ)- α_0_dot.approx(result.τ), label="finite", ax=ax)
fig, ax = plt.subplots()
for params, flow, ex_flow, keys in zip(multi_params, flow_hops, exact_flows, model_keys):
consistency = (flow).consistency(ex_flow)
pu.plot_with_σ(
params.IntP.t,
-1 * flow,
bath=0,
ax=ax,
label=rf"$α(0)={params.SysP.g[0].sum().real:.2f}$ $ω_c={keys['wc']}$ ${consistency}\%$",
)
ax.plot(params.IntP.t, ex_flow, linestyle="dotted", color="black")
plt.xlabel("$τ$")
plt.ylabel("$-J$")
# ax.plot(multi_params[-1].IntP.t, flow_τ, label="Analytic", linestyle="dotted", color="black")
ax.legend()
# fs.export_fig("convergence_zero")
pu.plot_convergence(
multi_params[-1].IntP.t,
-1 * flow_hops[-1].for_bath(0)
-flow_τ,
reference=np.zeros_like(flow_τ)
)
plt.axhline(0)
plt.legend()
# plt.plot(multi_params[-1].IntP.t, -1 * flow_hops[-1].for_bath(0).σ)

View file

@ -0,0 +1,324 @@
#+PROPERTY: header-args :session comparison_hops_gauss_new_thermal :kernel python :pandoc yes :async yes :tangle finite_t.py
#+begin_src jupyter-python :results none
import figsaver as fs
import plot_utils as pu
#+end_src
#+begin_src jupyter-python
import hops.util.bcf
import hops.util.bcf_fits
from hops.core.hierarchy_parameters import *
from stocproc import StocProc_FFT, StocProc_TanhSinh
from hops.util.abstract_truncation_scheme import TruncationScheme_Simplex
from hops.core.integration import HOPSSupervisor
import qutip
def ho_zero(max_HO_level=15, wc=2, s=1, bcf_terms=7, k_max=5, bcf_scale=1, T=1):
hops_bcf = hops.util.bcf.OhmicBCF_zeroTemp(
s,
1,
wc,
normed=True
)
g, w = hops_bcf.exponential_coefficients(bcf_terms)
integration = IntP(t_steps=(40, 1000))
q = qutip.operators.create(max_HO_level) + qutip.operators.destroy(max_HO_level)
p = (
qutip.operators.destroy(max_HO_level) - qutip.operators.create(max_HO_level)
) / 1j
system = SysP(
H_sys=0.25 * (p ** 2 + q ** 2).full(),
L=[0.5 * q.full()],
psi0=qutip.states.fock(max_HO_level, n=0).full().flatten(),
g=[g],
w=[w],
bcf_scale=[bcf_scale],
T=[.5],
description="A single HO at zero Temperature",
)
params = HIParams(
SysP=system,
IntP=integration,
HiP=HiP(
seed=0,
nonlinear=True,
result_type=ResultType.ZEROTH_AND_FIRST_ORDER,
truncation_scheme=TruncationScheme_Simplex(k_max),
save_therm_rng_seed=True,
auto_normalize=True,
accum_only=False,
rand_skip=None,
),
Eta=[
StocProc_FFT(
spectral_density=hops.util.bcf.OhmicSD_zeroTemp(
s,
1,
wc,
normed=True
),
alpha=hops_bcf,
t_max=integration.t_max,
intgr_tol=1e-5,
intpl_tol=1e-5,
negative_frequencies=False,
)
],
EtaTherm=[StocProc_TanhSinh(
spectral_density=hops.util.bcf.Ohmic_StochasticPotentialDensity(
s,
1,
wc,
beta=1/system.__non_key__["T"][0],
normed=True
),
alpha=hops.util.bcf.Ohmic_StochasticPotentialCorrelations(
s,
1,
wc,
beta=1/system.__non_key__["T"][0],
normed=True
),
t_max=integration.t_max,
intgr_tol=1e-5,
intpl_tol=1e-5,
negative_frequencies=False,
)],
)
return params
model_keys = [dict(bcf_scale=1, wc=ω_c) for ω_c in [1, 2]] + [dict(bcf_scale=s, wc=2) for s in [.5, 1, 2]]
multi_params = [ho_zero(**keys) for keys in model_keys]
multi_params = multi_params
#+end_src
#+RESULTS:
#+begin_example
stocproc.stocproc - INFO - Loaded instance from cache.
INFO:stocproc.stocproc:Loaded instance from cache.
stocproc.stocproc - INFO - non neg freq only
INFO:stocproc.stocproc:non neg freq only
stocproc.stocproc - INFO - get_dt_for_accurate_interpolation, please wait ...
INFO:stocproc.stocproc:get_dt_for_accurate_interpolation, please wait ...
stocproc.stocproc - INFO - requires dt < 1.953e-02
INFO:stocproc.stocproc:requires dt < 1.953e-02
stocproc.stocproc - INFO - yields N = 2049 (time domain)
INFO:stocproc.stocproc:yields N = 2049 (time domain)
stocproc.stocproc - INFO - find accurate discretisation in frequency domain
INFO:stocproc.stocproc:find accurate discretisation in frequency domain
stocproc.stocproc - INFO - wmax:5.151618328961643
INFO:stocproc.stocproc:wmax:5.151618328961643
stocproc.stocproc - INFO - Loaded instance from cache.
INFO:stocproc.stocproc:Loaded instance from cache.
stocproc.stocproc - INFO - non neg freq only
INFO:stocproc.stocproc:non neg freq only
stocproc.stocproc - INFO - get_dt_for_accurate_interpolation, please wait ...
INFO:stocproc.stocproc:get_dt_for_accurate_interpolation, please wait ...
stocproc.stocproc - INFO - requires dt < 9.766e-03
INFO:stocproc.stocproc:requires dt < 9.766e-03
stocproc.stocproc - INFO - yields N = 4097 (time domain)
INFO:stocproc.stocproc:yields N = 4097 (time domain)
stocproc.stocproc - INFO - find accurate discretisation in frequency domain
INFO:stocproc.stocproc:find accurate discretisation in frequency domain
stocproc.stocproc - INFO - wmax:5.665445489514328
INFO:stocproc.stocproc:wmax:5.665445489514328
stocproc.stocproc - INFO - Loaded instance from cache.
INFO:stocproc.stocproc:Loaded instance from cache.
stocproc.stocproc - INFO - Loaded instance from cache.
INFO:stocproc.stocproc:Loaded instance from cache.
stocproc.stocproc - INFO - Loaded instance from cache.
INFO:stocproc.stocproc:Loaded instance from cache.
stocproc.stocproc - INFO - Loaded instance from cache.
INFO:stocproc.stocproc:Loaded instance from cache.
stocproc.stocproc - INFO - Loaded instance from cache.
INFO:stocproc.stocproc:Loaded instance from cache.
stocproc.stocproc - INFO - Loaded instance from cache.
INFO:stocproc.stocproc:Loaded instance from cache.
#+end_example
#+begin_src jupyter-python
import ray
ray.shutdown()
ray.init()
#+end_src
#+RESULTS:
: RayContext(dashboard_url='', python_version='3.9.13', ray_version='1.13.0', ray_commit='e4ce38d001dbbe09cd21c497fedd03d692b2be3e', address_info={'node_ip_address': '192.168.100.160', 'raylet_ip_address': '192.168.100.160', 'redis_address': None, 'object_store_address': '/tmp/ray/session_2022-07-22_18-47-59_113889_118691/sockets/plasma_store', 'raylet_socket_name': '/tmp/ray/session_2022-07-22_18-47-59_113889_118691/sockets/raylet', 'webui_url': '', 'session_dir': '/tmp/ray/session_2022-07-22_18-47-59_113889_118691', 'metrics_export_port': 59610, 'gcs_address': '192.168.100.160:46332', 'address': '192.168.100.160:46332', 'node_id': '4145ec03d17533dce23caf24ec33bc176b6ba84fd6e6969da8d61160'})
#+begin_src jupyter-python
supervisors = []
for params in multi_params:
supervisor = HOPSSupervisor(
params,
10_000,
data_path="ho_data",
data_name="finite_t",
)
supervisor.integrate()
supervisors.append(supervisor)
#+end_src
#+RESULTS:
: 0it [00:00, ?it/s]
* Flow
#+begin_src jupyter-python :results none
from hopsflow import hopsflow
#+end_src
#+begin_src jupyter-python
flow_hops = []
for supervisor in supervisors:
hf_sys = hopsflow.SystemParams.from_hi_params(supervisor.params)
hf_therm = hopsflow.ThermalParams.from_hi_params(supervisor.params)
flow_hops.append(
hopsflow.heat_flow_from_data(
supervisor.get_data(True),
hf_sys,
thermal_params=hf_therm,
every=1000,
save=f"flow_finite",
overwrite_cache=True
)
)
#+end_src
#+RESULTS:
: Loading: 100% 22/22 [00:08<00:00, 2.54it/s]
#+begin_src jupyter-python
fig, ax = plt.subplots()
for params, flow, keys in zip(multi_params, flow_hops, model_keys):
pu.plot_with_σ(
params.IntP.t,
-1 * flow.for_bath(0),
ax=ax,
label=f"{len(params.SysP.g[0])}, {params.HiP.truncation_scheme.kmax}",
)
ax.legend()
#+end_src
#+RESULTS:
:RESULTS:
: <matplotlib.legend.Legend at 0x7f37bc65ceb0>
[[file:./.ob-jupyter/819c4a63463ef7a5ad06547f2423aa8920c8a021.svg]]
:END:
* Analytic
#+begin_src jupyter-python :results none
from hopsflow import gaussflow as gf
#+end_src
#+begin_src jupyter-python :results none
exact_flows = []
for params, keys in zip(multi_params, model_keys):
α = gf.BCF(
params.IntP.t_max,
hops.util.bcf.OhmicBCF_nonZeroTemp(
s=1,
eta=params.SysP.bcf_scale[0],
w_c=keys["wc"],
beta=1 / params.SysP.__non_key__["T"][0],
),
num_terms=6,
resolution=0.01,
)
α_0 = gf.BCF(
params.IntP.t_max,
# factors=params.SysP.g[0],
# exponents=params.SysP.w[0],
hops.util.bcf.OhmicBCF_zeroTemp(
s=1,
eta=params.SysP.bcf_scale[0],
w_c=keys["wc"],
normed=True,
),
# hops_bcf,
num_terms=6,
resolution=0.001,
)
α_0_dot = gf.BCF(
params.IntP.t_max,
factors=-α_0.factors * α_0.exponents,
exponents=α_0.exponents,
# lambda t: 2 / (1j * np.pi) * (wc / (1 + 1j * wc * t)) ** 3,
# num_terms=6,
# resolution=0.001,
)
sys = gf.SystemParams(Ω=1, η=1, α_0=α_0)
flow = gf.Flow(sys, α, α_0_dot, n=0)
flow_τ = flow(params.IntP.t)
exact_flows.append(flow_τ)
#+end_src
#+begin_src jupyter-python
_tt = np.linspace(0, 4, 1000)
fig, ax = pu.plot_complex(_tt, α(_tt)- α.approx(_tt), label="finite")
pu.plot_complex(_tt, α_0(_tt)- α_0.approx(_tt), label="finite", ax=ax)
# ut.plot_complex(result.τ, α(result.τ)- α.approx(result.τ), label="finite", ax=ax)
# ut.plot_complex(result.τ, α_0_dot(result.τ)- α_0_dot.approx(result.τ), label="finite", ax=ax)
#+end_src
#+RESULTS:
:RESULTS:
| hline | <AxesSubplot:> |
[[file:./.ob-jupyter/add7bd9b1175ad78784e79b3e3f6f12950eb9e3e.svg]]
:END:
#+begin_src jupyter-python
fig, ax = plt.subplots()
for params, flow, ex_flow, keys in zip(multi_params, flow_hops, exact_flows, model_keys):
consistency = (flow).consistency(ex_flow)
pu.plot_with_σ(
params.IntP.t,
-1 * flow,
bath=0,
ax=ax,
label=rf"$α(0)={params.SysP.g[0].sum().real:.2f}$ $ω_c={keys['wc']}$ ${consistency}\%$",
)
ax.plot(params.IntP.t, ex_flow, linestyle="dotted", color="black")
plt.xlabel("$τ$")
plt.ylabel("$-J$")
# ax.plot(multi_params[-1].IntP.t, flow_τ, label="Analytic", linestyle="dotted", color="black")
ax.legend()
# fs.export_fig("convergence_zero")
#+end_src
#+RESULTS:
:RESULTS:
: <matplotlib.legend.Legend at 0x7f37bc65c340>
[[file:./.ob-jupyter/5b1a80f6d49139558c8e638123e7da40d2b40504.svg]]
:END:
#+begin_src jupyter-python
pu.plot_convergence(
multi_params[-1].IntP.t,
-1 * flow_hops[-1].for_bath(0)
-flow_τ,
reference=np.zeros_like(flow_τ)
)
plt.axhline(0)
plt.legend()
# plt.plot(multi_params[-1].IntP.t, -1 * flow_hops[-1].for_bath(0).σ)
#+end_src
#+RESULTS:
:RESULTS:
: <matplotlib.legend.Legend at 0x7f37bc4c7a60>
[[file:./.ob-jupyter/b143bfbf41586019c77aa024349ade7ff6cc3b81.svg]]
:END: