Transformation theory.

Needs["Utilities`CleanSlate`"]
CleanSlate[]
Remove["Global`*"]

Cell 1:  Show that the matrix converting {(Subscript[Y, 1]^1), (Subscript[Y, 1]^0), (Subscript[Y, 1]^(-1 ))} to the orbitals {(Subscript[p, x]), (Subscript[p, y]), (Subscript[p, z])} is unitary.

We have the following definitions of the orbitals (Subscript[p, x]), (Subscript[p, y]), and (Subscript[p, z]) (see Metiu, Quantum Mechanics, p. 321) in terms of the Spherical Harmonics (Subscript[Y, 1]^(-1)), (Subscript[Y, 1]^0), (Subscript[Y, 1]^1)

	(Subscript[p, x]) = (1/(Sqrt[2]))((Subscript[Y, 1]^(-1))- (Subscript[Y, 1]^1))					(1)
	
	(Subscript[p, y]) = (i/( (Sqrt[2])))((Subscript[Y, 1]^(-1)) + (Subscript[Y, 1]^1))					(2)
	
	(Subscript[p, z]) =  (Subscript[Y, 1]^0) 							(3)
	
Do not confuse (Subscript[p, x]), (Subscript[p, y]), and (Subscript[p, z]) with momentum components. They are the molecular orbitals defined by Eqs. 1-3. 

We want to show that the matrix U that converts the basis set {(Subscript[Y, 1]^1),  (Subscript[Y, 1]^0), (Subscript[Y, 1]^(-1))} into the basis set {(Subscript[p, x]), (Subscript[p, y]), (Subscript[p, z])}  is a unitary matrix.  This means that a calculation using one of these basis sets will give the same results for any observable as a calculation using the other basis set.  The set {(Subscript[p, x]), (Subscript[p, y]), (Subscript[p, z])} is more frequently used in the theory of chemical bonds because the results are easier to interpret. The set (Subscript[Y, 1]^(-1)), (Subscript[Y, 1]^0), (Subscript[Y, 1]^1) is more convenient for spectroscopy.

Define a vector pv with the components p[1], p[2], p[3]

Clear[p]
pv=Table[p[i],{i,1,3}]

Define the vector Yv with the components Y[1], Y[0] and Y[-1] which are symbols for the spherical harmonics (Subscript[Y, ℓ]^m)(θ, ϕ), with ℓ=1 and m=1, 0, -1.

Yv=Table[Y[i],{i,-1,1,1}]

Define the matrix converting the Y kets into the px kets. I found this matrix by inspecting Eqs. 1-3

U={{1/Sqrt[2],0,-1/Sqrt[2]},{I/Sqrt[2],0,I/Sqrt[2]},{0,1,0}}

Display as a matrix

U//MatrixForm

Test that the matrix U transforms the basis set {Y[-1], Y[0], Y[1]} into {px, py, pz}

vec=U.Yv

If you compare this output to Eqs. 1-3 you see that U transforms the basis set {(Subscript[Y, 1]^1),  (Subscript[Y, 1]^0), (Subscript[Y, 1]^(-1))} into the basis set {(Subscript[p, x]), (Subscript[p, y]), (Subscript[p, z])} 

We need to determine now that U is unitary. If it is, the absolute value of its determinant should be equal to 1. Test that:

Det[U]//Abs

Just for fun let us convert our symbols Y[i] into the spherical harmonics (Subscript[Y, ℓ]^m)(θ, ϕ) to see that we get the correct results for (Subscript[p, x]), (Subscript[p, y]), and (Subscript[p, z]).

vec=vec/.Y[m_]->SphericalHarmonicY[1,m,θ,ϕ]

Simplify the expressions:

vec//ExpToTrig

