A versatile unitary transformation framework for an optimal bath construction in density-matrix based quantum embedding approaches

The performance of embedding methods is directly tied to the quality of the bath orbitals construction. In this paper, we develop a versatile framework, enabling the investigation of the optimal construction of the orbitals of the bath. As of today, in state-of-the-art embedding methods, the orbitals of the bath are constructed by performing a Singular Value Decomposition (SVD) on the impurity-environment part of the 1RDM, as originally presented in Density Matrix Embedding Theory (DMET). Recently, the equivalence between the SVD protocol and the use of unitary transformation, the so-called Block-Householder transformation, has been established. We present a generalization of the Block-Householder transformation by introducing additional flexible parameters. The additional parameters are optimized such that the bath-orbitals fulfill physically motivated constrains. The efficiency of the approach is discussed and exemplified in the context of the half-filled Hubbard model in one-dimension.


I. INTRODUCTION
Quantum chemistry is built on the foundation of solving the Schrödinger equation, an endeavor that quickly becomes computationally prohibitive for systems of practical interest.The challenge comes from the exponentially scaling cost with system size.In that context, quantum embedding methods [1,2] have emerged as a powerful tool, particularly for the study of strongly correlated systems [3][4][5], where traditional methods often fall short [6,7].In short, an embedding protocol consists in partitioning the original extended system into fragments.Each fragment, often referred to as an impurity, is complemented with bath orbitals to define an effective reduced system described by an effective Hamiltonian.The reduced system, comprising the fragment+bath (denoted the cluster ), contains N i orbitals of the fragment, which interact with the additional N b bath orbitals.Ideally, the fragment and bath reduced system is entirely decoupled from the rest of the system, which is referred to as the environment of the cluster.The cornerstone of these embedding methods is the construction of bath orbitals associated with a fragment in order to derive the effective Hamiltonian of the reduced system.This effective Hamiltonian captures the physics of the larger original environment while integrating out most of its degrees of freedom.Consequently, it provides a highly efficient approach to treat localized electron correlation without the need to explicitly consider the entire system, thus balancing computational feasibility with accuracy.Different embedding strategies have been proposed in the literature to construct such a reduced effective Hamiltonian.
Among these embedding strategies, some rely on the Green's function formalism, which focuses on singleparticle excitations and can naturally incorporate electronic correlation effects.In particular, within the dy- * quentin.marécat@umontpellier.frnamical mean-field theory (DMFT) formalism [8][9][10], the effective Hamiltonian is an effective Anderson impurity model (AIM) [11], and is derived self-consistently using the Green's function of the system.Note that various alternative embedding approaches using the Green's function have been developed to derive effective AIM [12][13][14] or other effective Hamiltonians [15,16].
In parallel, embedding approaches have been proposed within the one-body Reduced Density-Matrix formalism (1RDM).The aim is to design a protocol to construct the bath orbitals, and the corresponding effective Hamiltonian of the cluster, that become functional of the 1RDM.Among the different approaches recently proposed [17][18][19][20], the pioneering work of Knizia et al. [21] proposes to define the effective Hamiltonians by means of the Schmidt decomposition of a single Slater Determinant (SD) |Φ that is univocally associated with an idempotent 1RDM γ 0 , such that the Schmidt decomposition of |Φ defines a projector via a Singular Value Decomposition (SVD) of the fragment-environment part of γ 0 [22][23][24].Indeed, SVD is known to yield the most compact representation of all informations contained within the fragment-environment 1RDM.This approach imposes that the bath constructed through the SVD contains as many orbitals as the impurity fragment, and also that only idempotent 1RDM can be used to describe the extended system.Over the past decade, different variations and improvements of the original DMET have been investigated and benchmarked [25][26][27][28][29]. Recently, Sekaran et al. [23,30] proposed to use a specific unitary transformation, defined as a functional of the 1RDM, to construct the bath orbitals [30].The aforementioned unitary transformation is known as the Block-Householder transformation.They demonstrated that the SVD of the fragment-environment 1RDM is equivalent to the Block-Householder transformation, even for a non-idempotent 1RDM [24].Despite the appealing compactness of the resulting sub-space of the bath, it should be noted that there is no inherent physical reason to believe that such decomposition yields the "best" effective Hamiltonian for reproducing the interactions between the fragment and arXiv:2307.13446v1[cond-mat.str-el]25 Jul 2023 its environment.
In this contribution, we propose a versatile framework that generalizes the block-Householder transformation and introduces additional degrees of freedom for the construction of the orbitals of the bath and the derivation of the effective Hamiltonian.As a result, additional constraints are required for the bath orbitals such as maximally disentangling the cluster from the environment, or matching density matrices.The effects of the additional constraints to optimally construct the orbitals of the bath are benchmarked on the well-known but non-trivial halfband filled one-dimensional Hubbard model [31], following the divide and conquer algorithm proposed in a previous work by the authors [32].

