Yb$_2$Ti$_2$O$_7$ - Excitations in a quantum spin ice

Recreating experimental data in Phys. Rev. X 1, 021002
Simon Ward - Scientific Software Developer - ESS
Overview
paper data of Yb$_2$Ti$_2$O$_7$
Published data and simulation of Yb$_2$Ti$_2$O$_7$
Phys. Rev. X 1, 021002

Overview

We will be trying to recreate PhysRevX.1.021002

Yb$_2$Ti$_2$O$_7$ is a quantum spin ice with a pyrochlore structure

The Hamiltonian given by:


$$ \mathcal{H} = \frac{1}{2} \sum_{i, j} J_{i, j}^{\mu , \nu} S_i^\mu \cdot S_j^\nu - \mu_B H^\mu \sum_i g_i^{\mu, \nu} S_i^\nu $$

The Yb$_2$Ti$_2$O$_7$ model

Download the script here: sw_tutorial_06.m
Creating the lattice

The Yb$_2$Ti$_2$O$_7$ lattice

Download the paper from here we will be referencing it quite a bit.

SpinW does not contain the space group we need, so we have to manually enter it.

symStr = '-z, y+3/4, x+3/4; z+3/4, -y, x+3/4; ';
symStr = strcat(symStr, 'z+3/4, y+3/4, -x; y+3/4, x+3/4, -z; ');
symStr = strcat(symStr, 'x+3/4, -z, y+3/4; -z, x+3/4, y+3/4');

ybti = spinw;
a = 10.0307;
ybti.genlattice('lat_const', [a a a], 'angled', [90 90 90], 'spgr', symStr, 'label', 'F d -3 m Z');
ybti.addatom('label', 'Yb3+', 'r', [1/2 1/2 1/2], 'S', 1/2)
ybti.addatom('label', 'Ti4+', 'r', [0 0 0])
ybti.addatom('label', 'O2-', 'r', [0.3318 1/8 1/8])
ybti.addatom('label', 'O2-', 'r', [3/8 3/8 3/8])
plot(ybti, 'nMesh', 3)
swplot.legend('none')
A diversion into plotting options....

Making plots look better

We setup the number of triangular faces to be produced by the plot command using swpref. Each atom is an icosahedron, where each face is subdivided into triangles nmesh-times (see swplot.icomesh function). Thus nmnesh=3 will make all sphere to have 1280 faces. npatch determines the number of subdivision of the circle that is used to generate cylinders an arrows (swplot.arrow and swplot.cylinder).

pref = swpref;
pref.set({'nmesh', 'npatch'}, {3, 50})

Questions:

Take a note of the options in swpref. We will be using it later

A diversion into plotting options....

Drawing polyhedra

To the draw oxygen polyhedra around the Yb ions, we use the swplot.plotchem() function, that can draw polyhedra around arbitrary atoms on an existing crystal structure plot. Now we use center atom 'Yb' and polyhedra atoms 'O' for oxygen. Since the oxygen environment of Yb is octahedron, we set the limits to the 8 closes oxygen atom. The same can be achieved within the spinw.plot() function, adding the options below plus a string 'chem', such as 'chemAtom1', 'chemAtom2', 'chemLimit' and 'chemRange' and setting 'chemMode' to 'poly'.

swplot.plotchem('atom1','Yb','atom2','O','limit',8,'range',[0.1 0.9;0.1 0.9;0.1 0.9]);
Creating the Hamiltonian

Creating the Hamiltonian

We can remove the non-magnetic atoms from the spinw object with a single command using the unitcell() function (not to mix with the unit_cell property of the spinw object). The unitcell() function can return selected atoms from the list of symmetry inequivalent atoms in the unit cell. In our case the magnetic Yb ions are the first atom.

ybti.unit_cell = ybti.unitcell(1);
ybti.gencoupling
ybti.addmatrix('label','J1','value',1)
ybti.addmatrix('label','g0','value',1);

ybti.addcoupling('mat','J1','bond',1)
ybti.addg('g0')
Creating the Hamiltonian

Creating the Hamiltonian - Part 2

In the paper the anisotropic g-tensor is defined in the local coordinate system of the magnetic ions. Where the g_z component is along the local [1 1 1] direction, while the two perpendicular components are $g_xy$. In the lattice corrdinate system the g-tensor has the matrix form: [A B B;B A B;B B A]. One can check the eigenvalues of this matrix, that has to match with the published values: $g_{xy}=4.32$ and $g_z=1.8$. From the eigenvalue calculation we get: $g_{xy}=A-B$; $g_z = A + 2*B$. We store the calculated g-tensor in the sw object. When calculating the spin wave intensities, the code takes care the rotation of the g-tensor according the symmetry operators for every magnetic ion.

ybti.matrix.mat(:,:,2) =  -0.84*ones(3)+4.32*eye(3);
J1 = -0.09; J2 = -0.22; J3 = -0.29; J4 = 0.01;
ybti.setmatrix('mat','J1','pref',[J1 J3 J2 -J4]);

Follow the tutorial to make sure you understand what's going on.

Questions:

Check that the correct values have been entered into the exchange martix

Plotting the Hamiltonian

Plotting the Hamiltonian

With the plot() command, we can plot the magnetic bonds of Yb$_2$Ti$_2$O$_7$. The arrow pointing from an atom to another denotes the direction of the bond, while the thicker arrow at the middle of the bond denotes the direction of the Dzyaloshinskii-Moriya vector. Is is also possible to visualize the g-tensor by setting the 'ionMode' to 'g'.

plot(ybti,'ionMode','g')
swplot.zoom(1.3)
swplot.legend('none')
A diversion into plotting options....

Delving into plotting

The high level spinw.plot command calls lower level commands (swplot.plotatom, swplot.plotion, swplot.plotbond, etc). For details check the documentation of spinw.plot or any of the lower level functions. Any option of the lower level functions can be controlled by making new options as: lowlevelfunname + lowlevelfunoption, for example to set the color option of the swplot.plotatom function set 'atomColor' option in the spinw.plot method as in the example below. The low level functions can be also called separately.

swplot.plotchem('atom1','Yb','atom2','O','limit',8,'range',[0.1 0.9;0.1 0.9;0.1 0.9]);
A diversion into plotting options....

Delving into plotting - Manual polyhedra

To manually plot a polyhedra, you need to supply the vertices which describe the shape. And have them in a shape which MATLAB understands.

swplot.plot('type', 'polyhedron', 'position', permute(R,[2 3 1]));

Questions

Choose 4 atoms and make a tetrahedra or two. Make a note of the atoms x, y, z positions (in fractional lattice co-ordinates) to generate R. Hint, you might want to turn on the tooltip


We can also move the atoms by shifting positions by 10 Angstrom to the right.

swplot.plotatom('shift', [10 0 0]', 'mode', 'mag', 'replace', true)
Recreate the papers figure

Recreate the papers figure

Follow and understand what's happening in the script

  • Make Q, E slice
  • Apply a magnetic field and optimise th structure
  • Apply an instrumental resolution
  • Create a plot for each applied magnetic field