Identifying quantum spin liquids, magnon breakdown, or fractionalized excitations in quantum magnets is an ongoing challenge due to the ambiguity in interpreting excitation continua occurring in linear-response probes. Recently, it was proposed that techniques measuring higher-order response, such as two-dimensional coherent spectroscopy (2DCS), could resolve such ambiguities. Numerically simulating nonlinear response functions can, however, be computationally very demanding. We present an efficient Lanczos-based method to compute second-order susceptibilities directly in the frequency domain. Applying this to extended Kitaev models describing , we find qualitatively different nonlinear responses between intermediate magnetic field strengths and the high-field regime. To put these results into context, we derive the general 2DCS response of partially polarized magnets within the linear spin-wave approximation, establishing that is restricted to a distinct universal form if the excitations are conventional magnons. Deviations from this form, as predicted in our Lanczos-based simulations for , can hence serve in 2DCS experiments as direct criteria to determine whether an observed excitation continuum is of conventional two-magnon type or of different nature.