Accurate ground state potential of Cu 2 up to the dissociation limit by perturbation assisted double-resonant four-wave mixing

Perturbation facilitated double-resonant four-wave mixing is applied to access high-lying vibrational levels of the X 1 Σ +g ( 0 +g ) ground state of Cu 2 . Rotationally resolved transitions up to v ′′ = 102 are measured. The highest observed level is at 98% of the dissociation energy. The range and accuracy of previous measurements are significantly extended. By applying the near dissociation equation developed by Le Roy [R. J. Le Roy, J. Quant. Spectrosc. Radiat. Transfer 186 , 197 (2017)], a dissociation energy of D e = 16 270(7) hc cm − 1 is determined, and an accurate potential energy function for the X 1 Σ +g ( 0 +g ) ground state is obtained. Molecular constants are determined from the measured transitions and by solving the radial Schrödinger equation using this function and are compared with results from earlier measurements. In addition, benchmark multi-reference configuration interaction computations are performed using the Douglas–Kroll–Hess Hamiltonian and the appropriate basis of augmented valence quadruple ζ type. Coupled-cluster single, double, and perturbative triple calculations were performed for comparison.


I. INTRODUCTION
The bond energy of a molecule is a fundamental thermochemical quantity in chemistry: it tells us how much energy is released when a bond is made or, alternatively, how much energy is needed to break that bond. About one decade after the seminal paper of Gilbert Newton Lewis in 1916, 1 who suggested that a pair of electrons shared by two atoms is responsible for the formation of a chemical bond, quantum mechanics laid the ground for our understanding of the creation and destruction of molecules, 2-6 which has been continuously refined since then. A recent comprehensive review by Frenking and Shaik contains the relevant advances up to our present understanding of bonding. 7 Bond energies can be determined experimentally via traditional calorimetry. 8 However, the quantum mechanical understanding of the chemical bond has allowed us ever since its beginning to assess bond energies spectroscopically, such as with the more recently developed stimulated emission pumping experiments [9][10][11][12][13][14][15][16][17] or zero electron kinetic energy (ZEKE) experiments. [18][19][20][21] Quantum chemical models in combination with the rapid growth of computing power provide an alternative to experimental determination. For molecules composed of first and second row elements, accurate bond dissociation energies are accessible using the coupled-cluster methods with single and double excitations, while treating triple excitations perturbatively [CCSD(T) [22][23][24] ]. The highest possible accuracy is achieved by considering correlation consistent 25 or geminal basis sets 26 in addition to essential corrections from spin-orbit interactions, relativistic effects, and quantum electron dynamics. Coupled cluster theory combined with geminal basis sets can indeed be considered to be among the most accurate and powerful method to obtain accurate information on the chemical structure and properties of medium sized molecules, when the ground state is essentially single reference, even when late-transition metal atoms are included. 27 In a study on transition metal containing dimers, it is emphasized, however, that coupled cluster computations are usually impractically expensive for most problems in transition metal chemistry and suggest that the method yields similar but not better results than the much less expensive density functional theory (DFT). 28 When compared with experimentally determined dissociation energies, the mean unsigned deviations for the transition metal containing dimers were ≈4 kcal/mol (1700 hc cm −1 ) from that study (h is the Planck constant and c is the speed of light in vacuum).
To warrant our quantitative understanding of the chemical bond, the poor accuracy reported in Ref. 28 must be improved, both theoretically and experimentally. Such an understanding will reflect our knowledge of the nature of the chemical bond, which depends on the computational approaches used and the experimental verification of the calculated observables. In this context and in addition to covalent and ionic bonding, a third class of electron pair bond, called charge-shift (CS) bonding, was introduced. 29 CS bonding is characterized within valence bond theory by the occurrence of large resonance energies associated with the mixing of covalent and ionic components of the bonding wavefunctions, which is expected to take place with compact electronegative or lone-pair rich species, and we refer the reader to the original literature for further explanations on the concept of CS bonding. 29 Recent computational investigations on the nature of single bond transition-metal dimers 30 have shown that, for the isoelectronic coinage metal dimers, Cu 2 , Ag 2 and Au 2 , a significant contribution between 40% and 50% to the total bonding energy is attributed to CS bonding. This assertion can hardly be verified experimentally, as the aforementioned resonance energy is not observable. One can nevertheless address calculated and measured values for equilibrium bond lengths and bond energies.
Currently, perhaps the most reliable measurement for the dissociation energy for dicopper, one of the best studied transition metal dimers, dates back to 1986. 31 In that work, Rohlfing and Valentini applied UV excitation of Cu 2 produced by laser vaporization in a molecular beam and observed in emission long progressions to the ground state. Vibrational levels up to v ≤ 72 were measured by dispersed fluorescence. The determined vibrational origins were fitted to the near-dissociation equation reported by Le Roy and Lam. 32 The resulting value for the dissociation energy, De = (16 760 ± 200) hc cm −1 , at the vibrational quanta at the dissociation limit of v d = 128 ± 5 is considerably the most reliable value of De in the literature.
We report in the following on an experimental investigation, which yields a significantly increased accuracy for the bond dissociation energy by applying a non-linear spectroscopic method and by taking into account results from a recent deperturbation study of high lying energy levels of the copper dimer. 33 In addition to the determination of De, an accurate potential energy function is obtained to assess the structural and dynamical properties of this important transition metal species and to provide a further stringent assessment for high level computations. Multi-reference configuration interaction (MRCI) calculations were reported in a previous study. 34 A variation of those calculations as well as new CCSD(T) calculations will be discussed in the present work. Results from the present calculations as well as from previous theoretical work will be critically reviewed in light of the new experimental results.
Experimental investigations are based on the characterization of highly excited states in the spectral region between 37 400 and 38 050 cm −1 , which have opened ways to access spectroscopically vibrational levels of the ground state close to the dissociation limit. 33 Accurate molecular constants and the term symbols for the I and J states 35 have been determined. Furthermore, approximately 1000 rovibronic transitions in the vicinity of the J state of the two main isotopologues 63 Cu 2 and 65 Cu 63 Cu were assigned. The deperturbation of strongly interacting rotational levels of the J (v = 0-2) state revealed dark states that were preliminarily assigned to high vibrational levels of the G 0 + u state. As discussed in our previous work, 34 the G 0 + u state emerges from the avoided crossing of the B 1 Σ + u (0 + u ) state and an ion-pair state. Close to the equilibrium bond length, the B state has a stronger ionic character, while the G state is strongly covalent. At larger distances, however, in the region of the avoided crossing, the latter becomes ionic. Transitions from neutral into ionpair states are typically strong. Therefore, strong emission from the B 1 Σ + u (0 + u ) state to low vibronic levels of the ground state X 1 Σ + g (0 + g ) is observed. By implication, the shallow potential of the G 0 + u state should expose most ion-pair character at the outer wall at even larger internuclear distances beyond the avoided crossing. By the Franck-Condon principle, vertical transitions to low vibrational levels are not expected to occur for the G state. Alternatively, as will be shown in this work, significant transition strength to high vibrational levels in the ground state emerges and allows for spectroscopic detection of ground state ro-vibrational levels up to the asymptotic limit.
In the theoretical section (Sec. II), we briefly recall the key elements used in the MRCI and CCSD(T) calculations. In the experimental section (Sec. III), we outline the two-color resonant four-wave-mixing method to perform Perturbation-Facilitated Optical-Optical Double-Resonance (PFOODR) spectroscopy. The nonlinear spectroscopic technique is favorably suited to perform background-free and highly sensitive stimulated emission pumping (SEP) type investigations in a molecular beam environment. PFOODR has been applied initially on alkali metal dimers to access triplet states from the electronic singlet ground state. [36][37][38] Triplet ↔ singlet transitions are spin-forbidden, but both singlet and triplet states become accessible from the ground state by an intermediate singlet-triplet mixed state. Pulsed or continuous-wave lasers have been used to perform PFOODR spectroscopy by fluorescence excitation, 39 dispersed fluorescence, 38,40 ion detection, 41 or continuous-wave optical triple resonance spectroscopy. 42 Applying non-linear four-wave mixing to perform PFOODR has been demonstrated by some of us in an investigation of the dark triplet manifold of C 3 that exhibits a singlet ground state 43 and was also used to characterize the first quintet-quintet band of C 2 . 44