II. THEORY
Let us consider the paramagnetic Hubbard ring given by where ĉ † iσ (ĉ iσ ) corresponds to the creation (annihilation) of an electron of spin σ on the i-th orbital, niσ is the counting operator equal to ĉ † iσ c iσ .Indexes < ij > refer at nearest neighbor orbitals.−t corresponds to the hopping integral, while U stands for the Coulomb integral.

A. Quantum bath from the Block-Householder transformation
In this section we recall the generic construction of the unitary Block-Householder transformation R σ , following the work of Sekaran et al. [30].The unitary transformation R σ is defined as a functional of the spin 1RDM γ σ .It performs the rotation of the mono-electronic basis set such that the system is divided into a compact subset of N i + N b orbitals, the cluster, interacting as few as possible with a large number of N e orbitals.More precisely, the unitary transformation is design such that (i) the fragment in the full system is the same as in the cluster where stand for the creation and annihilation operators of an electron in the k-th orbital with spin σ expressed in the R σ representation, and (ii) the fragment is fully disconnected from the environment at the one-body level, Numerous unitary transformations satisfy conditions ( 2) and ( 4) resulting in identical bath orbitals subspace [33,34].Among all these transformations, we focus explicitly on the Block-Householder transformation R σ defined using an auxiliary matrix By construction, R σ is a normal involution (i.e.R σ −1 =R σ , R σ T = R σ and det(R σ ) = (−1) Ni ) with eigenvalues {1, −1}, and dim(Ker(R σ + Id.))= r = N i .The matrix V is given as following The superscript r (underscript Ni ) stands for th number of columns (lines) of the associated matrix, and the notation A c i:j corresponds to consider only the elements from i-th to the j-th lines for first c-th columns of the matrix A. For example, the square matrix γ r Ni:2Ni refers to the first r = N i -th columns of γ σ from the N i -th element to the 2N i -th.The null matrix 0 r Ni allows to fix the identity in equation ( 2) on the subspace corresponding to the fragment, and X r Ni is a square matrix to be determined in order to satisfy condition (4).To that aim, Rotella et al. [35] proposed a systematic construction with and , where the nonnegative scalar d i and the orthogonal matrix P are defined by and D (P) is eigenvalues (eigenvectors) matrix of the left sided matrix in equation (8), respectively.The construction of X r Ni in equation (7) holds only if the square matrix γ r Ni:2Ni is invertible.Then we obtain Equation ( 9) is strictly equivalent as the conditions defined in equations ( 2) and (4), i.e. the transformation preserves the orbitals of the fragment in the cluster and disconnects the fragment from the environment at the one-body level.
Recent study expounded on the mathematical equivalence between the bath orbitals constructed using the Block-Householder transformation of the 1RDM and the SVD of the fragment-environment 1RDM [24].As a result of this equivalence, the Block-Householder transformation provides the most compact representation, encapsulating information pertaining to the one-body level interactions between the fragment and its environment.Once the Householder transformation has been presented to construct the orbitals of the bath, we briefly recall the procedure based on R σ in order to define an effective, yet approximated, embedded Hamiltonian.Details and discussion on the straightforward generalization to the ab-initio Hamiltonian can be found in Ref. [32].Following equation (3), we perform the transformation of the full Hamiltonian (1) to derive the Hamiltonian H, where single-and two-bodies integrals in the block-Householder representation tij and Ũijkl , respectively, are calculated in the same manner, Considering that R σ maximally uncouples the fragment from the environment though the orbitals of the bath, an approximation for the effective Hamiltonian on the cluster Hc is obtained by projecting the Hamiltonian (10) onto the cluster orbitals.Part of the clusterenvironment two-body interactions are taken into account at a mean-field level in Hce MF .Finally, an effective chemical potential µ emb is added to preserve the total number of electrons in the fragment.All together, the cluster effective Hamiltonian reads and where the notation (kl) refers to pairs of orbitals in R σ representation such that at least k or l belongs to the environment.Likely influenced by the non-interacting character of the bath in DMFT, the non-locals interaction integrals Ũijkl that naturally arise in the bath have been initially neglected in DMET, leading to a Anderson impurity model in the cluster [21].This approximation is called non-interacting bath (NIB) in contrast to the interacting bath (IB) version of DMET where Ũijkl are explicitly taken into account [22].The similar distinction between NIB an IB can also be considered using the Block-Householder transformation to construct the orbitals of the bath [30,32].Interestingly for non-interacting Hamiltonian (U = 0), for which the associated 1RDM is idempotent (γ σ2 = γ σ ), the transformation R σ leads to a perfect decoupling of 1RDMs of the clusters (Eq.( 4) being fulfilled ∀i ∈ the cluster instead of the fragment), respectively.It corresponds to an exact factorization of the underlying associated wave-function |Ψ = |φ , where |φ refers to a single-Slater determinant.The factorization of the wavefunction gives where | φc (| φe ) is an anti-symmetrized product of orbitals that belong solely to the cluster (environment).Consequently, in this specific case, the projected cluster Hamiltonian proposed in equation ( 13) can be used to extract exact local properties of the impurity site, where the ground-state wave function | φ is equal to | φc [30].
In a more general context, particularly with interacting cases, the effective cluster Hamiltonian, as defined in equation ( 13), can serve as a useful approximation for local properties of the fragment.This Hamiltonian is de facto a functional of the 1RDM via the definition of R σ [γ σ ].Several studies have developed self-consistent schemes predicated on the local cluster 1RDM calculated with equation ( 13), as seen in DMET [21] and more recent divide and conquer algorithms [32].Regardless of the self-consistent matching or conquer strategy employed, the efficiency of the method hinges crucially on the ability of the effective cluster Hamiltonian (13) to locally mimic the full Hamiltonian on the fragment.On this subject, there is no unambiguously definitive approach to determine what exactly constitutes "mimicry".To address this challenge, we propose in the following sections to generalize equation ( 6), introducing variational parameters to explore different flavours of embedding.