These are the well known expressions for (Subscript[p, x]), (Subscript[p, y]) and (Subscript[p, z]) orbitals (see Metiu, Quantum Mechanics, p. 321). Clearly (Subscript[p, x ])(the first component of the vector) is proportional to x, (Subscript[p, z]) is proportional to z, and (Subscript[p, y]) is proportional to y.  Hence, the names (Subscript[p, x]), (Subscript[p, y]), and (Subscript[p, z]).

Cell 2: Show that the eigenvalues of (Subscript[(OverHat[L]), z]) are the same whether we use the basis set  {(Subscript[Y, 1]^1), (Subscript[Y, 1]^0), (Subscript[Y, 1]^(-1 ))}  or the basis set {(Subscript[p, x]), (Subscript[p, y]), (Subscript[p, z])}

To find the eigenstates of the operator (OverHat[L])  I will calculate the matrix of the operator in the basis set {(Subscript[p, x]), (Subscript[p, y]),(Subscript[p, z])} and in the basis set {(Subscript[Y, 1]^(-1)),(Subscript[Y, 1]^0), (Subscript[Y, 1 ]^1)} .  I will obtain two different matrices. Then I calculate the eigenvalues of these two matrices. I want to test that the two matrices have the same eigenvalues!

To calculate matrices I must calculate scalar products of the form 

	⟨(Subscript[ϕ, i])|(OverHat[(Subscript[L, z])|(Subscript[ϕ, j])⟩]) 
	
where |(Subscript[ϕ, i])⟩ are elements of one of the two basis sets we are working with. Therefore the first thing I do is define a function that calculate the scalar product and also one that defined how (OverHat[(Subscript[L, z])]) acts on the elements of the basis set.

Define a function scalar[f[n],f[m]] that calculates the scalar product of f[n] with f[m].

The function 

	scalar[f[i], f[j]]  
	
is defined to assume that the arguments f[i] (which are kets) form a orthonormal basis set:

Clear[scalar]
scalar[Y_[m_],Y_[mp_]]:=KroneckerDelta[m,mp]

Test it

{scalar[a[2],a[2]],scalar[a[x],a[y]],scalar[a[1],a[2]]}

The scalar product must satisfy
	⟨a f[i] | f[j]⟩ =(SuperStar[a]) ⟨f[i] | f[j]⟩ 
Force the function scalar to implement this equation

scalar[a_*Y_[m_],  Y_[mp_]]:=(a/.Complex[x_,y_]:>Complex[x,-y])* scalar[Y[mp],Y[m]]

Test that it works

scalar[Sqrt[q] f[x],f[y]]

This works as we want.

Implement 

	⟨f[i] | b f[j]⟩ =b ⟨f[i] | f[j]⟩

scalar[ Y_[m_], b_*Y_[mp_]]:= b*scalar[Y[ m],Y[mp]]

Test that it works as intended

scalar[ b[x], c b[y]]

Implement 

	⟨a f[i] | b f[j]⟩ =(SuperStar[a])b ⟨f[i] | f[j]⟩

scalar[ a_*Y_[m_], b_*Y_[mp_]]:= (a/.Complex[x_,y_]:>Complex[x,-y])*b*scalar[Y[m],Y[mp]]

scalar[b a[x],a[y] f]

scalar[ a (b+I 3) c[x], (d+I 2) m c[y]]

Implement 

 	⟨a+b|c⟩ = ⟨a|c⟩ + ⟨b|c⟩

scalar[a_+b_,c_]:=scalar[a,c]+scalar[b,c]

Implement 

 	⟨a|b+c⟩ = ⟨a|b⟩ + ⟨a|c⟩

scalar[a_,b_+c_]:=scalar[a,b]+scalar[a,c]

Test them:

scalar[a[x]+b[y],c[y]]

scalar[ a Sin[g] b[x]+Cos[k] b[y], Sqrt[2] b[z]]

Implement 
	
	⟨ a f[x] | 0⟩ =0   and ⟨ 0 | a f[x]⟩ =0

scalar[a_.*Y_[x_],0]:=0
scalar[0,a_.*Y_[x_]]:=0

Test this

