Angular momentum is classically defined as:

$$\mathbf L=\mathbf r\times\mathbf p$$

The quantum-mechanical operator is defined the same way, by just replacing the two quantities with their respective operators.

$$ \hat L_i =\hat x_j\hat p_k-\hat x_k\hat p_j =-i\hbar\left(x_j\frac{\partial}{\partial x_k}-x_k\frac{\partial}{\partial x_j}\right) $$

The three components of angular momentum are not compatible quantities with each other, since they don't commute.

$$\left[\hat L_i,\hat L_j\right]=i\hbar\hat L_k\Leftrightarrow\mathbf L\times\mathbf L=i\hbar\mathbf L$$

In order to fully describe the angular momentum of a quantum-mechanical system we need another quantity, which would be compatible. This quantity is the magnitude of the total angular momentum squared, $\hat L^2=\left|\hat{\mathbf{L}}\right|^2$.

Indeed, this quantity commutes with all the components.

$$\left[\hat L^2,\hat L_i\right]=0$$

Completely describing the angular momentum of a quantum system we use two operators, $\hat L^2$ for the magnitude and $\hat L_z$ for its projection on the $z$ axis, which is easier to play with in spherical coordinates.

Also, both of those operators commute with the Hamiltonian, and thus have common eigenfunctions. In the end, we have two eigenvalue problems.

$$\hat L^2\left|l,m\right>=\hbar^2\lambda\left|l,m\right>=\hbar^2l(l+1)\left|l,m\right>$$ $$\hat L_z\left|l,m\right>=\hbar\left|l,m\right>$$

The eigenvalue $\lambda=l(l+1)$ is written this way because it is easier to handle mathematically, as will be shown later on, without affecting the problem in any way. The quantum numbers $l,m$ are integers and define the so called Spherical Harmonics $Y_{lm}(\theta,\phi)$, which are the angular part of the wavefunctions $\psi_{lm}(r,\theta,\phi)=R(r)\cdot Y_{lm}(\theta,\phi)$.

The existence of a central potential $V(\mathbf r)=V(-\mathbf r)$ causes the quantization of the radial part of the wavefunction. Together with spin they form a quartet of quantum numbers ($n,l,m_l,m_s$) which fully defines the state of a particle around the central potential (e.g. a single electron around a nucleus).

The eigenvalue problem of $\hat L_z$ is similar to that of the Hamiltonian for the Quantum Harmonic Oscillator. We are looking for the raising and lowering operators $\hat L_\pm$ that will move us from the state $\left|l,m\right>$ to the state $\left|l,m\pm1\right>$.

The only difference is that $l$ should not be affected, meaning that the eigenvalue equation for $\hat L^2$ will still hold. $l$ represents the magnitude of the angular momentum, while $m$ its projection on an axis. It is a natural requirement that any action on $m$, which corresponds to axes rotation, won't affect $l$ in any way.

We start by factoring the magnitude:

$$L^2=L_x^2+L_y^2+L_z^2=L_z^2+\left(L_x+iL_y\right)\left(L_x-iL_y\right)$$

By comparing with the algebra of the harmonic oscillator (factoring $\hat a^\dagger\hat a$) we practically guess the existence of the operators:

$$\hat L_+=\hat L_x+i\hat L_y$$ $$\hat L_-=\hat L_x-i\hat L_y$$

It's obvious that $\hat L_\pm^\dagger=\hat L_\mp$. We now need to show that these two operators have the required properties:

• $\left[\hat L_\pm,\hat L^2\right]=0$
• $\left[\hat L_z,\hat L_\pm\right]=\pm\hbar\hat L_\pm$, which leads to $\hat L_\pm\left|l,m\right>\propto\left|l,m\pm1\right>$

The $\hat L_\pm$ operator is a linear combination of $\hat L_x, \hat L_y$ which commute with $\hat L^2$. So $\hat L_\pm$ commutes too.

$$ \left[\hat L_\pm,\hat L^2\right] =\left[\hat L_x,\hat L^2\right] \pm i\left[\hat L_y,\hat L^2\right] =0\pm i0=0 $$