B. Quantum bath from a versatile unitary transformation framework
Following the philosophy of exact diagonalization solver in DMFT [36], we would like to control the number of bath orbitals, independently of the number of impurity orbitals in the fragment, in order to systematically have a better description of the interactions of the fragment with the environment.As shown schematically in Fig. 1, the Block-Householder transformation leads the fragment (dark blue square) unchanged in the cluster (gray square), and interacts solely with the bath (orange square), where the number of bath orbitals is the same as in the fragment.In the following, we give a general and flexible definition of the unitary matrix R σ [γ σ ] to obtain an optimized cluster, composed of N i impurities coupled to N b ≥ N i bath orbitals.In what follows, spin index σ is omitted for clarity.Similarly as the Block-Householder transformation, we use the auxiliary matrix V with The auxiliary matrix V ∈ R N r , where N i ≤ r ≤ N refers to the rank of the auxiliary matrix, is constructed as following and indexes N i , N b , N e refer to the number of impurity orbitals in the fragment, the number of bath orbitals and the number of orbitals in the environment, respectively.γ r N b (γ r Ne ) corresponds to the r first columns of the spin 1RDM γ σ and the N b lines (N e ) after the N i -th (N i + N b -th) first one.At this point, condition (2) (preservation of the fragment in the cluster) is already satisfied.The matrix X r N b must be determined in order to preserve condition (4), i.e. it must satisfy the following equation, (X where the superscript N i refers to the N i first column of the matrix.In pratice, the nonlinear equation ( 18) can be solved numerically.In the following, we discuss only the single impurity case N i = 1, where equation ( 18) can be solved analytically.In this case, the matrix is reduced to a vector, and equation ( 18) becomes (X Equation ( 19) is nonlinear and constrains the norm of the vector X 1 N b .More precisely, this equation is fulfilled for any vectors X 1 where ||v|| = v T v refers to the norm of the vector v. Equation ( 20) exists only for r > N i and is linear for X j N b , corresponding to a scalar-product preservation, and can be written as follows, We introduce a spherical representation which allows to express the vectors X j N b with lengths l j and a complete set of angles {θ}.In this representation, length l 1 (l j ) is used to fulfill norm preservation (21) (scalar product ( 22)) for X  17) satisfying the conditions ( 2) and ( 4).The number of free parameters is equal to r ×(N b −N i ), with r ≥ N i .In Fig. 2, we illustrate all vectors that can be obtain for N i = 1, N b = 3 and r = 2. Dark blue vectors correspond to the special Block-Householder solution, where N b = 1 corresponds to the second impurity in fragment.The rank two Block-Householder transformation preserves the identity over the second impurity.As a result, the projection of the first (second) Block-Householder vector over the second impurity i = 2 axis gives γ σ 12 (γ σ 22 ) , as shown with the dashed dark blue line in Fig. 2. Beyond the special Block-Householder transformation, all vectors belonging to the orange sphere are solution of equation ( 21) such as presented with X 1 N b (orange vector) for example.The second vector (light blue vector) norm depends on its scalar product with X 1 N b .Thus the norm (and the direction) is fixed using equation ( 22) and all other angles left are free.More generally, when r > 1, every vectors X j N b , 1 < j ≤ r are independent from each other and are correlated to X 1 N b solely.In the case of r = N i , only one solution is available (with 0 degrees of freedom) and leads to the unique Block-Householder transformation.
Note that by considering r, N b and i.e. any rank r auxiliary matrices expressed using equation (17) with N b bath cannot be expressed with an auxiliary matrix with a greater number of bath N b .Similarly, if we consider N b , r and r > r, we get , which means that different ranks r lead to a a specific definition of the cluster.
At this stage we have proposed a generic construction of unitary transformations following equations ( 16) and ( 17) that generalize the block-Householder construction.Indeed, given an arbitrary number of bath orbitals N b ]: Eq. ( 22) (light blue and orange vectors) in spherical representation for N i = 1, N b = 3 and r = 2. Dark blue vectors correspond to exact Block-Householder vectors from a density-matrix of a non-interacting 1D-Hubbard model.The length l 1 is fixed according to equation ( 21), while lengths l j , r ≥ j > 1 depends on the vector X 1 N b following the equation (22).All angles parameters {θ j i }, 1 ≤ i ≤ r − N i are free.For example, all vectors belonging to the orange sphere are available as a choice of X we show that we can construct many unitary transformations that fulfill conditions (2) and (4) up to r ×(N b −N i ) free parameters.In the next section we propose to use these additional parameters in order to add physically motivated criteria to design bath orbitals.

