Abstract : In this paper, a class of linear parabolic systems of singularly perturbed second-order differential equations of reaction–diffusion type with initial and Robin boundary conditions is considered. The components of the solution → u 𝑢→ of this system are smooth, whereas the components of ∂ → u ∂ x 𝜕𝑢→𝜕𝑥 exhibit parabolic boundary layers. A numerical method composed of a classical finite difference scheme on a piecewise uniform Shishkin mesh is suggested. This method is proved to be first-order convergent in time and essentially first-order convergent in the space variable in the maximum norm uniformly in the perturbation parameters.