$$\begin{align} \left[\hat L_z,\hat L_\pm\right] &= \left[\hat L_z,\hat L_x\pm i\hat L_y\right] \\ &= \left[\hat L_z,\hat L_x\right]\pm i\left[\hat L_z,\hat L_y\right] \\ &= \left[\hat L_z,\hat L_x\right]\mp i\left[\hat L_y,\hat L_z\right] \\ &= i\hbar\left(\hat L_y\mp i\hat L_x\right) \\ &= \hbar\left(i\hat L_y\pm\hat L_x\right) \\ &= \pm\hbar\left(\hat L_x\pm i\hat L_y\right) \\ &= \pm\hbar\hat L_\pm \end{align}$$ The quantum number $m$ is the eigenvalue of $\hat L_z$, so in order to see if $\hat L_\pm$ indeed takes us from state $\left|l,m\right>$ to state $\left|l,m\pm1\right>$ we can act with it on state $\hat L_\pm\left|l,m\right>$.

$$\begin{align} \hat L_z\left(\hat L_\pm\left|l,m\right>\right) &= \left(\hat L_z\hat L_\pm\right)\left|l,m\right> \\ &= \left(\hat L_\pm\hat L_z\pm\hbar\hat L_\pm\right)\left|l,m\right> \\ &= \hat L_\pm\left(\left(\hat L_z\pm\hbar\right)\left|l,m\right>\right) \\ &= \hat L_\pm\left(\hbar m\pm\hbar\right)\left|l,m\right> \\ &= \hbar\left(m\pm1\right)\left(\hat L_\pm\left|l,m\right>\right) \end{align}$$ The eigenvalue is $\hbar\mu=\hbar\left(m\pm1\right)$, so the operator $\hat L_\pm$ indeed takes us from state $\left|l,m\right>$ to state $\left|,m\pm1\right>$.

We have shown that there is a relation between the two states, $\hat L_\pm\left|l,m\right>\propto\left|l,m\pm1\right>$. In order to find an exact equation we can use the fact that wavefunctions are normalized.

But first we need to calculate two operators that show up later on during the process.

$$\begin{align} \hat L_-\hat L_+ &= \left(\hat L_x-i\hat L_y\right)\left(\hat L_x+i\hat L_y\right) = \hat L_x^2+\hat L_y^2+i\left(\hat L_x\hat L_y-\hat L_y\hat L_x\right) \\ &= \hat L^2-\hat L_z^2+i\left[\hat L_x,\hat L_y\right] = \hat L^2-\hat L_z^2+i\left(i\hbar\hat L_z\right) \\ &= \hat L^2-\hat L_z^2-\hbar\hat L_z = \hat L^2-\hat L_z\left(\hat L_z+\hbar\right) \end{align}$$ $$\begin{align} \hat L_+\hat L_- &= \left(\hat L_x+i\hat L_y\right)\left(\hat L_x-i\hat L_y\right) = \hat L_x^2+\hat L_y^2-i\left(\hat L_x\hat L_y-\hat L_y\hat L_x\right) \\ &= \hat L^2-\hat L_z^2-i\left[\hat L_x,\hat L_y\right] = \hat L^2-\hat L_z^2-i\left(i\hbar\hat L_z\right) \\ &= \hat L^2-\hat L_z^2+\hbar\hat L_z = \hat L^2-\hat L_z\left(\hat L_z-\hbar\right) \end{align}$$ If we write $\hat L_\pm\left|l,m\right>=c_{lm}^\pm\left|l,m\pm1\right>$, we then can find the constant by taking the square magnitudes, since $\left\|\left|l,m\pm1\right>\right\|^2=1$ by definition.

