Equations and boundary conditions describing tidal deformations and the gravitational potential of prestressed elastic spherically symmetric bodies create series of boundary value problems that can be converted into a set of linear algebraic equations by using of pseudospectral schemes on Chebyshev grid. The method is used for an evaluation of the tidal Love numbers for models of Mars and Venus