scalar[3 Y[1],0]

scalar[(a +I b)f[1],b f[1]]

Calculate the matrix elements of (Subscript[L, z]) in the basis set {(Subscript[p, x]), (Subscript[p, y]), (Subscript[p, z])} and calculate the eigenvalues and the eigenvectors of the matrix

Define (Subscript[p, x]), (Subscript[p, y]), (Subscript[p, z]) as combinations of Y[1],Y[0] and Y[-1] where Y[1] is (Subscript[Y, 1]^1)(θ, ϕ) , Y[0] is (Subscript[Y, 1]^0)(θ, ϕ),  etc.

p[1]=(1/Sqrt[2]) Y[-1] -(1/Sqrt[2]) Y[1]
p[2]=(I/Sqrt[2]) Y[-1] +(I/Sqrt[2]) Y[1]
p[3]=Y[0]

These equations define the vectors  (Subscript[p, x]), (Subscript[p, y]) and (Subscript[p, z]) in the space defined by the basis set Y[1],Y[0] and Y[-1]. Show that the kets p[1], p[2], and p[3] are orthonormal. To test orthonormality we use the scalar product defined earlier.

Table[scalar[p[i],p[j]],{i,1,3},{j,1,3}]

Mathematica comment:  the first thing Mathematica does is evaluate the argument of a function. Therefore it transforms scalar[p[1],p[1]] into scalar[((Y[-1])/(Sqrt[2]))-((Y[1])/(Sqrt[2])),((Y[-1])/(Sqrt[2]))-((Y[1])/(Sqrt[2]))] and then uses the rules defining the function scalar to evaluate this expression.  If it did not do that, it would automatically give  scalar[p[i],p[j]] =(Subscript[δ, ij]) .

Define the operator Lz acting on the functions Y[m].  Lz assumes that its argument is an eigenstate of the operator (OverHat[(Subscript[L, z])]).

Clear[Lz]
Lz[Y_[m_]]:=ℏ m Y[m]

