27 Rotational symmetry and selection rules
Although it is a statement that doesn’t need repeating at this point, I will anyways: symmetry is a very powerful concept in physics. Rotational symmetry in particular has been crucial in understanding the physics of spin, solutions for the hydrogen atom, and more. However, not all systems in nature are “spherical cows”; an enormous amount of interesting physics occurs in systems which are very much not rotationally invariant.
A simple example occurs when we put a (hydrogen, for simplicity) atom in a background electric field \vec{E}. The presence of the E-field breaks rotational symmetry, since it gives a preferred direction in space. As we saw before, a weak E-field that we can treat as a perturbation gives a Hamiltonian of the form \hat{W} = e\vec{E} \cdot \hat{\vec{r}}. This gives matrix elements between two atomic eigenstates of the form \bra{n'l'm'} \hat{\vec{r}} \ket{nlm}. If we want to compute higher-order energy perturbations for an arbitrary E-field direction, these are what we need to compute. These matrix elements are also crucial ingredients for time-dependent perturbation theory, which we haven’t studied yet: in particular, they appear in the calculation of transition probabilities for the radiative transition between \ket{n'l'm'} and \ket{nlm} involving emission or absorption of light.
If we want to study such a radiative transition for excited states, the number of matrix elements proliferates quickly: for example, if we are interested in 3d \rightarrow 2p in hydrogen, there are 5 possible initial states, 3 final states, and 3 independent components of \hat{\vec{r}}, leading to a total of 45 matrix elements to calculate! We know everything at this point to do the calculation: we already know our wavefunctions for the energy eigenstates take the form \left\langle \vec{r} | nlm \right\rangle = R_{nl}(r) Y_l^m(\theta, \phi). Going to position space, each matrix element is given by an integral of the form \bra{n'l'm'} \hat{\vec{r}} \ket{nlm} = \int_0^\infty dr\ r^3 R_{n'l'}^\star(r) R_{nl}(r) \int d\Omega (Y_{l'}^{m'})^\star(\theta, \phi) Y_l^m(\theta, \phi) \left(\frac{\vec{r}}{r}\right). We already know enough to just open up Mathematica and compute all 45 matrix elements by brute force. However, if you do this, you will find that 27 of the matrix elements are immediately zero from the angular integral; this is a consequence of rotational symmetry, and if we’re smart about it, we can predict which ones are zero and save ourselves a lot of work. If we’re really smart about exploiting rotational symmetry, as we’re about to learn, then we only need to compute one matrix element, and we can obtain all the others from rotations.
You may have noticed that we’ve already computed diagonal matrix elements of the form \bra{nlm} \hat{\vec{r}} \ket{nlm}, when looking at splitting of hydrogen eigenstates under an applied electric field. In those previous examples (Section 21.4), we just waved our hands and said “hydrogen is rotationally invariant, so let’s just spin our coordinates around to put the E-field in the z direction.”
In a sense, the whole idea of this section is exactly the same: we’re going to study how matrix elements transform when we rotate our coordinates around. But we can’t always just rotate our coordinate system and be done with it. For the example of atomic transitions, the emitted light will be in a random direction. We could take the z-axis to be in the direction of each photon emitted, but then to collect our results together we have to figure out how the different coordinate systems are all related - which is just the hard way of re-doing the systematic treatment of rotation we’re about to do.
More generally, even in systems that aren’t rotationally invariant at all, we can make use of rotational symmetry. To see how, let’s remind ourselves of one of our most useful approaches to time evolution in quantum mechanics: we can take an arbitrary state \ket{\psi} and expand it in energy eigenstates, \ket{\psi} = \sum_{n} c_n \ket{E_n}. Since the states \ket{E_n} evolve simply in time, the time-evolved version of this state is easily written down as \ket{\psi(t)} = \sum_{n} c_n e^{-iE_n t/\hbar} \ket{E_n}. Mathematically, this is simple because the \ket{E_n} are eigenstates of the Hamiltonian operator \hat{H} which gives us time evolution. Physically, we’re exploiting a symmetry: the states \ket{E_n} are invariant (up to a phase) under time translation, which is the symmetry generated by \hat{H}.
In exactly the same way, we can take any state and expand it in angular momentum eigenstates, \ket{\psi} = \sum_{j,m} d_{jm} \ket{jm}, and then the effect of rotating the state is captured by how rotation acts on the \ket{jm} states: \hat{\mathcal{D}}(R) \ket{\psi} = \sum_{j,m} d_{jm} \sum_{m'} {\mathcal{D}}^{(j)}_{m'm}(R) \ket{jm'} where the object appearing instead of a phase is the Wigner D-matrix, which we’ve met before, {\mathcal{D}}^{(j)}_{m'm}(R) \equiv \bra{jm'} \exp \left( \frac{-i(\hat{\vec{J}} \cdot \vec{n}) \phi}{\hbar} \right) \ket{jm}. Although the basic idea is the same as for time evolution, the structure is more complicated because rotation symmetry is more complicated; states \ket{jm} with different m mix with each other when we rotate. Still, different j values don’t mix (easily seen from the definition of the D-matrix). Whereas time evolution is a diagonal matrix in the energy eigenbasis, rotation is only block diagonal in the angular momentum eigenbasis. But block diagonal is better than nothing!
27.1 Rotating vector operators
Comparing to how we use energy eigenstates to deal with time translation, there is one other important complication for rotations. In most of what we’ve studied, the Hamiltonian itself and most of the other physical observables are time-independent, so that all of the time evolution is in the state vectors. On the other hand, trying to only study rotation-independent operators is very restrictive; in particular, any vector like the position vector \hat{\vec{r}} transforms in a non-trivial way under rotations.
Let’s start with an arbitrary vector-valued operator \hat{\vec{V}}. We know that if we take an expectation value \left\langle \hat{\vec{V}} \right\rangle = \bra{\psi} \hat{\vec{V}} \ket{\psi} in some quantum state, the result is a classical vector, which means that it transforms under rotation of our system as (writing it out in Cartesian components) \left\langle \hat{V}_i \right\rangle \rightarrow \sum_j R_{ij} \left\langle \hat{V}_j \right\rangle, where R_{ij} is a standard 3x3 rotation matrix. However, we also have a quantum description of rotations: rotating our system maps our state \ket{\psi} \rightarrow \hat{\mathcal{D}}(R) \ket{\psi}. So we also have the relation \left\langle \hat{V}_i \right\rangle \rightarrow \left\langle \hat{\mathcal{D}}^{-1}(R) \hat{V}_i \hat{\mathcal{D}}(R) \right\rangle. Since all of the above is true for an arbitrary quantum state, we have an identity between sets of operators, \hat{\mathcal{D}}^{-1}(R) \hat{V}_i \hat{\mathcal{D}}(R) = \sum_j R_{ij} \hat{V}_j.
This holds for any vector which is an operator. (Vectors which are not operators, like background fields \vec{E} or \vec{B}, will either just rotate classically or not rotate at all, depending on how we define our rotations to act - be careful with this!)
While this is an accurate definition of how an arbitrary rotation should act on a quantum vector operator, a more practically useful version of the formula comes from considering an infinitesimal rotation. We can plug in explicitly to the above equation for such a rotation and then solve, but this requires some fiddling with tensor indices or cross products, so let’s skip to the result:
A vector operator \hat{\vec{V}} with Cartesian components \hat{V}_i transforms under a rotation R according to \hat{\mathcal{D}}^{-1}(R) \hat{V}_i \hat{\mathcal{D}}(R) = \sum_j R_{ij} \hat{V}_j.
Equivalently, the components of the vector operator satisfy the commutation relations [\hat{V}_i, \hat{J}_j] = i\hbar \epsilon_{ijk} \hat{V}_k, where \hat{\vec{J}} is the total angular momentum operator (i.e. the generator of rotations.)
Prove that the first equation in the box above implies the commutation relation, studying an infinitesimal rotation by angle \delta \phi about an arbitrary axis \vec{n}.
Answer:
Taking the definition of \hat{\mathcal{D}}(R) from above, if we consider the infinitesimal version of this relation, rotating by \delta \phi on the left-hand side we have (1 + \frac{i (\hat{\vec{J}} \cdot \vec{n}) (\delta \phi)}{\hbar}) \hat{V}_i (1 - \frac{i (\hat{\vec{J}} \cdot \vec{n}) (\delta \phi)}{\hbar} ) = \sum_j R_{ij}(\delta \phi) \hat{V}_j \\ \Rightarrow \hat{V}_i - \frac{i(\delta \phi)}{\hbar} [\hat{V}_i, \hat{\vec{J}} \cdot \vec{n} ] = \sum_j R_{ij}(\delta \phi) \hat{V}_j. We can also expand out the right-hand side, since the matrix R_{ij} simplifies for rotation by an infinitesimal angle about a given axis; the result is, referring to Merzbacher, \sum_j R_{ij}(\delta \phi) \hat{V}_j \approx \hat{V}_i - (\delta \phi) (\vec{n} \times \hat{\vec{V}})_i = \hat{V}_i - (\delta \phi) \epsilon_{ijk} n_j \hat{V}_k. The \hat{V}_i term cancels on both sides, and by comparing the \mathcal{O}(\delta \phi) terms, we find that the commutator with \hat{\vec{J}} gives a cross product with the normal vector, which in terms of the Cartesian components implies the commutation identity [\hat{V}_i, \hat{J}_j] = i\hbar \epsilon_{ijk} \hat{V}_k.
(If you don’t want to mess around with cross products, you can arrive at the same result just by setting \vec{n} equal to the coordinate axis directions in turn and plugging in the corresponding rotation matrices.)
The first thing to note about this commutation relation is that it’s the same as the one we derived previously for [\hat{J}_i, \hat{J}_j], which makes sense; \vec{J} itself is an example of a vector operator. The commutation relation itself may look a bit obscure, but the physics is clear if you think about breaking it into components. For example, the relation [\hat{V}_z, \hat{J}_z] = 0 is just a statement that if we apply a rotation about the z axis, the z component of our vector operator isn’t affected by it. While this is a pretty simple observation, it is also immediately useful. For example, suppose that we’re studying a system and we want matrix elements of \hat{V}_z in terms of \ket{jm} eigenstates. We can use the commutator as follows: \bra{j'm'} [\hat{V}_z, \hat{J}_z] \ket{jm} = 0 \\ \hbar (m - m') \bra{j'm'} \hat{V}_z \ket{jm} = 0 applying \hat{J}_z to the eigenstates on the left and right. So we see immediately that the matrix elements of \hat{V}_z vanish unless m = m'. This is an example of a selection rule; a condition that matrix elements have to meet in order to be non-zero. (It’s also simple physics: since \hat{V}_z is invariant under z rotations, it can only have non-vanishing matrix elements if the bra and ket involved don’t transform differently under z rotations, which is in turn determined by the \hat{J}_z eigenvalue.)
Singling out the z direction for our rotation makes sense (since it’s what we’ve chosen as part of our basis-defining CSCO), but we should also look at the other components of \hat{\vec{V}}. Writing the remaining two commutators out explicitly, we have [\hat{V}_x, \hat{J}_z] = -i\hbar \hat{V}_y, \\ [\hat{V}_y, \hat{J}_z] = i\hbar \hat{V}_x. These relations don’t immediately lead to simple selection rules; rotation about z, as we would expect, mixes up the x and y components of our operator. However, we can find linear combinations of the two components that at least map back into themselves in the commutator: it’s easy to see that [\hat{V}_x \pm i \hat{V}_y, \hat{J}_z] = i \hbar (-\hat{V}_y \pm i \hat{V}_x) = \mp \hbar (\hat{V}_x \pm i \hat{V}_y). The commutator being proportional to the operator we started with means that these combinations don’t mix at all under z-axis rotations. They also look familiar; these are exactly the combinations that, for angular momentum \hat{\vec{J}}, defined the ladder operators \hat{J}_{\pm}.
27.1.1 Spherical basis
The three combinations of components of \hat{\vec{V}} which map to themselves under commutation with \hat{J}_z (and thus, under rotations around the z-axis) represent a change of basis from Cartesian basis to what is known as the spherical basis, written \hat{V}_q^{(1)}: \begin{aligned} \hat{V}_1^{(1)} &= -\frac{1}{\sqrt{2}} (\hat{V}_x + i \hat{V}_y), \\ \hat{V}_0^{(1)} &= \hat{V}_z, \\ \hat{V}_{-1}^{(1)} &= \frac{1}{\sqrt{2}} (\hat{V}_x - i\hat{V}_y). \end{aligned} The factors of 1/\sqrt{2} and the minus sign on the first vector ensure that this is an orthonormal basis. By convention, the subscript runs through q=\{-1, 0, +1\} in this basis; the superscript (1) is part of a larger notation for “spherical tensors” that we’ll explain fully later, but for now it helps us keep track of whether we’re working in spherical basis or Cartesian basis.
We can write out the new basis that we’re changing to more explicitly: the basis vectors are \begin{aligned} \hat{e}_1 &= -\frac{1}{\sqrt{2}} \left( \hat{x} + i\hat{y}\right) \\ \hat{e}_0 &= \hat{z} \\ \hat{e}_{-1} &= \frac{1}{\sqrt{2}} \left( \hat{x} - i\hat{y} \right). \end{aligned} (I’m using the hats now to denote unit vectors, although you can take these as operator definitions too.)
The complex nature of this basis makes it obvious that this is designed to work in quantum mechanics, although if you’ve studied the formulation of classical electromagnetism using complex numbers, you may recognize the \hat{e}_{\pm 1} vectors as the unit vectors describing left and right circular polarization of a light wave.
We can check that these vectors indeed define an orthonormal basis, as long as we’re careful to keep track of complex conjugation: \hat{e}_q^\star \cdot \hat{e}_{q'} = \delta_{qq'} We can, of course, expand any arbitrary spatial vector in terms of the spherical basis components, \vec{X} = \sum_{q} \hat{e}_q^\star X_q^{(1)} where X_q^{(1)} = \hat{e}_q \cdot \vec{X}. We could have also expanded in the vectors \hat{e}_q instead of the conjugates \hat{e}_q^\star; for our present purposes this won’t change anything as long as we’re consistent, so we’ll use the convention above.
Since the motivation for spherical basis was how rotations act (i.e. commutators with \hat{\vec{J}}), let’s write those out explicitly. From what we already found above, matching the notation, the commutation relation with \hat{J}_z can be written compactly in the form [\hat{J}_z, \hat{V}_q^{(1)}] = \hbar q \hat{V}_q^{(1)}, flipping the commutator to avoid having to write the extra minus sign on the right. We haven’t really looked at the other rotational axes yet, but we have given all of the necessary definitions: a bit more algebra leads to a result which is compact when written in terms of the ladder operators, [\hat{J}_{\pm}, \hat{V}_q^{(1)}] = \hbar \sqrt{ (1 \mp q) (2 \pm q)} \hat{V}_{q \pm 1}^{(1)}. (Effectively, both \hat{\vec{V}} and \hat{\vec{J}} are in spherical basis here, where the rotational symmetry is most strongly manifest.)
If we write out the components of the position vector \vec{r} in spherical basis, the components turn out to be simply related to l=1 spherical harmonics: r_q^{(1)} = \sqrt{\frac{4\pi}{3}} r Y_1^q(\theta, \phi).
If we plug in explicitly and then think about converting between Cartesian and spherical coordinates, we can easily verify each of the three components by hand: \sqrt{\frac{4\pi}{3}} r Y_1^1(\theta, \phi) = -r \sqrt{\frac{1}{2}} \sin \theta e^{i\phi} = - \frac{x + iy}{\sqrt{2}} \\ \sqrt{\frac{4\pi}{3}} r Y_1^0(\theta, \phi) = r \cos \theta = z, \\ \sqrt{\frac{4\pi}{3}} r Y_1^{-1}(\theta, \phi) = r \sqrt{\frac{1}{2}} \sin \theta e^{-i\phi} = \frac{x-iy}{\sqrt{2}}.
In fact, it makes perfect sense that specifically the l=1 spherical harmonics are what appear here. These spherical harmonics transform under the 3-dimensional representation of the rotation group SO(3). Vectors also have three components, and they also transform in the same 3-dimensional representation. As a reminder, “transforms in a 3-dimensional representation” means explicitly that \vec{r} \rightarrow R \vec{r}, where R is a 3-by-3 matrix which is a representation of SO(3) (i.e. has the right algebraic structure to be a rotation.) Spherical basis is just a choice which makes the matrix R diagonal for rotations about the z-axis; this was already the case by construction for spherical harmonics, since they’re built as eigenstates \ket{jm} of the \hat{J}_z operator. Explicitly, for a rotation by \alpha, \ket{jm} \rightarrow \hat{\mathcal{D}}(R) \ket{jm} = \exp \left( -\frac{i \hat{J}_z \alpha}{\hbar} \right) \ket{jm} which is manifestly a rotation by a diagonal 3x3 matrix if j=1. This is the abstract Hilbert-space version, but we can also make this explicitly about spherical harmonics, which can be written as Y_l^m(\theta, \phi) = \left\langle \theta, \phi | lm \right\rangle. where \ket{\theta, \phi} is a position eigenstate in spherical coordinates with arbitrary radius r. (We can also write the position bra as \bra{\vec{n}}, where \vec{n} is a unit vector pointing in the direction given by the spherical angles.) Substituting in above for the case l=1 gives us the result Y_1^m(\theta, \phi) \rightarrow e^{-im \alpha} Y_1^m(\theta, \phi) i.e. for a rotation about the z axis the m=0 harmonic doesn’t transform at all, while the m=\pm 1 harmonics pick up a phase. This is, in fact, obvious if you look up the form of the harmonics; Y_1^0 is independent of \phi, while Y_1^{\pm 1} go as e^{\pm i \phi}.
27.2 Selection rules for vector operators
One of the main applications of this rotational machinery is the derivation of selection rules, which are rules for which matrix elements are automatically zero due to symmetry. (We’ve made use of similar rules for e.g. parity before.) We found one of these above specifically for \hat{V}_z, but let’s do the whole vector in spherical basis this time. We’ll start with the commutator with \hat{J}_z, \bra{j'm'} [\hat{J}_z, \hat{V}^{(1)}_q] \ket{jm} = \hbar q \bra{j'm'} \hat{V}^{(1)}_q \ket{jm} \\ = \hbar (m' - m) \bra{j'm'} \hat{V}^{(1)}_q \ket{jm} with the first equality using the commutator formula from above, and the second from expanding the commutator out and letting \hat{J}_z act on the eigenkets. Subtracting the two sides, we see that the matrix element of \hat{V}^{(1)}_q must be zero, unless we satisfy the relation m' = m + q known as the m-selection rule. Note that since for a vector |q| \leq 1, we can only connect two states whose m values differ by at most 1 unit with a vector operator.
This is only part of the story, since it’s only using a subset of all possible rotations; we also have rotations about the x and y axes. We can use the commutation relation in terms of \hat{J}_{\pm} written above, but in general this leads to a much more complicated analysis for selection rules, since \hat{J}_{\pm} sometimes changes the state on the left or right, and sometimes just annihilates the state and gives zero. To illustrate the basic idea, let’s just ask the question: which vector matrix elements with \ket{00} can be non-zero? Starting from the commutation relation again, we have \bra{j'm'} [\hat{J}_{\pm}, \hat{V}_q^{(1)}] \ket{00} = \hbar \sqrt{ (1 \mp q) (2 \pm q)} \bra{j'm'} \hat{V}_{q \pm 1}^{(1)} \ket{00} \\ = \hbar \sqrt{(j' \pm m') (j' \mp m' + 1)} \bra{j', m' \mp 1} \hat{V}_q^{(1)} \ket{00}. where in the second case, we get the opposite signs from normal since the \hat{J}_{\pm} operators are acting backwards on the \bra{jm} bra. Let’s try to plug in some numbers for q to make this more explicit. If we let q = 0 and take the + branch, then we get the equation \sqrt{2} \bra{j'm'} \hat{V}_1^{(1)} \ket{00} = \sqrt{(j'+m')(j'-m'+1)} \bra{j', m'-1} \hat{V}_0^{(1)} \ket{00} while taking q=+1 and the - branch gives the equation \sqrt{2} \bra{j'm'} \hat{V}_0^{(1)} \ket{00} = \sqrt{(j'-m')(j'+m'+1)} \bra{j', m'+1} \hat{V}_1^{(1)} \ket{00} If we apply the m-selection rule to both equations, the value of m' is set: \sqrt{2} \bra{j'1} \hat{V}_1^{(1)} \ket{00} = \sqrt{(j'+1)j'} \bra{j'0} \hat{V}_0^{(1)} \ket{00} \\ \sqrt{2} \bra{j'0} \hat{V}_0^{(1)} \ket{00} = \sqrt{j' (j'+1)} \bra{j'1} \hat{V}_1^{(1)} \ket{00} Plugging one equation into the other gives us the result 2 = j'(j'+1) which is solved only by j'=1. The other possible case to worry about is j'=0, but it’s straightforward to see that if we try j'=0 all the matrix elements we’re trying to compute are just zero automatically.
This is a special case of the more general selection rule, known as the triangular selection rule, which for a vector operator specifically states that \bra{j'm'} \hat{V}_q^{(1)} \ket{jm} is zero unless |j-1| \leq j' \leq j+1. (so indeed, we see that j=0 requires j'=1.)
One way to think of this physically is that with an operator in the middle, we can have non-zero overlap between states with different j, since the rotation of the operator can compensate the difference in rotation between the states. For example, \ket{00} doesn’t rotate at all, while \ket{10} does (it’s a vector component.) But if we have \bra{10} \hat{\vec{r}} \ket{00}, the part of the overall matrix element for which the rotation cancels between the vector operator and the j=1 state will still preserve rotational symmetry overall, and so can be non-zero.
There’s a much more mechanical approach to finding selection rules, which involves going to position space and then using algebraic identities for spherical harmonics. As an explicit example, for the dipole transition matrix elements putting \hat{\vec{r}} into spherical basis gives \bra{n'l'm'} \hat{r}_q \ket{nlm} = \int_0^{\infty} dr\ r^3 R_{n'l'}^\star(r) R_{nl}(r) \sqrt{\frac{4\pi}{3}} \\ \times \int d\Omega (Y_{l'}^{m'})^\star(\theta, \phi) Y_1^q (\theta, \phi) Y_l^m(\theta, \phi). The rotational selection rules are encoded in whether the integral over three spherical harmonics vanishes or not, but this doesn’t seem like a very clear condition. It is made much clearer by a powerful algebraic identity for decomposing products of two spherical harmonics as a sum over single harmonics, which leads to the integral formula \int d\Omega (Y_l^m)^\star (\theta, \phi) Y_{l_1}^{m_1}(\theta, \phi) Y_{l_2}^{m_2}(\theta, \phi) \\ = \sqrt{\frac{(2l_1+1)(2l_2+1)}{4\pi (2l+1)}} \left\langle l_1 l_2; 0 0 | l_1 l_2; l 0 \right\rangle \left\langle l_1 l_2; m_1 m_2 | l_1 l_2; l m \right\rangle. Here the inner products are Clebsch-Gordan coefficients. This is a beautiful result of getting around calculus using algebraic relations, and in fact it can be derived from our Hilbert space formalism quite nicely; however, the derivation is fairly technical, so I’ve decided to relocate it to an appendix (Appendix C) in these notes.
One nice thing about the Clebsch-Gordan based approach is since we already know the properties for C-G coefficients, we can just use them to read off the selection rules - in particular, the triangular selection rule is easy to see from this formula, whereas we were only able to easily get through the j=m=0 special case above.
27.3 Example: selection rules for atomic electric dipole transitions
Let’s make use of the results we’ve found so far to return to the study of electric dipole transitions in atoms. We haven’t studied the detailed derivation in time-dependent perturbation theory, but I’ll borrow one of the key results, which is that for spontaneous emission of light from an atom, the transition rate from initial state i to final state f is proportional to the square of the matrix element \bra{\alpha_f l_f m_f} \hat{\vec{r}} \ket{\alpha_i l_i m_i} where \alpha is a common notation standing in for “other quantum numbers”; we’ll look at this case first and then elaborate on what \alpha is exactly, but it doesn’t matter as long as they aren’t affected by rotations.
Rewriting the operator in spherical tensor basis, \bra{\alpha_f l_f m_f} \hat{r}_q^{(1)} \ket{\alpha_i l_i m_i}, we can immediately apply the two selection rules from above: |l_i-1| \leq l_f \leq l_i + 1; \\ m_f = m_i + q.
There is one more important symmetry constraint to consider, which is parity. So long as parity is a good symmetry of the underlying Hamiltonian (true for hydrogenic atoms), then we can use the fact that the position vector \hat{\vec{r}} is parity odd (meaning \hat{P} \hat{\vec{r}} \hat{P} = -\hat{\vec{r}}) to derive the relation \bra{\alpha_f l_f m_f} \hat{\vec{r}} \ket{\alpha_i l_i m_i} = - \bra{\alpha_f l_f m_f} \hat{P} \hat{\vec{r}} \hat{P} \ket{\alpha_i l_i m_i} \\ = (-1)^{l_f + l_i + 1} \bra{\alpha_f l_f m_f} \hat{\vec{r}} \ket{\alpha_i l_i m_i} using the fact that \ket{lm} eigenstates transform under parity as \hat{P} \ket{lm} = (-1)^l \ket{lm}. Comparing the two sides of this equation, we see that if l_i = l_f then the matrix element must be zero. (This proves as a side effect that if i and f are the same state, the dipole matrix element is zero, which tells us that no energy eigenstate of hydrogen can have a permanent electric dipole moment!) So putting angular momentum and parity together, the overall combined selection rule is l_f = l_i \pm 1.
Now let’s talk about the other quantum numbers for the specific example of atomic transitions in hydrogen. For labels like the principal quantum number n, there is no rotational or parity dependence, so the argument above doesn’t change; the values of n_i and n_f will affect the value of the matrix element, but not the selection rules. However, we should be more careful with the electron spin, which we’ve also ignored. We know that in our initial hydrogen Hamiltonian, the wavefunctions should really be written as \ket{nlmm_s}. Since the perturbation doesn’t couple to spin at all, we have the further selection rule \Delta m_s = 0. But this is a choice we’ve made for what the unperturbed Hamiltonian \hat{H}_0 is! If we want to study the hydrogen atom including fine-structure effects, or hydrogen in a strong magnetic field, we know a better choice of basis is total angular momentum \vec{J} = \vec{L} + \vec{S}, in which case the basis states are \ket{nljm}.
In the latter case, we need to be a bit more careful about what parity implies; the short version of the story is that parity still requires \Delta l = \pm 1. But it is possible to have cases where \Delta j = 0 without having \Delta l = 0. An example is the 2P_{1/2} \rightarrow 1S_{1/2} transition, with j=1/2 in both cases. Using Clebsch-Gordan coefficients, the initial state includes for example \ket{j=\tfrac{1}{2}, m=\tfrac{1}{2}} = \sqrt{\frac{2}{3}} \ket{m_l=1, m_s=-\tfrac{1}{2}} - \sqrt{\frac{1}{3}} \ket{m_l=0, m_s=+\tfrac{1}{2}}. So in the transition matrix, there is a piece with \Delta l \neq 0 that respects parity, even with \Delta j = 0. The only exception is the case where the transition is j=0 \rightarrow j=0, for which there is no mixing. (This case doesn’t show up in hydrogen, since the lone electron’s spin guarantees j \geq \tfrac{1}{2}.)
Putting these statements together, even before we compute anything quantitative, we have an important physics result:
Working in the dipole approximation, spontaneous emission transitions in hydrogen are only allowed between states for which \Delta l = \pm 1 and \Delta m_s = 0 (in the \ket{nlmm_s} basis), or |\Delta j| \leq 1 (in the \ket{nljm} basis.) In either case, we must also have |\Delta m| \leq 1.
(“Dipole approximation” means the leading order in a multipole expansion; more discussion can be found in Section 31.4 if you want to jump ahead.) Not a bad result, especially since we haven’t even set up the details of how to actually calculate the spontaneous emission rate for the cases where it’s allowed!
For a concrete test of this result, we can consider transitions down to the 1s ground state from the n=3 energy level of hydrogen. Our selection rule says that only the 3p \rightarrow 1s transition is allowed, while 3d \rightarrow 1s is forbidden (through the dipole interaction) due to rotational invariance. (So is 3s \rightarrow 1s by parity.) If we look up hydrogen transition lines in the NIST Atomic Spectra Database, taking a ratio so we don’t have to worry too much about units or conventions, we find that the ratio of decay rates from experiment is \frac{\Gamma(3d \rightarrow 1s)}{\Gamma(3p \rightarrow 1s)} \approx 9.5 \times 10^{-7}. In other words, the 3p \rightarrow 1s transition is about 1 million times more likely than 3d \rightarrow 1s! It isn’t exactly zero because the selection rule only forbids the electric dipole transition; higher multipoles can mediate the 3d \rightarrow 1s decay, so it still happens, just at a very suppressed rate.