C. Optimization of free parameters
The free parameters {θ} are variationaly optimized to adjust the bath orbitals using physical insights, unlike the systematic Block-Householder transformation.In what follows, the transformation R σ is a functional of the spin 1RDM, but also a function of a full set of free parameters {θ} and is denoted as R σ [γ σ ]({θ}), or in a more compact notation R σ (θ).The transformation of γ σ with the transformation R σ (θ) is denoted as γσ .
According to the construction of R σ (θ) fulfilling constraint (4), the fragment is disconnected from the environment at the one-body level.In the non-interacting case, we show analytically that the Block-Householder (r = N i ) disconnects bath orbitals from the environment.However, this is not the case for the correlated 1RDM.Therefore, the set of variational parameters is used to minimize the value of the so-called buffer-zone ∆, which gives a quantitative insight into the disentanglement of the environment cluster at the one body level and leads to the saddle-point equation Following a similar philosophy, one could minimize the single-particle Von-Neumann entropy of the truncated 1RDM of the cluster.
From the medium to the strong correlated regime, the buffer-zone ∆ is not able to give a quantitative value of the entanglement between the cluster and the environment, where two-body interactions dominate at this regime.In the following, we propose to minimize the square of the Hartree (i.e mean-field) contribution energy between the cluster and the environment where c (e) belongs to the cluster (environment), respectively, and the notation (jkl) refers to pairs of orbitals in R σ (θ) representation such that at least j, k or l belongs to the environment, and tce ( Ũ (θ) cjkl ) defined in equation (11) (equation ( 12)), respectively.In a practical way, the evaluation of E 2 H using equation ( 25) is numerically more expensive than the evaluation of ∆.
Finally, inspired by the DMET matching, we propose to design the transformation in order to enforce the matching between density matrix elements connected to the fragment in the transformation space, and in the cluster min with γσc refers as the ground-state spin σ 1RDM of the cluster obtained by solving Hc ef f defined in equation ( 13).As illustrated in Fig. 3, we test the cost functions proposed in equations ( 23), ( 25) and ( 28) with respect to two free parameters {θ} for N i = 1, N b = 3 and r = 2.A non-interacting one body reduced density matrix (idempotent) is used as a 1RDM test to obtain the transformation R σ (θ).The first rank r = 1 vector is fixed using Block-Householder solution.In that case, two free parameters are optimized, and cost functions can be represented in spherical coordinates.The minimization of the buffer zone ∆ (Fig. 3a) strictly cancels the value of ∆ and corresponds to the Block-Householder solution, as discussed in Sec.II A. This result is attributed to the fact that the trial density matrix is idempotent.However, this particular value is enclaved between regions of higher ∆ values, which can make optimization difficult depending on the starting point of the numerical minimization.Subsequently, the minimization of the meanfield term between the cluster and the environment (Fig. 3b) also yields the same solution in this case.However, the landscape is distinct from the previously studied cost function.Finally, the matching of the density matrices (Fig. 3c) presents a very specific landscape.Indeed, such a landscape is numerically very challenging to explore in order to obtain the global minimum.Contrary to the cost functions studied previously, the Block-Householder solution (dark blue triangle) is not the global minimum.It appears that there exists a continuous set of minimums, which further complicates the exploration of the landscape.Altogether, it highlights the non-trivial character of the resulting landscape that might display many local minima and quasi-flat regions.It results that the numerical optimization of the {θ} parameters might become challenging for large amount of parameters.