$$\begin{align} \left|c_{lm}^\pm\right|^2 &= \left\|\hat L_\pm\left|l,m\right>\right\|^2 \\ &= \left(\hat L_\pm\left|l,m\right>\right)^\dagger\left(\hat L_\pm\left|l,m\right>\right) \\ &= \left< l,m\left|\hat L_\mp\hat L_\pm\right|l,m\right> \\ &= \left< l,m\left|\hat L^2-\hat L_z\left(\hat L_z\pm\hbar\right)\right|l,m\right> \\ &= \left< l,m\left|\hbar^2\left[l(l+1)\right]-\hbar m\left(\hbar m\pm\hbar\right)\right|l,m\right> \\ &= \hbar^2\left[l(l+1)-m(m\pm1)\right]\left< l,m|l,m\right> \\ &= \hbar^2\left[l(l+1)-m(m\pm1)\right] \end{align}$$ We now have fully defined how the angular momentum raising and lower operators act on a state. $$\boxed{\hat L_+\left|l,m\right> = \hbar\sqrt{l(l+1)-m(m+1)}\left|l,m+1\right>}$$ $$\boxed{\hat L_-\left|l,m\right> = \hbar\sqrt{l(l+1)-m(m-1)}\left|l,m-1\right>}$$ ==== Accepted Values of $m$ and $l$ ====

There's a natural limit for the eigenvalue $m$. It isn't possible for the projection of the angular momentum on the $z$ axis to be of greater magnitude than the total angular momentum. Thus $m$ is bounded, and its minimum and maximum values are given by requiring the normalization term to be zero in the two formulas above.

$$\hat L_+\left|l,m_{max}\right> = 0\Rightarrow m_{max}=l$$ $$\hat L_-\left|l,m_{min}\right> = 0\Rightarrow m_{min}=-l$$

By proving that $\left[\hat L_z,\hat L_\pm\right]=\hbar\hat L_\pm$ we showed that the eigenvalues $\hbar m$ of $\hat L_z$ are equally spaced multiples of $\hbar$ and that the number $m$ is integer. The allowed values are all integer numbers in the interval $\left[-l,l\right]$. This interval has length $2l$, which also integer. So $l$ can be any integer divided by $2$, meaning an integer or half-integer.

The two paragraphs above is the reason we write the eigenvalue $\lambda$ from the problem $\hat L^2\left|l,m\right>=\lambda\left|l,m\right>$ as $\lambda=l(l+1)$. Without this formulation, finding a relation between the limits of $m$ and $\lambda$ would be much more difficult and less usable.

How is it possible for $l$ to theoretically be able to take half-integer values, while it's forced to be strictly integer since $m$ is strictly integer? During the analysis for $l$ and $\hat L^2$ above we silently assumed that by angular momentum we mean orbital angular momentum, as an equivalent to classical angular momentum. We can follow the exact same procedure for a more general angular momentum $j$, which include both the orbital angular momentum with $m=m_l$ and spin $m_s$.

The quantum number of spin $m_s$ can accept half-integer values, and ensures that all theoretically possible values of $j$ can be reached. Another impressive property is that it splits up all elementary particles in two categories, fermions and bosons, depending on if they have half-integer or strictly integer values of spin.

The momentum operator $\hat p^2=-\hbar^2\nabla^2$ represents the total momentum. When written in spherical coordinates, it is easy to spot a radial and an angular part.

$$\hat p^2 =-\hbar^2\nabla^2 =\underbrace{ -\hbar^2\left(\frac{\partial ^2}{\partial r^2}+\frac2r\frac{\partial}{\partial r}\right) }_{\large\hat p_r^2} + \frac{1}{r^2} \underbrace{ \left[-\hbar^2\left(\frac{1}{\sin\theta}\frac{\partial}{\partial\theta}\left(\sin\theta\frac{\partial}{\partial\theta}\right) +\frac{1}{\sin^2\!\theta}\frac{\partial^2}{\partial\phi^2}\right)\right] }_{\large\hat L^2}$$

This is how we can guess the form of $\hat L^2$. We can also try to guess how $\hat L_z$ looks like. Rotation around the $z$ axis is described exclusively from the azimuthal angle $\phi$. Taking the $\phi$ dependence from the formula above and ignoring the $\sin\theta$ scaling factor of the spherical Laplacian, we have a guess.

$$\hat L_z=-i\hbar\frac{\partial}{\partial\phi}$$

