  Stable and reliable numerical integration of second-order radial equations in
the fields of classical black holes can be performed by using the Pr\"ufer
transformation and by passing to the use of phase functions, which allows
uniquely selecting the solutions with physical asymptotics in numerical
calculations.
