The workflow for a polarizable embedding KS-DFT/AMOEBA wavefunction optimization using MiMiC includes modifications in the standard KS-DFT iterative cycle as shown in When running KS-DFT/AMOEBA calculations, there is a risk of electron spill-out. Since there are no electrons in the MM subsystem, there is no Pauli repulsion and the electron density of the QM subsystem can get close to the MM charges or even overlap with them, in particular when using extended basis sets such as PWs in combinations with real-space grids as illustrated in Fig. . The close proximity of electron density and (positive) MM charges leads to possible singularities of 1/r ij , when r ij → 0. To avoid those singularities, we use the modified Coulomb interaction tensor T mod , in which the covalent radii r cov of MM atoms determine the damping function. In the limit r cov → 0, T (n) mod and its derivatives recover the shape of the standard Coulomb interaction tensor. In the limit |r -r i | → 0 on the other hand, the modified Coulomb interaction tensor T (n) mod has a finite value while the standard Coulomb interaction becomes infinity, i.e., causes a singularity. In fact, the short-range limit for the 0-th order modified interaction tensor T (0) mod , is inversely proportional to the covalent radius (Eq. 21). As anticipated, the limit of higher derivatives of this tensor go to 0. Thus, we obtain a finite value that indirectly mitigates the electron spill-out atoms through the screened electrostatic potential for the multipolar expansion interactions mediated through the modified interaction tensor. 24