One can prove that this is indeed the operator if one starts from $\hat L_z=-i\hbar(x\partial_y-y\partial_x)$ and applies the chain rule on $\partial_x,\partial_y$ to end up with derivatives with respect to $r,\theta$ and $\phi$.

It was mentioned before that the two eigenvalue problems, $\hat L^2\left|l,m\right> =\hbar^2l(l+1)\left|l,m\right>$ and $\hat L_z\left|l,m\right> =\hbar m\left|l,m\right>$ must be solved at the same time. Thus, they have common eigenfunctions $Y_{lm}(\theta,\phi)$, which are called Spherical Harmonics.

In order to start solving the differential equations we initially try separating the variables and writing $Y_{lm}(\theta,\phi)=\Theta(\theta)\Phi(\phi)$ and start from the second eigenvalue problem.

$$\begin{align} \hat L_zY_{lm}(\theta,\phi) &= \hbar mY_{lm}(\theta,\phi)\\ -i\hbar\frac{\partial}{\partial\phi}\left(\Theta(\theta)\Phi(\phi)\right) &= \hbar m\Theta(\theta)\Phi(\phi) \\ \Theta(\theta)\frac{\partial\Phi}{\partial\phi} &= im\Theta(\theta)\Phi(\phi) \\ \frac{\partial\Phi}{\partial\phi} &=im\Phi(\phi) \end{align}$$ $$\boxed{\Phi(\phi)=e^{im\phi}}$$ Now, for the first eigenvalue problem: $$ \hat L^2Y_{lm}(\theta,\phi) = \hbar^2l(l+1)Y_{lm}(\theta,\phi) $$ $$\begin{align} -\hbar^2\left(\frac{1}{\sin\theta}\frac{\partial}{\partial\theta}\left(\sin\theta\frac{\partial}{\partial\theta}\right) +\frac{1}{\sin^2\!\theta}\frac{\partial^2}{\partial\phi^2}\right) \left(\Theta(\theta)e^{im\phi}\right) &= \hbar^2l(l+1)\Theta(\theta)e^{im\phi} \\ -\left(\frac{e^{im\phi}}{\sin\theta}\frac{d}{d\theta}\left(\sin\theta\frac{d\Theta}{d\theta}\right) +\frac{\Theta(\theta)}{\sin^2\!\theta}\frac{d^2\left(e^{im\phi}\right)}{d\phi^2}\right) &= l(l+1)\Theta(\theta)e^{im\phi} \\ \frac{m^2}{\sin^2\!\theta}\Theta(\theta)-\frac{1}{\sin\theta}\frac{d}{d\theta} \left(\sin\theta\frac{d\Theta}{d\theta}\right) &= l(l+1)\Theta(\theta) \end{align}$$ By doing a variable change of $\xi=\cos\theta, d\xi=-\sin\theta d\theta$ we can simplify the equation.

$$ \frac{m^2}{1-\xi^2}\Theta(\xi) -\frac{d}{d\xi}\left(\left(1-\xi^2\right)\frac{d\Theta}{d\xi}\right) = l(l+1)\Theta(\xi) $$

$$ \frac{d}{d\xi}\left(\left(1-\xi^2\right)\frac{d\Theta}{d\xi}\right) +\left(l(l+1)-\frac{m^2}{1-\xi^2}\right)\Theta(\xi) = 0 $$

This is the General (or Associated) Legendre Equation, and some of its solutions are the Associated Legendre Polynomials $P_l^m(x)$. From all the solutions, only the Associated Legendre Polynomials can describe this physical system. They are defined with respect to the normal, run of the mill Legendre Polynomials $P_l(x)$¹⁾. This is their historical definition, and only later on they were characterized as solutions to the differential equation above.

$$P_l^m(x)=(-1)^m\left(1-x^2\right)^{m/2}\frac{d^m}{dx^m}P_l(x)$$

The Legendre Polynomials are the solutions to the Legendre differential equation,

$$\frac{d}{dx}\left(\left(1-x^2\right)\frac{dP_l}{dx}\right) + l(l+1)P_l(x) = 0$$

and are generated by the Rodrigues formula.

