Cathinone and its synthetic derivatives belong to organic compounds with narcotic properties. Their structural diversity and massive illegal use create the need to develop new analytical methods for their identification in different matrices. NMR spectroscopy is one of the most versatile methods for identifying the structure of organic substances. However, its use could sometimes be very difficult and time-consuming due to the complexity of NMR spectra, as well as the technical limitations of measurements. In such cases, molecular modeling serves as a good supporting technique for interpreting ambiguous spectral data. Theoretical prediction of NMR spectra includes calculation of nuclear magnetic shieldings and sometimes also indirect spin–spin coupling constants (SSCC). The quality of theoretical prediction is strongly dependent on the choice of the theory level. In the current study, cathinone and its 12 fluorinated derivatives were selected for gauge-including atomic orbital (GIAO) NMR calculations using Hartree–Fock (HF) and 28 density functionals combined with 6-311++G** basis set to find the optimal level of theory for 1H, 13C, and 19F chemical shifts modeling. All calculations were performed in the gas phase, and solutions were modeled with a polarized-continuum model (PCM) and solvation model based on density (SMD). The results were critically compared with available experimental data.