  We develop the Fourier-Laplace Inversion of the Perturbation Theory (FLIPT),
a novel numerically exact "black box" method to compute perturbative expansions
of the density matrix with rigorous convergence conditions. Specifically, the
FLIPT method is extremely well-suited to simulate multiphoton pulsed laser
experiments with complex pulse shapes. The $n$-dimensional frequency integrals
of the $n$-th order perturbative expansion are evaluated numerically using
tensor products. The $N$ points discretized integrals are computed in $O(N^2)$
operations, a significant improvement over the $O(N^n)$ scaling of standard
quadrature methods.
