Skip to content

Commit

Permalink
Add SpinW tutorials 12, 13, 14 (#296)
Browse files Browse the repository at this point in the history
Co-authored-by: BhushanThipe <[email protected]>
  • Loading branch information
kbarros and 0ParadoX1 authored Aug 23, 2024
1 parent 2ff76ab commit 376872e
Show file tree
Hide file tree
Showing 10 changed files with 243 additions and 12 deletions.
2 changes: 1 addition & 1 deletion examples/02_LSWT_CoRh2O4.jl
Original file line number Diff line number Diff line change
Expand Up @@ -74,7 +74,7 @@ swt = SpinWaveTheory(sys_prim; measure=ssf_perp(sys_prim))
# path that connects high-symmetry points in reciprocal space.

qs = [[0, 0, 0], [1/2, 0, 0], [1/2, 1/2, 0], [0, 0, 0]]
path = q_space_path(cryst, qs, 400)
path = q_space_path(cryst, qs, 500)

# Select [`lorentzian`](@ref) broadening with a full-width at half-maximum
# (FWHM) of 0.8 meV. Use [`ssf_perp`](@ref) to calculate unpolarized scattering
Expand Down
2 changes: 1 addition & 1 deletion examples/03_LLD_CoRh2O4.jl
Original file line number Diff line number Diff line change
Expand Up @@ -138,7 +138,7 @@ qs = [[3/4, 3/4, 0],
[1/4, 1, 1/4],
[ 0, 1, 0],
[ 0, -4, 0]]
qpts = q_space_path(cryst, qs, 1000)
qpts = q_space_path(cryst, qs, 500)

# Calculate ``I(𝐪, ω)`` intensities along this path and plot.

Expand Down
2 changes: 1 addition & 1 deletion examples/07_Dipole_Dipole.jl
Original file line number Diff line number Diff line change
Expand Up @@ -38,7 +38,7 @@ plot_spins(sys_prim; ghost_radius=8, color=[:red, :blue, :yellow, :purple])

qs = [[0,0,0], [0,1,0], [1,1/2,0], [1/2,1/2,1/2], [3/4,3/4,0], [0,0,0]]
labels = ["Γ", "X", "W", "L", "K", "Γ"]
path = q_space_path(cryst, qs, 400; labels)
path = q_space_path(cryst, qs, 500; labels)

measure = ssf_trace(sys_prim)
swt = SpinWaveTheory(sys_prim; measure)
Expand Down
3 changes: 1 addition & 2 deletions examples/08_Momentum_Conventions.jl
Original file line number Diff line number Diff line change
Expand Up @@ -66,7 +66,7 @@ for _ in 1:nsamples
end
add_sample!(sc, sys)
end
path = q_space_path(cryst, [[0,0,-1/2], [0,0,+1/2]], 300)
path = q_space_path(cryst, [[0,0,-1/2], [0,0,+1/2]], 400)
res1 = intensities(sc, path; energies=:available, kT)

# Calculate the same quantity with linear spin wave theory at ``T = 0``. Because
Expand All @@ -76,7 +76,6 @@ res1 = intensities(sc, path; energies=:available, kT)
sys_small = resize_supercell(sys, (1,1,1))
minimize_energy!(sys_small)
swt = SpinWaveTheory(sys_small; measure=ssf_trace(sys_small))
path = q_space_path(cryst, [[0,0,-1/2], [0,0,+1/2]], 400)
res2 = intensities_bands(swt, path)

# This model system has a single magnon band with dispersion ``ϵ(𝐪) = 1 - D/B
Expand Down
10 changes: 6 additions & 4 deletions examples/spinw_tutorials/SW11_La2CuO4.jl
Original file line number Diff line number Diff line change
Expand Up @@ -2,7 +2,7 @@
# # SW11 - La₂CuO₄
#
# This is a Sunny port of [SpinW Tutorial
# 11](https://spinw.org/tutorials/15tutorial), originally authored by Sandor
# 11](https://spinw.org/tutorials/11tutorial), originally authored by Sandor
# Toth. It calculates the spin wave spectrum of La₂CuO₄.

# Load packages
Expand All @@ -14,7 +14,9 @@ using Sunny, GLMakie

units = Units(:meV, :angstrom)
latvecs = lattice_vectors(1, 1, 10, 90, 90, 90)
cryst = Crystal(latvecs, [[0, 0, 0]])
positions = [[0, 0, 0]]
types = ["Cu"]
cryst = Crystal(latvecs, positions; types)
view_crystal(cryst; dims=2)

# Build a spin system using the exchange parameters from [R. Coldea, Phys. Rev.
Expand Down Expand Up @@ -45,9 +47,9 @@ energies = range(0, 320, 400)
swt = SpinWaveTheory(sys; measure=ssf_perp(sys))
res = intensities(swt, path; energies, kernel=gaussian(fwhm=35))
res.energies .*= 1.18
plot_intensities(res)
plot_intensities(res; units)

# Plot instantaneous itensities, integrated over ω.

res = intensities_instant(swt, path)
plot_intensities(res; colorrange=(0,20))
plot_intensities(res; colorrange=(0,20), units)
67 changes: 67 additions & 0 deletions examples/spinw_tutorials/SW12_Triangular_easy_plane.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,67 @@
# # SW12 - Triangular lattice with easy plane
#
# This is a Sunny port of [SpinW Tutorial
# 12](https://spinw.org/tutorials/12tutorial), originally authored by Sandor
# Toth. It calculates the spin wave dispersion of a triangular lattice model
# with antiferromagnetic interactions and easy-plane single-ion anisotropy.

# Load packages

using Sunny, GLMakie

# Build a triangular lattice with arbitrary lattice constant of 3 Å.

latvecs = lattice_vectors(3, 3, 4, 90, 90, 120)
cryst = Crystal(latvecs, [[0, 0, 0]])

# Build a system with exchange +1 meV along nearest neighbor bonds.

S = 3/2
J1 = +1.0
sys = System(cryst, (3,3,1), [SpinInfo(1; S, g=2)], :dipole)
set_exchange!(sys, J1, Bond(1, 1, [1, 0, 0]))

# Set an easy-axis anisotropy operator ``+D S_z^2`` using
# [`set_onsite_coupling!`](@ref). **Important note**: When introducing a
# single-ion anisotropy in `:dipole` mode, Sunny will automatically include a
# classical-to-quantum correction factor, as described in the document
# [Interaction Strength Renormalization](@ref). The present single-ion
# anisotropy is quadratic in the spin operators, so the prescription is to
# rescale the interaction strength as ``D → (1 - 1/2S) D``. For purposes of this
# tutorial, our aim is to reproduce the SpinW result, which requires to "undo"
# Sunny's classical-to-quantum rescaling factor. Another way to achieve the same
# thing is to select system mode `:dipole_large_S` instead of `:dipole`.

undo_classical_to_quantum_rescaling = 1 / (1 - 1/2S)
D = 0.2 * undo_classical_to_quantum_rescaling
set_onsite_coupling!(sys, S -> D*S[3]^2, 1)
randomize_spins!(sys)
minimize_energy!(sys)
plot_spins(sys; dims=2)

# Plot the spin wave spectrum for a path through ``𝐪``-space.

qs = [[0, 0, 0], [1, 1, 0]]
path = q_space_path(cryst, qs, 400)
swt = SpinWaveTheory(sys; measure=ssf_perp(sys))
res = intensities_bands(swt, path)
plot_intensities(res)

# To select a specific linear combination of spin structure factor (SSF)
# components in global Cartesian coordinates, one can use [`ssf_custom`](@ref).
# Here we calculate and plot the real part of ``\mathcal{S}^{zz}(𝐪, ω)``.

measure = ssf_custom(sys) do q, ssf
return real(ssf[3, 3])
end
swt = SpinWaveTheory(sys; measure)
res = intensities_bands(swt, path)
plot_intensities(res)

# It's also possible to get data for the full 3×3 SSF. For example, this is the
# SSF for the 7th energy band, at the 10th ``𝐪``-point along the path.

measure = ssf_custom((q, ssf) -> ssf, sys)
swt = SpinWaveTheory(sys; measure)
res = intensities_bands(swt, path)
res.data[7, 10]
94 changes: 94 additions & 0 deletions examples/spinw_tutorials/SW13_LiNiPO4.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,94 @@
# # SW13 - LiNiPO₄
#
# This is a Sunny port of [SpinW Tutorial
# 13](https://spinw.org/tutorials/13tutorial), originally authored by Sandor
# Toth. It calculates the spin wave spectrum of LiNiPO₄.

# Load packages

using Sunny, GLMakie

# Build an orthorhombic lattice and populate the Ni atoms according to
# spacegroup 62 (Pnma).

units = Units(:meV, :angstrom)
a = 10.02
b = 5.86
c = 4.68
latvecs = lattice_vectors(a, b, c, 90, 90, 90)
positions = [[1/4, 1/4, 0]]
types = ["Ni"]
cryst = Crystal(latvecs, positions, 62, setting=""; types)
view_crystal(cryst)

# Create a system with exchange parameters taken from [T. Jensen, et al., PRB
# **79**, 092413 (2009)](https://doi.org/10.1103/PhysRevB.79.092413). The
# corrected anisotropy values are taken from the thesis of T. Jensen. The mode
# `:dipole_large_S` avoids a [classical-to-quantum rescaling factor](@ref
# "Interaction Strength Renormalization") of anisotropy strengths, as needed for
# consistency with the original fits.

S = 3/2
sys = System(cryst, (1,1,1), [SpinInfo(1, S=1, g=2)], :dipole_large_S)
Jbc = 1.036
Jb = 0.6701
Jc = -0.0469
Jac = -0.1121
Jab = 0.2977
Da = 0.1969
Db = 0.9097
set_exchange!(sys, Jbc, Bond(2, 3, [0, 0, 0]))
set_exchange!(sys, Jc, Bond(1, 1, [0, 0, -1]))
set_exchange!(sys, Jb, Bond(1, 1, [0, 1, 0]))
set_exchange!(sys, Jab, Bond(1, 2, [0, 0, 0]))
set_exchange!(sys, Jab, Bond(3, 4, [0, 0, 0]))
set_exchange!(sys, Jac, Bond(3, 1, [0, 0, 0]))
set_exchange!(sys, Jac, Bond(4, 2, [0, 0, 0]))
set_onsite_coupling!(sys, S -> Da*S[1]^2 + Db*S[2]^2, 1)

# Energy minimization yields a co-linear order along the ``c`` axis.

randomize_spins!(sys)
minimize_energy!(sys)
plot_spins(sys; color=[s[3] for s in sys.dipoles])

# Calculate the spectrum along path [ξ, 1, 0]

swt = SpinWaveTheory(sys; measure=ssf_perp(sys))
qs = [[0, 1, 0], [2, 1, 0]]
path = q_space_path(cryst, qs, 400)
res = intensities_bands(swt, path)
fig = Figure(size=(768, 300))
plot_intensities!(fig[1, 1], res; units);

# There are two physical bands with nonvanishing intensity. To extract these
# intensity curves, we must filter out the additional bands with zero intensity.
# One way is to sort the data along dimension 1 (the band index) with the
# comparison operator "is intensity less than ``10^{-12}``". Doing so moves the
# physical bands to the end of the array axis. Call the [Makie `lines!`
# function](https://docs.makie.org/stable/reference/plots/lines) to make a
# custom plot.

data_sorted = sort(res.data; dims=1, by= >(1e-12))
ax = Axis(fig[1, 2], xlabel="Momentum (r.l.u.)", ylabel="Intensity",
xticks=res.qpts.xticks, xticklabelrotation=π/6)
lines!(ax, data_sorted[end, :]; label="Lower band")
lines!(ax, data_sorted[end-1, :]; label="Upper band")
axislegend(ax)
fig

# Make the same plots along path [0, 1, ξ]

qs = [[0, 1, 0], [0, 1, 2]]
path = q_space_path(cryst, qs, 400)
res = intensities_bands(swt, path)
fig = Figure(size=(768, 300))
plot_intensities!(fig[1, 1], res; units)

data_sorted = sort(res.data; dims=1, by=x->abs(x)>1e-12)
ax = Axis(fig[1, 2], xlabel="Momentum (r.l.u.)", ylabel="Intensity",
xticks=res.qpts.xticks, xticklabelrotation=π/6)
lines!(ax, data_sorted[end, :]; label="Lower band")
lines!(ax, data_sorted[end-1, :]; label="Upper band")
axislegend(ax)
fig
69 changes: 69 additions & 0 deletions examples/spinw_tutorials/SW14_YVO3.jl
Original file line number Diff line number Diff line change
@@ -0,0 +1,69 @@
# # SW14 - YVO₃
#
# This is a Sunny port of [SpinW Tutorial
# 14](https://spinw.org/tutorials/14tutorial), originally authored by Sandor
# Toth. It calculates the spin wave spectrum of YVO₃.

# Load packages

using Sunny, GLMakie

# Build an orthorhombic lattice and populate the V atoms according to the
# pseudocubic unit cell, doubled along the c-axis. The listed spacegroup of the
# YVO₃ chemical cell is international number 62. It has been observed, however,
# that the exchange interactions break this symmetry. For this reason, disable
# all symmetry analysis by selecting spacegroup 1 (P1).

units = Units(:meV, :angstrom)
a = 5.2821 / sqrt(2)
b = 5.6144 / sqrt(2)
c = 7.5283
latvecs = lattice_vectors(a, b, c, 90, 90, 90)
positions = [[0, 0, 0], [0, 0, 1/2]]
types = ["V", "V"]
cryst = Crystal(latvecs, positions, 1; types)

# Create a system following the model of [C. Ulrich, et al. PRL **91**, 257202
# (2003)](https://doi.org/10.1103/PhysRevLett.91.257202). The mode
# `:dipole_large_S` avoids a [classical-to-quantum rescaling factor](@ref
# "Interaction Strength Renormalization") of anisotropy strengths, as needed for
# consistency with the original fits.

sys = System(cryst, (2,2,1), [SpinInfo(1, S=1/2, g=2), SpinInfo(2, S=1/2, g=2)], :dipole_large_S)
Jab = 2.6
Jc = 3.1
δ = 0.35
K1 = 0.90
K2 = 0.97
d = 1.15
Jc1 = [-Jc*(1+δ)+K2 0 -d; 0 -Jc*(1+δ) 0; +d 0 -Jc*(1+δ)]
Jc2 = [-Jc*(1-δ)+K2 0 +d; 0 -Jc*(1-δ) 0; -d 0 -Jc*(1-δ)]
set_exchange!(sys, Jab, Bond(1, 1, [1, 0, 0]))
set_exchange!(sys, Jab, Bond(2, 2, [1, 0, 0]))
set_exchange!(sys, Jc1, Bond(1, 2, [0, 0, 0]))
set_exchange!(sys, Jc2, Bond(2, 1, [0, 0, 1]))
set_exchange!(sys, Jab, Bond(1, 1, [0, 1, 0]))
set_exchange!(sys, Jab, Bond(2, 2, [0, 1, 0]))
set_onsite_coupling!(sys, S -> -K1*S[1]^2, 1)
set_onsite_coupling!(sys, S -> -K1*S[1]^2, 2)

# When using spacegroup P1, there is no symmetry-propagation of interactions
# because all bonds are considered inequivalent. One can visualize the
# interactions in the system by clicking the toggles in the
# [`view_crystal`](@ref) GUI.

view_crystal(sys)

# Energy minimization yields a Néel order with canting

randomize_spins!(sys)
minimize_energy!(sys)
plot_spins(sys)

# Plot the spin wave spectrum along a path

swt = SpinWaveTheory(sys; measure=ssf_perp(sys))
qs = [[0.75, 0.75, 0], [0.5, 0.5, 0], [0.5, 0.5, 1]]
path = q_space_path(cryst, qs, 400)
res = intensities_bands(swt, path)
plot_intensities(res; units)
4 changes: 2 additions & 2 deletions examples/spinw_tutorials/SW18_Distorted_kagome.jl
Original file line number Diff line number Diff line change
Expand Up @@ -77,7 +77,7 @@ energy_per_site(sys2) # < -0.7834 meV
# Define a path in q-space

qs = [[0,0,0], [1,0,0]]
path = q_space_path(cryst, qs, 512)
path = q_space_path(cryst, qs, 400)

# Calculate intensities for the incommensurate spiral phase using
# [`SpiralSpinWaveTheory`](@ref). It is necessary to provide the original `sys`,
Expand All @@ -93,7 +93,7 @@ plot_intensities(res; units)
radii = range(0, 2, 100) # (1/Å)
energies = range(0, 6, 200)
kernel = gaussian(fwhm=0.05)
res = powder_average(cryst, radii, 200) do qs
res = powder_average(cryst, radii, 400) do qs
intensities(swt, qs; energies, kernel)
end
plot_intensities(res; units)
2 changes: 1 addition & 1 deletion examples/spinw_tutorials/SW19_Different_Ions.jl
Original file line number Diff line number Diff line change
Expand Up @@ -39,7 +39,7 @@ plot_spins(sys)

swt = SpinWaveTheory(sys; measure=ssf_perp(sys))
qs = [[0,0,0], [1,0,0]]
path = q_space_path(cryst, qs, 512)
path = q_space_path(cryst, qs, 400)

# Plot three types of pair correlation intensities

Expand Down

0 comments on commit 376872e

Please sign in to comment.