Thermal and solvent effects on δ(57Fe) of [Fe(CN)5(NO)]2− (1) and [Fe(CN)6]4− (2) were simulated with two molecular-dynamics (MD)-based approaches. Employing the BP86 functional, Car–Parrinello (CP) and Born–Oppenheimer (BO) MD simulations were performed with and without periodic boundary conditions, respectively, for both free anions and aqueous solutions. In the BOMD simulations a QM/MM scheme was adopted, in which the solvent water was described by a force field with electrostatic embedding into the quantum-mechanical part. Magnetic shieldings were computed at the GIAO-B3LYP level and were averaged for snapshots along the respective trajectories. For 1, both approaches afford very similar thermal and solvent effects on δ(57Fe), and the simulated values in water agree well with experiment. For 2, very large effects are obtained, up to |Δδ| > 1200. While the BOMD-derived δ value in aqueous solution is close to the experimental value, the CPMD-based one is significantly underestimated. The computed trends for the chemical shifts can be rationalized by a very large sensitivity of the magnetic shieldings on the Fe-ligand distances, affording bond-length/shielding derivatives up to ∂σFe/∂rFeC ≈ −35 000 ppm A−1 for 2. There is only a small direct effect of the solvent on the δ(57Fe) values, which are governed by the geometric parameters that change upon solvation. Due to the limited box sizes in the periodic CPMD simulations, spurious interactions between molecules in neighboring boxes are especially large for the tetraanion 2, resulting in too short bond lengths and, thus, in an increased magnetic shielding. Nevertheless, both CP- and BO-MD approaches afford similar pictures of the respective hydration spheres around 1 and 2, as assessed by the average number of hydrogen-bonded water molecules. Both approaches differ, however, in details of the solvent structure around the anions, as revealed by suitable pair correlation functions.