At steady state, electroplating processes are governed by the dimensionless equations partial derivative/partial derivative x(d(i) partial derivative u(i)/partial derivative x + d(i)e(i) partial derivative phi/partial derivative x u(i)) = 0 (i = 1,..., n; 0 less-than-or-equal-to x less-than-or-equal-to 1), [GRAPHICS] where d(i), e(i), and u(i) are respectively the diffusion coefficient, charge, and concentration of the ith species. The extra electroneutrality condition SIGMA(i=1)u e(i)u(i) = 0 will determine the electric potential phi. This system of nonlinear differential equations is subject to the nonlinear boundary conditions modelling the actual electrode kinetics. The authors prove the existence of the solution and construct a computational algorithm. Numerical experiments are performed on practical data.