II. THEORETICAL METHODS
Technical details of the MRCI method used to obtain the theoretical results presented in this report on the X 1 Σ + g ground state of Cu 2 are given in Ref. 34. In essence, they are based on a restricted active space calculation, composed of 22 electrons in 18 orbitals (RAS-22,18) and using the Douglas-Kroll-Hess Hamiltonian 45,46 and the appropriate basis of augmented valence quadruple ζ type (AVQZ-DK) as defined in the MOLPRO program package. 47 These methods led to unprecedented accuracy of the theoretical description of the 0 + u and 1u manifolds in excited Cu 2 and a solid theoretical justification of the assignment of the A ′ state. We refer to Ref. 34 for further explanations on details of the calculations.
Such an accuracy was achieved because the active space was optimized at the MC-SCF level upon inclusion of a fourth state with 1 Σ + g and a third state with 1 Σ + u symmetry label. These are charge transfer states in the vicinity of the equilibrium structure (Ref. 34; Table I therein). We term this calculation "NpC." As the set of active orbitals include the 4p subshell, when the Cu-Cu distance is elongated, the fourth 1 Σ + g and third 1 Σ + u states change character and become the gerade and ungerade states of the "sp" asymptote, i.e., yielding the dissociation products Cu(3d 10 4s) + Cu(3d 10 4p). This asymptote is about 3.8 eV higher than the ground state asymptote leading to Cu(3d 10 4s) + Cu(3d 10 4s) 48 and 2.3 eV higher than the manifold of "sd" asymptotic states leading to Cu(3d 10 4s) + Cu(3d 9 4s 2 ), to which the several 0 + u and 1u states belong to, which were calculated in Ref. 34.
It turns out that the (22,18) active space is very likely not sufficient for a correct description of the "sp" asymptote. Higher lying orbitals would be needed to render this calculation more accurate. This would make the active space even larger and is beyond the scope of the present work. As a side effect, in the absence of higher orbitals, the quality of the ground state itself suffers when one attempts to optimize all states in the asymptotic region. To guarantee a sufficiently good description, in the present work, we decided to remove the fourth 1 Σ + g and third 1 Σ + u states from the active space during orbital optimization at the MC-SCF level. We term such a calculation "N." The active space contains effectively only the ground state and all singly excited states leading to neutral "ss" and "sd" asymptotes (see Table I in Ref. 34).
CCSD(T) calculations were performed on the basis of single reference Hartree-Fock functions for Cu 2 at the equilibrium bond length of about 225 pm and the separated fragments using the Douglas-Kroll-Hess Hamiltonian and the AVQZ-DK basis set compatible with it, as implemented in the MOLPRO program package. 47

