fix cuttoff dep

python/energy_flow_proper/08_dynamic_one_bath/cutoff_dep.py

@@ -15,6 +15,60 @@ from hops.util.logging_setup import logging_setup
 import logging
 logging_setup(logging.INFO, show_stocproc=False)
 
+import scipy
+import pathlib
+
+ωs = np.linspace(.5, 2.1, 15)
+ω_models = []
+Δ = 5
+
+for ω_c in ωs:
+    proto, strobe_t, strobe_indices = es.energy_shovel(
+        Δ,
+        periods=4,
+        modulate_system=False,
+        δ=0.1,
+        ω_c = ω_c,
+    )
+
+    proto.k_max = 5
+    proto.therm_method = "fft"
+
+    goal = 0.4
+
+
+    def cost(δ):
+        proto.δ = δ
+        aux.integrate(proto, 500)
+        with aux.get_data(proto) as data:
+            inter = proto.interaction_energy(data)
+
+        return (abs(inter.value).max() - goal) ** 2
+
+    cache_path = pathlib.Path(f"./.cache/{proto.hexhash}_delta_{goal}.npy")
+    if cache_path.exists():
+        with cache_path.open("rb") as cache:
+            δ = np.load(cache, allow_pickle=True)
+    else:
+        δ = scipy.optimize.minimize_scalar(cost,  bounds=(.1, 3), method="bounded", options=dict(xatol=1e-2, disp=3)).x
+        cache_path.parent.mkdir(parents=True, exist_ok=True)
+        with cache_path.open("wb") as cache:
+            np.save(cache, δ)
+
+    final, strobe_t, strobe_indices = es.energy_shovel(
+        Δ,
+        periods=50,
+        modulate_system=False,
+        δ=δ,
+        ω_c=ω_c
+    )
+    final.k_max = 5
+    final.therm_method = "fft"
+
+    ω_models.append(final.copy())
+
+aux.integrate_multi(ω_models, 10_000)
+
 print(es.models_table(ω_models))
 
 plt.plot(ωs, [model.full_thermal_bcf(0).real * model.bcf_scale for model in ω_models])