where all the local quantities $\argrebasis$ are obtained from the same wave function $\Psi$. In the limit of two non-overlapping and non-interacting dissociated fragments $\text{A}+\text{B}$, this integral can be rewritten as the sum of the integral over the region $\Omega_\text{A}$ and the integral over the region $\Omega_\text{B}$
Therefore, a sufficient condition to obtain size consistency is that all the local quantities $\argrebasis$ are \textit{intensive}, i.e. that they \textit{locally} coincide in the supersystem $\text{A}+\text{B}$ and in each isolated fragment $\text{X}=\text{A}$ or $\text{B}$. Hence, for $\br{}\in\Omega_\text{X}$, we should have
\begin{subequations}
\begin{equation}
n_\text{A+B}(\br{}) = n_\text{X}(\br{}),
\label{nAB}
\end{equation}
\begin{equation}
\zeta_\text{A+B}(\br{}) = \zeta_\text{X}(\br{}),
\label{zAB}
\end{equation}
\begin{equation}
s_\text{A+B}(\br{}) = s_\text{X}(\br{}),
\label{sAB}
\end{equation}
\begin{equation}
n_{2,\text{A+B}}(\br{}) = n_{2,\text{X}}(\br{}),
\label{n2AB}
\end{equation}
\begin{equation}
\mu_{\text{A+B}}(\br{}) = \mu_{\text{X}}(\br{}),
\label{muAB}
\end{equation}
\end{subequations}
where the left-hand-side quantities are for the supersystem and the right-hand-side quantities for an isolated fragment. Such conditions can be difficult to fulfil in the presence of degeneracies since the system X can be in a different mixed state (i.e. ensemble) in the supersystem $\text{A}+\text{B}$ and in the isolated fragment~\cite{Sav-CP-09}. Here, we will consider the simple case where the wave function of the supersystem is multiplicatively separable, i.e.
where $\otimes$ is the antisymmetric tensor product. In this case, it is easy to shown that Eqs.~(\ref{nAB})-(\ref{sAB}) are valid, as well known, and it remains to show that Eqs.~(\ref{n2AB}) and~(\ref{muAB}) are also valid. Before showing this, we note that even though we do not explicity consider the case of degeneracies, the lack of size consistency which could arise from spin-multiplet degeneracies can be avoided by the same strategy used for imposing the energy independence on $S_z$, i.e. by using the effective spin polarization $\tilde{\zeta}(n(\br{}),n_{2}(\br{}))$ or a zero spin polarization $\zeta(\br{})=0$. Moreover, the lack of size consistency which could arise from spatial degeneracies (e.g., coming from atomic p states) can also be avoided by selecting the same member of the ensemble in the supersystem and in the isolated fragement. This applies to the systems treated in this work.
with $\Gam{pq}{rs}=2\mel*{\wf{}{}}{\aic{r_\downarrow}\aic{s_\uparrow}\ai{q_\uparrow}\ai{p_\downarrow}}{\wf{}{}}$. As the summations run over all orbitals in the basis set $\Bas$, $n_{2{}}(\br{})$ is invariant to orbital rotations and can thus be expressed in terms of localized orbitals. For two non-overlapping fragments $\text{A}+\text{B}$, the basis set can then partitioned into orbitals localized on the fragment A and orbitals localized on B, i.e. $\Bas=\Bas_\text{A}\cup\Bas_\text{B}$. Therefore, we see that the on-top pair density of the supersystem $\text{A}+\text{B}$ is additively separable
Again, $f(\bfr{},\bfr{})$ is invariant to orbital rotations and can be expressed in terms of orbitals localized on the fragments A and B. In the limit of infinitely separated fragments, the Coulomb interaction vanishes between A and B and therefore any two-electron integral $\V{pq}{rs}$ involving orbitals on both $A$ and $B$ vanishes. We thus see that the quantity $f(\bfr{},\bfr{})$ of the supersystem $\text{A}+\text{B}$ is additively separable
with $\mu_{\text{X}}(\bfr{})=(\sqrt{\pi}/2) f_{\text{X}}(\bfr{},\bfr{})/n_{2,\text{X}}(\br{})$. The local range-separation parameter is thus a local intensive quantity.
In conclusion, if the wave function of the supersystem $\text{A}+\text{B}$ is multiplicative separable, all local quantities used in the basis-set correction functional are intensive and therefore the basis-set correction is size consistent.
The computational cost of the present approach is driven by two quantities: the computation of the on-top pair density and the $\murpsibas$ on the real-space grid. Within a blind approach, for each grid point the computational cost is of order $n_{\Bas}^4$ and $n_{\Bas}^6$ for the on-top pair density $n_{2,\wf{\Bas}{}}(\br{})$ and the local range separation parameter $\murpsibas$, respectively.
Nevertheless, using CASSCF wave functions to compute these quantities leads to significant simplifications which can substantially reduce the CPU time.
\subsection{Computation of the on-top pair density for a CASSCF wave function}
Given a generic wave function developed on a basis set with $n_{\Bas}$ basis functions, the evaluation of the on-top pair density is of order $\left(n_{\Bas}\right)^4$.
Nevertheless, assuming that the wave function $\Psi^{\Bas}$ is of CASSCF type, a lot of simplifications happen.
If the active space is referred as the set of spatial orbitals $\mathcal{A}$ which are labelled by the indices $t,u,v,w$, and the doubly occupied orbitals are the set of spatial orbitals $\mathcal{C}$ labeled by the indices $i,j$, one can write the on-top pair density of a CASSCF wave function as
Written as in eq. \eqref{def_n2_good}, the leading computational cost is the evaluation of $n_{2,\mathcal{A}}(\br{})$ which, according to eq. \eqref{def_n2_act}, scales as $\left( n_{\mathcal{A}}\right)^4$ where $n_{\mathcal{A}}$ is the number of active orbitals which is much smaller than the number of basis functions $n_{\Bas}$. Therefore, the final computational scaling of the on-top pair density for a CASSCF wave function over the whole real-space grid is of $\left( n_{\mathcal{A}}\right)^4 n_G$, where $n_G$ is the number of grid points.
At a given grid point, the computation of $\murpsibas$ needs the computation of $f_{\wf{}{}}(\bfr{},\bfr{})$ defined in eq. \eqref{eq:def_f} and the on-top pair density defined in eq. \eqref{eq:def_n2}. In the previous paragraph we gave an explicit form of the on-top pair density in the case of a CASSCF wave function with a computational scaling of $\left( n_{\mathcal{A}}\right)^4$. In the present paragraph we focus on simplifications that one can obtain for the computation of $f_{\wf{}{}}(\bfr{},\bfr{})$ in the case of a CASSCF wave function.
\textit{A priori}, for a given grid point, the computational scaling of eq. \eqref{eq:f_good} is of $\left(n_{\Bas}\right)^4$ and the total computational cost over the whole grid scales therefore as $\left(n_{\Bas}\right)^4 n_G$.
In the case of a CASSCF wave function, it is interesting to notice that $\Gam{pq}{rs}$ vanishes if one index $p,q,r,s$ does not belong
to the set of of doubly occupied or active orbitals $\mathcal{C}\cup\mathcal{A}$. Therefore, at a given grid point, the matrix $\mathcal{N}_{r}^s(\bfr{})$ has only at most $\left(n_{\mathcal{A}}+n_{\mathcal{C}}\right)^2$ non-zero elements, which is usually much smaller than $\left(n_{\Bas}\right)^2$.
As a consequence, in the case of a CASSCF wave function one can rewrite $f_{\wf{}{}}(\bfr{},\bfr{})$ as
Therefore the final computational cost of $f_{\wf{}{}}(\bfr{},\bfr{})$ is dominated by that of $\mathcal{V}_r^s(\bfr{})$, which scales as $\left(n_{\mathcal{A}}+n_{\mathcal{C}}\right)^2\left( n_{\Bas}\right)^2 n_G$, which is much weaker than $\left(n_{\Bas}\right)^4 n_G$.
