Why does the Method of Successive Approximations for a Differential Equation work? Time dependent perturbation theory in quantum mechanics is often derived using the Method of Successive Approximations for a Differential Equation.  I have not seen an explanation or  a more rigorous outline of proof as to why this should work?
Normally we start with this equation:
$$
i\hbar \dot d_f=\sum_n\langle f^0|H^1(t)|n^0\rangle e^{i\omega_{fn} t}d_n(t)
$$
and then we say that to zeroth order the RHS is 0. and then we plug that into the first order and so on, as outlined in Shankar:

and as what $d_f(t)$ is. To zeroth order, we ignore the right-hand side of Eq. (18.2.5) completely, becasue of the explicit $H^1$, and get 
  $$\dot d_f=0 \tag{18.2.7}$$
  in accordance with our expectations. To first order, we use the zeroth-order $d_n$ in the right-hand side because $H_1$is itself of first order. This gives us the first-order equation
  $$
\dot d_f(t)=\frac{-i}\hbar\langle f^0|H^1(t)|i^0\rangle e^{i\omega_{fi} t} \tag{18.2.8}
$$
  the solution to which is, with the right initial conditions, is
  $$
d_f(t) = \delta_{fi} - \frac{i}\hbar\int_0^t \langle f^0|H^1(t')|i^0\rangle e^{i\omega_{fi} t'}\text dt' \tag{18.2.9}
$$

But how could we explain why this should work more rigorously?

The series expansion, known as the Dyson series, is formally correct. However, to prove that it actually converges can be a major problem.