III. RESULTS AND DISCUSSION
In this section, the various cost functions outlined in Section II C are evaluated against the homogeneous paramagnetic Hubbard model at half-filling, and compared to exact Bethe-Ansatz (BA) results [37,38].The nonidempotent N-representable 1RDM space is spanned using the self-consistent protocol presented in [32].The various results pertain to the case of a single impurity.The generalization to multiple impurities, discussed in Section II B, will not be covered in this paper.We recall that the Block-Householder transformation corresponds to a particular case, where the set of parameters {θ} is set according to the equation (7).Additionally, the variational parameters are optimized using the numpy python library [39], particularly with the L-BFGS-B method, which is similar to the conjugate gradient optimization method.Finally, as explained in [32], a damping parameter is added to avoid drastic changes of the 1RDM and convergence issues, meaning that a fraction of the previous 1RDM obtained is kept in the new 1RDM.In the results presented here we used a damping of 60%.
In Fig. 4, we show the relative error of the kinetic en- where BA refers to the Bethe Ansatz solution, and the relative error of the double occupation ∆d = 100 × d BA − d /d BA as a function of the relative correlation strength U/(U + 4t), where 4t corresponds to the non-interacting band width.This figure presents the effect of the various cost functions outlined in Sec.II C (color-coded lines).In this case, the single impurity N i = 1 is embedded with three IB orbitals N b = 3 for a vector V of rank r = 2. Thus, there is a number of parameters θ equal to 2 × (3 − 1) = 4 to optimize.In the non-interacting limit U/t → 0, the kinetic energy and the double occupation are cor- Evaluation of the cost functions with respect to free parameters {θ} for Ni = 1, N b = 3 and r = 2 using a non-interacting one body reduced density matrix from a 1D Hubbard chain.The vector corresponding to the rank r = 1 is fixed and corresponds to the Block-Householder vector (see Fig. 2), and the second is presented with the dark blue triangle.In this case, two angles are available, and are represented using spherical coordinates.In (a) we evaluate the buffer-zone ∆ (see Eq. ( 23)).In (b), we evaluate the Hartree contribution E 2 H (see Eq. ( 25)) for U/t = 8, and in (c) we evaluate the density matrix matching (see Eq. ( 28)) for U/t = 8.

