I guess each physics simulator works with impulses implicitly because of the discrete nature of the simulations. The forces under consideration are used to keep bodies in resting contact.

I suppose you could integrate the equations symbolicly inside the time step. In classic RBD, you almost have to do this in order to conserve angular momentum properly, otherwise the system will leak energy. You can see this clearly when an object is nutating. A first order integrator will not work, or only work for a few seconds.

As far as stability, I haven't seen an LCP engine, so I can't say one way or another. I do know that adhoc methods have problems, because they have no concept that the solution depends on the entire system, rather work on pairs of objects separately from all the others.

This leads to impulses applied in opposite directions that grow in strength in order to componsate for the impulses working in the opposite directions. To use an analogy, think of people talking at a party, and having to talk louder and louder in order to be heard, meanwhile all the other people also start talking louder and louder so they can be heard, and so forth.