Define the rule

	(Subscript[(OverHat[L]), z])( |m⟩ +|m'⟩) = (OverHat[(Subscript[L, z])])|m⟩ +  (OverHat[(Subscript[L, z])])|m'⟩

Lz[Y_[m_]+Y_[n_]]:=Lz[Y[m]]+Lz[Y[n]]

Lz[a[m]+a[b]]

Define the rule 

 	(OverHat[(Subscript[L, z])]) a |m⟩ = a |m⟩

Lz[a_* Y_[m_]]:=a Lz[Y[m]]

Test it

Lz[ 3 Y[1]]

Define the rule

	(OverHat[(Subscript[L, z])]) (a |m⟩ +b |m'⟩ ) = (OverHat[a (Subscript[L, z])]) |m⟩ +(OverHat[b (Subscript[L, z])])  |m'⟩

Lz[a_*Y_[x_] +b_*Y_[y_]]:=a Lz[Y[x]]+b Lz[Y[y]]

See how Lz acts on p[1]

Lz[p[1]]

Test that Lz works when used in the argument of the function scalar

{scalar[Y[1],Lz[Y[1]]],scalar[Lz[Y[1]],Lz[Y[1]]],scalar[Y[0],Lz[Y[1]]]}

Those are the correct results.

Calculate the matrix elements 

	⟨p[i] | (OverHat[(Subscript[L, z])])|p[j]⟩

maLz=Table[scalar[p[i],Lz[p[j]]],{i,1,3},{j,1,3}]

This is the matrix of (OverHat[(Subscript[L, z])]) in the basis set |(Subscript[p, x])⟩, |(Subscript[p, y])⟩, |(Subscript[p, z])⟩.

Display the matrix

maLz//MatrixForm

Calculate the eigenvalues of the matrix (which are the eigenvalues of the operator (OverHat[(Subscript[L, z])])).

Eigenvalues[maLz]

The operator (Subscript[L, z])  has the same eigenvalues in the basis set {p[i]} and in the basis set {Y[i]}, as the theory says that it should.

Let us now calculate the eigenvectors of the matrix (Subscript[L, z]) in the basis set {p[1], p[2],p[3]} (this is {(Subscript[p, x]), (Subscript[p, y]), (Subscript[p, z])})

{eve[1],eve[2],eve[3]}=Eigenvectors[maLz]

Test whether they are normalized (they must be orthogonal since they are eigenvectors of a Hermitian matrix):

Make a function that returns the complex conjugate of the vector

Clear[conjV]
 conjV[eve_?VectorQ]:= Map[Conjugate,eve]

Check that it works properly

Table[conjV[eve[i]],{i,1,3}]

Calculate the overlap matrix of the eigenvectors. Remember that the scalar product for vectors with complex components is
{v1, v2, v3}.{w1, w2, w3} =  (∑_(i=1)^3)(Subscript[v, i]^( * ))×(Subscript[w, i])

Table[conjV[ eve[i]].eve[j],{i,1,3},{j,1,3}]//MatrixForm

The eigenvectors are orthogonal but are not normalized. Normalize them:

Table[eveN[i]=eve[i]/Sqrt[conjV[eve[i]].eve[i]],{i,1,3}]

Test orthonormality:

Table[conjV[ eveN[i]].eveN[j],{i,1,3},{j,1,3}]//MatrixForm

Since the matrix is written in the basis {p[1], p[2],p[3]} the eigenvectors above are the coordinates of the eigenvectors in this basis set.   Hence

eigenvector1= 0 p[1] + 0 p[2] + 1 p[3]
eigenvector2  = (i/(Sqrt[2]))p[1] + (1/(Sqrt[2]))p[2] + 0 p[3]
eigenvector3 = ((-i  )/(Sqrt[2]))p[1] +(1/(Sqrt[2]))p[2] + 0 p[3]
However since we defined p[1], p[2] and p[3] in terms of Y[1], Y[2] and y[3], when we write the above equations for the eigenvectors Mathematica will automatically write p[i] in terms of Y[i], giving us the eigenvectors in therms of Y[-1], Y[0] and Y[1].  Let's see what these eigenvectors look like.

Table[v[i]= eveN[i].{p[1],p[2],p[3]}//FullSimplify,{i,1,3}]

It may surprise you that we do not get the eigenvectors to be Y[0],Y[1] and Y[-1]. However the eigenvectors we got are correct and are orthonormal as you can see below. They are as good as Y[0], Y[-1] and Y[1] because i = exp[ I π/2] is a phase factor as you can see below

Exp[I   Pi/2]

Table[scalar[v[i],v[j]],{i,1,3},{j,1,3}]

Cell 3:  Hybridization: (sp^3)

You have learned in the introductory chemistry about the hybrid orbitals (sp^3) of the carbon atom. They are given by

h1=(1/2) (s+(Subscript[p, x]) +(Subscript[p, y])+(Subscript[p, z]))

h2=(1/2) (s+(Subscript[p, x]) -(Subscript[p, y])-(Subscript[p, z]))

h3=(1/2) (s-(Subscript[p, x]) +(Subscript[p, y])-(Subscript[p, z]))

h4=(1/2) (s-(Subscript[p, x]) -(Subscript[p, y])+(Subscript[p, z]))

Create a 4-dimensional vector ov with the atomic orbitals and hv with the hybrid orbitals.

ov={s, px,py,pz}
hv={h1,h2, h3, h4}

Define the matrix U which converts ov into hv.

U=(1/2){{1,1,1,1},{1,1,-1,-1},{1,-1,1,-1},{1,-1,-1,1}}

U/(1/2)//MatrixForm

Test that U gives the correct hybrid orbitals

hv=U.ov

If U is unitary the absolute value of its determinant should be 1 and the absolute values of its eigenvalues should be 1. See below that they are

Det[U]//Abs

Eigenvalues[U]//Abs