  We introduce a numerical method based on an integral equation formulation for
simulating drops in viscous fluids in the plane. It builds upon the method
introduced by Kropinski in 2001, but improves on it by adding an interpolatory
quadrature approach for handling near-singular integrals. Such integrals
typically arise when drop boundaries come close to one another, and are
difficult to compute accurately using standard quadrature rules. Adapting the
interpolatory quadrature method introduced by Helsing and Ojala in 2008 to the
current application, very general drop configurations can be handled while
still maintaining stability and high accuracy. The performance of the new
method is demonstrated by some challenging numerical examples.
