  We investigate the accuracy of the second-order time-convolutionless (TCL2)
quantum master equation for the calculation of linear and nonlinear
spectroscopies of multichromophore systems. We show that, even for systems with
non-adiabatic coupling, the TCL2 master equation predicts linear absorption
spectra that are accurate over an extremely broad range of parameters and well
beyond what would be expected based on the perturbative nature of the approach;
non-equilibrium population dynamics calculated with TCL2 for identical
parameters are significantly less accurate. For third-order (two-dimensional)
spectroscopy, the importance of population dynamics and the violation of the
so-called quantum regression theorem degrade the accuracy of TCL2 dynamics. To
correct these failures, we combine the TCL2 approach with a classical ensemble
sampling of slow microscopic bath degrees of freedom, leading to an efficient
hybrid quantum-classical scheme that displays excellent accuracy over a wide
range of parameters. In the spectroscopic setting, the success of such a hybrid
scheme can be understood through its separate treatment of homogeneous and
inhomogeneous broadening. Importantly, the presented approach has the
computational scaling of TCL2, with the modest addition of an embarrassingly
parallel prefactor associated with ensemble sampling. The presented approach
can be understood as a generalized inhomogeneous cumulant expansion technique,
capable of treating multilevel systems with non-adiabatic dynamics.
