  Multi-reference calculations along the lines of the Generator Coordinate
Method or the restoration of broken symmetries within the nuclear Energy
Density Functional (EDF) framework are becoming a standard tool in nuclear
structure physics. These calculations rely on the extension of a
single-reference energy functional, of the Gogny or the Skyrme types, to
non-diagonal energy kernels. There is no rigorous constructive framework for
this extension so far. The commonly accepted way proceeds by formal analogy
with the expressions obtained when applying the generalized Wick theorem to the
non-diagonal matrix element of a Hamilton operator between two product states.
It is pointed out that this procedure is ill-defined when extended to EDF
calculations as the generalized Wick theorem is taken outside of its range of
applicability. In particular, such a procedure is responsible for the
appearance of spurious divergences and steps in multi-reference EDF energies,
as was recently observed in calculations restoring particle number or angular
momentum. In the present work, we give a formal analysis of the origin of this
problem for calculations with and without pairing, i.e. constructing the
density matrices from either Slater determinants or quasi-particle vacua. We
propose a correction to energy kernels that removes the divergences and steps,
and which is applicable to calculations based on any symmetry restoration or
generator coordinate. The method is formally illustrated for particle number
restoration and is specified to configuration mixing calculations based on
Slater determinants.
