A FORTRAN-77 program is presented which solves the Sturm-Liouville problem for a system of coupled second-order differential equations by the finite difference method of the second order using the iterative Richardson extrapolation of the difference eigensolutions on a sequence of doubly condensed meshes. The same extrapolational procedure and error estimations are applied to the eigenvalues and eigenfunctions. Zero-value (Dirichlet) or zero-gradient (Neumann) boundary conditions are considered.