Within DFT the vicinity of the solution to the Kohn-Sham equations is explored treating small changes of nuclei positions as a perturbation. The restoring forces that act in the system upon such perturbations give us information about the force constant matrix, from which we can get the harmonic vibrational frequencies as explained elsewhere, e.g. in the book of Wilson and Decius. Basically, you need the second derivatives.

A key point is to do it efficiently. Doing it numerically is very costly. Since the basis functions are known and they have some nice properties one can do it analytically. Thus, researchers have devised a set of coupled equations to compute the effect of perturbations on the density matrix, which is the core of the calculation.

You will find a detailed explanation in this article by F. Neese et al.

Within DFT the vicinity of the solution to the Kohn-Sham equations is explored treating small changes of nuclei positions as a perturbation. The restoring forces that act in the system upon such perturbations give us information about the force constant matrix, from which we can get the harmonic vibrational frequencies as explained elsewhere, e.g. in the book of Wilson and Decius. Basically, you need the second derivatives.

A key point is to do it efficiently. Doing it numerically is very costly. Since the basis functions are known and they have some nice properties one can do it analytically. Thus, researchers have devised a set of coupled equations to compute the effect of perturbations on the density matrix, which is the core of the calculation.

You will find a detailed explanation in this article by F. Neese.