Infinitely dilute aqueous solutions of Nd3+, Gd3+, and Yb3+ have been studied with molecular dynamics simulations, using ab initio effective ion-water pair potentials based on the polarizable continuum model. Structural results, as first peak positions of g(r) are in good agreement with experimental data. We obtain a coordination number of 9 for all the cations. This value agrees with experimental measurements for Nd3+ and Gd3+ but overestimates them for Yb3+. Significant differences between Yb3+ and the other two ions have however been observed in a detailed analysis of the solvation shell structure, based on the diagonalization of the inertia tensor. A similarity index based on a linear combination of the inertia tensor eigenvalues is proposed. Going from Nd3+ to Yb3+, structures like a capped square antiprism (CSQA) are more favored than a tricapped trigonal prism. In contrast to the lanthanide contraction observed for the most probable ion-oxygen distances, in CSQA structures the distance from the center of mass of the polyhedron to the capped O increases from Nd3+ to Yb3+. Dynamics of internal rearrangements has been studied and residence times of waters in the first and second hydration shell have been calculated. Diffusion coefficients for the ions appear somewhat underestimated with respect to the experimental data, which agree with the value we obtain for hydration shell waters. Linear and angular velocity correlation functions and spectra of first shell waters show a large increase of intensity at high frequency compared to bulk solvent. The Yb3+ solution consistently presents these features in enhanced form.