An integral method is developed to solve the two-dimensional Stokes problem with Neumann boundary conditions for multiply connected domains in which the inside hole area can shrink and disappear. The method is applied to simulate viscous sintering. In particular the sintering of glasses can be modelled as such, i.e. a viscous incompressible Newtonian volume flow driven solely by surface tension. A Boundary Element Method is applied to solve the integral equations of Stokes flow involved, and the time integration is carried out by a variable-step, variable-order Backward Differences Formulae method. The derived numerical algorithm is demonstrated for several arbitrarily shaped multiply connected sintering domains. In particular some cylindrical packings are considered. The latter simulations provide a justification for the use of ‘unit problems’ in the theory of sintering.