Download this example as Julia file or Jupyter notebook.

SW15 - Ba₃NbFe₃Si₂O₁₄

This is a Sunny port of SpinW Tutorial 15, originally authored by Sandor Toth. It calculates the linear spin wave theory spectrum of Ba₃NbFe₃Si₂O₁₄.

Load packages

using Sunny, GLMakie

Specify the Ba₃NbFe₃Si₂O₁₄ Crystal cell following Marty et al., Phys. Rev. Lett. 101, 247201 (2008).

units = Units(:meV, :angstrom)
a = b = 8.539 # (Å)
c = 5.2414
latvecs = lattice_vectors(a, b, c, 90, 90, 120)
types = ["Fe", "Nb", "Ba", "Si", "O", "O", "O"]
positions = [[0.24964,0,0.5], [0,0,0], [0.56598,0,0], [2/3,1/3,0.5220],
             [2/3,1/3,0.2162], [0.5259,0.7024,0.3536], [0.7840,0.9002,0.7760]]
langasite = Crystal(latvecs, positions, 150; types)
cryst = subcrystal(langasite, "Fe")
view_crystal(cryst)
Example block output

Create a System and set exchange interactions as parametrized in Loire et al., Phys. Rev. Lett. 106, 207201 (2011).

sys = System(cryst, [1 => Moment(s=5/2, g=2)], :dipole; seed=0)
J₁ = 0.85
J₂ = 0.24
J₃ = 0.053
J₄ = 0.017
J₅ = 0.24
set_exchange!(sys, J₁, Bond(3, 2, [1,1,0]))
set_exchange!(sys, J₄, Bond(1, 1, [0,0,1]))
set_exchange!(sys, J₂, Bond(1, 3, [0,0,0]))

The final two exchanges are setting according to the desired chirality $ϵ_T$ of the magnetic structure.

ϵT = -1
if ϵT == -1
    set_exchange!(sys, J₃, Bond(2, 3, [-1,-1,1]))
    set_exchange!(sys, J₅, Bond(3, 2, [1,1,1]))
elseif ϵT == 1
    set_exchange!(sys, J₅, Bond(2, 3, [-1,-1,1]))
    set_exchange!(sys, J₃, Bond(3, 2, [1,1,1]))
else
    error("Chirality must be ±1")
end

This compound is known to have a spiral order with approximate propagation wavevector $𝐤 ≈ [0, 0, 1/7]$. Search for this magnetic order with minimize_spiral_energy!. Due to reflection symmetry, one of two possible propagation wavevectors may appear, $𝐤 = ± [0, 0, 0.1426…]$. Note that $k_z = 0.1426…$ is very close to $1/7 = 0.1428…$.

axis = [0, 0, 1]
randomize_spins!(sys)
k = minimize_spiral_energy!(sys, axis)
3-element StaticArraysCore.SVector{3, Float64} with indices SOneTo(3):
 -6.00541550923142e-11
  7.384862266239884e-11
  0.14264604651041948

We can visualize the full magnetic cell using repeat_periodically_as_spiral, which includes 7 rotated copies of the chemical cell.

sys_enlarged = repeat_periodically_as_spiral(sys, (1, 1, 7); k, axis)
plot_spins(sys_enlarged; color=[S[1] for S in sys_enlarged.dipoles])
Example block output

One could perform a spin wave calculation using either SpinWaveTheory on sys_enlarged, or SpinWaveTheorySpiral on the original sys. The latter has some restrictions on the interactions, but allows for our slightly incommensurate wavevector $𝐤$.

measure = ssf_perp(sys)
swt = SpinWaveTheorySpiral(sys; measure, k, axis)
SpinWaveTheorySpiral(SpinWaveTheory(System([Dipole mode], Supercell (1×1×1)×3, Energy per site 9.601), Sunny.SWTDataDipole(StaticArraysCore.SMatrix{3, 3, Float64, 9}[[0.8526317395268106 -0.40188677348270807 0.3339253480212965; -0.30205335680855694 0.1423724260235337 0.9425995236307179; -0.42636004314683157 -0.9045535438037 2.736855799247858e-10], [0.664920841999281 0.3690821829375078 0.6493524590798527; 0.5677515137203324 0.31514573578159544 -0.7604875961427303; -0.4853230805795652 0.8743348943378396 -6.159333581702433e-11], [-0.021234304141621734 -0.1808697279412665 -0.9832778070525509; 0.11465047861443019 0.9765707929257129 -0.1821119275553544; 0.9931789224913622 -0.11660029124704963 -3.57481316600529e-10]], StaticArraysCore.SVector{3, Float64}[[1.7052634790536207, -0.8037735469654159, 0.6678506960425931] [1.3298416839985618, 0.7381643658750154, 1.2987049181597055] [-0.04246860828324342, -0.3617394558825326, -1.9665556141051002]; [-0.6041067136171134, 0.2847448520470672, 1.8851990472614348] [1.1355030274406641, 0.6302914715631904, -1.5209751922854597] [0.22930095722886024, 1.9531415858514247, -0.3642238551107087]; [-0.8527200862936629, -1.8091070876073996, 5.473711598495715e-10] [-0.9706461611591302, 1.7486697886756788, -1.2318667163404864e-10] [1.9863578449827235, -0.23320058249409914, -7.149626332010577e-10]], Sunny.StevensExpansion[Sunny.StevensExpansion(0, [0.0], [0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]), Sunny.StevensExpansion(0, [0.0], [0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0]), Sunny.StevensExpansion(0, [0.0], [0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0])], [1.5811388300841898, 1.5811388300841898, 1.5811388300841898]), MeasureSpec, 1.0e-8), [-6.00541550923142e-11, 7.384862266239884e-11, 0.14264604651041948], [0.0, 0.0, 1.0])

Calculate broadened intensities for a path $[0, 1, L]$ through reciprocal space

qs = [[0, 1, -1], [0, 1, -1+1], [0, 1, -1+2], [0, 1, -1+3]]
path = q_space_path(cryst, qs, 400)
energies = range(0, 6, 400)
res = intensities(swt, path; energies, kernel=gaussian(fwhm=0.25))
plot_intensities(res; units, saturation=0.7, colormap=:jet, title="Scattering intensities")
Example block output

Use ssf_custom_bm to calculate the imaginary part of $\mathcal{S}^{2, 3}(𝐪, ω) - \mathcal{S}^{3, 2}(𝐪, ω)$. In polarized neutron scattering, it is conventional to express the 3×3 structure factor matrix $\mathcal{S}^{α, β}(𝐪, ω)$ in the Blume-Maleev polarization axis system. Specify the scattering plane $[0, K, L]$ via the spanning vectors $𝐮 = [0, 1, 0]$ and $𝐯 = [0, 0, 1]$.

measure = ssf_custom_bm(sys; u=[0, 1, 0], v=[0, 0, 1]) do q, ssf
    imag(ssf[2,3] - ssf[3,2])
end
swt = SpinWaveTheorySpiral(sys; measure, k, axis)
res = intensities(swt, path; energies, kernel=gaussian(fwhm=0.25))
plot_intensities(res; units, saturation=0.8, allpositive=false,
                 title="Im[S²³(q, ω) - S³²(q, ω)]")
Example block output