(* 
          Cartesian and spherical tensors in NMR Hamiltonian 
  P .  P .  Man
  www . pascal - 
            man . com
  April 2013
*) 

 Table[a[k, q],  {k,  0,  2},  {q,  -k,  k}];  
        Table[t[k, q],  {k,  0,  2},  {q,  -k,  k}];

  
h1c = "bx sx txx + by sx txy + bz sx txz + \
bx sy tyx + by sy tyy + bz sy tyz + bx sz tzx + by sz tzy + bz sz tzz";
  
h2c = "bx sx txx + bx sy txy + bx sz txz + \
by sx tyx + by sy tyy + by sz tyz + bz sx tzx + bz sy tzy + bz sz tzz";

  
t[0, 
              0 ] :=  -(txx  +  tyy  +  tzz)/Sqrt[3];
t[1, 1] := 
            plusOrMinusT*(tzx  -  txz  +  I (tzy  -  tyz))/2;
t[1, 0]  := 
            plusOrMinusT*I (txy  -  tyx)/Sqrt[2];
t[1, -1] := 
            plusOrMinusT*(tzx  -  txz  -  I (tzy  -  tyz))/2;
t[2, 
              2] := (txx  -  tyy  +  I (txy  +  tyx))/2;
t[
              2, -2] := (txx  -  tyy  -  I (txy  +  tyx))/
              2;
t[2, 
              1] := -(txz  +  tzx  +  I (tyz  +  tzy))/2;
t[
              2, -1] := (txz  +  tzx  -  I (tyz  +  tzy))/
              2;
t[2, 
              0] := (3 tzz  -  (txx  +  tyy  +  tzz))/
              Sqrt[6];

  
a[0, 
              0 ] :=  -(bx*sx  +  by*sy  +  bz*sz)/Sqrt[3];
a[1, 1] := 
            plusOrMinusA*(sz*bx  -  sx*bz  +  
                    I (sz*by  -  sy*bz))/2;
a[1, 0]  :=  
            plusOrMinusA*I (sx*by  -  sy*bx)/Sqrt[2];
a[1, -1] := 
            plusOrMinusA*(sz*bx  -  sx*bz  -  
                    I (sz*by  -  sy*bz))/2;
a[2, 
              2] := (bx*sx  -  by*sy  +  I (bx*sy  +  by*sx))/
              2;
a[
              2, -2] := (bx*sx  -  by*sy  -  
                  I (bx*sy  +  by*sx))/2;
a[2, 
            1] := -(bx*sz  +  bz*sx  +  
                  I (by*sz  +  bz*sy))/2; 
        a[2, -1] := (bx*sz  +  bz*sx  -  
                I (by*sz  +  bz*sy))/2;
a[2, 
              0] := (3 bz*sz  - (bx*sx  +  by*sy  +  bz*sz))/
              Sqrt[6];

Print["========================================================"]

plusOrMinusT = 1;
plusOrMinusA = 1;
Print["Tensor components T{1}-A{1}"]

Print["in H{1,s},"]
Print["Result:"]
H  =  
          Expand[Sum[((-1)^(k + q))*a[k, q]*
                t[k, -q],  {k,  0,  2},  {q,  -k,  k}]]
HH  =  ToString[H];
Print["corresponding to "]
If [Equal[HH,  h1c], 
          Print["H{1,c}"],  Print["H{2,c}"]]

Print["---------"]
Print["in H{2,s},"]
Print["Result:"]
G = 
          Expand[Sum[((-1)^q)*a[k, q]*t[k, -q],  {k,  0,  
                2},  {q,  -k,  k}]]
GG  =  ToString[G];
Print["corresponding to "]
If [Equal[GG,  h2c], Print["H{2,c}"], 
          Print["H{1,c}"] ]

Print["========================================================"]

plusOrMinusT = 1;
plusOrMinusA = -1;
Print["Tensor components T{1}-A{2}"]

Print["in H{1,s},"]
Print["Result:"]
H  =  
          Expand[Sum[((-1)^(k + q))*a[k, q]*
                t[k, -q],  {k,  0,  2},  {q,  -k,  k}]]
HH  =  ToString[H];
Print["corresponding to "]
If [Equal[HH,  h1c], 
          Print["H{1,c}"],  Print["H{2,c}"]]

Print["---------"]
Print["in H{2,s},"]
Print["Result:"]
G = 
          Expand[Sum[((-1)^q)*a[k, q]*t[k, -q],  {k,  0,  
                2},  {q,  -k,  k}]]
GG  =  ToString[G];
Print["corresponding to "]
If [Equal[GG,  h2c], Print["H{2,c}"], 
          Print["H{1,c}"] ]

Print["========================================================"]

plusOrMinusT = -1;
plusOrMinusA = 1;
Print["Tensor components T{2}-A{1}"]

Print["in H{1,s},"]
Print["Result:"]
H  =  
          Expand[Sum[((-1)^(k + q))*a[k, q]*
                t[k, -q],  {k,  0,  2},  {q,  -k,  k}]]
HH  =  ToString[H];
Print["corresponding to "]
If [Equal[HH,  h1c], 
          Print["H{1,c}"],  Print["H{2,c}"]]

Print["---------"]
Print["in H{2,s},"]
Print["Result:"]
G = 
          Expand[Sum[((-1)^q)*a[k, q]*t[k, -q],  {k,  0,  
                2},  {q,  -k,  k}]]
GG  =  ToString[G];
Print["corresponding to "]
If [Equal[GG,  h2c], Print["H{2,c}"], 
          Print["H{1,c}"] ]

Print["========================================================"]

plusOrMinusT = -1;
plusOrMinusA = -1;
Print["Tensor components T{2}-A{2}"]

Print["in H{1,s},"]
Print["Result:"]
H  =  
          Expand[Sum[((-1)^(k + q))*a[k, q]*
                t[k, -q],  {k,  0,  2},  {q,  -k,  k}]]
HH  =  ToString[H];
Print["corresponding to "]
If [Equal[HH,  h1c], 
          Print["H{1,c}"],  Print["H{2,c}"]]

Print["---------"]
Print["in H{2,s},"]
Print["Result:"]
G = 
          Expand[Sum[((-1)^q)*a[k, q]*t[k, -q],  {k,  0,  
                2},  {q,  -k,  k}]]
GG  =  ToString[G];
Print["corresponding to "]
If [Equal[GG,  h2c], Print["H{2,c}"], 
          Print["H{1,c}"] ]
Null