Electroless plating in micro-channels is a rising technology in industry. In many electroless plating systems, hydrogen gas is generated during the process. A numerical simulation method is proposed and analyzed. At a micrometer scale, the motion of the gaseous phase must be addressed so that the plating works smoothly. Since the bubbles are generated randomly and everywhere, a volume-averaged, two-phase, two-velocity, one pressure-flow model is applied. This fluid system is coupled with a set of convection–diffusion equations for the chemicals subject to flux boundary conditions for electron balance. The moving boundary due to plating is considered. The Galerkin-characteristic finite element method is used for temporal and spatial discretizations; the well-posedness of the numerical scheme is proved. Numerical studies in two dimensions are performed to validate the model against earlier one-dimensional models and a dedicated experiment that has been set up to visualize the distribution of bubbles.