Figure 3
rectly reproduced for all cost functions.In regards to the atomic limit U/t → ∞, this is also the case with an asymptotic limit of −8ln(2)t 2 /U for the kinetic energy, and 4ln(2)t 2 /U 2 for the double occupation [40].
For intermediate regimes U/4t 1, there is a strong competition between electronic delocalization which increases kinetic energy, and the electron-electron repulsion strength which penalises the number of double occupation.For state-of-the-art embedding methods, such as DMET [21] or the projected site-occupation embedding theory (PSOET) [18], describing accuratly this regime is very challenging.For weakly to intermediate correlated regimes (U/(U +4t) < 0.6, U/t < 6), the results from the minimization of constrain (25) (dark blue line) are consistently better than those obtained by minimizing the constrain (23) (orange line), which are in turn better than the Block-Householder solutions (yellow line).However, the numerical cost associated with the minimization of constrain ( 25) scale as O N c N 4 e , which is significantly higher than the cost associated with constrain (23) scaling as O N c N 2 e , where N e (N c ) represents the number of orbitals belonging to the environment of the cluster (cluster), respectively.Concerning the one body reduced density matrix matching proposed in equation ( 28) (blue line), the results obtained are similar to the constrain (25) from low to middle correlated regime, and deviate for values of U/(U + 4t) 0.6 for kinetic energy and double occupancy.For strongly correlated regimes (U/(U + 4t) ≥ 0.6), the interaction energy E int = U × d dominates.Although the double occupancy is similar for the Block-Householder method and the buffer-zone minimization at this regime, the minimization of the meanfield term improves the results.The double occupancy is badly described by the density matrix matching and exhibits nonphysical numerical instabilities for values of U/(U + 4t) 1.
In Fig. 5, we show the evolution of the different variational parameters θ for different cost functions as a function of relative correlation strengh U/(U + 4t) for a two rank r and a three bath orbitals calculation.In the left panel, we show the two variational parameters associated with rank r = 1, i.e., the first column of the matrix V in equation (17), while in the right panel, the parameters are associated with rank r = 2. Solid lines correspond to the first angles for all ranks, and dotted lines to the second ones.We found that all parameters except for the first of rank two (solid line in the right panel) do not change significantly with respect to U/t.Interpreting this result is challenging.It might be due to specific symmetries of the system, such as translation invariance and electron-hole symmetry at half-filling.Regarding the first parameters of rank two, they all follow the same trend, with the Hartree cost function slightly lower than the others.In light of these results, it is conceivable to simplify the variational optimization of the different parameters by considering only a reduced set of parameters (in this case the first of the second rank) varying significantly in the process.This simplification could greatly improve the numerical optimization of the parameters, and therefore the numerical efficiency of the method in general.
As originally presented in DMET [21], the effective Hamiltonian was an AIM, which implies that interaction terms within the bath are neglected, a scenario referred to as the NIB approximation.In this context, we present in Fig. 6 the relative error of the ground-state energy per site E gs = E k + U d, with respect to the correlation strength for a system with a single impurity in the fragment, three bath orbitals, and a rank-two vector.In the low correlation regime (U/(U +4t) < 0.3), both NIB (represented by dashed colored lines) and IB (represented by full colored lines) yield similar results across all presented cost functions.However, for the intermediate to strong correlation regime (U/(U + 4t) > 0.6), the NIB approximation fails to provide an accurate description of the ground-state energy of the system.It should be noted that the Block IB (full yellow line) appears to offer the most accurate representation of the ground-state, contradicting the observations made in Fig. 4.However, this apparent accuracy is misleading and arises from a larger error compensation between the kinetic and interaction energies.Additionally, we found the same error compensation for the NIB approximation, leading to inadequate  results for both kinetic and interaction energies.It is also noteworthy to optimize variational parameters to minimize the strength of Coulomb repulsion in the bath.This method could potentially close the gap between finite interacting baths in DMET or DaC approaches, and the infinite but non-interacting baths present in DMFT.We present in Fig. 7 relevant optimized parameter (the first parameters of rank two, see results presented in Fig. 5) with respect to the relative correlation strengh for different cost functions (colored lines).The two sce- narios, NIB (dashed lines) and the IB (full lines), are presented.For the NIB approximation, optimized parameters follow the same trend as the IB approximation, but are overestimated for all relative correlation strengh.As explained in Fig. 6, both NIB and IB yield similar results across all presented cost functions in the low correlated regime and equivalent parameters for both noninteracting (U/t → 0) and atomic (U/t → ∞) limits.Interestingly, the interaction terms in the bath (IB), that emerge naturally with an embedding scheme based on a unitary transformation, are mandatory to obtain a correct description, in particular at large U/t.
In the following, we focus on the IB case and the minimization of the buffer zone (see Eq. ( 23)) to explore the influence of the number of bath orbitals and the rank of the matrix V used to define the unitary transformation.In Fig. 8, we present the per site kinetic energy scaled with the non-interacting kinetic energy per site E 0 K = −4/π (upper panel) and the per site double occupation (lower panel) as a function of relative repulsion strength U/(U + 4t).Results correspond to the cases N i = 1, N b = 3 and are given for different rank r of the vector V , ranging from one up to three.For instance, the rank two (blue line) corresponds to the results previously shown in Fig. 4 (orange line).The number of variational parameters is equal to 2 for r = 1, equal to 4 for r = 2, and equal to 6 for r = 3.For a rank 1 vector, we are in the special case where the rank equals the number of impurities N i .In this case, the construction of vector V needs to satisfy only the norm preservation in equation (21).As demonstrated previously in Section II B, increasing the rank of vector V does not system- atically improve the solutions, as the accessible solution spaces are disjoint.Indeed, we show that the rank 2 results are the closest to the exact results for both kinetic energy and double occupation, and this applies to all correlation regimes U/t.For strongly correlated regimes, rank 3 slightly improves the results of rank 1, but it exhibits numerical instabilities due to the optimization of a larger amount of variational parameters.According to this figure, it is not necessary to increase the rank of vector V in order to systematically improve the results.Moreover, the best results are obtained for r = 2, corresponding to the number of singular values in DMET [21], or the number of columns of the vector V = V for the Block-Householder method presented here [32].In Fig. 9, we present the scaled kinetic energy per site (upper panel) and the per site double occupation (lower panel) as a function of relative repulsion strength U/(U + 4t).The results are shown for a rank r = 3 vector V , for different number of orbitals in the bath, ranging from one up to five (colored lines).For spin symmetry reasons, we only consider cases where the number of orbitals in the cluster is even.The case with a single orbital in the bath (yellow line) is very particular, as there are no variational parameters to be optimized in this case.In the other cases, the number of variational parameters corresponds to 6 for N b = 3 orbitals, and 12 for N b = 5 orbitals.Importantly, increasing the number of orbitals in the bath systematically improves the kinetic energy and double occupation for all correlation regimes U/t.However, for strongly correlated regimes U/t >> 1, we observe oscillations of the solutions for five bath orbitals.In the latter case, there are a large number of variational parameters, and their optimization is numerically challenging with the proposed iterative process.

