A methodology is developed to fit semi-empirical interatomic potentials aimed at obtaining a consistent thermodynamic behavior. The procedure is based on the cluster variation method theory which is seamlessly integrated to the other more standard equations of the fitting technique. A new interatomic potential for the Fe-Cu system is thus built within the framework of the embedded atom method, to be used in studies of the microstructure evolution of Fe-Cu alloys under irradiation. The potential is shown to reproduce very reasonably the Cu solubility curve in the Fe matrix as well as to lead to enhanced description of the point defect kinetics behavior with respect to previous interaction models. Limitations of the fitting technique and possible ways of improvement are discussed.