Benchmark Study of the Electronic States of the LiRb Molecule: Ab Initio Calculations with the Fock Space Coupled Cluster Approach

Accurate potential energy curves (PECs) are determined for the twenty-two electronic states of LiRb. In contrast to previous studies, the applied approach relies on the first principle calculations involving correlation among all electrons. The current methodology is founded on the multireference coupled cluster (CC) scheme constructed within the Fock space (FS) formalism, specifically for the (2,0) sector. The FS methodology is established within the framework of the intermediate Hamiltonian formalism and offers an intruder-free, efficient computational scheme. This method has a distinctive feature that, when applied to the doubly ionized system, provides the characteristics of the neutral case. This proves especially beneficial when investigating PECs in situations where a closed-shell molecule dissociates into open-shell fragments, yet its double positive ion forms closed-shell species. In every instance, we successfully computed continuous PECs spanning the entire range of interatomic distances, from the equilibrium to the dissociation limit. Moreover, the spectroscopic characteristic of various electronic states is presented, including relativistic effects. Relativistic corrections included at the third-order Douglas-Kroll level have a non-negligible effect on the accuracy of the determined spectroscopic constants.


Introduction
The chemistry and the physics of polar ultracold molecules have been extensively studied in the last decade [1][2][3][4][5].These types of molecules attract scientists mostly due to their high dipole moment, which induces such phenomena as long-range and anisotropic interactions [6].An example of a diatomic polar molecule is LiRb, the object of this study.The research finds that the LiRb molecule possesses a dipole moment of around ∼4.2 debye [7] and it has vast potential applications.In 2010, Herrera et al. proposed that the mixture of cold LiRb and LiCs molecules trapped on an optical lattice and controlled using an electric field produces a perfect system for potential applications in quantum simulation, quantum computation, studies of negative refraction of microwave fields, or studies of the formation of wave-vector space crystals of excitons [8].In 2011, Kuznetsova and her coworkers presented the application of the LiRb molecule as an atom-molecule platform for quantum computing [9].They discovered that using LiRb, it is possible to encode a qubit in atomic states that allowed an easy initialization, readout, one-qubit operations, and mapping of the qubit state onto a photon for quantum communication.Besides that, LiRb, being a polar molecule composed of alkali metals, can possibly be applied in many fields, e.g., few-body [10] and many-body physics [11], studies of chemical reactions [12], or to expand our knowledge of quantum mechanics itself [13].These studies have driven us to explore the properties of this molecule more accurately than ever before using the state-of-the-art method [14] described in the Section 3.
The most important experimental works for our study are those concerning potential energy curves (PECs) and/or spectroscopic constants of the LiRb molecule.This consists of the following: two papers by Ivanova et al. in which they used Fourier-transform spectroscopy to evaluate PECs and spectroscopic constants of various states of LiRb [15,16]; the study of Dutta et al. in which they used similar technique [17]; and the studies which were conducted using resonantly enhanced multiphoton ionization spectroscopy [18][19][20][21].Values obtained in these papers were cited in the Section 2 to evaluate the accuracy of our theoretically obtained data.
The current work is focused on the quantum chemical description of the bond-breaking process for LiRb.We have in mind an approach that would be able to provide correct energy values at the dissociation limit as well as an accurate description of the potential energy curves (PECs).
Fundamental information about LiRb can be obtained out of PECs, which allows the evaluation of selected spectroscopic constants of the molecule, e.g., equilibrium distance R e , adiabatic excitation energy T e , well depth D e , etc.The detailed knowledge of PECs in the ground and excited states is crucial in the study of ultracold molecules [52].
A theoretical description is not trivial for alkali metals diatomics using most of the well-established methods.The difficulty comes from the fact that the molecule with a closed-shell character at the equilibrium distance alters into an open-shell system at the dissociation limit.Therefore, the restricted Hartree-Fock scheme (RHF) should not be used for long-range distances in such cases.Usually, the unrestricted Hartree-Fock (UHF) or the restricted open-shell Hartree-Fock (ROHF) methods are used; however, such calculations are often faced with convergence problems of the HF and post-HF solutions.Theoretical chemistry came up with several solutions to this problem, but none of them is ideal.For instance, researchers often use the multireference configuration interaction (MRCI) method, which, in its approximate variants, is not size-extensive.Therefore, they decide to carry out the calculations for valence electrons only, freezing the inner-shells' electrons, which also affects the accuracy of the results.To overcome this drawback, these electrons are replaced in some studies using the ECP (Effective Core Potential) method, which makes it possible to consider the valence electrons as moving in the potential of atomic cores [36,38,53].This approximation neglects the subtle effects connected with the correlation of all electrons and, on the other hand, introduces new parameters specific to the considered system to mimic the potential of electrons originating from the inner shells.When discussing interatomic potentials, we should also mention the methods based on the algebraic approaches using the well-known models for interatomic potentials like Morse, Lennard-Jones, or Pöschl-Teller and their usefulness of the characteristic of vibrational spectra [54][55][56][57].On the other hand, the "method of choice", frequently applied in the studies of excited states, i.e., EE-EOM-CC (Excitation Energy Equation-of-Motion Coupled Cluster) [58][59][60][61][62], is not size-extensive.This is the major reason that, in the studies of the potential energy curves where the size extensivity plays a prominent role, the latter method is less frequently applied in the calculations of alkali metal dimers [63][64][65].
In our study, we use the Fock space multireference coupled cluster (FS-MRCC) method, which is-on the one hand-size-extensive (an important feature in the PECs calculations) and-on the other-has a built-in capability to provide-by selecting a proper FS sectorcorrelation corrections for systems with an altered (with respect to the Hartree-Fock (HF) reference) number of electrons.For example, the FS(m, 0) sector produces results pertaining to the system with m electrons added to the HF function.Assuming a neutral HF reference, these results would correspond to the m-tuply negative anion.Smooth and well-behaved PECs are obtained, e.g., in the case of the van der Waals molecules, when the closed shell structure dissociates into the closed shell fragments.This example has been a sort of inspiration when dealing with the problem of the dissociation of alkali metal diatomics.Instead of dissociating the dimer Me 2 (Me 2 →Me + Me ), we can dissociate the doubly ionized molecule Me 2+  2 according to the following equation: and we see here that the closed shell structure Me 2+ 2 dissociates into closed shell fragments, isoelectronic with noble gas atoms, specifically discussing an alkali metal dimer.So, this scheme allows us to obtain continuous and smooth curves for the reference system.In the next step, the FS-MRCC scheme is employed to solve the FS equations for the valence (2,0) sector, recovering the original neutral structure.Consequently, at the correlated level, we obtain results for the Me 2 molecule.
The state-of-the-art research is documented with many results obtained for the twovalence sector of the Fock space [14,[66][67][68].Such a strategy, applied to the study of the PECs of a large number of electronic states for a series of alkali metal diatomics [66][67][68], provided the results of incomparable accuracy.The results obtained for the Rb 2 dimer (on the basis of the (2,0) sector of Fock space) have already been utilized in the experimental works carried out in ultra-low temperatures [69].
In this study, we present accurately calculated potential energy curves and spectroscopic constants for the 22 lowest-lying electronic states of the LiRb molecule using the FS-CCSD (2,0) (S-Singles, D-Doubles) method in the Intermediate Hamiltonian (IH) formulation [14] to eliminate so-called intruder-state problems [70,71].Although the structure of the LiRb molecule is similar to that of NaLi, studied in [68], the presence of the heavy Rb atom requires special attention in the choice of the model space and the calculation of the relativistic corrections.The latter effects were included by combining the DK3 (third order Douglas-Kroll) [72] formalism and the IH-FS-CCSD (2,0) scheme [14].The local version of the IH-FS-CCSD (2,0) program was connected to the GAMESS (General Atomic and Molecular Electronic Structure System) package.This required the development of the necessary software [68].
Obtained results are presented in the Section 2, followed by the more detailed method description in the Section 3, and the study is summarized in the Section 4.

Results and Discussion
All calculations were performed using ACES II [73] ver.2.7.0 and GAMESS [74] ver.2021 R2 Patch 1 packages, both supplemented with our own modules in the local version for IH-FS-CCSD(2,0) calculations [14,68].The spectroscopic constants were computed using 8.0 version of Robert J. LeRoy's LEVEL program [75].For the main calculations, we used the uncontracted ANO-RCC [76] basis set with additional diffuse functions, which we called unANO-RCC+.The exponents for the six additional diffuse functions for the lithium are as follows: 0.0027497, 0.0009619 for the s shell; 0.0017173, 0.0006010 for the p shell; and 0.0067528, 0.0023635 for the d shell.Analogously, for the rubidium atom, the exponents for the new diffuse functions are as follows: 0.0012643, 0.0005057, 0.0002023, 0.0000809 for the s shell; 0.0034222, 0.0013689, 0.0005476, 0.0002190 for the p shell; and 0.0727457, 0.0327356, 0.0147310, 0.0066289 for the d shell.We have obtained the exponents using the even-tempered scheme to provide the correct ordering of atomic states.The final size of unANO-RCC+ basis set was 326 basis functions.The spherical harmonic polarization functions were used in all of our calculations and all electrons were correlated.
To evaluate the contribution of relativistic effects in spectroscopic constants of the LiRb molecule, we used the following scheme: we computed potential energy curves and spectroscopic constants of LiRb using Sapporo-QZP-2012-diffuse basis set [77] (167 basis functions), and, in the next step, our calculations were again completed using the scalar relativistic corrections of Douglass-Kroll third order (DK3) [72] and applying the Sapporo-DKH3-QZP-2012-diffuse basis set (also the same number of basis functions), which is the relativistic equivalent of the former basis set.Finally, we calculated the differences between the relativistic and non-relativistic results and these values were added to the constants obtained using the unANO-RCC+ basis set since our previous works show that spectroscopic constants derived from PECs calculated with the family of Sapporo basis sets are not as accurate as those derived from PECs obtained using augmented ANO-RCC basis sets [68,78].
As a reference system for all of the double electron attachment calculations, we used the LiRb 2+ cation.The reference function was always RHF.The size of the active space was set to 91 (i.e., 91 lowest virtual orbitals have been selected as active) for IH-FS-CCSD(2,0)/unANO-RCC+ calculations, which resulted in the model space size equal to 8281.For IH-FS-CCSD(2,0)/Sapporo-QZP-2012-diffuse, as well as for IH-FS-CCSD(2,0) DK3/Sapporo-DKH3-QZP-2012-diffuse calculations, the active space size was set to 50 (model space size = 2500).
The dipole moment µ of the ground state of LiRb was calculated using the numerical Finite Field approach [79].We used the following equation: where E(F) is the energy of the system in the weak uniform external electric field F. The value of F was equal to 0.001 a.u. and its vector was aligned with the C ∞ axis of the molecule.The external electric field perturbation in IH-FS-CCSD(2,0) computations was introduced using ACES II.In order to estimate the contribution of relativistic effects for the dipole moment, we used a similar procedure as in the case of spectroscopic constants, this time using the difference between CCSD/Sapporo-QZP-2012-diffuse and CCSD/Sapporo-DKH3-QZP-2012-diffuse + DK3 results.We used OpenMolcas quantum chemistry software package [80], which allowed us to combine both an external electric field and DK3 scalar relativistic correction.

Atomic Energies at the Dissociation Limit
The important feature of the IH-FS-CCSD(2,0) method is its size-extensivity.This property enforces that the energies of electronic states of the LiRb molecule should converge at an infinite distance to the sum of atomic values.Table 1 presents computed values in the unANO-RCC+ basis set.The first and the second columns present values of energies of Li and Rb-these were calculated using EA-EOM-CCSD-which is equivalent to the IH-FS-CCSD(1,0) method.The third column presents the sum of those energies, and the fourth column shows the energy in the respective dissociation limit obtained using the IH-FS-CCSD(2,0) method.We can clearly see that the results are equal and our method is strictly size-extensive.Thus, the FS-CC method shows correct separability.This feature is particularly useful in the studies of the dissociation process.

Dipole Moment
The calculated dipole moment of the ground state of the LiRb molecule is presented in Table 2 alongside experimental and theoretical reference values.The LiRb molecule has one of the largest dipole moments among heteronuclear alkali metal dimers, being second only to the species containing cesium [7], and our results comply with this fact.Our value computed using IH-FS-CCSD(2,0) DK3/unANO-RCC+ is very close to the experimental values of Tarnovsky et al. [81] and Stevenson et al. [21].We can clearly see that the methods based on pseudopotentials usually overestimate the values of the dipole moment of LiRb, and our results are close to the quality of CCSDT [27,40] with significantly smaller computational effort needed to solve the equations.a Method used in this work, which is IH-FS-CCSD(2,0) DK3/unANO-RCC+ as described in the text.The contribution of the relativistic correction is given in the parenthesis.
The computed PECs are shown in Figure 1-we chose separate point types for each symmetry and multiplicity, and states assigned to different dissociation limits have distinct colors.The applied computational method, which uses the closed shell configuration as the reference, has the advantage that the calculations can be completed equally easily for any internuclear distance.The total molecular energies of each state for the various interatomic distances are available in the Supplementary Materials.All of the states are bound, and four of them have double minima, i.e., 4 3 Σ + , 5 1 Σ + , 6 1 Σ + , 6 3 Σ + (see Figure 2), which is in accordance with the results obtained in Ref. [38] via pseudopotential calculations.The energy barriers/positions are 300 cm −1 /4.75 Å, 2314 cm −1 /8.05 Å, 2214 cm −1 / 6.20 Å, and 3057 cm −1 /8.00 Å, respectively (see Table 3).The higher-lying states have a specific form due to avoided crossing (see Table 4 and Figures 1 and 2).  4. Avoided crossing positions obtained using IH-FS-CCSD(2,0)/un-ANORCC+ compared with Ref. [38].

State Position in This Work (Å)
Position in [38] (Å) We also present Figure 3 in which the comparison with experimental curves of Refs.[15,16] is shown.We were able to compare our theoretical PECs with the experimental ones in two cases: 1 1 Σ + and 2 1 Π.We can see the exceptional alignment of our IH-FS-CCSD(2,0) values with the experimental ones.Although the experimental data are accessible only in a limited range, we can see that the respective parts overlap each other very well.This accuracy is also revealed in spectroscopic constants, which are described in the next subsection.In addition, we computed the C 6 coefficient for the ground state and obtained a value of 2457 a.u.comparing well with the experimental one equal to 2550 a.u.[82].
Finally, just for comparison purposes, we also present PECs calculated using IH-FS-CCSD(2,0) DK3/Sapporo-DKH3-QZP-2012-diffuse (see Figure 4).As we can see, the shapes of the curves are similar to those obtained using a significantly larger unANO-RCC+ basis set (Figure 1).
We observe, and this is an important feature of the results, that the curves collected in Figures 1 and 4 dissociate for each group into a common limit, and the shapes of our PECs are consistent with those available from the literature [38].

Spectroscopic Constants
The spectroscopic constants of the LiRb molecule are shown in Table 5.These include equilibrium distance R e , well depth D e , adiabatic excitation energy T e , harmonic frequency ω e , anharmonicity constant ω e x e , and equilibrium rotational constant B e .The data obtained in the current work include the relativistic contributions, while the pure relativistic corrections are shown in parentheses.In addition to our IH-FS-CCSD(2,0) results, we present values obtained in previous papers for comparison purposes.In order to undertake reliable analysis of our results, we decided to limit the number of cited references to fairly recent and the most relevant theoretical works, i.e., the paper of Berriche et al. in which they use a well-established method based on pseudopotential [38], the two publications of You et al., which apply multireference configuration interaction method (MRCI) [41,42], and the recent article of Kozlov et al. in which they combine effective core potentials with the relativistic variant of FS-CCSD method [43].Also, we present the experimental values of Refs.[15,16,18,19].[41] was based on the multireference configuration interaction method (MRCI).Basis set for lithium: aug-cc-pwCV5Z, for rubidium: ANO-type basis set.All-electron DKH relativistic formalism was used.d Method used in [42] was similar as in [41].e Method used in [43] was Fock-space relativistic coupled-cluster with singles and doubles.Basis set for lithium: cc-pwCVQZ.For rubidium semilocal shape-consistent effective core potentials and [7s7p5d3f2g] basis set was used.f Mean absolute error calculated for the states where experimental results are available.For the 3 1 Σ + state the [16] is treated as a true value while for the ω e x e constant the [19] is treated as a true value.
Out of the results cited in the Table 5, we were able to perform a statistically significant mean absolute error (MAE) analysis of experimental values only with the paper by Jendoubi et al. [38].For T e , D e and ω e x e , our results agree better with experimental values, while for R e , the results in the cited paper are better.The value of ω e is similar (the difference is about 1.7%).This is the expected behavior-pseudopotential-based methods produce quite accurate results near the equilibrium distance but are not size-extensive, therefore cannot properly reproduce the PECs for the distances far from the minimum, in contrast to our IH-FS-CCSD(2,0) method, which behaves properly in the whole spectrum of interatomic distances.Also, for the case of the ground state of LiRb, our results are closer to the experimental values of [15], with the perfect alignment of the R e value.
For the ground state, our method provides results more consistent with the multireference configuration interaction (MRCI) method [41,42] than with the method based on pseudopotential [38].For the other states where MRCI results are available, it is not always the case, i.e., a part of our results is closer to the pseudopotential method.We were also able to compare our results to the relativistic variant of FS-CCSD method by Kozlov et al. [43] for 2 1 Σ + and 1 3 Π states.The alignment is exceptionally good with the difference of less than 0.2%, except of the R e value for the 1 3 Π state where the difference is slightly higher but still less than 2%.
We have correctly identified 4 3 Σ + , 5 1 Σ + , 6 1 Σ + , 6 3 Σ + states as having double minima, which we confirmed by comparison with the study of Jendoubi et al. [38].Generally, the agreement of R e , D e , and T e values with the cited paper is good.The best agreement of R e can be seen for the singlet Π and ∆ states, in the case of D e for the 4 1 Π state (the difference of 0.4%), and the absolute difference in T e for the second minimum of 5 1 Σ + state is only 10 cm −1 .The good match in R e implies the same for the derived B e values.However, ω e have some exceptions where the differences are larger than expected, i.e., 4 1 Σ + , 4 3 Σ + , and 4 3 Π.
Relativistic corrections given in parentheses in the Table 5 are not regular-they fluctuate between states and spectroscopic constants.In some cases, they contribute less than 1% of the value of the given spectroscopic constant (e.g., T e for 1 3 Σ + or R e for 1 3 ∆), while in others, they contribute more than 50% of the computed value (e.g., 506 vs. 915 cm −1 for D e of 3 3 Σ + state).Rubidium is an atom where relativistic effects cannot be neglected.Thus, their contribution to the energies of the LiRb dimer should always be included.

Methods
Within the single reference (SR) coupled cluster (CC) theory, the wave function is defined via exponential Ansatz [83,84] |Ψ = e T |Φ 0 where T is a cluster operator defined at the CCSD level as follows: Operators T 1 , T 2 are responsible for single, double excitations from the reference function |Φ 0 .It should be stressed that in our calculations, the reference system is understood as a doubly ionized structure LiRb +2 , and both RHF and CCSD solutions are referred to this reference.The a...i... are the second-quantized operators, removing an electron from the occupied level i and placing it on the virtual one a.The t a• i• amplitudes are the solution of the CC equations obtained by a projection of the H(≡ e −T He T ) operator against excited configurations Φ a... i... |: The multireference [85][86][87][88][89][90][91][92][93][94][95][96][97][98][99][100][101][102] formalism assumes that the configurational space is divided into two subspaces: the model space M with the projector operator P and its orthogonal complement M ⊥ defined by the projector Q.The advantage of the FS-MRCC approach lies in the fact that the diagonalization of the Hamiltonian H within full configurational space (size: millions, billions,...) is replaced by the diagonalization of the effective Hamiltonian H e f f operator in the model space.The size of the model space is drastically smaller and the diagonalization of the H e f f can be carried out in most cases with the standard diagonalization techniques.The effective Hamiltonian is defined as follows: where Ω is a wave operator, which, on acting on the model function Ψ 0 , generates an exact wave function: In order to use particle-hole second-quantized operators, the configuration adopted as a Fermi vacuum must be selected.Then, the model determinants can be generated by the action of the particle-hole creation operators within the valence one-particle levels.The K-valence space is defined by indication of the N v valence electrons and n v valence unoccupied levels The MR-CC considered here, i.e., formulated in the Fock space (also known as a valence universal), admits the configuration with a variable number of electrons.The FS model is obtained by a distribution of 0,1,2,...,K valence electrons among K valence orbitals.It is obvious that it includes configurations containing 0 valence electrons, i.e., with N v electrons removed from the system, as well as K valence electrons, i.e., n v -tuply ionized anion.The valence universality means that i) all model configurations are defined with respect to the same Fermi vacuum; ii) the wave operator Ω is defined identically for all reference determinants.An important step in the FS approach relates to the selection of the Fermi vacuum, which-on the one hand-defines the reference system for which the one-particle states are obtained (in our approach, the Hartree-Fock solutions) and-on the other-determines the sector structure of the model space.The sector is defined by the number of valence particles and valence holes created with respect to the vacuum, e.g., (k, l) sector indicates that the reference configuration contains k valence particles and l valence holes.The (0,0) sector indicates the Fermi vacuum.
The flexibility of the FS approach-indicated already in the Section 1-can be used in all situations where we want to replace the open-shell reference (requiring UHF function) with the closed shell one, conveniently described by the RHF method.
The basic general formula for the MR approach is a Bloch equation of the form Operating from the left with the P operator, we obtain the expression for the H e f f (H e f f = PHΩP) (the identity PΩ = P has been used) while the projection against the configurations belonging to the orthogonal subspace (represented by the Q operator) generates MR-CC equations.The Ω operator can be expressed as and the braces {} indicate that the normal ordering is imposed on the product of S operators.
A specificity of the Fock space method lies in the hierarchical structure of the CC solutions.So, when we want to solve the Fock space equations for the given (k, l) sector, the lower rank sectors (i, j) solutions must be known.That is why in the current work, since we use the FS-CCSD (2,0) variant, the solutions for the (0,0) and (1,0) sectors are needed.The (0,0) sector corresponds to the single reference solution for the reference system (LiRb +2 ion), and the (1,0) sector is responsible for the electron affinity calculations, and the (2,0) sector is used for double electron affinity ones.
The particular sectors are defined by the configurations belonging to them.So, if the model space contains configurations Φ α with a single electron placed on the virtual level, α, we have the (1,0) sector; for configurations Φ αβ with two additional electrons placed on the virtual levels, α, β, we have the (2,0) sector, etc.Note that the convention we adopted here indicates that the indices a, b, ... refer to virtual one-particle levels, and their Greek equivalence α, β, ... refer to active particles, while the ā, b, ... refer to inactive particles.Thus, the DEA (double electron attached) states refer to the (2,0) sector.Consequently, only particle active levels belong to the active space and their number determines the size of the active space.
The energy values of the DEA states are obtained by the diagonalization of the effective Hamiltonian within the Φ αβ configurational space: e f f = P (2,0) He S (0,0) +S (1,0) +S (2,0) P (2,0) (10) where the projection operator P (2,0) is defined as The cluster operators S (0,0) (≡ T), S (1,0) , and S (2,0) correspond to the sector indicated by the superscripts.As we mentioned before, in order to construct the effective Hamiltonian for the (2,0) sector H (2,0) e f f , the amplitudes from the sectors (1,0) and (0,0) are needed, and the diagonalization of the H (2,0) e f f is performed only within the (2,0) sector.The S (1,0) operator can be constructed by solving the FS-MRCC equations for the (1,0) sector, but the identical result can be obtained by solving the EOM-CC problem set up for the electron affinity (EA).It is known that the Fock space approach at the one-valence level, (1,0) and (0,1), is equivalent to the EOM-CC scheme applied to the EA and IP (Ionization Potentials) cases, respectively.This equivalence means that the eigenvalues in both approaches are identical while the eigenvectors can be obtained from each other by a simple transformation.
A well-known problem that complicates the MR formulations of the coupled cluster theory is the presence of intruder states.They arise when the excited determinants from the orthogonal space (usually low in energy) are close to those belonging to the model space.This may cause numerical instabilities and, consequently, may be a reason for making the MRCC equations divergent.
In order to eliminate convergence problems in the (2,0) sector due to intruder states [70,71], we applied the intermediate Hamiltonian (IH) strategy [98][99][100].Within IH formalism, we select a part of the orthogonal space, M ⊥ as an intermediate space, M I , connected with the P I (projector onto the subspace defined by the operators S (2,0) ).Thus, in this formulation, we have the following subspaces: I⊥ and projectors: and wave operators, with X (2,0) = Z (2,0) + Y (2,0) .Z (2,0) , operating on the model space P (2,0) , generates determinants belonging to the M (2,0) I space, while the Y (2,0) operator connects the P (2,0) and M (2,0) I⊥ subspaces.For the IH-FS-CCSD case, they take the following form: The eigenvalues of the effective Hamiltonian [99], defined within the FS-CC, are the same as those of the simpler intermediate Hamiltonian operator, i.e., H (2,0) I [14]: This is a matrix representation of H in M (2,0) 0 with the P (2,0) 0 − P (2,0) part modified by dressing, which is constructed from H and the cluster operators S known from the lower sector, i.e., (1,0) (for more details see the paper devoted to the IH-FS-CCSD (2,0) method [14]).The respective matrix to be diagonalized is obtained in the electron-affinity EOM-CC [103][104][105][106] calculations in order to omit convergence problems in the FS onevalence sector.Diagonalization of the H (2,0) I operator provides a subset of eigenvalues that are identical to those obtained by diagonalization of the H (2,0) e f f operator.It follows from the above that the diagonalization of the IH matrix can replace the iterative solution of the H e f f -based equations.The final advantage of this is the elimination of the intruder state problem from the FS approach.In addition, the (2,0) step scales as n 5 (neglecting the limited number of active orbitals).The scaling of the CCSD solution for the ground state is n 6 , and the (1,0) part is n 5 .

