The correlation theory of the chemical bond

The quantum mechanical description of the chemical bond is generally given in terms of delocalized bonding orbitals, or, alternatively, in terms of correlations of occupations of localised orbitals. However, in the latter case, multiorbital correlations were treated only in terms of two-orbital correlations, although the structure of multiorbital correlations is far richer; and, in the case of bonds established by more than two electrons, multiorbital correlations represent a more natural point of view. Here, for the first time, we introduce the true multiorbital correlation theory, consisting of a framework for handling the structure of multiorbital correlations, a toolbox of true multiorbital correlation measures, and the formulation of the multiorbital correlation clustering, together with an algorithm for obtaining that. These make it possible to characterise quantitatively, how well a bonding picture describes the chemical system. As proof of concept, we apply the theory for the investigation of the bond structures of several molecules. We show that the non-existence of well-defined multiorbital correlation clustering provides a reason for debated bonding picture.

In this work we provide the true multiorbital correlation theory, consisting of a framework for handling the structure of multiorbital correlations, a toolbox of true multiorbital correlation measures, and the definition together with an algorithm for the multiorbital correlation clustering. The presented theory significantly outgrows the multipartite entanglement theory 27 , on which it is based, namely, in the last three items mentioned just above. (The detailed construction is presented in the Supplementary). We adopt the principle that bonds are where the electrons can freely move among atoms, and this is reflected in the correlations of occupations of localised orbitals. Then we show illustrative results by investigating the multipartite correlations in several molecules showing multiorbital bonds. We give quantitative characterisation how well a bonding picture describes the chemical systems. We also illustrate that in the debated case of the dicarbon molecule, there is no well-defined multiorbital correlation clustering, which provides a reason for the ambiguous bonding picture [29][30][31][32][33] . This is not only the first true multiorbital correlation based study of the chemical bond, but also the first application of true multipartite correlation based techniques in physics.
We emphasise that the notions of correlations are basis dependent. We employed two basis sets in this study, standard STO-3G and STO-6G with optimised exponents, which latter provided HF energy close to HF/cc-pVTZ level of theory (see Methods section), together with the localisation procedure of Pipek and Mezey 34 to produce atomic-like orbitals. All results discussed about the correlation structure of localised orbitals are understood with respect to this localisation. For the prototypical molecules which were considered for illustrating our theory, we have found that employing the minimal unchanged STO-3G basis set is sufficient for the description of bonding, and using a basis set closer to the Complete Basis Set limit have not changed the bonding picture. The results using optimised STO-6G basis set are presented in the main text, while the results using unchanged STO-3G basis set are also presented in the Supplementary for comparison.
We note that our work is not connected to previous works of de Giambiagi, Giambiagi and Jorge 35 regarding generalised bond indices based on density-density correlation functions.

