Matrix representations of spin operators.

S = sop(Spins,Comps)
[S1,S2,...] = sop(Spins,Comps1,Comps2,...)
... = sop(...,'sparse');

sop returns the matrix representation of spin operators, defined in Comps (or Comps1, Comps2, etc), for a set of spins given in Spins.

The representation of the operators is in terms of product Zeeman basis states |m₁,m₂,m₃,..>, where m₁ is for the first spin, m₂ is for the second spin, etc. The basis functions are in descending order of m_i, left to right. E.g. for a system with a spin-1/2 as the first spin and a spin-1 as the second spin, the states are in the following order: |1/2,+1>, |1/2,0>, |1/2,-1>, |-1/2,+1>, |-1/2,0>, |-1/2,-1>.

Spins is a vector containing the spin quantum numbers (for example [1/2 1]) or, alternatively, a spin system structure.

Comps is a character or numeric array specifying the particular spin operator. There are several different ways to do this:

1. One component for each spin: Give one operator component for each spin in Spins. Elements of Comps can be 'e' (identity), 'x' (S_x), 'y' (S_y), 'z' (S_z), '+' (S₊) or '-' (S_-). For spin-1/2, 'a' and 'b' correspond to S_α and S_β. The order in Comps corresponds to the order in Spins.

   For example, for a spin system with one spin-1/2 and one spin-1 (Spins = [1/2 1]), Comps = 'xe' would give S_1x, Comps = 'ex' gives S_2x, and Comps = 'xy' gives S_1xS_2y.

2. Component plus spin index: Give a list of components, each followed by the index of the spin.

   For example, for a 3-spin system with Spins = [1/2 1/2 1/2], Comps = 'x2' gives S_2x, Comps = 'z1y3' gives S_1zS_3y, and Comps = 'z1+2-3' gives S_1zS₂₊S_3-. If indices are present, the order of the entries is irrelevant: 'z1y3' and 'y3z1' are identical. Optionally, you can use commas to separate the components: 'z1y3' and 'z1,y3' give the same result.

3. Component and transition (and index): For transition-selective operators, specify the pair of levels after the component in the form (level1|level2), optionally followed by the spin index. The level indices refer to the standard descending m order, so level 1 corresponds to m = S, level 2 to m = S-1, etc.

   For example, for a single spin, 'x(1|3)' indicates an Sx operator connecting levels 1 and 3. 'e(2|2)' (or 'e(2)') indicates a matrix with a single non-zero element at (2,2), corresponding to the population of level 2. For more than one spin, commas are required, and spin indices are necessary if you don't give a component for each spin. For example, '+(1|3),z' for a 2-spin system gives an S₊ operator that connects levels 1 and 3 of spin 1, times an S_z operator for spin two. 'y(1|4)1,x(2|3)3' for a 3-spin system indicates an S_y operator connecting levels 1 and 4 of the first spin, times an S_x operator connecting levels 2 and 3 of the third spin.

4. Spin index and component code: Give an Nx2 array, where each row contains a spin index and a letter code for the component, where 0 corresponds to 'e', 1 is 'x', 2 is 'y', 3 is 'z', 4 is '+', and 5 is '-'. For example [1 3; 2 4] indicates 'z' for the first spin and '+' for the second spin (and thus equivalent to 'x+' if this is a two-spin system, or 'x1+2'.) This numeric syntax is a bit bulkier than the character-based ones above, but can be useful in loops.

If more than one operator is provided (Comps1, Comps1, etc), sop computes all operators at the same time. Of course, there must be a corresponding number of output arguments.

If 'sparse' is given as last argument, spin operator matrices are returned in sparse instead of full form.

Here are a few basic spin operators:

Splus = sop(1/2,'+')

Splus =
     0     1
     0     0

S2z = sop([1/2 1/2],'z2')

S2z =
    0.5000         0         0         0
         0   -0.5000         0         0
         0         0    0.5000         0
         0         0         0   -0.5000

SzIx = sop([1/2 1/2],'zx')

SzIx =
         0    0.2500         0         0
    0.2500         0         0         0
         0         0         0   -0.2500
         0         0   -0.2500         0

Ip = sop([1/2 1],[2 4])

     Ip =
     0       1.4142            0            0            0            0
     0            0       1.4142            0            0            0
     0            0            0            0            0            0
     0            0            0            0       1.4142            0
     0            0            0            0            0       1.4142
     0            0            0            0            0            0

S1zI3p = sop([1/2 1/2 1/2],'z1+3')

S1zI3p =
         0    0.5000         0         0         0         0         0         0
         0         0         0         0         0         0         0         0
         0         0         0    0.5000         0         0         0         0
         0         0         0         0         0         0         0         0
         0         0         0         0         0   -0.5000         0         0
         0         0         0         0         0         0         0         0
         0         0         0         0         0         0         0   -0.5000
         0         0         0         0         0         0         0         0

To compute multiple operators with one call, use for example

[Sx,Sy,Sz] = sop(1,'x','y','z');

Sx =
         0    0.7071         0
    0.7071         0    0.7071
         0    0.7071         0
Sy =
        0                  0 - 0.7071i        0
        0 + 0.7071i        0                  0 - 0.7071i
        0                  0 + 0.7071i        0
Sz =
     1     0     0
     0     0     0
     0     0    -1

Instead of the spin quantum numbers, a spin system structure can be provided.

System.S = 1;
System.g = 2;
System.D = [1000 0];
sop(System,'-')

ans =
         0         0         0
    1.4142         0         0
         0    1.4142         0

Here are a few examples of transition-selective operators:

Sx13 = sop(1,'x(1|3)')      % Sx operator connecting levels 1 and 3

Sx13 =
         0         0    0.5000
         0         0         0
    0.5000         0         0

SpIb = sop([1/2 1/2],'+,e(2)')    % S+ operator for first spin, only level 2 (beta) of second spin

SpIb =
     0     0     0     0
     0     0     0     1
     0     0     0     0
     0     0     0     0

Sxc = sop(7/2,'x(4|5)')    % Sx operator for central -1/2<->+1/2 transition of an S=7/2 (e.g Gd3+)

Sxc =
         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.5000         0         0         0
         0         0         0    0.5000         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

spin operators, spin system, isto spinvec, stev