III. EXPERIMENT
The experiments are performed in a molecular beam apparatus designed for simultaneous linear and non-linear spectroscopic measurements of stable and transient species and have been described in detail previously. 34,43,44,49,50 Briefly, Cu 2 is prepared in a homebuilt metal cluster source by laser-vaporization. 50 The second harmonic of a Nd:YAG laser (Continuum, NY81; ≈100 mJ/pulse) is focused through a 500 mm focal lens onto a copper disk target (99% purity). The target is translated and rotated by electric motors (Maxon) to control carefully the rate at which fresh surface is sampled. A pulsed valve (Series 9 general valve, Parker-Hannifin) is used to introduce helium (6.0, Messer Schweiz AG) carrier gas that is expanded with the copper plume in a near supersonic expansion through a 1 mm diameter nozzle into the vacuum. The helium backing pressure behind the pulsed valve is 50 bars. Copper dimers are probed with the two-color resonant four-wave mixing (TC-RFWM) technique ≈5 mm downstream from the nozzle.
The optical setup has been described recently 43 and is only summarized here. TC-RFWM is performed by using two dye lasers (NarrowScan, Radiant Dyes), which are simultaneously pumped by a single Nd:YAG laser (QuantaRay Pro 270-10, Spectra-Physics). The bandwidth of the dye lasers is specified to ≈0.04 cm −1 . A combination of optical components establishes a forward BOXCARS 51 configuration. 52 Three parallel propagating laser beams pass along the three main diagonals of a parallelepiped and cross at a small angle of ∼1 ○ . Doppler broadening is minimized by arranging these beams orthogonally to the propagation direction of the molecular beam. The two laser beams of equal frequencies are referred to as PUMP beams and the third beam as PROBE beam. Phase matching conditions govern the direction of the SIGNAL beam that propagates roughly along the fourth (dark) diagonal of the parallelepiped. The SIGNAL beam is detected by a photomultiplier after removing scattered light and unwanted fluorescence by spatial and spectral filters on its nearly 5 m path. Fluorescence signals are observed perpendicular to both the molecular beam and laser propagation through a 1 m spectrometer (SPEX) and detected by an additional photomultiplier tube.

