SunnySuite / Sunny.jl

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

Fix sign of momentum transfer #271

Closed kbarros closed 1 month ago

kbarros commented 1 month ago

This PR fixes the direction of momentum transfer for structure factor intensities calculated using SampledCorrelations. The issue is documented in #270.

The PR additionally unifies notation, adds test, and adds a new tutorial to illustrate the Sunny sign convention for momentum transfer.

(Early versions of this PR also attempted to swap the sign of $q$ in LSWT calculations, but that was erroneous and has been reverted.)

Lazersmoke commented 1 month ago

An example script to exercise the sublattice part of this:

Example script ```julia using Sunny, GLMakie fig = Figure() for type = [:bravais,:sublattice] D = 1 B = 1 if type == :bravais latvecs = lattice_vectors(2, 2, 1, 90, 90, 90) cryst = Crystal(latvecs, [[0,0,0]], "P1") sys = System(cryst, (1,1,32), [SpinInfo(1,S=1,g=1)], :dipole; units=Units.theory) set_exchange!(sys, dmvec([0, 0, D]), Bond(1, 1, [0, 0, 1])) elseif type == :sublattice latvecs = lattice_vectors(2, 2, 32, 90, 90, 90) cryst = Crystal(latvecs, [[0,0,x/32] for x = 0:31], "P1") sys = System(cryst, (1,1,1), [SpinInfo(i,S=1,g=1) for i = 1:32], :dipole; units=Units.theory) for i = 1:32 set_exchange!(sys, dmvec([0, 0, D]), Bond(i, mod1(i+1,32), [0, 0, i == 32 ? 1 : 0])) end end set_external_field!(sys, [0, 0, B]) randomize_spins!(sys) minimize_energy!(sys) display(plot_spins(sys)) @assert energy_per_site(sys) ≈ -5/4 density = 400 Qs = (type == :bravais ? 1 : 32) * [[0,0,-1/2], [0,0,+1/2]] path, xticks = reciprocal_space_path(cryst, Qs, density) swt = SpinWaveTheory(sys) formula = intensity_formula(swt, :trace; kernel=gaussian(; fwhm=0.1)) energy_max = 4 energies_swt = 0:0.02:energy_max intens_swt = intensities_broadened(swt, path, energies_swt, formula) dt = 0.05 kT = 0.02 damping = 0.1 langevin = Langevin(dt; kT, damping) suggest_timestep(sys, langevin; tol=1e-2) for _ in 1:10_000 step!(sys, langevin) end sc = dynamical_correlations(sys; dt, nω=400, ωmax=energy_max) add_sample!(sc, sys) nsamples = 20 for _ in 1:nsamples for _ in 1:1_000 step!(sys, langevin) end add_sample!(sc, sys) end intens_classical = intensities_interpolated(sc, path, intensity_formula(sc, :trace; kT)) energies_classical = available_energies(sc); row = type == :bravais ? 1 : 2 scale = type == :bravais ? 1 : 32 ax = Axis(fig[row,1]; aspect=1.4, ylabel="ω (meV)", xlabel="𝐪 (r.l.u.)", xticks, xticklabelrotation=π/10, title = type == :bravais ? "Bravais SWT" : "Sublattice SWT") heatmap!(ax, 1:size(intens_swt, 1), energies_swt, intens_swt, colormap=:gnuplot2, colorrange=(0, 10 * scale)) ax = Axis(fig[row,2]; aspect=1.4, ylabel="ω (meV)", xlabel="𝐪 (r.l.u.)", xticks, xticklabelrotation=π/10, title = type == :bravais ? "Bravais Classical" : "Sublattice Classical") heatmap!(ax, 1:size(intens_classical, 1), energies_classical, intens_classical, colormap=:gnuplot2, colorrange=(0, 0.2 * scale)) end fig ```

(N.B. it fails the @assert stochastically, sorry about this, just need to rerun until it works)

kbarros commented 1 month ago

In addition to the new tutorial 8 and Sam's example above, here is a script that shows agreement between LSWT and SampledCorrelations. Passes on both this branch, and (for accidental reasons) on the release Sunny 0.5.11.

Example script ```jl using Sunny, GLMakie latvecs = lattice_vectors(2, 2, 1, 90, 90, 90) cryst = Crystal(latvecs, [[0,0,0]], "P1") sys = System(cryst, (1,1,4), [SpinInfo(1,S=1,g=1)], :dipole; units=Units.theory) D = 1 B = 1 set_exchange!(sys, dmvec([0, 0, D]), Bond(1, 1, [0, 0, 1])) set_external_field!(sys, [0, 0, B]) randomize_spins!(sys) minimize_energy!(sys) @assert energy_per_site(sys) ≈ -5/4 density = 400 path, xticks = reciprocal_space_path(cryst, [[0,0,-1/2], [0,0,+1/2]], density) swt = SpinWaveTheory(sys) formula = intensity_formula(swt, :trace; kernel=gaussian(; fwhm=0.1)) energy_max = 4 energies = 0:0.02:energy_max intens = intensities_broadened(swt, path, energies, formula) sys2 = resize_supercell(sys, (1,1,64)) @assert energy_per_site(sys2) ≈ -5/4 dt = 0.05 kT = 0.02 damping = 0.1 langevin = Langevin(dt; kT, damping) suggest_timestep(sys2, langevin; tol=1e-2) for _ in 1:10_000 step!(sys, langevin) end sc = dynamical_correlations(sys2; dt, nω=400, ωmax=energy_max) add_sample!(sc, sys2) nsamples = 100 for _ in 1:nsamples for _ in 1:1_000 step!(sys2, langevin) end add_sample!(sc, sys2) end intens2 = intensities_interpolated(sc, path, intensity_formula(sc, :trace; kT)) energies2 = available_energies(sc); fig = Figure() ax = Axis(fig[1,1]; aspect=1.4, ylabel="ω (meV)", xlabel="𝐪 (r.l.u.)", xticks, xticklabelrotation=π/10) heatmap!(ax, 1:size(intens, 1), energies, intens, colormap=:gnuplot2, colorrange=(0, 10)) ax = Axis(fig[1,2]; aspect=1.4, ylabel="ω (meV)", xlabel="𝐪 (r.l.u.)", xticks, xticklabelrotation=π/10) heatmap!(ax, 1:size(intens2, 1), energies2, intens2, colormap=:gnuplot2, colorrange=(0, 0.2)) fig ```
image