Skip to main content

Symbolic programming

Quantum spins​

The beauty of Wolfram Language is its flexibillity. One can shape the language making a dialect of it for a specific problem to solve

Download original notebook

Wavefunctions​

Let us use built-in syntax for Ket and Bra vectors. Apart from the decorations, there is no predefined rules for those expressions

{Ket[1], Bra[1]}
{(*BB[*)(Ket[1])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFaYwgCV4gEZaZWu6SmpxflFiSXxTMDBTxTi2B6GMBEkGlOanBIIZHamIKTBoAFg4VsQ=="*)(*]BB*),(*BB[*)(Bra[1])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFaYwgCV4gEZaZWu6SmpxflFiSXxTMDBRxKkqE6GMBEkGlOanBIIZHamIKTBoAFBUVkw=="*)(*]BB*)}

Let us define an ortogonallity rule on Dot operation aka

⟨Ms∣Ms′⟩=δMs,Ms′\langle Ms | Ms^\prime \rangle = \delta_{Ms, Ms^\prime} 
Ket /: (*BB[*)(Bra[Ms1_])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMaeiRJfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABBAEzo="*)(*]BB*).(*BB[*)(Ket[Ms2_])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*) := KroneckerDelta[Ms1, Ms2]
Ket /: Conjugate[(*BB[*)(Ket[Ms2_])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*)] :=  (*BB[*)(Bra[Ms2])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMaeiRJfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABBAEzo="*)(*]BB*)

now test it in the basis of a spin S=2S=2

∣S=2,Ms⟩| S=2, Ms \rangle 
basis = Table[(*BB[*)(Ket[Ms])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*), {Ms, -2,2}];

Table[Conjugate[i].j, {i, basis}, {j, basis}] 
{{1,0,0,0,0},{0,1,0,0,0},{0,0,1,0,0},{0,0,0,1,0},{0,0,0,0,1}}

or more appealing view

% // MatrixForm 
((*GB[*){{1(*|*),(*|*)0(*|*),(*|*)0(*|*),(*|*)0(*|*),(*|*)0}(*||*),(*||*){0(*|*),(*|*)1(*|*),(*|*)0(*|*),(*|*)0(*|*),(*|*)0}(*||*),(*||*){0(*|*),(*|*)0(*|*),(*|*)1(*|*),(*|*)0(*|*),(*|*)0}(*||*),(*||*){0(*|*),(*|*)0(*|*),(*|*)0(*|*),(*|*)1(*|*),(*|*)0}(*||*),(*||*){0(*|*),(*|*)0(*|*),(*|*)0(*|*),(*|*)0(*|*),(*|*)1}}(*]GB*))

Ladder operators​

Let us define a typical lowering and raising operators

Ket /: Sp.(*BB[*)(Ket[Ms_])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*) := (*SqB[*)Sqrt[2(2+1) - (Ms+1)Ms](*]SqB*) (*BB[*)(Ket[Ms+1])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*)
Ket /: Sm.(*BB[*)(Ket[Ms_])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*) := (*SqB[*)Sqrt[2(2+1) - (Ms-1)Ms](*]SqB*) (*BB[*)(Ket[Ms-1])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*)
Ket /: Sz.(*BB[*)(Ket[Ms_])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*) := Ms (*BB[*)(Ket[Ms])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*)

Sx = (*FB[*)((1)(*,*)/(*,*)(2))(*]FB*)(Sp + Sm);
Sy =  (*FB[*)((1)(*,*)/(*,*)(2I ))(*]FB*)(Sp - Sm);

Let's test

Sp.(*BB[*)(Ket[1])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFxTxAMe/UEpfU5PyixJL8IohKFiARVJqTGgxieKQmpgQzQ9QBABLAE1g="*)(*]BB*)
2 ((*BB[*)(Ket[2])(*,*)(*"1:eJxTTMoPSmNiYGAo5gcSAUX5ZZkpqSn+BSWZ+XnFaYwgCV4gEZaZWu6SmpxflFiSXxTMDBTxTi2B6GMBEkGlOanBIIZHamIKTBoAFg4VsQ=="*)(*]BB*))

Linear algebra​

In Woflram Language a linear properties of Dot operator are not predefined, i.e.

a . (2 b) === 2 a.b 
a . (b + c) === a.b + a.c 
False
False

But it is relatively easy to define them

Unprotect[Dot];
Unprotect[Conjugate];

Conjugate[a_ b_] := Conjugate[a] Conjugate[b]
Conjugate[a_ + b_] := Conjugate[a] + Conjugate[b]

Dot[a_, b_?NumericQ] := b a
Dot[b_?NumericQ, a_] := b a
Dot[a_ + b_, c_] := a.c + b.c 
Dot[c_, a_ + b_] := c.a + c.b  
Dot[a_, b_?NumericQ c_] := b a.c
Dot[b_?NumericQ c_, a_] := b c.a

Protect[Dot];
Protect[Conjugate];

Now it is better

a . (2 b) === 2 a.b 
a . (b + c) === a.b + a.c 
True
True

Spin Hamiltonian​

For example, one can try to solve a well-known spin-hamiltonian for a spin system with an uniaxial anisotropy in magnetic field perpendicular to it

H^=S^xh−dS^z2\hat{H} = \hat{S}_x h - d \hat{S}_z^2

If d is negative, the spin system has easy-plane anisotropy, i.e. it is more preferable for a spin to lie in xy plane.

In a code it looks like this

H[h_, d_] := Sx h - d Sz.Sz

Let us find an average projection of a spin to the field direction

SxMatrix = Table[Conjugate[i].Sx.j, {i, basis}, {j,basis}];

Compute for each field the ground state and find an average of Sx

estimate[d_] := Table[
  With[{
    matrix = Table[Conjugate[i].H[h, d].j, {i, basis}, {j,basis}]
  },
    With[{
      system = Eigensystem[matrix]
    },
      With[{
        gState = {SortBy[system // Transpose, First] // First // Last}
      },
        {h, Conjugate[gState].SxMatrix.Transpose[gState] // First // First} 
      ]
    ]
  ]
, {h, 0, 5, 0.1}];


ListLinePlot[{estimate[-1], estimate[1]}, PlotLegends->{"easy-plane", "easy-axis"}, AxesLabel->{"h", "S_{X}"}]
(*VB[*)(Legended[ToExpression[FrontEndRef["6295e4b2-f7d7-4184-b83c-691d6b88cdc2"], InputForm], Placed[LineLegend[{Directive[PointSize[1/72], RGBColor[0.24, 0.6, 0.8], AbsoluteThickness[2]], Directive[PointSize[1/72], RGBColor[0.95, 0.627, 0.1425], AbsoluteThickness[2]]}, {"easy-plane", "easy-axis"}, LegendMarkers -> {{False, Automatic}, {False, Automatic}}, Joined -> {True, True}, LabelStyle -> {}, LegendLayout -> "Column"], After, Identity]])(*,*)(*"1:eJytUstKAzEUrW8tCuoHCILbWfRhrYsyVFtfVKwdcZ9J7mgwTUqS0Y57P0K37vyCrv0AF278AN2IIPgHJhmtVAQRzOKQ3NzHuefexVC0opFMJqNmDBxSOKsBFhJpIYMJY2nAEXCSj4atS9ZAnVDzZx2jIWubN7AhBdd1TupdwLFGIYNgyZhL+dVlKIZ5L1ohK14xVy56YbmAvdJqjpTCchkT/JF41EArNmGT9gKI7HGWOOuBjCHlN26gyRAGEo1/kmlQDinDrzwNqnQaMWWgRiVgTU8hZWtNTUG5Dug5pDGuJNJUcMSodaJbBtIM7m9zbV0wIWVv4eJlv3fny4I7j768urTn2U9zzxmohkqwWMPBMcUnHJSiw/1k/0wncufVlzdv97vh7JMvK9mH607l9nc6A1oFVkhAKvE6DHEIpj6fqEvV4HjcjqR67yJ5AlJ9k33gpcbsaiCmwPVZjbVom8bwX71cZTvyHaMUfJt0f0m+tmUw0K0JCoEFOmEQZX7m61yn+901UCJiHdiiRuy4zR3LaqRNy3YI2wS4pjp5BzruzsA="*)(*]VB*)

One can modify it and estimate the projection for different anisotropy term d magnitudes

Taking it from 0 up to ~h

Table[
  With[{
    matrix = Table[Conjugate[i].H[h, d].j, {i, basis}, {j,basis}]
  },
    With[{
      system = Eigensystem[matrix]
    },
      With[{
        gState = {SortBy[system // Transpose, First] // First // Last}
      },
        {h, Conjugate[gState].SxMatrix.Transpose[gState] // First // First} 
      ]
    ]
  ]
, {d, 0, 10, 1}, {h, 0, 10, 0.2}];

ListLinePlot[%, ImageSize->400, Frame->True, PlotStyle->"DarkRainbow", PlotLegends->Table["H="<>ToString[SetPrecision[h, 1]], {h, 0, 10, 0.2}]]
(*VB[*)(Legended[ToExpression[FrontEndRef["97462b5c-bc3b-4580-a58a-dcbbbc3dce0e"], InputForm], Placed[LineLegend[{Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.237736, 0.340215, 0.575113]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.253651, 0.344893, 0.558151]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.264425, 0.423024, 0.3849]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.291469, 0.47717, 0.271411]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.416394, 0.555345, 0.24182]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.624866, 0.673302, 0.264296]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.813033, 0.766292, 0.303458]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.877875, 0.731045, 0.326896]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.812807, 0.518694, 0.303459]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.72987, 0.239399, 0.230961]], Directive[PointSize[1/120], AbsoluteThickness[2], RGBColor[0.72987, 0.239399, 0.230961]]}, {"H=0", "H=0.2", "H=0.4", "H=0.6", "H=0.8", "H=1.", "H=1.", "H=1.", "H=2.", "H=2.", "H=2."}, LegendMarkers -> {{False, Automatic}, {False, Automatic}, {False, Automatic}, {False, Automatic}, {False, Automatic}, {False, Automatic}, {False, Automatic}, {False, Automatic}, {False, Automatic}, {False, Automatic}}, Joined -> {True, True, True, True, True, True, True, True, True, True, True}, LabelStyle -> {}, LegendLayout -> "Column"], After, Identity]])(*,*)(*"1:eJxTTMoPSmNmYGAo5gUSYZmp5S6pyflFiSX5RcHsQBGf1PTUvBSjNCaQEi4g4ZqSCZQDKUxjBIkJAQm3ovy8Ete8FNeK1OTSksSknNRgFaCwpbmJmVGSabJuUrJxkq6JqYWBbqKpRaJuSnJSElAoJTnVIBViMAuQCCoFauMAMVITU/zzcirBoiFFpakQ97EBiYCcxOTUlDQ2mGN8MvNSIS5M44aZ45NZXALRwQkkXDKLUpNLMstSIa4FCQXkZ+aVBGdWQe0GW5lYkpmfl5iTCVKUWQEkIMoFgYRjUnF+TmlJakhGZnJ2XmpxcSZIG8QGsF53J+f8nPyiIqeqLhmlvHP2RV8/L1MSPXbVvqhA5UJCcPoj+wFxznKHtPSL1hfsi4QWPxfbJXTNvmh1ZHpN/J2HA+Occ6vOXw17A3SOski042Wh2/ZFD/bxzTFedGNgnNPDv3l67vJL9kXvyjW4vnTesy/6PtnQ7kzcxYFxjiWXAKvxslv2RRkv5j9IOvrQvig79ZD+18/nBsY5kjY+v5/9e2xfdHjWqUUbup/aFz1aLXvU8tWFgXFOpLHRiViWV/ZFDt+fnilrfWFfNCEzV+p29uWBcU6K9X3/XvE39kUGH+5U70h5bl+0PNbH+O7rKwPjHNftcZebmYChs8DkiobKzAf2RXsDtX3uD1ToXJjdFSIRDQwU7Z83YhcsAZaHL5MueylOODvonINSZQSDtHrYGgSzQmg9IzjLBM4yg7MsglnALEM9nAwjLAzU6g5c50LqL99EYI4vKk7jQnYTQjWIVwyy2y0xpzgVHGCOpSX5ucAQSh5VRVtV4JgCNTm8gEk0Fa2lAW+k0I6B5hJwuycxKTUnuKQyJzWNAXuCASvlgScvn8TK/NKSYJAvgNmiNDcP7G3HtBJgmgNlF8+U1LySzJJKAOIInDE="*)(*]VB*)
Last updated on