IV. CONCLUSION
In this study, we have thoroughly investigated the performance of various 1RDM based embedding methods to construct the orbitals of the bath.A particular emphasis is placed on the characteristics of the resulting reduced and effective Hamiltonian.Indeed this Hamiltonian is tasked with accurately reproducing the interactions between the fragment of the system and its environment within a downscaled cluster.
While DMET employs the SVD of the fragmentenvironment 1RDM to define the effective Hamiltonian, we have demonstrated that the compact subspace is not the optimal setting for deriving the effective Hamiltonian.By generalizing the Block-Householder equations, we introduce a significant amount of additional flexible parameters, notably by adding bath orbitals that are nearly independent of the number of fragment orbitals or by exploring different transformation domains via rank augmentation.To efficiently leverage these additional degrees of freedom, we have proposed cost functions that, in most cases, effectively disconnect the cluster, containing an integer number of electrons, from the environment.
These cost functions were tested on the half-filled Hubbard model Hamiltonian, with a single impurity orbital in the fragment for which the equations are simplified.The results showed significant improvements over the Block-Householder outcomes.
Nevertheless, these improvements imply numerical optimization, which often proves challenging due to the complex landscape of cost functions.The complexity of these landscapes likely contributes to the fact that we can currently only achieve half-filled results.Moreover, this complexity occasionally makes it difficult to obtain continuous solutions for all relative correlation strengths, resulting in certain non-physical instabilities.Therefore, we encourage further research into the development of efficient cost functions that can derive an optimized effective Hamiltonian to describe the fragment, offering a smoother landscape than its counterparts.A compelling challenge for future research would be to test the method on multiple-impurity fragments and, importantly, to derive linearized equations to define the unitary transformation more effectively.