Multiorbital correlations
For the quantum mechanical description of the molecule, we use the second quantized picture, that is, the Hilbert space of the electronic system is built up by the one-orbital Hilbert spaces, describing the occupation of the orbitals 36 . An orbital can be unoccupied, occupied with one electron of spin up or down, or doubly occupied by two electrons of spin up and down, resulting in four-dimensional one-orbital Hilbert spaces. In the Hilbert space formalism of quantum mechanics, any linear combination of orbitals is an orbital, however, the interpretation, or physical properties single out some of them. Correlations among orbitals are not invariant under such nonlocal (among-orbital) operations. In order that the correlations express some connection among local objects (atoms), it is necessary that the orbitals are localised on the atoms.
So, for the description of (the electronic system of) the molecule, we consider m localised, atomic-like orbitals. Let M, stands for "molecule", denote the set of (the labels of) these orbitals. We aim at describing the correlations in an L ⊆ M set of orbitals (cluster). (If L = M then the correlations in the whole molecule is considered). In general, the state of the full electronic system of the cluster L can be described by the density operator 1, 37, 38 ρ L . The reduced state of a (sub)cluster X ⊆ L can be described by the reduced density operator 1, 37, 38 ρ X . If the cluster of orbitals L can be described by a state vector |ψ L 〉 (for example, when a given eigenstate of the whole molecule is considered), then its density operator is of rank one, ρ L = |ψ L 〉〈ψ L |, called a pure state. Its reduced density operator is mixed (not of rank one) in general, which is the manifestation of the entanglement 4 of (sub)cluster X and the rest of the cluster L\X.
The correlation can be defined with respect to a split of the M set of the orbitals 27,39 . Let ξ denote such a split, that is, a partition 40 ξ = {X 1 , X 2 , …, X |ξ| } ≡ X 1 |X 2 |…|X |ξ| , where the clusters X ∈ ξ, called parts, are disjoint subsets of the cluster L, and their union is the full cluster L. A natural comparison of partitions is the "refinement": we say that partition υ is finer than partition ξ, if the parts of υ are contained in the parts of ξ. The set of the partitions of L is denoted with Π(L). (For illustrations, see Supplementary Fig. S1). The measure of correlation with respect to this split is the ξ-correlation 27 , Here S(ρ) = −tr ρ ln ρ is the von Neumann entropy 1,37 . As a special case, the i|j-correlation being the well-known (two-orbital) mutual information 1,37,41 , has already been used 7, 8, 10, 12-20, 26, 28 . The ξ-correlation is zero if the state is uncorrelated with respect to ξ (it can be written in a product form of reduced states of clusters X ∈ ξ); and nonzero otherwise, characterising the strength of the correlation among the parts X ∈ ξ. This comes from the information-geometrical meaning of this quantity: it characterises how "far" the state is from the states uncorrelated with respect to ξ. (For more details of the construction, see the Supplementary). Note that C ξ is larger for finer partitions, (this is called multipartite monotonicity 27 ), it is zero for the trivial split ξ = = L { }  , and it takes its maximum, C ⊥ , for the finest split ξ = ⊥ = {{i} | i ∈ L}. The latter quantity is also called total correlation [42][43][44][45][46] , (Note that if cluster L is described by a pure state, e.g., if L = M, then S(ρ L ) = 0, and the correlation is entirely quantum entanglement 4,5,27 ). It is easy to check the following sum rule 46 , valid for any partition ξ, X X X L L , that is, the total correlation is the sum of the total correlations inside the parts plus the correlation with respect to the partition.
We would also like to characterise the correlations in an overall sense, that is, without respect to a given partition. There are several ways of this 27 , here we consider two of them. Let us introduce the k-partitionability correlation and the k-producibility correlation, respectively, for 1 ≤ k ≤ |L|. These characterise two different (one-parameter-) notions of multiorbital correlations. The k-partitionability correlation is zero if the cluster can be split into at least k parts which are uncorrelated with one another, and the correlations are restricted only inside those parts; and nonzero otherwise, characterising the strength of this kind of correlation. In general, C k-part increases with k, and it jumps after the number k of parts into which L can approximately be split. The k-producibility correlation is zero if the cluster L contains correlated (sub)clusters of size not larger than k; and nonzero otherwise, characterising the strength of this kind of correlation. In general, C k-prod decreases with k, and it jumps at the size k of the largest part in the partition into which L can approximately be split. As special cases, C |L|-part = C 1-prod = C ⊥ grabs all the correlations, it is zero if there is no correlation at all in the cluster L, that is, its state is a product of the states of orbitals; and nonzero otherwise. On the other hand, C 2-part = C (|L| − 1)-prod is sensitive only for the strongest correlations, it is nonzero if the cluster L is globally correlated, that is, its state is not a product of states of two (or more) clusters; and zero otherwise.
For other values of k there are no such coincidences among the partitionability and producibility correlations, however, the relation C k-part ≥ C (|L|−k+1)-prod holds. Also, the bounds The tools (1) and (5), despite being so simple, are proven very useful in a wide range of applications for the characterisation of multiorbital correlations in the electronic system of molecules. In the sequel, we show four of these applications. Illustrating these, we present numerical results for several prototypical molecules, namely benzene, pyrrole, borole, cyclobutadiene, furan, thiophene, and the sequence C 2 H 2x for x = 1, 2, 3 and C 2 .