A. Dispersed fluorescence and four-wave mixing
As mentioned in the introduction, high-lying vibrational levels of the G 0 + u state gain transition strength by perturbation with the J 0 + u state and are, therefore, accessible from the ground state of Cu 2 . 34 Fig. 1 depicts a dispersed fluorescence (DLIF) spectrum upon excitation of the P(18) G (62) − X 1 Σ + g (0 + g ) (v ′′ = 0) transition of 65 Cu 63 Cu at 37 445.03 cm −1 (lower trace, blue). The notation G (62) is used to refer to the level that is presumably the v ′ = 62 vibrational level of the G state. The assignments and characterization of the high lying levels of the G state are subject of a separate publication. 33 The DLIF spectrum shown in Fig. 1 demonstrates the large Franck-Condon overlap between G (62) and high lying vibrational levels of the ground state X 1 Σ + g (0 + g ). Upon excitation, substantial emission up to v ′′ ≈ 102 for J ′′ = 18 and 16 is observed. Above the asymptotic atom limit, broad diffuse features are present that originate from bound-free transitions. 53 The Condon internal diffraction produces an interference pattern that extends substantially beyond the dissociation limit. A detailed report of the phenomenon in Cu 2 is subject of a forthcoming publication.
Due to the small rotational constant (B ′′ ≈ 0.1 cm −1 ), the DLIF spectrum is not rotationally resolved. However, the long progression with large emission to high v ′′ levels opens ways to perform rotationally resolved and isotopologue-specific four-wave mixing experiments to characterize precisely the ground state potential up to ≈98% of the dissociation limit. Individual TC-RFWM spectra are shown in Fig. 1 (upper traces, the ordinate is shifted to ease visual comparison). For the applied SEP scheme, the PUMP laser is operated at a fixed-wavelength resonant with a particular rotational of Chemical Physics transition within the target species, while the PROBE laser is scanned over the spectral region of interest. The SEP scheme makes use of the same excited intermediate state, and therefore, the possible final states are governed by strict two-photon selection rules from the initial state. 50,54 As for the DLIF spectrum in the figure, the P(18) G (62) − X 1 Σ + g (0 + g )(v ′′ = 0) transition is selected to populate J ′ = 17 of the upper state with the PUMP laser. Since the excited G (62) is an 0 + g state in the appropriate Hund's case (c) notation, only one R and one P branch transition to J ′′ = 16 and 18 of each ground state vibration is optically allowed. To exemplify a typical measurement to determine accurate molecular constants, the spectral region for transitions to v ′′ = 101 is reproduced at full rotational resolution in the inset of Fig. 1. The P(16) and R(18) lines are unambiguously assigned in a straightforward manner, and their spectral positions yield the origin and the rotational constant for the X 1 Σ + g (0 + g ), v ′′ = 101 of the 65 Cu 63 Cu isotopologue. In general, several low lying J ′′ -levels are measured for each vibration v ′′ by pumping selected rotational levels in the G (62) state. The molecular constants are then determined by performing non-linear least squares fits to the line positions using the PGOPHER 55 software from Colin Western.
To characterize the ground state potential, 347 SEP type transitions have been measured for the two main isotopologues of Cu 2 spanning the range of vibrational levels from v ′′ = 36 to 98 and v ′′ = 36 to 102 for 63 Cu 2 and 65 Cu 63 Cu, respectively. Levels with v ′′ ≧ 83 are obtained by pumping selected J ′ of the G (62) intermediate level. In some cases, the v ′ = 0 level of the J 0 + u state is accessed as the common intermediate level. The level mixing in the J, v ′ = 0 ∼ G (62) system provides sufficient ion-pair state character for specific rotational levels such that they are equally suitable to access high lying vibrational states of X 1 Σ + g (0 + g ). Lower ground state vibrational levels (v ′′ = 36-61) are obtained by exciting rotational levels of v ′ = 8, 20, and 28 of the G 0 + u state. These levels have sufficient FC-overlap with the ground state X 1 Σ + g (0 + g ), v ′′ = 0 level and can be excited directly. In fact, Rohlfing and Valentini 31 have reported excitations of the G-states up to v ′ = 39 and vibrationally resolved fluorescence emission to levels as high as v ′′ = 72.
In addition to the experimentally determined molecular constants in this work, high-resolution data from Ram et al. 56 are included for the analysis. The authors applied Fourier transform emission spectroscopy and obtained precise constants in the range of v = 0-3 for both isotopologues. The complete line list is provided in the supplementary material to this publication. The resulting vibronic origins Gv and rotational constants Bv from the leastsquares fit 55 for the two main isotopologues are listed in Tables IV and V. The root mean square (rms) value of the fit is 0.02 cm −1 , which is less than the linewidths of the dye lasers. Statistical uncertainties are given in units of the last significant figure. No allowance is given for systematic errors in the calibration. The constants include the centrifugal distortion constants (up to Ov), which are obtained by the analysis of the effective Gv and Bv values, as described in Sec. IV B. Figure 2 displays a Birge-Sponer plot. The measured ΔG v ′′ + 1 vs v ′′ data are plotted for 63 Cu 2 (red circles) and 65 Cu 63 Cu (blue squares). The heavier isotope values are shifted by applying isotope relations. The green line represents a Dunham extrapolation for v < 60 assuming a Morse potential with the parameters ωe and ωexe. To account for the negative curvature for v > 60, the Dunham parameters ωeye and ωeze are included for the extrapolation procedure (red line). Birge-Sponer extrapolation by using Eq. (1) of order 6 (light blue). The physically correct 1/R 6 behavior at the asymptotic limit is included for two NDE curves: NDE(LIF) is obtained from vibrationally resolved fluorescence measurements (Ref. 31) for v ′′ ≤ 72. The black solid trace depicts the NDE approach to the rotationally resolved data up to v = 102 in this work. The inset shows an enlarged view of the dissociation region (see text for details).
isotopologue are shifted in accord with the usual isotope relationships. Up to v ≈ 60 ΔG v ′′ + 1 2 is, in a first approximation, linear. In this region, the potential is described by a Morse potential with Dunham parameters 57 ωe = 266.30(5) cm −1 and ωexe = 0.9939(9) cm −1 . A fit with these parameters gives the straight line [Morse (v < 60), solid green line] in Fig. 2. The dissociation energy De for a Morse potential is hc ω 2 e /4ωexe = 17 838(18) hc cm −1 at the vibrational level ωe/2ωexe ≈ 134. It is obvious, that the extrapolation strongly overestimates the dissociation energy and fails to describe ΔG v ′′ +1/2 values for higher ground state vibrations. A relatively large rms value of the residuals (1.3 cm −1 ) indicates that the linear approximation is inaccurate by neglecting the negative curvature of the higher lying vibrational levels. For an improved result, it is necessary to extend the Dunham parameters to ωeye and ωeze (DH, red line). The fit results are listed in Tables I and II. The fit parameters ωe, ωexe, ωeye and ωeze are in agreement with the results obtained from the vibrationally resolved fluorescence measurements from Rohlfing and Valentini. 31 The first anharmonicity parameter, ωexe, is close to the Fourier transform value from Ram, 56 while the second anharmonicity parameter ωeye is higher. The parameters are expected to be more accurate by considering the larger range of vibrational levels. The parameter pertaining to the rotational constant, Be, and hence the equilibrium distance re is in accordance with the high resolution study of Ram et al. The vibration-rotation interaction parameter, αe is about 2% smaller than the reported value of Ram. Again, the difference is most likely due to the larger dataset measured in this work. The parameters listed in the last column ["theory (this work)"] were obtained by using the analytical potential Vr defined in Eq. (19) below. The agreement with the parameters obtained from the analysis of the experimental data is satisfactory.
An estimation of the dissociation energy is obtained by a Birge-Sponer extrapolation to the expression Applying the equation with polynomials of the order n = 4, 5, and 6 and integrating under the fitted ΔG(v + 1 2 ) curves yields an average dissociation limit, De, of 16 118(25) hc cm −1 . The average (v + 1 2 ) intercept is evaluated to vD = 108.2(3). The Birge-Sponer fit of order 6 is shown in Fig. 2 (solid, light blue). However, it has been argued that the Birge-Sponer method is limited because the correct longrange behavior of the potential is not adequately described. 32,[58][59][60][61][62] In Sec. IV B, the dissociation energy is evaluated more accurately by taking into account the appropriate asymptotically dominant inverse-power contribution to the potential.

