  Recent progress in experimental techniques has made the quantum regime in
plasmonics accessible. Since plasmons correspond to collective electron
excitations, the electron-electron interaction plays an essential role in their
theoretical description. Within the Random Phase Approximation, this
interaction is incorporated through a system of equations of motion, which has
to be solved self-consistently. For homogeneous media, an analytical solution
can be found using the Fourier transform, giving rise to Lindhard theory. When
the medium is spatially inhomogeneous, this is no longer possible and one often
uses numerical approaches, which are however limited to smaller systems. In
this paper, we present a novel semi-analytical approach for bulk plasmons in
inhomogeneous media based on the semiclassical (or WKB) approximation, which is
applicable when the charge density varies smoothly. By solving the equations of
motion self-consistently, we obtain the expressions of Lindhard theory with a
spatially varying Fermi wavevector. The derivation involves passing from the
operators to their symbols, which can be thought of as classical observables on
phase space. In this way we obtain effective (Hamiltonian) equations of motion
for plasmons. We then find the quantized energy levels and the plasmon spectrum
using Einstein-Brilllouin-Keller quantization. Our results provide a
theoretical basis to describe different setups in quantum plasmonics, such as
nanoparticles, quantum dots and waveguides.
