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₁₄. The ground state is an incommensurate spiral, which can be directly studied using the functions minimize_spiral_energy! and SpinWaveTheorySpiral.

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)
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 set 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):
 -5.696178589183529e-13
 -1.0464674673019045e-12
  0.8573539534382282

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.424665049576274 -0.5445830267803297 -0.7232488663047842; -0.4447496101337116 0.5703391156996577 -0.6905874871360212; 0.7885793428402068 0.6149330207800746 -2.1049216223493393e-12], [0.13100592856269436 0.24865865999063697 0.9596907405477492; 0.44732799078563457 0.8490606492782332 -0.2810581479105578; -0.8847231857146696 0.4661168143919359 -3.2079322177671286e-13], [0.9698207691111254 -0.059522398347763884 -0.23644187424319718; 0.23599780831390804 -0.01448427998444506 0.971645635046216; -0.06125936884867228 -0.998121881199016 1.1242730670753908e-12]], StaticArraysCore.SVector{3, Float64}[[0.8493300991525476, -1.0891660535606589, -1.4464977326095683] [0.2620118571253886, 0.49731731998127376, 1.919381481095498] [1.9396415382222492, -0.11904479669552767, -0.47288374848639386]; [-0.8894992202674227, 1.1406782313993147, -1.381174974272042] [0.8946559815712687, 1.6981212985564655, -0.5621162958211153] [0.4719956166278159, -0.028968559968890112, 1.943291270092431]; [1.5771586856804132, 1.229866041560149, -4.209843244698678e-12] [-1.7694463714293387, 0.9322336287838716, -6.415864435534256e-13] [-0.1225187376973445, -1.9962437623980311, 2.2485461341507808e-12]], 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), [-5.696178589183529e-13, -1.0464674673019045e-12, 0.8573539534382282], [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