B. Determination of the ground state potential function of Cu 2
To analyze the vibrational energies and rotational constants (Gv and Bv, respectively) listed in Tables IV and V, the procedures reported by Ji and co-workers 63 and Liu et al. 64 are adapted and briefly outlined in the following. The applied "near-dissociation expansion" (NDE) has been introduced by Le Roy (Ref. 65 and references therein) and contains the theoretically known dissociation behavior. Compared to the conventional Dunham expansion, 66 the NDE expressions are more reliable for an extrapolation to high vibrational levels beyond the observed data.
The X 1 Σ + g (0 + g ) state dissociates to the Cu( 2 S) + Cu( 2 S) atomic limit. Thus, the vibrational spacing lying near dissociation is mainly governed by a long-range van der Waals dispersion (or induceddipole-induced-dipole) potential of the form where De is the dissociation energy and C 6 is the leading long-range potential coefficient. A reliable value estimate for the latter is available from a recent ab initio study. 67 The dispersion interaction is computed at the coupled cluster level [CCSD(T)] and amounts to   Tables IV and V and J ≤ 50. Both isotopologues are fitted simultaneously. The values for the heavier isotopologue are shifted by using the usual isotope relationships. The rms uncertainty of the residuals is 0.03 cm −1 , which is smaller than the specified laser bandwidth. All values are in cm −1 except re in Å (see text for details).

Parameter
Rohlfing 1.154 × 10 6 hc cm −1 Å 6 . An accuracy of this theoretical value was not established by the authors. However, as has been observed by Le Roy and Lam 32 and verified in this work, the results are very insensitive to this value (vide infra).
The equation for the vibrational energies of isotopologue-α is given by and for the rotational constant, D is a non-integer effective vibrational index at dissociation and μα is the reduced mass. For the C 6 -potential of Cu 2 adapted in this work, K 0 and K 1 are given by the equation where and Xm are known numerical factors (7931.949 and 546.64 for m = 0 and 1, respectively). 65 The rational polynomial with orders L and M, is applied for the vibrational energies. For the rotational constants, an exponential expansion is used. The expansion parameters p 0 i , p 1 i , q 0 i , the effective vibrational index, vD, and the energy at the dissociation limit De are fitted by taking into account the constants of both main isotopologues (Tables IV and V) simultaneously. The expansion parameters for the isotopologues are related by p The exponent power S is S = 1 for an "outer" expansion or S = 2n/(n − 2) for an inner expansion n = 6, the power of the asymptotically dominant leading inverse power term appropriate for the ground-state dissociation of Cu 2 . Together with the order of the polynomials L and M, the chosen model determines how fast the transition to the limiting behavior takes place. In order to determine the best model that reproduces the term energies to within the experimental accuracy, all possible combinations of the expansion   Tables IV and V by applying the lmfit software. 69 192 models for the rational polynomial orders L = 0, 1, . . ., 7 and M = 0, 1, . . ., 5 in Eq. (7) with a leading expansion coefficient of τ = 1 62 were computed for inner and outer Padé expansions in order to obtain model-averaged estimates of De and vD. Only four models reproduced the term energies of both isotopologues simultaneously to within an accuracy of rms <0.1 cm −1 . The best fit by applying an inner Padé approximant with nine parameters (L = 7 and M = 2) is listed in Table III and reproduces the data within the experimental accuracy (rms = 0.05 cm −1 ). The complete dataset of the model calculations is provided in the supplementary material to this report. The averaged dissociation energy of the four models is 16 270 (7) hc cm −1 and the averaged limiting vibrational level at dissociation is vD = 113.6(10). For the rotational constants, Bv, an exponential NDE expansion with eight parameters is applied [Eq. (8)]. An excellent fit is obtained with an uncertainty for the residuals of rms = 2 × 10 −4 cm −1 .
The NDE coefficients are then used to compute a Rydberg-Klein-Rees (RKR) representation of the potential energy function The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp up to the dissociation limit employing the software RKR1 70 from Le Roy. Since the measured transitions contain levels with high rotational quantum numbers and large internuclear distances (see the line list in the supplementary material), higher-order centrifugal distortion constants are included in the analysis. Therefore, (i) the radial Schrödinger equation is solved for the evaluated potential by using the program LEVEL. 71 (ii) The computed centrifugal distortion constants (up to Ov) are held fixed in a subsequent least-squares fit procedure, as described in Sec. IV A, in order to obtain refined Gv and Bv values. (iii) The RKR potential and centrifugal distortion constants are re-evaluated by applying the optimized NDE parameters. The procedures (i) to (iii) are repeated until conversion is achieved, which occurs after two iterations. The important advantage of the NDE approach is its capability to extrapolate Gv and Bv values beyond the measured values and construct a physically reliable RKR potential up to the dissociation asymptote. The ΔG v+ 1 2 curve vs v ′′ shown in Fig. 2 (black solid line) illustrates the effect of the inverse-power expansion in the long range. The point of inflection at v ≈ 110 yields a higher number of vibrational levels below the dissociation. Note that the earlier application of the NDE equations based solely on vibrational term values up to v ≤ 72 31 and lacking rotational characterization significantly overestimates the dissociation limit (violet solid line in Fig. 2).
The RKR potential obtained with the converged Gv and Bv values in Tables IV and V is depicted in Fig. 3 (NDE, red solid line). The turning points from the RKR inversion procedure are marked with red dots in the figure. In addition, rotationless vibrational wave functions are shown for a few vibrational levels (v = 0, 20, 50, and 102). The wave function for the highest observed vibrational level at v = 102 shows clearly that the molecule remains predominantly close to the outer turning point. The probability of residing close to the inner wall is much smaller. Obviously, Franck-Condon factors for vertical transitions from low-lying vibrational levels are strongly impeded. The optical-optical double-resonance method applied  in this study circumvents this hindrance by taking advantage of gateway states of mixed J ∼ G character, which enable optical pumping from initially low (v = 0-2) to high (v ≥ 62) lying vibrational states.

