Molecular dynamics (MD) is a cornerstone of computational materials science and biology, but it suffers from a severe timescale problem: standard MD is unable to simulate the long timescales required for macroscopic phenomena to occur. This is the algorithmic manifestation of metastability, wherein a system remains trapped in a sub-region of the configuration space for very long times before transitioning to another state. While accelerated MD methods exploit mathematical properties of local equilibria to reach macroscopic timescales, their practical efficiency depends heavily on the user's ability to provide adequate definitions of these metastable states. In this presentation, we will explore the definition of optimal states. We view this is a spectral shape optimization problem in the high-dimensional configuration space, associated with the underlying Langevin dynamics. We tackle this problem through two complementary approaches. First, we explore the problem analytically in the low-temperature regime. Using semiclassical analysis, we establish sharp quantitative asymptotics for the eigenvalues in temperature-dependent domains. This yields a geometry-sensitive generalization of the celebrated Eyring--Kramers formula, and allows us to derive asymptotically optimal, temperature-dependent definitions of metastable states. Second, we attack the shape optimization problem numerically. A difficulty here is that degenerate eigenvalues make the shape-functional non-smooth, and dedicated numerical methods have to be constructed. To handle the high dimensionality of realistic molecular systems, we propose a dimensionality-reduction scheme using collective variables, and apply the method to a biomolecular system.