This paper is the first attempt to build CGC/saturation model based on the next-to-leading order corrections to linear and non-linear evolution in QCD. We assume that the renormalization scale is the saturation momentum and found that the scattering amplitude has geometric scaling behaviour deep in the saturation domain with the explicit formula of this behaviour at large $tau = r^2 Q^2_s$. We built a model that include this behaviour, as well as the ingredients that has been known: (i) the behaviour of the scattering amplitude in the vicinity of the saturation momentum, using the NLO BFKL kernel, (ii) the pre-asymptotic behaviour of $lnLb Q^2_sLb Y RbRb$, as function of $Y$ and (iii) the impact parameter behaviour of the saturation momentum, which has exponential behaviour $propto expLb -, m, bRb$ at large $b$.We demonstrated that the model is able to describe the experimental data for the deep inelastic structure function. Despite this, our model has difficulties that are related to the small value of the QCD coupling at $Q_sLb Y_0Rb$ and the large values of the saturation momentum, which indicate the theoretical inconsistency of our description.