$$P_l(x)=\frac{1}{2^ll!}\frac{d^l}{dx^l}\left[\left(x^2-1\right)^l\right]$$

By replacing the Rodrigues formula on the definition of the associated polynomials we get:

$$P_l^m(x)=(-1)^m\frac{1}{2^ll!}\left(1-x^2\right)^{m/2}\frac{d^{l+m}}{dx^{l+m}}\left[\left(x^2-1\right)^l\right]$$

The formula above stands only for positive values of $m$. It can expanded to negative values through the formula below:

$$P_l^{-m}(x)=(-1)^m\frac{(l-m)!}{(l+m)!}P_l^m,\quad m>0$$

After taking all of the above into account, and some calculations, one can find the exact form of the spherical harmonics.

$$\boxed{ Y_{lm}(\theta,\phi) = (-1)^m\frac{1}{2^ll!} \sqrt{\frac{2l+1}{4\pi}\frac{(l-m)!}{(l+m)!}} \left(1-\cos^2\!\theta\right)^{\left|m\right|/2} e^{im\phi} \left(\frac{d}{d(\cos\theta)}\right)^{l+\left|m\right|} \left[\left(\cos^2\!\theta-1\right)^l\right] }$$

The maximum degree of the polynomial $\left(\cos^2\!\theta-1\right)^l$ is $2l$, so the constraint regarding the values of $m$ is also obeyed. For any $\left|m\right|>l$ there will be the differentiation of a constant which will immediately zero the spherical harmonic.

The raising and lower operators $\hat L_\pm$ obey all physical constraints for the spherical harmonics as they have been defined. So an easier way to find said harmonics is to solve the differential equation for one of the termination conditions and find all other states one by one.

$$\hat L_+Y_{l,l}(\theta,\phi)=0,\quad \hat L_-Y_{l,-l}(\theta,\phi)=0$$

First, we have to find formulas of $\hat L_x, \hat L_y$ with respect to $\theta, \phi$ the same way we found $\hat L_z=-i\hbar\partial_\phi$.

$$\hat L_x = i\hbar\left(\sin\phi\frac{\partial}{\partial\theta} + \frac{\cos\phi}{\tan\theta}\frac{\partial}{\partial\phi}\right), \quad \hat L_y = i\hbar\left(-\cos\phi\frac{\partial}{\partial\theta}+\frac{\sin\phi}{\tan\theta}\frac{\partial}{\partial\phi}\right)$$

Calculating the raising and lowering operators through $\hat L_\pm=\hat L_x\pm i\hat L_y$ we find:

$$\hat L_+ =\hbar e^{i\phi}\left(\frac{\partial}{\partial\theta}+i\frac{1}{\tan\theta}\frac{\partial}{\partial\phi}\right), \quad \hat L_-=\hbar e^{-i\phi}\left(-\frac{\partial}{\partial\theta}+i\frac{1}{\tan\theta}\frac{\partial}{\partial\phi}\right)$$

The lowest state will have the form $Y_{l,-l}=\Theta_l(\theta)e^{-il\phi}$. Solving the termination condition we have:

$$\hat L_-Y_{l,-l}(\theta,\phi)=0$$ $$\hbar e^{-i\phi}\left(-e^{-il\phi}\frac{d\Theta}{d\theta}+i\frac{\Theta(\theta)}{\tan\theta}(-il)e^{-il\phi}\right)=0$$ $$-e^{il\phi}\left(\frac{d\Theta}{d\theta}-\frac{l}{\tan\theta}\Theta(\theta)\right)=0$$

We then perform a $u$ substitution, $u=\sin\theta, du=\cos\theta d\theta$.

$$\cos\theta\frac{d\Theta}{du}-\frac{l\cos\theta}{u}\Theta(u)=0$$ $$\frac{d\Theta}{du}-\frac{l}{u}\Theta(u)=0$$ $$\Theta(u)=N_lu^l$$ $$\Theta(\theta)=N_l\sin^l\!\theta$$

Due to normalization we require that:

