  The interplay of nuclear and electronic dynamics characterizes the
multi-dimensional electronic spectra of various molecular and solid-state
systems. Theoretically, the observable effect of such interplay can be
accounted for by response functions. Here, we report analytical expressions for
the response functions corresponding to a class of model systems. These are
characterized by the coupling between the diabatic electronic states and the
vibrational degrees of freedom resulting in linear displacements of the
corresponding harmonic oscillators, and by nonadiabatic couplings between pairs
of diabatic states. In order to derive the linear response functions, we first
perform the Dyson expansion of the relevant propagators with respect to the
nonadiabatic component of the Hamiltonian, then derive and expand with respect
to the displacements the propagators at given interaction times, and finally
provide analytical expressions for the time integrals that lead to the
different contributions to the linear response function. The approach is then
applied to the derivation of third-order response functions describing
different physical processes: ground state bleaching, stimulated emission,
excited state absorption and double quantum coherence. Comparisons between the
results obtained up to sixth order in the Dyson expansion and independent
numerical calculation of the response functions provide an evidence of the
series convergence in a few representative cases.
