In order to employ continuum models in the analysis of the flow behaviour of a viscous Newtonian fluid at micro scale devices, it is necessary to consider at the wall surfaces appropriate slip boundary conditions instead of the classical non-slip condition. To account for the slip condition at the nano-scale, we used the Navier's type boundary condition that relates the tangential fluid velocity at the boundaries to the tangential shear rate. In this work a boundary integral equation formulation for Stokes slip flow, based on the normal and tangential projection of the Green's integral representational formulae for the Stokes velocity field, which directly incorporates into the integral equations the local tangential shear rate at the wall surfaces, is presented. This formulation is used to numerically simulate concentric and eccentric rotating Couette mixers and a Single rotor mixer, including the effect of thermal creep in cases of rarefied gases. The numerical results obtained for the Couette mixers, concentric and eccentric, are validated again the corresponding analytical solutions, showing excellent agreements.