Semiclassical trajectory-based methods are computational techniques that approximate solutions to the time-dependent Schrödinger equation by using information from classical trajectories equipped with evolving phases or width matrices. These additional structures to the classical trajectory enable semiclassical methods to go beyond purely classical descriptions and to capture quantum effects, at least approximately. While single-trajectory semiclassical methods, such as the thawed Gaussian approximation, are efficient and well suited for studying weakly anharmonic systems, more computationally demanding multi-trajectory methods are required in highly anharmonic systems.

This thesis deals with the efficient and accurate simulation of the multi-trajectory Herman-Kluk approximation and its application to spectroscopy. The Herman-Kluk approximation is a continuous superposition of frozen Gaussian wavepackets and is one of the most accurate semiclassical methods. However, its practical application requires the numerical evaluation of a high-dimensional integral over the phase space. Conventional grid-based approaches are generally not feasible for this task, but the structure of the frozen Gaussian wavepackets makes probabilistic techniques well suited. In the first part of this thesis, we investigate Monte Carlo approaches for the Herman-Kluk wavefunction, where the sampling of initial conditions is either based on the Husimi density or on its square root.

Quantum-mechanical expectation values, which require the knowledge of the wavefunction, can give further information about the quantum system. The special structure of the Herman-Kluk approximation allows to evaluate such expectation values without evaluating the wavefunction by working on a double phase space. We extend the two wavefunction sampling approaches from single to double phase spaces, and we prove sufficient conditions for the convergence of the estimator and provide error estimates. Moreover, we introduce an observable-dependent sampling that minimizes the variance of the Monte Carlo integrand at initial time.

We next introduce a refined cellularization scheme for the Herman-Kluk wavepacket autocorrelation function, a quantity that plays a central role in the simulation of spectra. This new scheme is based on cellularization (Filinov filtering), a technique that enhances convergence of estimators by suppressing unimportant contributions. Particularly, in the new methodology, the sampling density depends on the cell size and number of cells, allowing for an effective sampling of the phase space. The main advantage of the presented approach is that, unlike the standard cellularization, it converges to the original Herman-Kluk result in the limit of an infinite number of trajectories and to the thawed Gaussian approximation when a single trajectory is used.

To have more accurate spectra simulations, theoretical models must account for temperature effects. We therefore combine the Herman-Kluk approximation with coherence thermofield dynamics for the simulation of nonzero-temperature spectra. In coherence thermofield dynamics, nonzero-temperature effects are naturally incorporated through an augmented initial pure state. Particularly, spectra can be obtained directly from a single wavepacket autocorrelation function.