$$\oint\left|Y_{lm}(\theta,\phi)\right|^2\ d\Omega=1$$ $$\require{cancel}\oint N_l^2\sin^{2l}\!\theta\cancelto{1}{\left|e^{im\phi}\right|^2}\ d\Omega = 1$$ $$N_l=\frac{1}{2^ll!}\sqrt{\frac{(2l+1)!}{4\pi}}$$

The lowest state is:

$$\boxed{Y_{l,-l}(\theta,\phi)=(-1)^l\frac{1}{2^ll!}\sqrt{\frac{(2l+1)!}{4\pi}}\ \sin^l\!\theta\ e^{-il\phi}}$$

and the rest can be found through the raising operator.

$$\boxed{Y_{l,m+1}(\theta,\phi)=\frac{1}{\hbar\sqrt{l(l+1)-m(m+1)}}\hat L_+Y_{l,m}(\theta,\phi)}$$

Below are listed the first few spherical harmonics:

$$Y_{00}(\theta,\phi)=\frac{1}{\sqrt{4\pi}}$$ $$Y_{10}(\theta,\phi)=\sqrt{\frac{3}{4\pi}}\ \cos\theta$$ $$Y_{1,\pm1}(\theta,\phi)=\mp\sqrt{\frac{3}{8\pi}}\ \sin\theta\ e^{\pm i\phi}$$ $$Y_{20}(\theta,\phi)=\sqrt{\frac{5}{16\pi}}\ (3\cos^2\!\theta-1)$$ $$Y_{2,\pm1}(\theta,\phi)=\mp\sqrt{\frac{15}{8\pi}}\ \cos\theta\sin\theta\ e^{\pm i\phi}$$ $$Y_{2,\pm2}(\theta,\phi)=\sqrt{\frac{15}{32\pi}}\ \sin^2\!\theta\ e^{\pm i2\phi}$$

Calculating quantities is trivial if their operator can be described as a function of the raising and lowering operators. Take for example the $x$ and $y$ components of the angular momentum.

$$\hat L_x=\frac{\hat L_++\hat L_-}{2},\quad\hat L_y=\frac{\hat L_+-\hat L_-}{2}$$

For example below are the calculations for $\left<\hat L_x\right>,\left<\hat L_y\right>,\left<\hat L_x^2\right>, \left<\hat L_y^2\right>$.

$$\begin{align} \left<\hat L_x\right> &= \left< l,m\left|\hat L_x\right|l,m\right> \\ &\propto \left< l,m\left|\hat L_+\right|l,m\right> + \left< l,m\left|\hat L_-\right|l,m\right> \\ &\propto \left< l,m|l,m+1\right> + \left< l,m|l,m-1\right> \\ &= 0 \end{align}$$ The system is $\phi$-independent, since it has a central potential. This means symmetry around the $z$ axis and thus $\left<\hat L_x\right>=\left<\hat L_y\right> =0$. Similarly $\left<\hat L_x^2\right>=\left<\hat L_y^2\right>$. Then:

$$ \left<\hat L^2\right> = \left<\hat L_x^2+\hat L_y^2+\hat L_z^2\right> = \left<\hat L_x^2\right> + \left<\hat L_y^2\right> + \left<\hat L_z^2\right> = 2\left<\hat L_x^2\right> + \left<\hat L_z^2\right> $$

$$\begin{align} \left<\hat L_x^2\right> = \left<\hat L_y^2\right> &= \frac12\left(\left<\hat L^2\right> - \left<\hat L_z^2\right> \right) \\ &= \frac12\left(\left< l,m\left|\hat L^2\right|l,m\right> -\left< l,m\left|\hat L_z^2\right|l,m\right>\right) \\ &= \frac12\Big(l(l+1)\left< l,m|l,m\right> - m^2\left< l,m|l,m\right>\Big) \\ &= \frac{l(l+1)-m^2}{2} \end{align}$$ This way we can find the dispersion around the expected value. $$\require{cancel} \Delta L_x = \Delta L_y = \sqrt{\left<\hat L_x^2\right> - \cancelto{0}{\left<\hat L_x\right>^2}} = \sqrt{\frac{l(l+1)-m^2}{2}} $$