During the development and application of a scattering algorithm, its accuracy is normally validated by comparing with results of spherical particles given by the exact Mie theory. Being the simplest shape, sphere supports morphology-dependent resonances (MDRs), which cause sharp variations of the scattering properties in narrow size ranges. We show that MDRs may mislead the validation of any volume- or surface-discretization methods, including the discrete dipole approximation (DDA) and, thus, should be explicitly avoided. However, the brute-force DDA simulations can actually capture the narrow peaks in the extinction efficiency over the size parameter, but only if a dipole size parameter is smaller than twice the MDR width. That is much more computationally intensive than typical DDA simulations. We find that a single Lorentzian MDR peak may be split into two due to the symmetry breaking by the DDA discretization. Furthermore, instead of time-consuming high-resolution DDA simulations for reproducing MDR, we developed and validated a significantly more computationally efficient method. It is based, first, on fitting simulated data with one or two Lorentzian peaks combined with a cubic baseline. Second, we use Richardson extrapolation of peak parameters to zero dipole size, exploiting the smooth convergence of these parameters towards the reference Mie values. When applied to two MDRs with relative widths 2 × 10−3 and 9 × 10−4, the developed workflow, powered by intensive simulations, reproduces the peak positions with unprecedented accuracy – errors less than 0.07% and 0.4% of their widths, respectively. This extends the way for studying the evolution of the MDR under non-axisymmetric deformations of a sphere or a spheroid.