  A model is developed, based on the density functional perturbation theory and
the inverse Kohn-Sham method, that can be used to improve relativistic nuclear
energy density functionals towards an exact but unknown Kohn-Sham
exchange-correlation functional. The improved functional is determined by
empirical exact ground-state densities of finite systems. A test of the model
and an illustrative calculation are performed, starting from two different
approximate functionals, to reproduce the parameters and density dependence of
a target functional, using exact ground-state densities of symmetric N=Z
systems.