Figure 1 :
Figure 1: Schematic representation of the Block-Householder transformation and the generalization, both functionals of the 1RDM.Small dark blue squares correspond to the fragment 1RDM in the system, identical as in the cluster, and orange squares correspond to bath 1RDM in the cluster.The gray squares correspond to the cluster.

1 N
b (X j N b ) respectively.The set of angles {θ} are thus completely free and can take any value between [−π, π[.Consequently, we have a complete set of {θ} parameters to construct the corresponding set of different auxiliary matrices defined in equation (

1 N
b .The number of free parameters is equal to r × (N b − N i ) = 4 (orange and light blue dashed curves).

Figure 4 :
Figure 4: Relative error for the kinetic energy ∆E k (top panel) and per site double occupation ∆d (bottom panel) in percent with respect to correlation strength U/(U + 4t) for one impurity orbital in the fragment, three bath orbitals and a rank two vector.Colored lines correspond to different cost functions.Black dashed line corresponds to Bethe Ansatz.

Figure 5 :
Figure 5: Optimized free angle parameters with respect to the relative correlation strengh U/(U + 4t) for a rank r = 2 and a three bath orbitals calculation, for different cost functions (colored dashed lines).The left panel corresponds to the first rank, while the right panel to the second.Solid lines correspond to the first parameter of each rank, and dashed lines to the second.

Figure 6 :
Figure 6: Relative error of the ground state energy per site ∆E gs = 100 × (1 − E gs /E BA ) with respect to correlation strength U/(U + 4t) for one impurity and three bath orbitals and a rank two vector.Colored lines correspond to different cost functions for the IB case, and dashed lines for the NIB case.Black solid line corresponds to Bethe Ansatz.

Figure 7 :
Figure 7: Optimized free angle parameters with respect to the relative correlation strengh U/(U + 4t) for the first angle of the rank two vector of the set of parameters {θ}, for different cost functions (colored lines).Solid lines correspond to the IB case, and dashed lines to the NIB case.

Figure 8 :
Figure 8: Renormalized kinetic energy per site E k /E 0 k (top panel) and per site double occupation d (bottom panel) with respect to correlation strength U/(U + 4t) for one impurity and three bath orbitals and the minimization of the buffer zone as a cost function.Colored lines correspond to the rank of the vector up to three.Black solid line correspond to Bethe Ansatz.

Figure 9 :
Figure 9: Renormalized kinetic energy per site E k /E 0 k (top panel) and per site double occupation d (bottom panel) with respect to correlation strength U/(U + 4t) for a rank three vector and the minimization of the buffer zone as a cost function.Colored lines correspond to the number of bath orbitals up to five.Black solid line corresponds to Bethe Ansatz.