The distribution of eigenvalues of the wave equation in a bounded domain is known as Weyl's problem. We describe several computational projects related to the cumulative state number, defined as the number of states having wavenumber up to a maximum value. This quantity and its derivative, the density of states, have important applications in nuclear physics, degenerate Fermi gases, blackbody radiation, Bose–Einstein condensation, and the Casimir effect. Weyl's theorem states that, in the limit of large wavenumbers, the cumulative state number depends only on the volume of the bounding domain and not on its shape. Corrections to this behavior are well known and depend on the surface area of the bounding domain, its curvature, and other features. We describe several projects that allow readers to investigate this dependence for three bounding domains—a rectangular box, a sphere, and a circular cylinder. Quasi-one- and two-dimensional systems can be analyzed by considering various limits. The projects have applications in statistical mechanics, but can also be integrated into quantum mechanics, nuclear physics, or computational physics courses.

## REFERENCES

*Spectra of Finite Systems: A Review of Weyl's Problem, the Eigenvalue Distribution of the Wave Equation for Finite Domains and its Applications on the Physics of Small Systems*

*Über die partielle differentialgleichung u+k*

^{2}y=0 und deren auftreten in der mathematischen physik*An Introduction to Nuclear Physics*

*Quantum Wells, Wires and Dots: Theoretical and Computational Physics of Semiconductor Nanostructures*

*An Introduction to Thermal Physics*

This expression is for bosons; an extra factor of two is usually included when considering a degenerate Fermi gas (Ref. 20).

See, for example, Refs. 8, 22, and 23. Results related to Eq. (8) are sometimes given in terms of the density of states (rather than the cumulative state number). Also, many of the formal results related to the distribution of eigenvalues of Eq. (1) are derived by solving related problems, such as the heat diffusion equation. See Ref. 22 for a discussion of this topic.

See Ref. 1 for a similar expression for the cubic case.

For a complementary description, we refer interested readers to Ref. 12.

*Data Reduction and Error Analysis for the Physical Sciences*

We can also use embedded Do loops; it is tricky to determine the upper limits for each of the loops.

It is possible to miss degeneracies due to rounding or other issues, which can be checked after the fact by comparing subsequent *κ* values to determine if they are very close to each other.

*Mathematical Methods in the Physical Sciences*

More mathematically sophisticated approaches, which correctly account for the averaging procedure, are described in Refs. 1, 8, 22, and 23. We agree with Lambert that, “As an alternative to the rather inaccessible general proof, it seems worthwhile to give explicit proofs for the sphere and the cylinder since these shapes…will tend to make the general result more plausible (Ref. 30).” Reference 22 gives a beautiful introduction to a method involving the diffusion equation, Green's functions and a Laplace transform.

See also Ref. 3 for an early calculation of the volume term for the spherical and cylindrical cases in the context of blackbody radiation.

*A Treatise on the Theory of Bessel Functions*

*Handbook of Mathematical Functions with Formulas, Graphs and Mathematical Tables*

One might wonder whether a correction is required for the points with $\u2113=0$. These points form a line in *n*-*m*-$\u2113$ space, so the corresponding correction would be of order *k*. Note that we have inserted the degeneracy factor of $2\u2113+1$ as a sum in Eq. (C6), rather than as an integral. If we had done the integral from $\u2212\u2113$ to $\u2113$ instead, we would have obtained $2\u2113$ instead of $2\u2113+1$, which would have yielded the same result at order *k*^{3}, but would have differed at order *k*^{2}. The reason for the difference is that the bounding surfaces suspended above the *n*-*m* plane (see Fig. 4) contain the points with $m=\xb1\u2113$, and these points all lie precisely on these surfaces. If we had done the integral in the *m* direction, we would only have captured half of the points on these surfaces (see the discussion between Eqs. (14) and (15) in Sec. II A). To correctly capture all of these points, we can extend the range of the integration to be between $\u2212\u2113\u22121/2$ to $\u2113+1/2$, which gives the same result as the sum.

*American Journal of Physics*and

*The Physics Teacher*as a member benefit. To learn more about this member benefit and becoming an AAPT member, visit the Joining AAPT page.