SunnySuite / Sunny.jl

Spin dynamics and generalization to SU(N) coherent states
Other
86 stars 19 forks source link

Estimate variance on classical SF calculations #154

Closed ddahlbom closed 1 year ago

ddahlbom commented 1 year ago

The only new thing that this PR publicly exposes is the keyword calculate_errors to dynamical_correlations and instant_correlations. However, I have not documented it because I don't want to make it official yet.

If calculate_errors is set to true, Sunny will use a modified Welford algorithm to iteratively estimate the complex variance on each matrix element of S(q,w). The exact approach is documented in a lyx file. This variance is recorded into an array called variance in SampledCorrelations. Its type is Union{Nothing, Array{Float64, 7}} and it is only allocated if the keyword is set to true.

The provisional (non-exposed) way to retrieve error estimates is by creating a error_formula, which is simply an intensity_formula that sets a flag when creating a ClassicalIntensityFormula. This flag in turn ensures that the calc_intensity closure uses the correct buffer and correct basis reduction function. Error information is then retrieved with errors_interpolated(sc, qs, err_formula).

Adding full error propagation through intensities calculation opens up the possibility of more complexity (e.g., potentially parallel contraction functions), and I don't really want to confront that at this point. There are other higher priority issues. But I'd like to get this merged so I can save structure factor data with variance estimates. This will get me unstuck in another project. We can reexamine the interface and error propagation with more care at a later point.

@Lazersmoke I don't think this affects you very much, but I did remove the copybuf from SampledCorrelations. I believe you added this in one of your initial optimization passes. Writing the relevant loop explicitly and storing intermediate values on the stack seems to be slightly faster, and it saves memory. accum_sample! still has some small allocations due to the iteration on the SortedDict (this hasn't changed), but this doesn't seem to be performance sensitive. (I hardcoded a particular case that removes the iteration on the sorteddict, and it wasn't noticeably faster even though it eliminated the allocations.)

I'm letting you know in case I overlooked some aspect of this.

kbarros commented 1 year ago

Looks good to me. The new features are isolated and not exposed to the user, so this seems safe.

ddahlbom commented 1 year ago

Closes #148