Calculates energy levels of a spin system in an external field.

E = levels(SpinSystem,Ori,B)
E = levels(SpinSystem,phi,theta,B)
[E,V] = levels(...)

This function calculates the energy levels of the spin system specified in the structure SpinSystem by diagonalising its Hamiltonian. The direction and strength of the external magnetic field is specified in the other input arguments.

theta is the angle between the z axis of the molecular frame and the magnetic field vector, phi is the angle between the molecular x axis and the projection of the magnetic field vector onto the molecular xy plane. Both angles are in radians.

Ori is a 2xn or nx2 array that contains both phi and theta, equivalent to [phi; theta]. Alternatively, the three principal directions x, y, and z and intermediate orientations can be specified by shorthands: In Ori, give 'x' for [0 pi/2], 'y' for [pi/2 pi/2], 'z' for [0 0], 'xy' for [pi/4 pi/2], 'xz' for [0 pi/4], 'yz' for [pi/2 pi/4], and 'xyz' for [pi/4 acos(1/sqrt(3))].

B is an array containing values of the external magnetic field, in mT. If using a list of orientations, this can also be a single number.

E is an array of the energies of the eigenstates of SpinSystem, in ascending order. The energies are in MHz.

phi, theta and B can be vectors, in which case energy levels for all possible combinations of their values are calculated. Depending on the dimensions of the three parameters, the output E can be an array with up to 4 dimensions. The energies are always along the last dimension.

A scalar theta together with vectors phi and B yields a 3D array E. A set of energies for a single magnetic field orientation is then in E(k,l,:) with k and l indexing phi and B, respectively.

If the second output Vis requested, the eigenvectors corresponding to the eigenvalues in E are returned.

To generate a plot of the energy levels of a spin-5/2 system as functions of the external magnetic field, use

clear
Sys.S = 5/2;
Sys.g = 2;
Sys.D = [-250 -1750 2000];  % MHz

B = 0:500;  % field, mT
Ori = [10 15]*pi/180;  % phi and theta, rad

E = levels(Sys,Ori,B);
plot(B,E)

Note that several level anticrossings are visible.

levelsplot, resfields_eig, resfields, resfields_perturb, ham