A multigrid (MG) numerical solution procedure for a mathematical model for epitaxial growth in metalorganic chemical vapor deposition reactors is developed. The model is governed by the full compressible Navier–Stokes equations for two-dimensional (plane and axisymmetric) laminar flows extended by transport equations for the chemical species, the energy equations and the equation of state, together with boundary conditions providing information on the reactor geometry and experimental conditions. The MG numerical solution is implemented in a finite element procedure that uses the same linear finite elements to discretize both convection and diffusion fluxes by adding to the model a stabilizing term, resulting in easy implementation of the procedure. A deposition process for a model of GaAs growth is calculated with our procedure showing consistency with experimental data.