IV. CONCLUSIONS In the present paper we have numerically calculated the evolution of both homogeneous background and density fluctuations in a viable f ( R ) DE model based on the specific func- tional form proposed in Ref. [14]. We have found that viable f ( R ) gravity models of present DE and accelerated expansion of the Universe generically exhibit phantom behaviour during the matter-dominated stage with crossing of the phantom boundary w$_{DE}$ = − 1 at redshifts z ≲ 1. The predicted time evolution of w$_{DE}$ has qualitatively the same behaviour as that was recently obtained from observational data in Ref. [3] . However, it is important that the condition of stability, or even metastability, of the future de Sitter epoch strongly re- stricts possible deviation of w$_{DE}$ from − 1 by several percents in these models. Thus, the DE phantomness should be small, if exists at all, that agrees with present observational data. Still for the models considered, it is not so hopelessly small as in the case of the similar model[12] with n = 0 . 5 recently considered in Ref. [40] using data on cluster abundance. Note also, that in contrast to Ref. [41], we do not impose the so called thin-shell condition | Δ( f ′ ( R ) − 1) | ≲ | Φ$_{N}$ | , where Φ$_{N}$ is the Newtonian potential of matter inhomogeneities and Δ means change in the quantity in question, for scales exceeding galatic ones where a background matter density approaches the cosmological one. On the other hand, this condition is satisfied automatically for matter overdensities more than 10 for the parameter range n ≥ 2 considered in our paper. As for the density fluctuations, we have numerically confirmed our previous analytic results of a shift in the power spectrum index for larger wavenumbers which exceed the scalaron mass during the matter-dominated epoch[15], while for smaller wavenumbers fluc- tuations have the same amplitude as in the ΛCDM model. Once more, the future de Sitter epoch stability condition bounds possible increase in density fluctuations for cluster scales (compared to the ΛCDM model) by ∼ 20% for n ≥ 2. On the contrary, if it is proven from observational data that this increase is less than 5%, then the background evolution should be practically indistinguishable from the ΛCDM one: | w$_{DE}$ + 1 | < 10 − 4 for n = 2. This shows that σ$_{8}$ and related density perturbations tests are the most critical ones for the f ( R ) DE models considered in the paper. We have obtained that the upper limit on | w$_{DE}$ + 1 | for n = 2 and λ = 8 is 4 . 4 × 10 − 5 when z = 0 . 16, which is of the same order as B (0). We have also investigated the growth index γ ( k, z ) of density fluctuations and have pre- 16