The Feynman path-integral variational approach to the polaron problem, along with the associated Feynman-Hellwarth-Iddings-Platzman (FHIP) linear-response mobility theory, provides a computationally amenable method to predict the frequency-resolved temperature-dependent charge-carrier mobility, and other experimental observables in polar semiconductors. We show that the FHIP mobility theory predicts non-Drude transport behavior, and shows remarkably good agreement with the recent diagrammatic Monte Carlo mobility simulations of Mishchenko et al. [Phys. Rev. Lett. 123, 076601 (2019)] for the abstract Fr\"ohlich Hamiltonian. We extend this method to multiple phonon modes in the Fr\"ohlich model action. This enables a slightly better variational solution, as inferred from the resulting energy. We carry forward this extra complexity into the mobility theory, which shows a richer structure in the frequency and temperature-dependent mobility, due to the different phonon modes activating at different energies. The method provides a computationally efficient and fully quantitative method of predicting polaron mobility and response in real materials.