We would like to thank the reviewers for their carefull reading of our manuscript, and we reply in the present document to their remarks and criticisms.
The basis-set correction can indeed be applied to any approximation to FCI such as multireference perturbation theory, similarly to what was done for weakly correlated systems in P.-F. Loos, B. Pradines, A. Scemama, J. Toulouse, and E. Giner, J. Phys. Chem. Lett. 10, 2931 (2019) in which the basis-set correction was applied to CCSD(T) calculations. As for weakly correlated systems, we expect a similar speed-up of the basis convergence and we will report on this in a forthcoming paper. We have added the following sentence in the Computational Details section:\\
{\it We note that, even though we use near-FCI energies in this work, the DFT-based basis-set correction could also be applied to any approximation to FCI such as multireference perturbation theory, similarly to what was done for weakly correlated systems for which the basis-set correction was applied to CCSD(T) calculations [51].}
This is a very high-quality manuscript, in which a DFT-based correction for finite basis-set errors in wavefunction theory is extended to the case of strongly correlated molecular systems, by looking at molecular dissociation curves. The beauty of the theory behind this DFT basis-set correction is that it can be applied to any wave function method. The authors use their very smart idea (published in previous works) of a suitable mapping of the interaction projected in the finite basis to a long-range only interaction (erf(mu*r)/r), determining the value of the local range separated parameter mu by matching the interaction at r=0.
Here the authors further extend the DFT part of their theory by using explictilyt the on-top (OT) pair density dependence of the short-range functionals of multideterminant range separated DFT to correct the basis set error. They also show that the OT dependence can fully eliminate the spin-polarization dependence.
The idea is very neat and well explained, and the results are excellent. I recommend publication with the minor revisions below.
1) The authors should comment more the results of Table 1 for atomization energies. Although their theory is guaranteed to converge to the exact result for the complete basis set limit, the results of table 1 show that for atomization energies this does not always occur in a monotinic way. There are cases (H$_{10}$ and N$_2$ with the OT functionals, O$_2$ with all the functionals and F$_2$ with augmented basis with all the functionals) in which the results are worse for a quadruple-zeta basis than for triple zeta. Could this be due to the fact that when augmenting the basis the OT pair density from the CASSCF wavefunction does not necessarily improve but might oscillate around the Coulomb cusp creating artefacts in the extrapolated OT? After all the basis sets are optimized on the energy and do not really care about the cusp. Or is the effect due to a different accuracy of the correction for the molecule and the atoms? It would be nice if the authors could say something on that.
2) Appendix A 1: eqs. A3, as discussed by the authors may break down for the case of degenearcy of the isolated fragments. The authors comment that the spin degeneracy can be avoided by considering the functionals without spin polarization and that the spatial degeneracy can be also cured by considering for the fragment the same element of the ensemble as for the suprasystem.
However I still wonder if this is ok. For the spin I agree with the authors, but for the spatial degeneracy I still have some doubts. Could the authors comment at least qualitatively on what would happen for C$_2$ and B$_2$ for example? Is the DFT-based basis set correction going to give the same result for the spherically-averaged C or B atom than for the atoms with the occupied p orbital oritented alsong the bond axis? Or does the DFT-based correction show on this the same problems of standard DFT?
\item[1] As pointed out by reviewer 2, the atomization energies computed with the DFT-based basis set correction does not converge in a monotonic way in some cases, as for instance the error is larger for N$_2$ in aug-cc-pVQZ (-1.5 mH) than for the aug-cc-pVTZ basis set (-0.3 mH). The reviewer suggested that this could come from a non monotonic behaviour of the on-top pair density itself with respect to the basis set. The authors have investigated in that direction and reported several levels of approximations for the on-top pair density for N$_2$ and the Nitrogen atom in increasing basis sets in order to quantify that. It was found that average on-top pair density at near FCI level for a given system (\textit{i.e.} its integral in real-space) decays in a monotonic way as one enlarges the basis set, which is the signature that the depth of the coulomb hole is enlarged when one increases the basis set in a monotonic way.
Nevertheless, the on-top pair density used in the SU-PBE-OT functional is not the near FCI one but the on-top pair density at the CASSCF level extrapolated towards the exact on-top thanks to the known large $\mu$ behaviour and the $\mu_{\mathcal{B}}({\bf r })$ computed for a given point in real space in a given basis set. It was found by the authors that, unlike the near-FCI one, such extrapolated CASSCF on-top increases with the basis set as the CASSCF on-top is almost basis set independent and the global $\mu_{\mathcal{B}}({\bf r})$ increases with the basis set. On the other hand, the same extrapolation scheme based on the near FCI on-top is almost constant and can therefore be considered as a reliable estimate of the exact (\textit{i.e.} near CBS) on-top.
Comparing with estimated exact values, it was observed that extrapolated CASSCF on-top overestimates the on-top pair density, and that this overestimation is more pronounced in the molecular system than in the atomic limit. Acknowledging that the correlation energy is a growing function of the integral of the on-top pair density, this can explain the overestimation of the atomization energy.
To further confirm this statement, we perform a SU-PBE-OT calculation based on the largest variational CIPSI wave functions and using the extrapolated on-top pair density at such level. It was found that the atomization energies were clearly improved, which confirms that the unbalanced overestimation of the on-top pair density between the molecule and the atomic limit was the origin of the overestimation of the atomization energies.
\item[2] Our DFT-based basis set correction does not generally preserve spatial degeneracy for arbitrary states or ensembles. So, in that regard, it has the same problem than standard DFT. Specifically, for the example of C$_2$ and B$_2$, it will not give the same result for the spherically-averaged C or B atom as for the atoms with the occupied p orbital oritented along the bond axis. However, what we write in the Appendix is that for the systems treated in this work the CASSCF wave function dissociates into fragments in simple identified pure states and we can thus choose to perform the calculation on the isolated atom with the same pure state in order to preserve size-consistency. Of course, this may not be always possible for other more complicated systems. We have clarified this point in the Appendix. The new paragraph is now:\\
{\it Moreover, for the systems treated in this work, the lack of size consistency which could arise from spatial degeneracies (coming from atomic $p$ states) can also be avoided by selecting the same state in the supersystem and in the isolated fragment. For example, for the F$_2$ molecule, the CASSCF wave function dissociates into the atomic configuration $\text{p}_\text{x}^2\text{p}_\text{y}^2\text{p}_\text{z}^1$ for each fragment, and we thus choose the same configuration for the calculation of the isolated atom. The same argument applies to the N$_2$ and O$_2$ molecules. For other systems, it may not be always possible to do so.}