Progress

scripts/01_renormalization_factors_for_sum_rule.jl

@@ -19,7 +19,48 @@ sim_params = (;
 )
 
 sys, cryst = FeI2_sys_and_cryst((4, 4, 4))
-
-κ = 1.070
+κ = 1.0
 kT = 0.1
-@time estimate_sum(sys, κ, kT, sim_params)
+
+@time estimate_sum(sys, κ, kT, sim_params; observables = observable_matrices())
+
+
+kTs = range(0.1, 100.0, 3)
+κs = zero(kTs)
+
+ref = 16/3
+thresh = 0.05
+global_bounds = (1.0, 2.0)
+
+for (n, kT) in enumerate(kTs)
+ println("kT = $kT")
+
+ κ0 = n == 1 ? 1.0 : κs[n-1]
+ bounds = (1.0, 5.0)
+
+ sum = estimate_sum(sys, κ0, kT, sim_params)
+
+ @time while abs(sum - ref) > thresh
+ sum = estimate_sum(sys, κ0, kT, sim_params)
+ println("\t κ=$κ0, sum=$sum")
+ if sum < ref
+ if κ > global_bounds[2]
+ κ0 = global_bounds[2]
+ continue
+ end
+ bounds = (κ0, bounds[2])
+ κ0 = (bounds[2] + κ0)/2
+ continue
+ elseif sum > ref
+ if κ < global_bounds[1]
+ κ0 = global_bounds[1]
+ continue
+ end
+ bounds = (bounds[1], κ0)
+ κ0 = (bounds[1] + κ0)/2
+ continue
+ end
+ end
+
+ κs[n] = κ0
+end

src/structure_factor_utils.jl

@@ -51,10 +51,10 @@ end
 
 # Estimate S(𝐪,ω) at temperature kT and evaluate the sum rule using both the
 # classical-to-quantum correspondence factor and moment renormalization.
-function estimate_sum(sys::System{N}, κ, kT, sim_params) where N
+function estimate_sum(sys::System{N}, κ, kT, sim_params; observables=nothing) where N
  (; nω, ωmax, Δt, Δt_therm, dur_therm, dur_decorr, nsamples, λ) = sim_params
 
- sc = dynamical_correlations(sys; Δt, nω, ωmax)
+ sc = dynamical_correlations(sys; Δt, nω, ωmax, observables)
  saved_coherents = copy(sys.coherents)
  ndecorr = round(Int, dur_decorr/Δt_therm)
  ntherm = round(Int, dur_therm/Δt_therm)