  We suggest a new procedure for extrapolating the parton distributions from
HERA to much higher energies. The procedure suggested consists of two steps.
First, we solve the non-linear evolution equation. Second, we introduce a
correcting function for which we write a DGLAP-type linear evolution equation.
The nonlinear equation is solved numerically and estimates for the saturation
scale, as well as for the gluon density at THERA and LHC energies are made.
