The hydrodynamic stability equations for natural convection flows adjacent to a vertical isothermal surface in cold or warm water (Boussinesq or non-Boussinesq situation for density relation), constitute a two-point-boundary-value (eigenvalue) problem, which was solved numerically using the simple shooting and the orthogonal collocation method. This is the first instance in which these stability equations have been solved using a computer code COLSYS, that is based on the orthogonal collocation method, designed to solve accurately two-point-boundary-value problem. Use of the orthogonal collocation method significantly reduces the error propagation which occurs in solving the initial value problem and avoids the inaccuracy of superposition of asymptotic solutions using the conventional technique of simple shooting.