A calculation scheme of the nuclear Fukui function via a coupled perturbed Hartree-Fock approach is proposed avoiding the finite difference approach in DFT-based descriptors. Nucleophilic and electrophilic nuclear Fukui functions are compared with the numerical approximation for the nuclear Fukui function (FF) as the negative derivative of the chemical potential with respect to the atomic coordinates and as the derivative of the Helman-Feynman force with respect to the total number of electrons. The results for a set of diatomic molecules are shown. Analytical and numerical techniques do show a high correlation. Overall, values from both numerical methods are larger than those from the analytical one. The analytical results can be interpreted in terms of the character of the orbital involved during ionization or adding of electrons; the change in the equilibrium bond length upon ionization, which is positive for bonding orbitals and negative for antibinding orbitals is connected with the negative or positive values of the left-hand-side nuclear FF, respectively. The nucleophilic nuclear FF is positive for all cases except CO indicating a systematical increase of the bond length after addition of an electron.