A. Ground state potential energy function
The black and green solid lines shown in Fig. 3 are spline interpolations V spline of the X 1 Σ + g ground state energies calculated ab initio within the "NpC" and "N" scheme, respectively. The reference energy is −3307.4167 E h for the "NpC" and −3307.4174 E h for the "N" calculation, both occurring at r = 220 pm. The dissociation energy in the X 1 Σ + g ground state is approximately De ≈ 16 700 hc cm −1 for the "NpC" and 16 600 hc cm −1 for the "N" calculation.
Close to equilibrium, the RKR and ab initio potential energy functions agree well on the scale of the graph. The agreement is also reasonably good asymptotically for the potential function obtained from the "N" calculation (green line), which leads to a dissociation energy that is about 300 hc cm −1 larger than the RKR result. From the "N" calculations, one obtains an interpolated equilibrium bond length of 222 pm, and the first three vibrational terms Gv for the 63 Cu 2 isotopologue are 272 cm −1 , 533 cm −1 , and 788 cm −1 for v = 1, 2, and 3, respectively, while the RKR potential yields the same value for the equilibrium bond length and term values that are 8 cm −1 , 6 cm −1 , and 1 cm −1 lower (see Table IV). The terms from the "NpC" calculations are 258 cm −1 , 519 cm −1 , and 779 cm −1 for v = 1, 2, and 3, respectively. Qualitatively, the potential function from the "N" calculation hence agrees fairly well with the RKR result. We understand the poor agreement of the "NpC" and RKR potential functions in the intermediate region between 300 pm and 500 pm as stemming from the lack of accuracy of the former, which is induced by a too strong compromise in attempting to optimize higher lying states leading to the "sp" asymptote, as discussed above.
A comparison with results from previous studies might be useful. In Ref. 72, the "experimental" RKR potential for the electronic ground state of K 2 is compared with a theoretical potential function obtained from a one-electron pseudo-potential treatment of the electronic structure. The result presented in Fig. 1 of that work indicates an apparent good agreement. However, the lines for the two potentials compared in that figure are scaled such that the asymptotes are identical, leaving a difference of about 200 hc cm −1 at the equilibrium. When the potentials are shifted so that the equilibrium energies are equal, as in Fig. 3, the ab initio potential turns out to be flatter than the RKR potential.
In Fig. 2 of Ref. 73, a potential function obtained from a very high level ab initio calculation for the ground state of F 2 is compared with a RKR potential. The agreement is perfect. Closer analysis of Table X of Ref. 74 shows that the excellent agreement is related to inclusion of correlation energy obtained from an extrapolation to the full configuration interaction (FCI) and complete basis set limit that the authors of those works dubbed "correlation energy extrapolation by intrinsic scaling (CEEIS)." This additional correlation energy that, in essence, mimics FCI has an important variation of the order of 10 mE h (10 millihartree, 1 mE h ≈ 219 hc cm −1 ) toward smaller interatomic distances, which effectively renders the theoretical potential energy function stiffer, leading to a good agreement with the "experimental" RKR potential similar to the potential function from the "N" calculations for Cu 2 in the present work. Currently, the implementation of CEEIS to the copper dimer case is not practicable, however.
We note that, while in F 2 , spin-orbit coupling also contributes to a variation on the order of mE h along the ground state potential energy function (see Table IV in Ref. 73), this coupling is much smaller in the ground state of Cu 2 from the present calculations, which is also confirmed in Ref. 75.

B. Dissociation energy
The dissociation energies are De ≈ 16 700 hc cm −1 for the "NpC" and 16 600 hc cm −1 for the "N" MRCI-calculation. These values are about 400 hc cm −1 and 300 hc cm −1 larger, respectively, than the highly accurate value obtained from the measurements presented in the present work and the corresponding value of the RKRpotential. Table VI presents 30 In addition to the MRCI treatment, all electron CCSD(T) calculations were carried out using the same AVQZ-DK basis set, which is comparable to the aug-cc-pVQZ-PP basis of Ref. 75 The pseudo-potential treatment is thus only slightly superior to the treatment adopted in the present work. The CCSD(T) calculations are inappropriate to render the potential energy functions depicted in Fig. 3, however, because of the multi-configurational character of the wave function for r → ∞, as discussed below.
We might conclude that, albeit rather costly and involved, calculations at the MRCI or CCSD(T) level are necessary to warrant an accuracy of at least 300 hc cm −1 to 400 hc cm −1 for the calculated value of the dissociation energy (1 kcal mol −1 ).