Applications
Application 1: Molecule, formed by bonds. Our fundamental principle is that, in the equilibrium, the bonds are almost uncorrelated with one another, and the orbitals involved in a bond are strongly (multiorbital) correlated. Using the tools introduced above, we formulate this principle, and we demonstrate it for the aforementioned molecules.
An ansatz for the bond structure is given by a partition β = B 1 |B 2 |…|B |β| ∈ Π(M) (bond split), representing the molecule as a set of bonds (represented by B ⊆ M sets of orbitals), together with some nonbonding orbitals (e.g., core orbitals or lone pairs, for those, |B| = 1). Then the β-correlation C β (ρ M ), given in (1), characterising the correlation with respect to the ansatz β, expresses how well this ansatz describes the physical situation: the lower the C β the better the ansatz from a purely information-theoretical point of view. The aim of this application is to find the bond split β (if exists) from ab initio data, without taking into account anything which can a priory be known about the bond structure in quantum chemistry. We call this multiorbital correlation clustering.
Since in a real electronic system one cannot expect such a simple ansatz to be exactly valid (C β = 0), we actually pose the question, which ansatz β is the best choice for the description of the bonds from a physical chemical point of view. Being the best, however, is a delicate question. Note that, on the one hand,  C = 0, for the trivial split  β = , which takes the whole molecule to be one big bond. On the other hand, C β grows with respect to the refinement of β, and takes its maximal value C ⊥ , the total correlation (3), for the finest split β = ⊥, which excludes nontrivial bonds. These extremal cases, obviously, do not give proper descriptions of the bond structure of a molecule, since, on the one hand, there can be clusters weakly correlated but not uncorrelated with the remaining part of the molecule, on the other hand, it is not allowed to neglect strong correlations. Instead of these, we have to be able to split the molecule into weakly correlated clusters consisting of strongly correlated orbitals. In order to grasp the meaning of the multipartite correlation clustering, we have to be able to decide about a given ξ, if it is a good ansatz, or it is worth considering a ξ′, which is "a bit" finer than ξ. That is, we have to investigate the difference C ξ′ − C ξ , where ξ′ is finer than ξ, and there is no other partition between them. We seek β, for which, while ξ is coarser than β, this difference is small if ξ′ is coarser than β, but large, if ξ′ is not coarser than β. If there exists such a β, then it is meaningful to consider the electronic system to be weakly correlated bonds consisting of strongly correlated orbitals, and this is described by β. (For the whole construction, see the Supplementary).
Here we face the problem that verifying that this definition holds for a given partition (calculating S(ρ X ) for all clusters X ⊆ M, needed for the calculation of C ξ for all partitions ξ) is numerically prohibitive. We can decrease the demands by successive refinement of the partitioning (bipartitioning of one cluster in each step) following the smallest increase in C ξ . One can show that if there exists a β satisfying the definition above, then the successive refinement goes through β, C ξ increases slowly until β, and rapidly after β. (For the whole construction, see the Supplementary).
We can also have a hint for the path of the successive bipartitioning. Consider the γ = G 1 |G 2 |…|G |γ| clustering based on the "connectivity" with respect to the two-orbital correlations (2). That is, the parts G ∈ γ are the sets of orbitals being the vertices of the connected components of the two-orbital correlation graph 8,28 . (This is the graph with vertices being the orbitals i ∈ M and with edges of weights C i|j (ρ {i,j} ) above a threshold T b . We call this two-orbital correlation clustering). It is proven to be a good strategy to do the successive bipartitioning with respect to γ, that is, not to split apart the parts of γ. Following this strategy, one reaches γ; C ξ increases rapidly after γ, but it is not sure that C ξ increases slowly before γ. This is because of the possibility of the existence of hidden correlations, which is an interesting feature of the multiorbital setting. For example, there exist states in which all two-orbital correlations (2) are zero, but the states are correlated as a whole, that is, they cannot be written in a nontrivial product form (see the Supplementary). This means that if we follow the above strategy, then C ξ may change rapidly before we reach γ. In this case, β does not equal to γ, but coarser than γ.
We have investigated the two-orbital and the multiorbital correlation clustering for the aforementioned molecules.
The two-orbital correlations are drawn by different shades of grey lines in subfigures (a) of Figs 1, 2 and 3. The two-orbital correlation clusterings γ are based on the appropriate threshold values T b . The distributions of two-orbital correlations, and the possible two-orbital correlation thresholds T b leading to the known bond structure in the given cases are shown in subfigures (b). For C 2 H 2 , there is a much wider range for ′ T b , leading to triple bond in γ′, than for T b , leading to double bond in γ, and for C 2 , there is a much wider range for ′ T b , leading to quadruple bond in γ′, than for T b , leading to triple bond in γ. A drawback of the two-orbital correlation clustering method is that, although the two-orbital correlation (2) is bounded by 0 ≤ C i|j ≤ 2(ln4) uniformly, a uniform threshold covering all the cases is contained in a quite narrow range 0.269(ln4) < T b ≤ 0.307(ln4). The reason for this is that an orbital seems to be forced to share its (two-orbital-)correlations among the ones strongly correlated with it, which may be a manifestation of the monogamy of entanglement 47 in correlations. (Different thresholds for the different cases may be obtained based on the separation of the correlation scales, however, this leads to a bond-interpretation rather arbitrary).
The multiorbital correlation clustering β, determined by the use of the method described above, is drawn by solid red lines in subfigures (a) of Figs 1, 2 and 3. The values of C ξ during the successive bipartitioning are shown in subfigures (c). In the cases of the cyclic molecules (Figs 1 and 2) we could find a well-defined β bond split, after which the value of C ξ jumps about at least twice as large as the maximal step before that. Note that in certain positions, some humps are observed in the tendencies C ξ (designated with red rectangles). These are the effects of correlated clusters of size more than two: When the successive bipartitioning reaches such a cluster, following the first large step, smaller steps become possible, leading to this concave behaviour. Such humps are coming from two origins. The more characteristic one is the cluster of the inner bonding 2p z orbitals (denoted with X in in the figures). In the cases of borole and cyclobutadiene, these humps can be found directly before β, they are not steep enough to keep these orbitals together during the successive bipartitioning, contrary to those in the cases of benzene, pyrrole, furan and thiophene, when these humps can be found after β. In this way we can distinguish between aromaticity and antiaromaticity. On the other hand, in the cases of furan and thiophene (Fig. 2), there is an additional weaker multiorbital correlated structure in each case besides the aromatic rings, due to the hyperconjugative interaction of the lone pair with the adjacent σ-bonds 19 (denoted with X hc in the figures). The correlations in them are not strong enough to keep these orbitals together during the successive bipartitioning. (The almost uniform increase after β comes from the bipartitioning of the two-orbital clusters). In the cases of the C 2 H 2x sequence (Fig. 3), it can be seen how the correlation picture becomes more and more fuzzy. Interestingly, for the case of C 2 H 2 , investigating the tendency C ξ during the successive bipartitioning in subfigure (c), one can see that it changes significantly at the partition leading to a double bond β (=γ), and there is a much less significant change at the partition leading to a triple bond β′ (=γ′). For bipartite correlation clustering, γ′ were more plausible than γ, however, here β seems to be more plausible than β′. This is indeed a very interesting observation, which might be to some extent an indication of hidden correlation (β is coarser than γ′). Note however, that despite not being divided in two in our multiorbital correlation point of view, the four-orbital bond does not contradict the classical picture of triple-bonded C 2 H 2 , as it contains four electrons. For the case of C 2 , investigating the tendencies C ξ during the successive bipartitioning in subfigure (c), one cannot give a well-defined bond split β besides the 1 + 8 + 1-orbital partition, because of the high multiorbital correlation of the eight orbitals. (The splits γ and γ′, given by the two-orbital correlation clustering, are drawn by solid orange lines in subfigure (a). In subfigure (c), we show two different paths of the successive refinement in the partitioning of the eight bonding orbitals, the upper one shows a more significant hump, while the lower one leads through the triple bond γ). That is, according to our observations, there exists no well-defined multiorbital correlation clustering. The bonding situation in the multireference C 2 is well known as a long-standing puzzle, and several bond orders have been suggested, including the extensively debated quadruple bond [29][30][31][32][33] . In spite of the fact that four strong two-orbital correlations have been found, from the reason mentioned above, we cannot give a decisive multiorbital correlation answer on the bond order of C 2 . Application 2: Bonds, formed by orbitals. The bonds, that is, the highly correlated clusters, given by the parts B ∈ β, are identified in the previous point. Now, we can investigate the correlations inside the bonds B ∈ β. For this purpose, we use the k-partitionability and k-producibility correlation C k-part,B and C k-prod,B , respectively, (see (5) is high in the four aromatic cases, that is, the orbitals in X in cannot be split even into two parts, and, accordingly, the greatest part is of size |X in |; on the other hand, − C X 2 prod, in and − C X 2 part, in , respectively, − C X 3 part, in are low in the two antiaromatic cases, that is, the orbitals in X in can be split into parts of size at most 2, their number are 2, respectively 3. So we can distinguish between aromaticity and antiaromaticity also in this way. The orbitals participated in the hyperconjugative interaction in furan and thiophene (contained in X hc ⊂ M) show weaker correlation than the ones in the aromatic ring. In the cases of the C 2 H 2x sequence, the orbitals participating in the bonds between the carbon atoms (contained in X in ) are getting more and more multiorbital-correlated. For C 2 H 4 , − C X 2 part, in and − C X 2 prod, in are near zero, that is, the two two-orbital bonds can be considered independent. For C 2 H 2 , − C X 2 part, in and − C X 4 prod, in are negligibly low, while − C X 3 part, in is significant, leading again to a double-bond, containing a four-orbital one. For C 2 , − C X 2 part, in , although being relatively low, does not seem to be completely negligible. In the latter two cases, we can see now from a local point of view, which was proposed in the previous Section from a global point of view, that X in is not divided completely into independent bonds.    Fig. 1. (d) The correlations for the orbitals participated in the hyperconjugative interaction, contained in X hc are also shown. (e) For the thiophene, the correlations among the valence orbitals, contained in ′ A S are shown. this can be represented by the split α = A 1 |A 2 |…|A |α| ∈ Π(M) (atomic split) of the molecule. Here an important quantity is the α-correlation C α (ρ M ), and the α-coarsened k-partitionability and k-producibility correlations C k-part,α (ρ M ), C k-prod,α (ρ M ). These characterise the different aspects of the correlations with respect to the atomic split α.
The atomic split α for the aforementioned molecules are drawn by dashed blue lines in subfigures (a) of Figs 1, 2 and 3. The values of C α are also shown. Calculating C k-part,α (ρ M ) and C k-prod,α (ρ M ) is infeasible, due to the large density matrices of high entropy, however, note that we already have the largest members of these hierarchies, since C |α|-prod,α (ρ M ) = C 1-prod,α (ρ M ) = C α (ρ M ). The value of this is near C ⊥ (ρ M ), that is, as can be expected, the atoms are strongly correlated with one another in the molecules. Application 4: Atoms, formed by orbitals. The orbitals localised on given atoms are collected in the parts A ∈ α in the previous point. Now, we can investigate the correlations in the atoms A ∈ α. For this purpose, we use the k-partitionability and k-producibility correlation C k-part,A and C k-prod,A , respectively, (see (5)), both of them are considered with respect to the splits Π(A).
We have investigated the nontrivial (non-H) atoms in the aforementioned molecules. Although not all the positions of the C atoms are equivalent in the molecules, the correlation measures take roughly the same values for those. The full increasing and decreasing tendencies C k-part,A (ρ A ) and C k-prod,A (ρ A ) are shown in subfigures (e). Note that the values of these are usually smaller than the correlations in the bonds, by about two orders of magnitude. In the sequence C 2 H 4 , C 2 H 2 and C 2 we can also see, how the increase of the multiorbital correlations leads to more and more strong correlations among the orbitals localised on the same C atom. The hyperconjugative interaction leads to the same results on the O and S atoms in furan and thiophene. Remarks on the applications. Having the results of all the four applications in hand, we can now observe how the sum rule (4) works. In the first two applications, when we considered the β bond split, ∑ B∈β C ⊥,B (ρ B ) was large and C β (ρ M ) small; while in the second two applications, when we considered the α atomic split, ∑ A∈α C ⊥,A (ρ A ) was small and C α (ρ M ) large, and these are connected by the sum rule (4) as Based on these, we can consider the molecule as the weakly correlated set of strongly correlated bonds, or the strongly correlated set of weakly correlated atoms. Note that this holds for the equilibrium structure, which is the only one considered here. If the internuclear distances are altered, which is a method for the investigation of bond-formation 12,[14][15][16]20 , we expect that the above picture is altered accordingly, however, the sum rule (6) holds with altered numerical values.
On the other hand, we may give a definition of the molecule from a correlation point of view: the orbitals M form a molecule, if there exists no nontrivial partition which is coarser than α, describing the atoms, and β, describing the bonds, that is,  α β ∨ = . (In these cases, intermolecular bonds do not appear in β).

Conclusions and outlook
We have presented a novel theory of the chemical bond which is inspired by quantum information theory and based on multiorbital correlations. Contrary to the literature, where only two-orbital notions were considered, we have invented and used true multiorbital notions. Illustrating the use of this theoretical toolbox, we have investigated several small prototypical molecules and showed how in a black-box manner the bonding picture of a molecule naturally comes out from the multiorbital correlations of occupations of localised atomic-like orbitals.
We have identified the bonds with the strongly correlated clusters, and characterised quantitatively how well a given bonding picture describes these molecules. Our tools are, e.g., able to distinguish between aromaticity and antiaromaticity in cyclic conjugated systems. On the example of the sequence C 2 H 2 , C 2 H 4 , and C 2 H 6 , we have seen that the increase of wide-range multiorbital correlations results in the decrease of the well-posedness of multiorbital correlation clustering. In the extreme case of C 2 , this leads to the nonexistence of a well-defined multiorbital correlation clustering, which provides a reason for the debated bonding picture.
We would like to emphasise again that the treatment in terms of true multiorbital correlations seems to be a very natural point of view in the investigation of bonding among more than two atoms. The multiorbital correlation based quantities have their statistical meaning on their own right, and we have already seen several of their applications. However, it would be interesting to relate them to other standard quantities in quantum chemistry, quantifying, e.g., bond strength or aromaticity. Besides aromaticity, this treatment may find applications also in multicenter transition metal cluster chemistry.
We have seen how the multiorbital correlations characterise the chemical bonds, if the orbitals are localised. We note, however, that the theory can also be applied to any (orthonormalized) sets of orbitals, then it characterises the correlations among those orbitals. These, of course, do not have to be related to the chemical bonds, but may be related to other chemical properties.
From the point of view of theoretical power and beauty, the multiorbital correlation theory provides a much more natural and flexible treatment for multiorbital situations than using only two-orbital correlations, done in preceding works. An example supporting this is given by the (6) application of the sum rule (4). Contrary to this, a treatment based only on two-orbital correlations is theoretically hard to grasp, due to monogamy-like issues of correlations in quantum systems. This is why the notion of hidden correlations is not well-defined, and to formulate a clear-cut (quantitative and/or operative) definition is an open question.

Methods
For the numerical results shown in this paper we have performed calculations using the quantum chemistry version of the density matrix renormalization group (QC-DMRG) algorithm 36,[48][49][50][51][52][53][54][55][56] . We have controlled the numerical accuracy using the dynamic block state selection (DBSS) approach 36,44,57 and the maximum number of block states varied in the range of 500-4000 for an a priory set quantum information loss threshold value χ = 10 −6 . The ordering of molecular orbitals along the one-dimensional topology of the DMRG was optimised using the Fiedler approach 10, 15 and the active space was extended dynamically based on the dynamically extended active space (DEAS) procedure 6,36 . We have used DMRG to obtain the optimised MPS wavefunction, which was then used to construct the reduced density matrices, from which the correlation measures (1) and (5) were calculated.
Geometries have been optimised at HF/cc-pVTZ level of theory which yielded sufficient geometries in accordance with higher level methods. To obtain the localised atomic orbitals for the DMRG procedure, we first optimised the exponents of the STO-6G basis set using the MRCC program 58-60 which approach resulted in similar HF energy to the cc-pVTZ basis set result within 10 −2 Hartree. Then we used the Pipek-Mezey procedure 34 implemented in MOLPRO 61 Version 2010.1, with tight threshold 10 −12 , and minimised the number of atomic orbitals contributed to each localised orbitals. All localised orbitals have been used in the DMRG procedure thus, as a result, we have carried out calculations at the FCI limit for all molecules. Then the results close to the FCI limit have been analysed in the paper 19 . We note that we also calculated all results using HF/STO-3G geometry and localised STO-3G orbitals, with literature value exponents, and found neglectable difference compared to the results presented in the manuscript. (These results are presented in the Supplementary for comparison). This suggests that our analysis is very robust in general. We note that this robustness is not entirely surprising. Mayer has shown 62-64 that extracting chemical information from molecular wavefunctions such as bond orders and valence indices could also be obtained using only STO-3G basis set.
The datasets generated during and/or analysed during the current study are available from Sz.Sz. on reasonable request.