Polynomial interpolation of the function of two variables with large gradients in the boundary layers

The problem of interpolation of the function of two variables with large gradients in the boundary layers is investigated. It is assumed that the function has large gradients near the boundaries of a rectangular domain. Such function corresponds to the solution of the elliptic equation with a small parameter in the highest derivatives. It is known that the error of polynomial interpolation on a uniform grid for the function can be of the order of $O(1)$. It is suggested to use the two-dimensional Lagrange interpolation on the piecewise uniform Shishkin mesh, which is dense in the boundary layers. The Lagrange polynomial with $k_1$ interpolation nodes on $x$ and with $k_2$ interpolation nodes on $y$ is used. The error estimate which is uniform in the small parameter is obtained. Results of the numerical experiments are discussed.