C. The nature of the Cu-Cu bond
With these results in mind, we may now turn to the discussion of the nature of the bond in the copper dimer. The MRCI wave functions are of the general form where ΦI (CSF) are spin and symmetry adapted configuration state functions (CSFs). Table VII gives the main coefficients for two wave functions of 1 Σ + g character having the lowest energies. CSF 1 describes essentially a σg[4s]-bond, CSF 2 is the σu[4s] anti-bonding configuration state, CSF 3 describes an attenuated σg[4s]-bond and an enhanced σg[3d z 2 ]-bonding configuration state, and CSF 4 is the anti-bonding state corresponding to CSF 3.
Close inspection of Table VII shows that the ground state (state 1) is essentially a single-configuration state with a doubly occupied σg[4s] orbital, a typical σ-bond case generated by two 4s frontier orbitals. Roughly, the Cu-Cu bond is to 90% a σg[4s]-bond and to 10% a σg[3d z 2 ]-bond.
The CSFs in Table VII are molecular orbital configuration state functions. They may be expressed in terms of valence bond covalent and ionic states as follows: where At equilibrium, state 1 is hence a nearly 50% mixture of the covalent and ionic valence bond states formed by the two 4s frontier orbitals with normalized weights w (cov) ≈ 0.686 2 + 0.066 2 ∝ 61% and w (ion) ≈ 0.550 2 + 0.049 2 ∝ 39%.
At the asymptote, combining with the coefficients from Table VII, one obtains At r → ∞, state 1 represents hence essentially the dissociation of a bond formed by two 4s orbitals ("ss"-channel), while state 2 describes the dissociation of a bond formed by the 4s and 3d z 2 orbitals ("sd"-channel).
Not shown in Table VII is that, at the asymptote, the ground state results from an important mixing of the reference states, which describe the "ss" and "sd"-channels at the multi-configurational self-consistent-field level of calculation. In a single reference multiconfigurational interaction calculation, this mixing would not appropriately be taken into account, and consequently, the energy at the asymptote would potentially be considerably higher than that calculated here, leading to a larger bond dissociation energy.

D. Discussion
The description of the electron pair bond as a superposition of covalent and ionic states was originally proposed by Pauling 76 within valence bond theory. When there is an overwhelming weight of one of the two forms in the total wave function, the bond may be termed either covalent or ionic, depending on whether the covalent or ionic configuration state functions are preponderant. When the weights of these forms become similar, the bond is better described as a chargeshift bond, which is a typical situation for electron rich bonding partners such as in F 2 . 29 Our results reported in Eqs. (15) and (16) indicate that the bonding in Cu 2 should be considered to be of the charge-shift type, which is in agreement with results from valence bond calculations. 30 In the latter work, however, the weight of the covalent state (structure 1 in Table II of Ref. 30) is 71%, whereas that of the ionic states (structures 2 and 3) is only 26%. The weight of the covalent part in a valence bond calculation, where frontier orbitals rather than molecular orbitals are optimized, is somewhat larger than that obtained within a molecular orbital calculation, such as in the present work, where molecular orbitals rather than frontier orbitals are optimized. The best value obtained in Ref. 30 for the binding energy in Cu 2 is 14 235 hc cm −1 (40.7 kcal mol −1 ; Table I of  Ref. 30).
From Table VII, one concludes that the configuration state functions involving singly occupied 4s orbitals contribute to (nonnormalized) ∼(0.874 2 + 0.096 2 ) ∼ 77% of the ground state wave function, whereas those involving each singly occupied 4s and 3d z 2 orbital contribute to less than 1%-the remaining 22% are distributed among many excited configuration state functions not shown in the table. The small contribution from the singly occupied 4s and 3d z 2 orbitals is only one part of the "orbital splitting" inferred in Ref. 30. By this procedure, the occupation of two different orbitals is meant, which are localized on the same center. It allows us to improve on seizing electron correlation. The other part stems from 4s and 4pz correlation, visible in the MOLPRO output, but not influential in the characterization of the wave function. Both the covalent and ionic portions of the ground state wave function of dicopper stem from the interaction between the 4s orbitals.
One should also note that a similar analysis of the MRCI ground state of dihydrogen yields w (cov) ≈ 58% and w (ion) ≈ 42%. In Ref. 29, the covalent part in the total wave function of H 2 is w (cov) ≈ 76%. Hence, the weight of the covalent part in the total wave function of H 2 turns out to be larger in a valence bond calculation than in a molecular orbital calculation, too, similar to the case of Cu 2 . However, while the bond in the latter is considered to be of charge-shift type, 30 that of the former is considered to be of covalent character. 29 From the analysis of the present molecular orbital calculations, both molecules would be considered to be of charge-shift type. One obvious consequence of this analysis is that one cannot rely on the relative weights of the covalent and ionic portions of the wave function alone to conclude about the nature of the bond. Analysis of bond critical points is necessary in addition, as well as the partitioning of the total energy into pure covalent, ionic, and charge-shift contributions.
Another consequence is that the charge-shift attribute of a bond depends on the calculation method, as much as does the analysis of the bond critical points and energy partitioning. These are not observable quantities, however, and so, the "disclosure of the true nature of the chemical bond," as suggested in Ref. 29, might not be possible from experimentally available information.

E. Alternative analytical potential energy function
The RKR potential energy function discussed in Sec. IV B is one possibility to represent the ground state potential of Cu 2 in terms of turning points at given energies. Such representations are very accurate but not easily available, in general, and deriving them for polyatomic molecules is practically not feasible. For polyatomic molecules, one often needs to resort to simplified, compact analytical representations of potential energy surfaces in one or more dimensions. One dimensional analytical potential energy functions along a bonding coordinate, such as the Morse potential, are valuable ingredients of representations of potential energy hyper-surfaces 77 and it is useful to test how accurately analytical representations can describe the interaction potential of molecules such as Cu 2 , where the spectroscopic structure can be very well assessed experimentally up to the dissociation energy.