Conclusions
The studies of the potential energy curves of diatomic molecules have a well-established place in the quantum chemical literature.The problem is challenging since upon a dissociation of the closed shell structure, usually the open shell fragments are formed, which complicates the calculations at the correlated level.
In this work, we focus on the Fock space realization of the multireference coupled cluster theory, specifically applied to investigate the electronic structure of the LiRb molecule.The central aspect of our approach lies in the Intermediate Hamiltonian formulation, providing a convenient way to circumvent the iterative solution inherent in the standard Bloch equation.Instead, we took advantage of the direct diagonalization of a suitably constructed matrix ('dressed' H, i.e., IH).This methodology eliminates complications associated with intruder states, a persistent hindrance that has restrained the broader adoption of this method over the years.
Promising findings were obtained in the studies of the dissociation of a single bond for LiRb with a solution based on the (2,0) sector of FS as an indirect application of the DEA calculations.Using this scheme, excellent results were obtained in the determination of the spectroscopic constants and PECs for the twenty-two lowest electronic states.Notably, calculations were performed for the first time across the entire range of interatomic distances using the RHF reference function, correlating all 40 electrons.This approach produced smooth and accurate PECs from the equilibrium distance to infinity.The mean absolute error from experimental D e is equal to 39 cm −1 , while a corresponding value for T e was 46 cm −1 .Additionally, good agreement with experimental data was observed for ω e , with an average deviation of 13 cm −1 .Thus, the molecular data obtained in this study are well justified and, in the absence of experimental data, may be treated as reference values.This work presents the capability of the IH-FS-CCSD (2,0) method to provide comprehensive characteristics of the investigated molecule without relying on additional approximations.
Author Contributions: G.S. carried out theoretical calculations, co-interpretation of the results, preparing all figures and tables, co-editing the manuscript; M.M. conceived the research plan, cointerpreted the results, co-editing the manuscript.All authors have read and agreed to the published version of the manuscript.

Funding:
The research activities co-financed by the funds granted under the Research Excellence Initiative of the University of Silesia in Katowice.

Table 1 .
Energies of electronic states at the dissociation limit of the LiRb molecule compared to the atomic energies in the unANO-RCC+ basis set.

Table 2 .
The dipole moment of the X 1 Σ + state of the LiRb molecule.

Table 5 .
Spectroscopic constants of LiRb.D e , T e , ω e , ω e x e and B e given in cm −1 ; R e given in Å. DK3 relativistic corrections are given in parentheses.
a Method used in this work, which is IH-FS-CCSD(2,0) DK3/unANO-RCC+ as described in the text.b Method c Method used in