The equations of the two-fluid model of low-pressure plasmas with warm ion gas are taken into consideration including collisions between charged particles and neutrals, the charge exchange, and the ionization. The basic equations contain a removable singularity at the ion sonic speed. These equations are ill-conditioned in the subsonic interval of the ion flux, but they are well-conditioned in the transsonic one. First, several transformations and auxiliary functions are introduced in order to eliminate the singularity at the ion sound speed. The resulting boundary value problem is numerically solved by a multi-shooting method for one of the versions of the transformed equations. Second, an improved one-fluid-model is well-conditioned wherein the space charge density is calculated additionally using the electric field and the Poisson equation. The numerical solution yields usable approximated results in the subsonic interval and suitable initial values for the solution of the two-fluid model in the transsonic interval. Third, the unknown functions are expanded as a power series in the relation of the ion temperature to the electron temperature. These equations can be numerically integrated throughout both intervals without serious difficulties. A set of parameters is given describing subsonic intervals extending over the whole plasma. Results obtained by means of the used methods confirm that Bohm's sheath criterion loses its meaning in collision-dominated plasmas. The scopes of application of the different methods are treated by means of examples.