ARTICLE scitation.org/journal/jcp
In the case of Cu 2 , the Morse potential was found to be insufficiently flexible to describe simultaneously even the lowest spectroscopic terms and the dissociation energy. Instead, the following analytical function, Vr(r), which describes a Morse-like potential, is more flexible: where The quantities V e , r e , A, B, E, and R 6 are adjustment parameters. They are varied in the definition domains V e ≥ 0, r e ≥ 0, A ≥ 0, B ≥ 0, 0 ≤ E ≤ 1, and R 6 ≥ 0. The function Vr(r) was originally proposed in Ref. 78, where graphical representations are discussed for varied forms that this function can take, depending on the choice of parameter values. The term described by the function e(r), Eq. (22), captures the correct asymptotic Vr ∼ −C 6 /r 6 behavior related to the dispersion or induction interaction between two neutral fragments upon dissociation of a diatomic molecule, as discussed in Refs. 78 and 79.
The adjustment to the ab initio energy points was carried out with a modified version of the Levenberg-Marquardt algorithm 80 in which additional non-linear constraints among adjustable parameters can be incorporated. These are the dissociation energy, the C 6 coefficient, and the Morse potential expression for the fundamental transition, which is only an approximation for the actual value that can be obtained for it from the solution of the Schrödinger equation [in Eq. (25), μ is the reduced mass of the dicopper molecule].
The function Vr is depicted in Fig. 3. It overlaps almost perfectly with the RKR potential in the vicinity of the equilibrium distance as well as in the asymptotic region. The agreement is somewhat poorer between 350 pm and 500 pm. Vr yields accurate fit parameters from a Dunham analysis, which are reported in the column "theory (this work)" in Table I. The theoretical parameter values deviate within 1% from the experimental ones with the exception of ωeye, ωeze, and γe, which could be related to the poor representation of the potential by the analytical function in the intermediate range of bond distances. All in all, the function defined by Eq. (19) can be considered to yield a satisfactory overall description of the potential.

VI. CONCLUSIONS
In this work, 347 double-resonant transitions to excited vibrational levels in the X 1 Σ + g (0 + g ) state have been measured with rotational resolution. By accessing high-lying perturbed intermediate states, SEP-type transitions to vibrational levels of the ground state X 1 Σ + g (0 + g ) up to v ′′ = 102 were accessible. The substantially increased accuracy of the measured band positions and the extended vibrational energy range up to ≈98% of the dissociation energy allow for a precise determination of the potential function and the dissociation energy.
The NDE analysis applied simultaneously for the 63 Cu 2 and 65 Cu 63 Cu isotopologues reproduces the measured level energies within the experimental accuracy and leads to the determination of a faithful Rydberg-Klein-Rees (RKR) representation of the potential function up to the dissociation. A compact analytical representation in the form of Eq. (19) is also proposed, which yields slightly less accurate ro-vibrational term values in the higher energy domain, while reproducing the experimental dissociation energy and thus being valuable for its simplicity and compactness. In parallel, high-level ab initio calculations at the internally contracted multi-reference configuration-interaction level of theory have been performed. Calculations based on an active space containing the ground and singly excited states that lead to the neutral "ss" and "sd" asymptotes yield a ground state potential curve in good agreement with the experiment. Attempts to account for the higher lying states leading to the "sp" limit without inclusion of higher lying orbitals, such as with the "NpC" calculations reported in Ref. 34, compromise the accuracy of the ground state potential in the asymptotic region. To obtain the missing correlation energy calculations within larger configuration interaction spaces, approaching the complete basis set limit, such as with the CEEIS method, 74 would be required. Such calculations, however, are not currently feasible for Cu 2 . We hope to be able to improve, in future work, the accuracy of The Journal of Chemical Physics ARTICLE scitation.org/journal/jcp both the ground and excited state calculations of dicopper within a more adaptive approach of the full configuration interaction space by a selective inclusion of 4p orbitals. The congruence of the theoretical and experimental results sheds some light on the nature of the Cu-Cu bond. We clearly showed that bond formation in Cu 2 is governed by the interaction of two 4s frontier orbitals leading to a doubly occupied σ-orbital in the ground state at equilibrium. Close to equilibrium, the wave function thus obtained is an almost 1:1 mixture of covalent and ionic configuration state functions, which is quite typical for homonuclear bonds of essentially monovalent species. Configurations involving pairs of 4s and 3d orbitals play a role in the assessment of the correlation energy but are of minor importance with respect to the bond type, which in valence bond theory would be best termed a mixed covalent-ionic bond.
The attribute "charge-shift" was given to the bond type in dicopper as a result of valence bond calculations. 30 It is compatible with the findings from the present calculations using molecular orbital theory. However, the bond in dihydrogen could be classified in the same way within molecular orbital theory, if solely the weights are considered by which the ionic and covalent configuration state functions intervene in the multi-reference configuration interaction wave function of the ground state at the equilibrium distance. Yet, dihydrogen is said to be covalent within valence bond theory. 29 Clearly, bond types are calculation method dependent attributes depending on the system. Observables must not depend on the calculation method, when calculations are sufficiently accurate. To the best of our knowledge, a direct experimental observation of the type of a homonucleus chemical bond has not yet succeeded.

SUPPLEMENTARY MATERIAL
See the supplementary material for observed and calculated lines for the two main isotopologues of Cu 2 in the range of vibrational levels from v ′′ = 36-98 and v ′′ = 36-102 for 63 Cu 2 and 65 Cu 63 Cu, respectively, a table that contains also the PUMP transitions to the intermediate excited state for the SEP-type double resonance measurements, and a table that lists the expansion parameters for the 192 model calculations, which were performed to determine NDE polynomials that reproduce the term energies for both isotopologues simultaneously to within an accuracy of rms <0.1 cm −1 .

DATA AVAILABILITY
The data that supports the findings of this study are available within the article and its supplementary material.