text:
# Atomistics of Tensile Failure in Fused Silica: Weakest Link Models Revisited

## 1 INTRODUCTION

Leonardo da Vinci knew that longer wires fail at lower tensile loads than shorter ones. For a simple chain this result is so obvious it hardly bears mention: a chain is only as strong as its weakest link; if any link fails the chain fails. Adding links to a chain will weaken it if one of the added links is weaker than any of the initial links, or will leave its strength unchanged if all added links are stronger than the weakest initial link. Adding links to a given chain will inevitably weaken it if enough links are added, if the strengths of the links are drawn randomly from a non-singular distribution. The weakest link found in a large sample (a long chain) will on average be weaker than the weakest link found in a small sample (a short chain).

These results are obvious to us, but may have been less obvious in Leonardo’s day, before the invention of the science of mechanics (and statistics), when the mechanics of a chain was imperfectly understood. In fact, it is not completely trivial, for it depends on the loading conditions, the properties of the material of which the links are made, and the definition of failure. For example, if the links are capable of work-hardening and failure is defined as extension of the chain by a suitable fixed length, then a long chain may actually be able (depending on its detailed properties) to sustain a greater load than a short one.

## 2 WEAKEST LINK MODELS

From the picture of a chain in tension a simple idealized model, known as the weakest link model, is derived. In this model a fiber (wire, filament or rod) loaded in tension is described as a series of independent mechanical elements ordered in a single file, each of which transmits a load (force, or momentum flux) from its immediate neighbor on one side to its immediate neighbor on the other side. The model is quasi-static: all elements are assumed to be subjected to the same tensile load, which will be the case if there are no accelerations and no forces on the links except those exerted by their immediate neighbors. It is also local, in the sense that force is transmitted only across the boundaries between two elements; this assumption may fail if the elements are coupled by long range forces, such as electromagnetism, as in a current-carrying or charged chain or one made of piezoelectric material, although these effects are unlikely to be significant in practice. Failure of the entire chain is defined as the failure of a single element. Given these assumptions, the strength of a chain is the strength of its weakest link.

The statistics of weakest link models were apparently first considered by Peirce, who worked for the British Cotton Industry Research Association and published in the Journal of the Textile Institute. This is curious because weakest link models are unlikely to be applicable to ropes and cables made up of many parallel strands. These strands divide the load amongst themselves, so that the failure of a single strand need not, and usually does not, imply the failure of the entire rope. We are all familiar with ropes that are partly frayed, but which still carry a load. Multi-stranded ropes are also non-local, because the shift of load from a failed strand to the remaining strands occurs by static friction over an extended length of the rope, and depends in a complex manner on the lay of the rope. In addition, load shifting occurs as a result of differential initial length, elastic stiffness and creep between strands, even when there has been no failure.

The statistics of weakest link models were reviewed by Epstein, whose paper remains useful today. He made a number of significant points:

1. If the links have a statistical distribution of strengths, then the mean strength of the entire chain (the strength of the weakest link) is a decreasing, albeit generally very slowly decreasing, function of the number of links in the chain.
2. The distribution of chain strengths is generally a skew function, with a longer “tail” in the direction of weaker chains; this is a consequence of the fact that a single unusually weak link makes a weak chain, while an unusually strong link has no effect on the strength of the chain.
3. For normally distributed link strengths, the dispersion in chain strengths is a slowly decreasing function of the number of links.
4. In many physically important problems the weakest link model is inapplicable; he cites notches on the surface of a bent beam as an example. A single notch serves as a site of stress concentration, while a large number of closely spaced notches (a roughened surface), if of equal depth, reduce the stress in the matter below them to low values, thus reducing the stress at the apex of any one notch.

Weakest link models are expected to be applicable to the tensile failure of brittle materials in many, but not all, circumstances. The reason for this, known since the work of Griffith, is that brittle materials usually fail in tension because of stress concentration at imperfections—surface scratches, voids, or heterogeneities in the bulk. Once local failure occurs, the size of the imperfection and its degree of stress concentration increase, and the object ruptures.

Suppose a macroscopic brittle object may be divided into small cells, each of which contains no more than one imperfection, and each imperfection is contained within a single cell. In most materials these cells must be very small because the density of imperfections is large, but in pristine glass fiber the imperfections may be so rare that a cell can encompass the entire diameter of the fiber and extend in length for many diameters. If such a division into cells is possible, with negligible interaction between the different cells and the imperfections within them, then we would expect a weakest link model to be applicable, because the failure of any one cell would lead to the failure of the entire object. This is true even though the cells need not be arrayed in a line like links in a chain; even cells side-by-side are independent potential sites of failure.

An essential assumption is that imperfections do not interact elastically. This is valid if the imperfections are small and well separated in all their dimensions, but may not hold if they are geometrically extended, such as grain boundaries in a ceramic or surface scratches. In such more complex geometries it may again be possible to use weakest link models if the cells are large enough to contain many imperfections, for then each cell may act as an independent link, even though the individual imperfections interact strongly. In such a case, the resulting statistics will be those of the cells, not those of the individual imperfections.

It is a curious feature of models of this kind that a renormalization group transformation is not possible (their statistics are not invariant under such a transformation) because the statistics of weakest link models depend on the total number of links, which is not invariant. This follows from the assumption that the flaws have a characteristic (small) size, and are not “fractal”. It also assumes that a single local failure is sufficient to produce general failure; in other words, that general failure is not a consequence of the interaction of many failures in the type of process now known as self-organized criticality.

The dependence of the statistics of weakest link models on the number of participating links means that it may be possible to learn something about the microphysics of failure from the statistics of the strength of macroscopic objects. That is the subject of this paper.

## 3 FUSED SILICA

Fused silica and related silica glasses are materials of enormous technological importance. They are also among the few materials (the others are carbon nanotubes and, when creep is negligible, metallic whiskers) which may approach their ideal limiting strengths because they can be produced apparently without volume or surface stress-concentrating flaws.

### 3.1 History

Griffith was apparently the first to realize the importance of flaws in reducing the strength of brittle materials. Griffith also found that finer glass fibers were systematically stronger than thicker ones, which may naturally be attributed (in a weakest link model) to the presence of more links in the thicker fibers. He compared this result to the similar results of Karmarsch from 1858 on metal wires, although we would now consider brittle (glass) and ductile (metallic) failure to be so different that a close comparison is not justified. For historical reviews see the work of Holloway and Kurkjian.

Griffith measured tensile strengths under ordinary room conditions of short lengths of very fine (3–4$`\mu `$ diameter) glass fibers of about 500,000 psi, which are comparable to those measured under similar conditions today. He also noted the effect of aging and of static fatigue in reducing these strengths, which remain important subjects of research. He speculated that the strength of these fine fibers could be extrapolated to infinitesimal diameters to determine the limiting strength of ideal glass. This extrapolation is inconsistent with his hypothesis that he had produced flaw-free samples and not very meaningful because it was based on data obtained under ambient conditions for aged samples; we now know that ideal strength requires an inert environment. Despite this, his numerical value ($`1.6\times 10^6`$ psi) for the ideal strength of glass is quite close to modern values.

Some of Griffith’s data are also surprising, given modern understanding. He found (Table 1) a fairly smooth dependence of strength on diameter, while we would expect fibers to be either without mechanical damage and very strong, or damaged and very weak; finer fibers would be expected to be more frequently found undamaged, but not systematically stronger than undamaged thicker fibers. We would also expect the strength to be essentially independent of diameter for the finest (rarely damaged) fibers, so that no extrapolation to zero diameter would be necessary or appropriate. After the lapse of 80 years it is impossible to interpret all his data, but some of his results may be a consequence of his testing of aged samples and his use of fiber diameter as an independent variable, in contrast to modern experiments which generally use communications fiber of a single diameter.

At the time of Griffith’s work no theory, reasonable in modern terms, of the ideal strength of materials existed. Approximate theories were soon developed. Polanyi used a simple energetic argument to predict the ideal strength of crystals from the energy required to form new surfaces by rupture. Although he did not refer to the (one year earlier) work of Griffith, his surface energy-based argument resembles Griffith’s argument for the growth of a crack, differing chiefly in that Griffith considered (using Inglis’s analytic solution) stress concentration at the crack tip while Polanyi made the naïve assumption of a cleavage with plane parallel faces. From this Polanyi correctly concluded that for most materials the measured strength is far below any reasonable theoretical estimate, but did not point out stress concentration as the explanation.

### 3.2 Properties

The development of quantum mechanics and the theory of interatomic forces made possible more quantitative models of the theoretical strength of materials, beginning with the work of Frenkel. Modern band structure theory makes reliable predictions of ideal strengths possible for crystalline materials, although these are generally unobservable because their much lower real strengths are determined by stress concentration at flaws for brittle materials and by dislocation motion for ductile materials. The modern methods are applicable to carbon nanotubes, which may perhaps be made free of defects because they are so small. They may also undergo plastic flow, which may prevent them from reaching their ideal strength, even if free of defects, and they are known to buckle reversibly on bending.

Fused silica and silica glasses have apparently not been the subject of modern theoretical strength calculations, probably because it is difficult to specify their non-crystalline structure. This structure can be calculated by molecular dynamics methods, but the result is dependent on the model for the interatomic forces and represents only a small subvolume of the bulk material; extension to the bulk requires a questionable extrapolation (such as periodic boundary conditions). $`\beta `$-cristobalite, a high-temperature crystalline phase of silica with cubic symmetry and a density close to that of fused silica, may be a useful model, but its properties do not appear to have been calculated from first principles and even its elastic constants have not been measured.

The development of fiber optic communications produced a renewal of research on fused silica. Essentially defect-free fiber is now available in enormous quantity, generally encapsulated in a protective coating which prevents mechanical damage and reduces the rate of environmental chemical attack. The tensile strength of this material, directly measured at low temperature (77$`^{}`$K) where chemical attack is believed to be negligible, may reasonably be taken to be its ideal strength, and is about 140 Kbar. Low-temperature measurements in which the stress is inferred from the strain in two-point bending give similar results, but are less direct because the reduction of the data depends on applying thin beam theory to a geometry in which it is not quantitatively applicable (because the radius of curvature is so small) and because the nonlinear stress-strain relation of fused silica is somewhat uncertain. The measured low-temperature strength is about 0.2 of the Young’s modulus, in accord with rough estimates of the strengths of the chemical bonds. Strengths measured under ambient conditions are about one third of those measured at low temperatures, and depend on the temperature, humidity, and duration of load, indicating complex processes of stress corrosion by water.

## 4 MICROSCOPIC PHYSICS

Now that the theory of strength-reducing cracks is well understood, and the strength of fused silica has been measured reasonably well, what more can be learned?

### 4.1 Mechanisms

It took many years before the importance of surface damage, as opposed to volume flaws, for the strength of glass was appreciated. However, the strength of an unflawed (pristine) sample of glass, its ideal strength, is likely determined in its bulk, by the strength of its network of covalent bonds. This does not answer the question of which component of the glass, or which physical process, actually determines that strength; it is possible that its strength is actually determined by the bonds at its surface.

A number of candidate strength-limiting components or processes need be considered. On the largest scale, a small heterogeneity may provide a slight degree of stress concentration and therefore determine the strength. The extraordinarily small measured dispersion ($`1`$%) of the strength of pristine glass and the at least semi-quantitative approach of the measured low-temperature strength to its estimated ideal value appear to make this implausible, unless the degree of stress concentration is very slight, as might be produced by heterogeneities of small amplitude (a small amplitude fluctuation in elastic constants). In optical communications fiber the germania-doped core provides such a fluctuation, but it is ordered rather than random and homogeneous (to good accuracy) along the length of the fiber.

Some heterogeneity is implicit in the fact that glasses are amorphous, and on the smallest spatial scale significant heterogeneity arises from the fact that Si—O bonds are not all equivalent; each has a slightly different length, and the O—Si—O and Si—O—Si bond angles differ from their ideal values. Thus, in a piece of glass subjected to tension in a specified direction different bonds will fail at different stresses.

In addition to SiO2, real glass, even fused silica, contains other components. There are a variety of impurities, reduced to very low levels in communications-grade optical fiber, but present nonetheless. These introduce bonds other than Si—O bonds, and these other bonds will generally be weaker (as is shown by the fact that mixed silicate glasses are less strong than fused silica). Some of these impurities (such as Na2O and CaO in soda-lime glass) are present at high abundances, while others (such as transition metal oxides in communications-grade fiber) are very rare, but all are sites of weakness. There may be other kinds of strength-limiting defects, such as non-stoichiometrically bonded atoms, interstitials, etc.. In the core of optical communications fiber several percent of the Si has been replaced by Ge. Because Ge atoms are slightly larger than Si atoms the bond network is further distorted there, and at the boundary between the pure SiO2 and the GeO2-substituted regions there is systematic strain, analogous to imperfect lattice parameter matching in epitaxial growth of crystals. These effects are distinct from heterogeneity of bulk elastic constants referred to above.

Finding which of these mechanisms actually determines the strength of glass, particularly fused silica, is of both scientific interest and practical importance, for some mechanisms may be affected by changes in composition or processing, potentially leading to stronger fibers. Analogous considerations may apply to the more complicated problem of room-temperature strength, where the resistance to stress-assisted chemical attack determines the strength, rather than ideal mechanical failure alone.

### 4.2 Can Failure Statistics Help?

For each of the possible modes of failure a macroscopic (or even mesoscopic) piece of glass contains a large number of well-separated elastically independent regions. If the failure of a single region leads to general failure, a plausible but unproven assumption, then a weakest link model will apply. The specific mechanisms of failure differ greatly in the number $`N`$ of independent potential failure sites present per unit volume, or in a given sample. For example, in a 1 m length of 125$`\mu `$ diameter fiber the number of elastic heterogeneities or geometrical irregularities may be very small ($`N1`$$`10^4`$), the number of Si—O bonds is $`N10^{21}`$, and the number of rare impurities may (depending on the impurity) be in the range $`0N<10^{15}`$.

One means of determining the sources of failure is to perform measurements on specimens in which the number of each type of potential failure site is varied. This is difficult in practice, partly because it may not be possible to control their number, and partly because each series of experiment would require the production of a batch of fiber of custom composition.

It is a general feature of weakest link models that as $`N`$ increases the strength of the specimen (as in the classical chain or wire) decreases. If the statistics of the strengths of the individual links were known this fact, by itself, would be used to determine $`N`$ from the dependence of strength on size, but these statistics are rarely known, and not in the present case. It is also often difficult to perform experiments with a sufficiently large range of specimen sizes. Similarly, the dispersion of strength among specimens depends on $`N`$ and could be used to determine it, but again only if the statistics of the individual links are known. These two effects are closely related, and do not give independent information.

Quantitative analysis depends on the assumed statistics of the links, and many possible distributions have been considered. One possible model, closely related to classical Weibull statistics, describes each link by a fracture readiness parameter $`\alpha `$, which may be considered a stress concentration factor, or the reciprocal of a link’s ideal strength. The values of $`\alpha `$ are distributed according to a distribution $`f(\alpha )`$, which is defined for $`\alpha >0`$ and normalized

$$_0^{\mathrm{}}f(\alpha ^{})𝑑\alpha ^{}=1.$$
(1)
The sample fails if the largest fracture readiness parameter found among the $`N`$ links, $`\alpha _{max}`$, exceeds a value $`\alpha _0`$. The probability that it does not fail is, to good approximation if $`N1`$,

$$P(\alpha _0)\mathrm{exp}\left(_{\alpha _0}^{\mathrm{}}Nf(\alpha ^{})𝑑\alpha ^{}\right).$$
(2)
The probability that failure occurs for a value of $`\alpha _0`$ between $`\alpha `$ and $`\alpha +d\alpha `$ is $`P(\alpha )d\alpha `$, where

$$P(\alpha )=Nf(\alpha )\mathrm{exp}\left(_\alpha ^{\mathrm{}}Nf(\alpha ^{})𝑑\alpha ^{}\right).$$
(3)
The most useful way to parametrize the results is as the ratio of the width $`w`$ of $`P(\alpha )`$ to the value $`\alpha _{max}`$ at which $`P(\alpha )`$ is a maximum; in terms of the Weibull modulus $`m`$

$$\frac{w}{\alpha _{max}}\frac{1}{[m(m1)]^{1/2}},$$
(4)
and the approximation is almost exact if $`w`$ is defined as the dispersion of $`P(\alpha )`$

$$w\left(\left|\frac{d^2\mathrm{ln}P(\alpha )}{d\alpha ^2}\right|_{\alpha =\alpha _{max}}\right)^{1/2}$$
(5)
and $`m1`$.

The stretched exponential function is defined

$$f(\alpha )=\frac{C(\nu )}{\alpha _0}\mathrm{exp}\left[\left(\frac{\alpha }{\alpha _0}\right)^\nu \right].$$
(6)
This is a general form widely used when the actual functional form is unknown, and includes the simple exponential and Gaussian as special cases. The normalizing constant $`C(\nu )\nu /\mathrm{\Gamma }(1/\nu )`$ and $`N^{}NC(\nu )/\nu `$. The distribution $`P(\alpha )`$ is shown in Figure 1 for $`\nu =2`$ (a Gaussian).

By successive approximations,

$$\frac{w}{\alpha _{max}}\frac{\left\{1\frac{1}{2\mathrm{ln}N^{}}\left[\frac{\nu 1}{\nu }\left(\frac{\nu 1}{\nu }\right)^2\mathrm{ln}\mathrm{ln}N^{}\right]\right\}}{\nu \mathrm{ln}N^{}\left(1\frac{\nu 1}{\nu }\frac{\mathrm{ln}\mathrm{ln}N^{}}{\mathrm{ln}N^{}}\right)^{1/\nu }}\frac{1}{\nu \mathrm{ln}N^{}}.$$
(7)
It is evident that if $`f(\alpha )`$ is, or can be fitted to, a stretched exponential or to one of its special cases useful and plausible estimates of $`N\mathrm{exp}[\alpha _{max}/(w\nu )]`$ can be obtained. However, the inferred value of $`N`$ depends very sensitively on $`\nu `$, and measurement of $`m`$ alone for a sample of test specimens of a single size does not determine $`\nu `$. Measurement of two or more populations of very different-sized test specimens of the same material (for which $`N`$ is proportional to the size) may determine both $`N`$ and $`\nu `$, and may be feasible; for example, in two-point bending experiments on optical fiber of 125$`\mu `$ diameter the number of atoms $`N_a`$ at significant risk of initiating fracture (those with stresses within about $`1/m`$ of the maximum, where the Weibull modulus $`m>100`$, which are found only close to the outside of the sharpest part of the bend) is $`N_a5\times 10^{14}`$, while in tensile loading of a 50 m gauge length of fiber $`N_a4\times 10^{22}`$ atoms are uniformly loaded. Two or more measurements of $`m`$ in which very different numbers of atoms are stressed permit simultaneous determination of both the ratio $`N/N_a`$ and $`\nu `$, although no extant data serve the purpose, in part because bending and tensile measurements are affected differently by variations in fiber diameter.

### 4.3 Skewness

In order to obtain useful results, even with data of unprecedented accuracy, it will probably be necessary to measure more than the dispersion of the strength. It has long been known, both theoretically within weakest link models and empirically, that distributions of failure strengths are strongly skewed, with outliers preferentially found in the direction of anomalously weak specimens. The reason for this is, of course, that an unusually weak link produces an unusually weak specimen, while an unusually strong link has no effect at all on the specimen strength. Quantitative measurements of the skewness of the strength distribution may constrain the statistics, and, more importantly, the number of individual links contributing.

The skewness of $`P(\alpha )`$ is defined:

$$s\frac{_0^{\mathrm{}}(\alpha \alpha _{max})^3P(\alpha )𝑑\alpha }{w^3_0^{\mathrm{}}P(\alpha )𝑑\alpha }.$$
(8)
The skewness is not small; see Fig. 1. It is also not readily estimated analytically because the Taylor expansion of $`P(\alpha )`$ does not converge sufficiently rapidly, but it may be calculated numerically. Values of the skewness, as a function of $`\nu `$ and $`N`$, are shown in Fig. 2.

The values plotted were computed using cutoffs on the integrals of $`\pm 5w`$ from $`\alpha _{max}`$, with $`w`$ defined self-consistently using the same cutoff. The reason for this is that $`P(\alpha )`$ has a long tail extending toward increasing $`\alpha `$ which contributes significantly to the skewness, but which is unlikely to be observed in a real experiment with a reasonable number of specimens because there will probably be no specimens that far out in the the tail of the distribution of strengths. The skewness computed without this cutoff is significantly larger, typically by $`O(10\%)`$.

When comparing experimental data to Fig. 2 a similar cutoff must be applied to the data. This will, in addition, exclude samples which are anomalously weak because of mechanical damage or other gross flaws, which otherwise must be excluded ad hoc. If the number of measurements is not large it may be necessary to choose a narrower cutoff, and to recompute Fig. 2 accordingly.

### 4.4 Beyond Skewness

Skewness is not the only statistical parameter beyond dispersion which may be useful. In principle, there are an infinite number of moments of the strength distribution, but the higher moments (and the skewness, too) are difficult to calculate reliably and accurately, even from data of high quality, because they are very sensitive to outliers. A data set of reasonable size will not adequately sample rare outliers, and, in addition, there will generally be an additional population of outliers resulting from damaged specimens or errors of measurement. It is difficult, and perhaps impossible, to decide in an unbiased manner if an outlier should be discarded as a probably consequence of a bad specimen or error, or retained in calculating the skewness and higher moments, and even a single outlier may have a substantial effect on their calculated values. Therefore, it is better to truncate the data distribution deliberately in defining the moments (as was done for the skewness in Fig. 2).

Rather than use moments, it may be more informative (and robust) to compare the complete distributions of strength to models. This avoids the problem of outliers, bur requires large data sets. While the truncated skewness can be determined to useful accuracy with comparatively few data (1000 measurements will typically determine it to $`\pm <0.1`$ if outliers are excluded), comparing distributions divided among many bins may require more data, simply because with few data per bin statistical fluctuations are more important. However, optical communications fiber is cheap, and automated testing machines may permit the measurement of thousands of specimens.

## 5 DISCUSSION

In order to use successfully the methods discussed here it will be necessary to obtain strength measurements of unprecedented accuracy and number. For example, it is believed that the dispersion in existing measurements of the strength of fused silica optical fiber is largely the result of dispersion in fiber diameter, rather than dispersion in the actual strengths of the fiber. This problem is likely to occur whenever small dispersions in strength are measured, and may have contributed to Leonardo’s observation that long wires are weaker than short ones: a long wire is more likely to contain a region of narrower diameter, just as it is more likely to contain an intrinsic flaw, and it is not straightforward to distinguish the two effects.

Equation (7) indicates that the dispersion in intrinsic strengths is unlikely to be much less than 1% unless the individual links have a fractional dispersion in strength $`1`$. Well controlled fused silica fibers have a dispersion in diameter of 0.4% and in cross-section of 0.8%, and an apparent dispersion in strength comparable to that in cross-section. In order to measure the real dispersion in strength much tighter control of fiber diameters will be required. These data already hint that the intrinsic fractional dispersion in strength of individual links may be small, which gives significant information about the nature of the weakest links—perhaps surprisingly, it points towards the SiO2 matrix, in which all chemical bonds are nearly the same, and away from rare impurities, which might be expected to involve a more heterogeneous range of bond strengths. It should be remembered, however, that these measurements were performed at room temperature, and reflect failure by stress corrosion rather than the ideal materials strength.

Weakest link models have been widely discussed for over half a century, but their applicability to real materials remains unproven. In part, this is because analysis of failure statistics rarely goes beyond calculation of the dispersion and Weibull modulus, and plotting on Weibull coordinates. There are few detailed quantitative studies of failure statistics, which would include measurement of the dependence of the dispersion and mean strength on specimen size, as well as quantitative study of the distribution of failure stresses. For specimens affected by heterogeneities and damage the data are unlikely to justify more careful analysis, but optical communications fiber offers both the opportunity to investigate the ultimate limit of ideal materials strength and the possibility of homogeneous data of sufficient quality to make such an investigation possible.

In fact, there are qualitative reasons for questioning the applicability of weakest link models to the ideal strength of glass. As discussed, extant data indicate a surprisingly small dispersion in measured strength. The surprise is increased by the fact that no sample can be completely free of chemical impurities and radiation damage, which would be expected to introduce anomalously weak links. The three-dimensional network of Si—O bonds may not be described by weakest link models, but may instead show a well-defined collective failure limit which leads to homogeneous macroscopic behavior from a microscopically heterogeneous material.

Understanding the statistics of the ideal strength of glass at low temperatures is of interest to the physicist. Understanding its strength at ambient conditions is of much greater practical interest, and nearly all data have been obtained under these conditions. It is also unclear if weakest link models are applicable to stress corrosion; the simple picture of the failure of a single small link leading to a crack propagating across the entire specimen is not obviously applicable. The kinetics of stress corrosion of glass are not understood quantitatively, and the collection of more accurate statistics of its strength distribution may be a new and useful means of attacking that problem.

## ACKNOWLEDGMENTS

I thank A. E. Carlsson for discussions and DARPA for support.


format: text

text:
# Higgs production at Hadron colliders as a probe of new physics

## Abstract

The coupling of the Higgs with gauge bosons $`g_{hVV}`$ ($`V=g,\gamma ,W,Z`$) can receive non-decoupling corrections due to heavy quanta. Deviations from the SM predictions are described by a set of parameters $`ϵ_V`$, which can be related to the parameters $`S,T`$ for specific models. The Higgs production by gluon fusion can be used to probe $`ϵ_g`$ at Tevatron, whereas the promising Higgs decay into $`WW^{}`$ can probe $`ϵ_W`$. We find that the resulting bounds may imply the exclusion of heavy particles that receive their mass directly from the SM Higgs, including additional standard or mirror families, chiral colored sextets and octet quarks. Within the MSSM, we also find that gluon fusion is a sensitive probe for the spectrum of squarks masses.

1.- Introduction. The discovery of the Higgs boson as the remmant of the mechanism of electro-weak symmetry breaking (EWSB), is one of the most cherised goals of present and future high-energy experiments. Within the minimal standard model (SM), the mass of the physical Higgs particle is a free-parameter, but present data seems to favor a moderate mass ($`105<m_h<220`$ GeV) , however this conclusion may be changed by the presence of new phyics with an scale $`\mathrm{\Lambda }O(1)`$ TeV . On the other hand, the minimal supersymmetric version of the SM (MSSM), which is one of its most appealing extensions, predicts a light Higgs boson, with an upper mass bound of about 130 GeV . Detection of the full spectrum of Higgs bosons in the SM and beyond, has been extensively studied, because it constitutes an important test of the possible realization of the Higgs mechanism in fundamental physics.

The characteristic Higgs boson couplings determine the strategies employed for its search at present and future colliders. For instance, at LEP the large couplings of the Higgs with massive gauge bosons allows to use the reaction $`e^+e^{}Z+h`$; whereas at next linear colliders (NLC) it will also be possible to study Higgs production by WW and ZZ fusion. On the other hand, Hadron colliders can also test these couplings, either through the reactions $`ppW+h,Z+h`$, or through the decays $`hWW,ZZ`$; vector fusion can also be used for heavy Higgs masses. The couplings of the Higgs with the heavier fermions, can be studied either by open production of $`t\overline{t}h,b\overline{b}h`$ or by the loop-induced coupling with gluon and photon pairs ($`hgg,h\gamma \gamma `$) . Any additional heavy particle that receive its mass from the SM Higgs mechanism, will couple to the Higgs with strength proportional to the particle mass itself, which means that such heavy quanta will induce non-decoupling contributions to the 1-loop vertices $`hgg,h\gamma \gamma `$. Moreover, since these effects will compete at the same order with the SM loop, it can be a significan effect that may be probed at future colliders, as has been explored in previous studies . However there is another important effect not considered previously, namely that the presence of these new particles could also induce non-decoupling corrections to the tree vertices $`hff`$, $`hWW`$ and $`hZZ`$ , which can affect the decay rate of detectable signatures, and thus must be included in the analysis.

In this papaer we study the effect of additional colored particles in the loop vertex $`hgg`$, and the bounds that could be obtained at Tevatron (RUN-II). We shall identify two cases, in the first (Scenario-I) we consider the situation when the mass of heavy colored particles comes from the SM Higgs boson, but which are weakly bounded by electroweak precision meassurements, i.e. we simply assume that there is some unspecified physics at the scale of the mass of the new particles that makes their contribution to the Peskin-Takeuchi parameters $`S,T`$ to be within experimental range . The second case (Scenario-II) will include a heavy fourth family, for which we find that the corrections induced by the heavy fermions on the vertices $`hgg`$ and $`hWW/hZZ`$ are correlated, and in fact can be written in terms of the parameter $`T`$. We also study gluon fusion for the MSSM higgs bosons, as a possible test for the mass spectrum of squarks. Large effects are obtained when the loop amplitude includes non-degenerate squark masses, which could be realized in many particular models . These results suggest that the discovery of a Higgs boson can be turned into a new tool for high precision electroweak physics studies.

2.- Parametrization of new physics in Higgs couplings. In order to probe heavy scales through their effect on the Higgs coupling, we shall consider extensions of the SM that include additional particles, but with a minimal Higgs sector consisting of one doublet. We are interested in the possible representations that can receive their mass from the SM Higgs mechanism. Because of the quantum numbers of the SM Higgs, we are restricted to consider only fermion doublets and singlets. Then the possibilities reduce to: a) additional families with SM quantum numbers b) additional mirror families and c) a combination of the above. Within cases a and b, new quarks must have the same color properties as in the SM, but since for case c it is possible to cancel anomalies among the quarks, they could lay in larger $`SU(3)_c`$ representations, like sextets or octets. Since these new states have electroweak charges, they will contribute in general to the parameters S,T,U, but not necessarily. Present global fits for electroweak data seem to exclude more than one additional SM family , however this conclusion relies on the assumption that no other physics occurs at the energy scale of the new fermion masses, which may not be the case in SUSY models. For instance, consider a model with a complete SUSY 4th family, in the limit $`m_A>>m_Z`$, under which the light Higgs boson $`h^0`$ behaves like the SM Higgs, and also assume a mass degeneracy for the components of the fermion doublets, as well as for the sfermions, this gives $`T=0`$; whereas fermions give $`S=2N_c/6\pi `$, the contribution of sfermions to $`S`$ is given by $`S=\frac{N_c}{36\pi }log(\frac{m_{\stackrel{~}{u}_{Li}}^2}{m_{\stackrel{~}{d}_{Li}}^2})`$ , thus by choosing the appropriate masses one could decrease the total value of S, and satisfy present experimental constraints. This first scenario, where we have new particles that receive their mass from the SM Higgs mechanism, but are weakly constrained by electroweak precision data, is the right place where the study of the Higgs signal can test the presence of such heavy particles. For the second case, when the corrections to the different vertices $`hVV`$ are correlated, one must include all the corrections to the relevant branching ratios for the analysis of bounds; here we can also obtain bounds on the heavy quanta, though they will be somewhat limited.

We shall describe the effects of heavy quanta on the higgs couplings to gauge bosons, by writing them as: $`g_{hVV}=g_{hVV}^{SM}(1+ϵ_V)`$, where $`V=W,Z,g,\gamma `$, and similarly for the fermion couplings $`g_{hff}=g_{hff}^{SM}(1+ϵ_f)`$. The parameters $`ϵ_{f,V}`$, encompass the effects of heavy quanta. In general, the expression for $`ϵ_{g,\gamma }`$ will be give by a ratio of complicated loop expresions, however since we are interested in the limit $`m_h<<M_{heavy}`$, one can use the low-energy theorems , to relate $`ϵ_{g,\gamma }`$ to the beta-function coefficients for the corresponding coupling constant ($`\beta _{I,X}`$). Thus, the higgs-gluon coupling is given by $`ϵ_g=\beta _{3,X}/\beta _{3,t}`$, where $`\beta _{3,t(X)}`$ denote the contribution of top and heavy quanta to strong beta function. The values of $`ϵ_g`$ arising from a pair of heavy color triplets, sextets and octet quarks are : $`ϵ_g=2,\mathrm{\hspace{0.17em}5}`$ and 6, respectively, whereas a new SM family plus its mirror partner gives $`ϵ_g=4`$. It is remarkable that coming stages of Tevatron will be able to test these values at significant levels.

The correction induced by the heavy quanta on the couplings $`hWW,hZZ`$ can also be derived using the low-energy theorems. However, for the case of a fourth heavy family, one finds an interesting relation between $`ϵ_{f,W,Z}`$ and the corresponding expression for the parameter $`T`$, namely,

$`ϵ_f`$ $`=`$ $`{\displaystyle \frac{\alpha T}{2}}+{\displaystyle \frac{N_cG_F}{12\sqrt{2}\pi ^2}}(m_1^2+m_2^2)`$ (1)
$`ϵ_Z`$ $`=`$ $`{\displaystyle \frac{N_cG_F}{6\sqrt{2}\pi ^2}}(m_1^2+m_2^2)`$ (2)
$`ϵ_W`$ $`=`$ $`ϵ_Z+{\displaystyle \frac{\alpha T}{2}}`$ (3)
where $`T`$ includes the contribution of the fourth family fermions with masses $`m_{1,2}`$.

We could also consider effects arising from heavy particles that receive their mass from the breaking of new symmetries characterized by a large scale $`V_{new}>>v=246`$ GeV, induced by another set of heavy Higgs bosons ($`\mathrm{\Phi }_{new}`$). Supposse this is communicated to the SM Higgs by including in the scalar potential a mixing term of the form: $`\lambda |\mathrm{\Phi }_{new}|^2|\varphi _{sm}|^2`$. Then one can also include the contribution of these ultra-heavy particles to the parameters $`ϵ_V`$. However, we found that in this case the mixing angle that transforms the Higgs weak to mass-eigenstate basis, is suppressed by the large scale ($`V_{new}`$), and it induces decoupling effects at low-energies, namely $`ϵ_gv/V_{new}`$.

3.- The mechanism of gluon fusion and bounds on $`ϵ_g`$ . The contribution coming from the new particles will modify the cross section for gluon fusion, as follows: $`\sigma (pph+X)=\sigma _{SM}(1+ϵ_g)^2`$ where the SM cross-section $`\sigma _{SM}`$ can be written in terms of the Higgs decay width $`\mathrm{\Gamma }(h_{SM}gg)`$; details of the equations can be found in the literature , Since the new particles modify the decay width into $`hgg`$, it will enhance the production by gluon fusion reaction, but this enhancement could work against the final signal rate, which usually involves some decay of the Higgs boson into a mode with clear signature, like $`hWW^{}`$ or $`\gamma \gamma `$ for the intermediate Higgs mass region. Both effects can be taken properly into account in the analysis by writing the product of the cross-section times the branching ratio of the signal as $`[\sigma \times B.R.(hVV)]_{new}=R_V\times [\sigma \times B.R.(hVV)]_{SM}`$, where $`R_V`$ is given by ,

$$R_V=\frac{(1+ϵ_g)^2(1+ϵ_V)^2}{[1+_Y(ϵ_Y^2+2ϵ_Y)B.R.(h_{SM}YY)]}$$
(4)
where the sum in the denominator runs over $`Y=t,b,g,A,W,Z`$.

We can obtain bounds on the parameters $`ϵ_{g,W}`$ using gluon fusion and the decay $`hWW^{}`$ at Tevatron, which was studied in detail in ref. , this work concluded that it is possible to detect a SM Higgs boson with an integrated luminosity of 30 fb<sup>-1</sup>, provided that an optimized selection of cuts is implemented; here we shall only use their first stage of cuts, namely: transverse lepton momentum $`p_t(e,\mu )>10`$ GeV, psedorapidity $`\eta _{e,\mu }<1.5`$, lepton invariant mass $`m_{ll}>10`$ GeV, jet resolution $`\mathrm{\Delta }R(lj)>0.4`$, missing transverse energy $`E_T>10`$ GeV. For the case when the corrections to $`ϵ_W`$ can be neglected (Scenario-I), the above cuts already allow us to get interesting bounds on the parameter $`ϵ_g`$ at 95 % c.l., as it is shown in figure 1 ( for 2 and 10 $`fb^1`$). These bounds will in turn constrain the presence of heavy colored particles coming from cases a, b and c, provided that the Higgs mass lays in the intermediate mass range. The interception of the exclusion lines with the straight dashed line (which corresponds to $`ϵ_g=2`$), shows the values of Higgs masses where it will be possible to limit the presence of a pair of heavy color triplets. Furthermore, since sextets and octet quarks give a larger value of $`ϵ_g`$, their presence could be excluded too.

However, these bounds need some clarfication for the case when the new physics also modifies the parameters $`ϵ_{W,Z}`$ , since this will affect the decay rates into $`hWW^{}`$. For the fourth family case (Scenario-II), the corrections to the vertex $`hVV`$ are negative (as compared with the positive tree-level value) and grow with the mass of the heavy quanta; for instance for fermion masses of order 500 (700) GeV the deviations from the SM tree-level couplings are of order - 27 (-52) % , which seems to imply that the signal will no longer be detectable i.e. one can only probe heavy but not ultraheavy quanta. With 10 $`fb^1`$ of integrated luminosity, it will be possible to exclude only a limited range of fermion masses, up to about 900 GeV at best, as it is shown in table 1. However since for these heavy masses one enters into the non-perturbative domain one can not draw a definite conclusion, unless it is found an scheme to treat the non-perturbative Higgs effects. On the other hand, the decay into photon pair can also be used to probe higgs couplings at LHC; but in this case the W-loop dominates the amplitude, and the contribution from addditional fermions has opposite sign, which will make it difficult to separate the effects, or even worst, may conspire to reduce the signal.

3.- Probing a Non-universal sfermion spectrum. The MSSM includes two Higgs doublets, and the Higgs spectrum consists of two neutral CP-even scalars $`h^0`$ and $`H^0`$, one CP-odd pseudoscalar $`A^0`$ and a charged pair $`H^\pm `$ . The Higgs sector of the model is completely determined at tree level by fixing two parameters, conventionally chosen to be $`\mathrm{tan}\beta `$ and the pseudoscalar mass $`m_A`$. At loop levels, the radiative corrections, mainly from top and stop, modify the tree-level mass bound ($`m_h<m_Z\mathrm{cos}2\beta `$), allowing to have $`m_h^{max}130`$ GeV . The corrections to Higgs couplings can also lead to modification of its production mechanisms at hadron colliders, as previously studied . Although conventional wisdom states that the contribution of top/bottom quarks (squarks) dominate the amplitude for the gluon fusion, and squarks from first and second generations can be neglected, it should be mentioned that this result only holds within the minimal SUSY breaking scenarios, where it is typically assumed that sfermions are mass degenerate, as required to respect FCNC constraints . However, mass non-degeneracy for sfermions within the same family but different isospin can be acceptable, since it will only be midly constrained by the parameter $`T`$. In fact, some non-degeneracy can arises even in the minimal SUSY-GUT, by imposing universal boundary conditions at the Planck scale, rather than the GUT scale. Other cases that require non-degeneracy are models with additional D-terms , and in Finite Grand Unification .

Thus, it is interesting to find some process that could test the sfermion mass spectrum. In this letter we consider the decay $`hgg`$ (which determines the gluon fusion mechanism) to probe the structure of soft-breaking masses predicted by models of SUSY breaking. The expression for the decay width of $`hgg`$ that includes a non-universal mass spectrum for squarks, can be written as,

$$\mathrm{\Gamma }(hgg)=\frac{G_F\alpha _s^2m_h^3}{4\sqrt{2}\pi }|G_h|^2$$
(5)
where

$$G_{h^0}=G_{t,\stackrel{~}{t}}+G_{b,\stackrel{~}{b}}+G_{LR}+G_{UD}$$
(6)
For $`h=h^0`$, $`G_h`$ is given by:

$$G_{t,\stackrel{~}{t}}=\frac{m_t^2}{m_h^2}\frac{\mathrm{cos}\alpha }{\mathrm{sin}\beta }[f_1(\lambda _t)+f_3(\lambda _{\stackrel{~}{t}_L})+f_3(\lambda _{\stackrel{~}{t}_R})]$$
(7)
$$G_{b,\stackrel{~}{b}}=\frac{m_b^2}{m_h^2}\frac{\mathrm{sin}\alpha }{\mathrm{cos}\beta }[f_1(\lambda _b)+f_3(\lambda _{\stackrel{~}{b}_L})+f_3(\lambda _{\stackrel{~}{b}_R})]$$
(8)
$`G_{LR}`$ $`=`$ $`{\displaystyle \frac{m_W^2s_w^2}{m_h^2c_W^2}}\mathrm{sin}(\beta +\alpha )[Q_u(f_3(\lambda _{\stackrel{~}{u}_R})f_3(\lambda _{\stackrel{~}{u}_L})`$ (10)
$`+Q_d(f_3(\lambda _{\stackrel{~}{d}_R})f_3(\lambda _{\stackrel{~}{d}_L})]`$
$$G_{UD}=\frac{m_W^2}{2m_h^2c_W^2}\mathrm{sin}(\beta +\alpha )\mathrm{\Sigma }_{u,d}[f_3(\lambda _{\stackrel{~}{u}_L})f_3(\lambda _{\stackrel{~}{d}_L})]$$
(11)
The expression for $`h=H^0`$ are obtained by making the replacements $`\mathrm{cos}\alpha \mathrm{sin}\alpha `$, $`\mathrm{sin}\alpha \mathrm{cos}\alpha `$, $`\mathrm{sin}(\alpha +\beta )\mathrm{cos}(\alpha +\beta )`$; the explicit form of $`f_1(z),f_3(z)`$ can be found in ref. . The contribution proportional to the fermion masses are kept only for the stop and sbottom ($`G_{t,\stackrel{~}{t}},G_{b,\stackrel{~}{b}}`$), whereas the remmaining SUSY-breaking effects are kept for all sfermions ($`G_{LR},G_{UD}`$). To reduce the number of parameters, we choose our $`\mu ,A_t`$ values in such a way that one can neglect L-R mixing, which was found previously to give small effects. We have used the previous equation to evaluate the contribution of squarks to the parameter $`ϵ_g`$, for $`m_A=200`$ GeV, $`\mathrm{tan}\beta =5,20`$, and squark masses covering the range from 200 to 700 GeV. For simplicity, we apply the 1-loop leading log formula to the Higgs masses and parameters. Results are shown in table 2, and for comparision we have also included the values obtained with universal squark masses. These results show that the effect of non-universal masses can modify the $`B.R.(hgg)`$, by values of order $`\pm 20\%`$.

4.- Conclusions. We have studied the Higgs interaction with fermions and gauge bosons, and found that deviations from the SM prediction, induced by heavy quanta, can be described by a set of parameters $`ϵ_{f,V}`$, which can be tested at future colliders. The production cross-section by gluon fusion can be used to probe $`ϵ_g`$ at Tevatron, whereas the promising Higgs decays into $`WW^{}`$ can probe $`ϵ_W`$. We find that the resulting bounds may imply the exclusion of heavy particles that receive their mass directly from the SM Higgs, including a 4th standard or mirror family or chiral colored sextets and octets. For the 4th family case we disply a relations between $`ϵ_{f,V}`$ and the Peskin-Takeuchi parameter $`T`$. Within the MSSM, we also find that gluon fusion is a sensitive probe for a non-universal spectrum of squarks masses, which allows for the possibility to test the sfermion masses predicted in the models of SUSY breaking. Whether the discovery of a Higgs signal will lead Particle Physics into a dark middle-age or into a renaissance remains to be seen, but the answer will depend crucially on our hability to meassure the Higgs properties, which could reveal us the path to the promised land of new physics or may be just the boring solitude of the SM up to very heavy scales.

We thank M. Chanowitz and C. Kolda for useful discussions. The support from “Fundacion Mexico-USA para la ciencia” and the hospitality of LBL are also acknwoledged. JLDC is also supported by CONACYT-SNI (Mexico).

FIGURE CAPTION

Fig 1. 95 % exclusion contours for the parameter $`ϵ_g`$ that can be obtained from Higgs search at Tevatron RUN-II, with 2 (dashes) and 10 (solid) $`fb^1`$ of integrated luminosity. The straight (dashed) lines corresponds to $`ϵ_g=2`$.

TABLE CAPTION

Table 1. Upper limits to 4th generation fermion masses that can be obtained from Higgs search at Tevatron RUN-II, with 10 $`fb^1`$.

Table 2. Values predicted for the parameter $`ϵ_g`$ for the MSSM with $`\mathrm{tan}\beta =5,20`$. $`M_{\stackrel{~}{Q}}`$ represents the range of values taken for $`M_{\stackrel{~}{U}_R}`$, $`M_{\stackrel{~}{D}_L}`$, $`M_{\stackrel{~}{D}_R}`$.

Table. 1

| $`m_h`$ \[GeV\] | $`m_{4th}^{max}`$ \[GeV\] |
| --- | --- |
| 120. | 350. |
| 130. | 520. |
| 140. | 650. |
| 150. | 730. |
| 160. | 870. |
| 170. | 710. |
| 180. | 670. |
| 190. | 600. |
| 200. | 510. |

Table. 2

| $`m_{\stackrel{~}{U}}`$ \[GeV\] | $`m_{\stackrel{~}{Q}}`$ \[GeV\] | $`ϵ_g`$ (non-univ) | $`ϵ_g`$ (univ) |
| --- | --- | --- | --- |
| 200. | $`200700`$ | $`0.10.4`$ | 0.4 |
| 500. | $`200700`$ | $`0.050.14`$ | 0.05 |
| 700. | $`200700`$ | $`0.050.18`$ | 0.1 |


text:
# Local Quasiparticle States around an Anderson Impurity in a 𝑑-Wave Superconductor: Kondo Effects

\[

## Abstract

The mutual influence between a Kondo impurity and a $`d`$-wave superconductor is numerically studied within the slave-boson mean-field approach. The Bogoliubov-de Gennes type equations are derived and solved using exact diagonalization. We show that a critical coupling strength, above which the Kondo effect takes place, exists regardless of whether the band particle-hole symmetry is present or not. In the Kondo regime the double resonant quasiparticle peaks are found in the local density of states (LDOS) both directly at the impurity and around its neighbors, which is in sharp contrast to the case of nonmagnetic unitary impurities, where the LDOS vanishes at the impurity site.

\]

The identification of $`d`$-wave pairing state in high-$`T_c`$ superconductors  has renewed much interest in low-energy quasiparticle states in unconventional superconductors. Different from conventional $`s`$-wave superconductors having a constant energy gap, the $`d`$-wave superconductors have an anisotropic energy gap which changes sign along the nodal directions on the essentially cylindrical Fermi surface. The anisotropy of the energy gap makes the $`d`$-wave superconductors very sensitive to impurity or other potential scattering, in contrast to the $`s`$-wave superconductors which, as dictated by the Anderson theorem , are almost unaffected by the presence of nonmagnetic impurities. As an important signature of the sign change of $`d`$-wave order parameter (OP), Balatsky, Salkola and co-workers , based on the T-matrix approximation, theoretically predicted that a strong atomic-like nonmagnetic impurity can induce a virtual bound state. Later studies  showed that the resonant state would show up in the local density of states (LDOS), which could be measured by the local differential tunneling conductance. Recent scanning tunneling spectroscopy measured on high-$`T_c`$ superconductor Bi2Sr2CaCu2O8+δ (BSCCO) with individual atomic-scale defects  or strong impurities , has provided strong evidence for the existence of low-energy quasi-bound states. Especially, the experiment by Pan et al.  on the observation of quasiparticle resonant states has a high spatial and energy resolution and the identity of impurity atoms is known in a controlled way. Very interestingly, the imaging of the effect of a single zinc (Zn) impurity atom showed the strongest intensity of the LDOS at resonance directly on the Zn atom, which is opposite to the earlier theories . In our previous work , we explained the novel pattern of the tunneling conductance in terms of the blocking effect of the BiO and SrO layers which exist between the tunneling tip and the CuO2 layer being probed. This explanation may not exclude other possibilities and needs further experimental justification . Although almost all existing theories are treating the Zn atom as a static potential scatterer, there exists a few experiments  indicating that Zn impurity may possess a magnetic moment. It is now an appropriate time to consider the effects of dynamic impurities on the quasiparticle properties of $`d`$-wave superconductors. The study of a dynamic impurity in a superconductor began with the works of Yu  and Shiba  on magnetic impurities in $`s`$-wave systems. The solving of this problem is quite challenging mainly due to the difficulty in treating the dynamical correlations of the coupled impurity-conduction electron system with pair correlation. The problem of Kondo impurities in a $`d`$-wave system has been studied by focusing on either the moment formation  or the suppression of superconducting transition temperature by them . Quite recently, within a Kondo spin model, Tsuchiura et al.  have studied the LDOS around a magnetic impurity by treating the impurity spin as an effective local magnetic field and the screening of the impurity spin by the quasiparticles was not considered there. However, other new theories  argued that the moment may be screened deep in a $`d`$-wave superconducting state. Experimentally, the specific heat measurement of Zn-doped YBa2Cu3O6.95  provided a possible evidence of Kondo screening in the superconducting state.

In this paper we study the Kondo effects of an Anderson impurity coupled with the conduction electrons of a $`d`$-wave superconductor. Within a slave-boson mean-field approach, we are able to exactly solve the composite system with the self-consistency for the $`d`$-wave OP included. It is shown that a critical coupling strength, above which the Kondo effect emerges, exists regardless of the band particle-hole (PH) symmetry being present or not. We also show in the Kondo regime that a double resonant-peak structure shows up in the LDOS both at the impurity and at its nearest neighboring sites. This result is different from those for the static impurities that the LDOS at the impurity site exhibits only a single peak below (or above) the Fermi energy for a weak repulsive (or attractive) scattering center while vanishes in the unitary scattering limit.

To describe the impurity and its interaction with the band electrons of the $`d`$-wave superconductor defined on a two-dimensional lattice, we start with an Anderson-like Hamiltonian, which consists of three parts: $`H=H_S+H_{imp}+H_{S\text{-}imp}`$. Here $`H_S=_{\mathrm{𝐢𝐣},\sigma }(t_{\mathrm{𝐢𝐣}}+\mu \delta _{\mathrm{𝐢𝐣}})c_{𝐢\sigma }^{}c_{𝐣\sigma }+_{\mathrm{𝐢𝐣}}(\mathrm{\Delta }_{\mathrm{𝐢𝐣}}c_𝐢^{}c_𝐣^{}+\mathrm{\Delta }_{\mathrm{𝐢𝐣}}^{}c_𝐣c_𝐢)+_{𝐢𝐣}\frac{|\mathrm{\Delta }_{\mathrm{𝐢𝐣}}|^2}{g_{\mathrm{𝐢𝐣}}}`$, is the Hamiltonian for the host $`d`$-wave superconductor, where $`c_{𝐢\sigma }^{}`$ ($`c_{𝐢\sigma }`$) are the creation (annihilation) operators of a conduction electron with spin $`\sigma `$ at the $`𝐢`$th site, $`t_{\mathrm{𝐢𝐣}}`$ is the hopping integral between two nearest neighbor sites, $`\mu `$ is the chemical potential, $`g_{\mathrm{𝐢𝐣}}`$ represents the strength of $`d`$-wave pairing interaction, and $`\mathrm{\Delta }_{\mathrm{𝐢𝐣}}`$ is the bond pairing amplitude which is in turn determined self-consistently as $`\mathrm{\Delta }_{\mathrm{𝐢𝐣}}=g_{\mathrm{𝐢𝐣}}c_𝐢c_𝐣`$; $`H_{imp}=_\sigma ϵ_dd_\sigma ^{}d_\sigma +Un_dn_d`$ is the Hamiltonian for the Anderson impurity, where $`d_\sigma ^{}`$ ($`d_\sigma `$) is the creation (annihilation) operator of the impurity electron with spin $`\sigma `$, $`n_{d\sigma }=d_\sigma ^{}d_\sigma `$, $`ϵ_d`$ is the bare impurity energy level with respect to the chemical potential of band electrons, and $`U`$ is the on-site repulsive interaction of the impurity electrons; and $`H_{S\text{-}imp}=𝒱_\sigma (c_{𝐢=\mathrm{𝟎}+𝜹,\sigma }^{}d_\sigma +d_\sigma ^{}c_{𝐢=\mathrm{𝟎}+𝜹,\sigma })`$ represents the hybridization of the impurity with the conduction electrons, where without loss of generality, we have chosen the impurity to be located at the origin $`𝐢=\mathrm{𝟎}`$, $`𝜹`$ are the unit vectors of the square lattice, and $`𝒱`$ is the coupling strength. For our purpose of studying the low temperature physics, we assume that the on-site interaction $`U`$ is infinite so that the double occupancy on the impurity is forbidden, which allows us to apply the slave-boson mean-field theory . This approach has been applied successfully to the low-temperature properties of a Kondo impurity in the presence of normal state conduction electrons , and to transport properties of a quantum dot connected to a normal and to a (conventional) superconducting lead , or embedded into a (conventional) Josephson junction . By writing the impurity electron operator as $`d_\sigma =b^{}f_\sigma `$, where $`f_\sigma `$ and $`b`$ are the spinon and holon operators and subject to the single occupancy constraint $`_\sigma f_\sigma ^{}f_\sigma +b^{}b=1`$, in the mean-filed approximation, we are able to obtain $`H_{imp}=_\sigma (ϵ_d+\lambda _0)f_\sigma ^{}f_\sigma +\lambda _0(b_0^21)`$, and $`H_{S\text{-}imp}=𝒱b_0_\sigma (c_{𝐢=\mathrm{𝟎}+𝜹,\sigma }^{}f_\sigma +f_\sigma ^{}c_{𝐢=\mathrm{𝟎}+𝜹,\sigma })`$. Here the holon operator $`b`$ has been replaced by a $`c`$-number $`b_0`$, $`\lambda _0`$ is the Lagrange multiplier introduced to enforce the single occupancy constraint. Then the composite Hamiltonian can be diagonalized by performing the canonical transformation for both the conduction and the impurity electrons, which leads to the Bogoliubov-de Gennes equations:

$$\underset{𝐣}{}\left(\begin{array}{cc}H_{\mathrm{𝐢𝐣}}& \mathrm{\Delta }_{\mathrm{𝐢𝐣}}\\ \mathrm{\Delta }_{\mathrm{𝐢𝐣}}^{}& H_{\mathrm{𝐢𝐣}}\end{array}\right)\left(\begin{array}{c}u_𝐣^n\\ v_𝐣^n\end{array}\right)=E_n\left(\begin{array}{c}u_𝐢^n\\ v_𝐢^n\end{array}\right),$$
(1)
where $`(u_𝐢^n,v_𝐢^n)^{\text{Transpose}}`$ is the quasiparticle wavefunction with eigenvalue $`E_n`$,

$$H_{\mathrm{𝐢𝐣}}=\stackrel{~}{t}_{\mathrm{𝐢𝐣}}\mu \delta _{𝐢\mathrm{𝟎},𝐣}+(ϵ_d+\lambda _0)\delta _{\mathrm{𝐢𝟎}}\delta _{\mathrm{𝐢𝐣}},$$
(2)
with $`\stackrel{~}{t}_{\mathrm{𝐢𝐣}}=𝒱b_0`$ for $`(\mathrm{𝐢𝐣})=(𝜹\mathrm{𝟎})`$ or $`(\mathrm{𝟎}𝜹)`$ and $`t`$ otherwise, $`\mathrm{\Delta }_{\mathrm{𝐢𝐣}}`$ and $`b_0`$ are subject to the self-consistent condition,

$$\mathrm{\Delta }_{\mathrm{𝐢𝐣}}=\frac{g_{\mathrm{𝐢𝐣}}}{2}\underset{n}{}(u_𝐢^nv_𝐣^n+u_𝐣^nv_𝐢^n)\mathrm{tanh}(E_n/2k_BT),$$
(3)
and

$$b_0^2=12\underset{n}{}\{|u_\mathrm{𝟎}^n|^2f(E_n)+|v_\mathrm{𝟎}^n|^2[1f(E_n)]\},$$
(4)
with $`f(E)=[1+\mathrm{exp}(E/k_BT)]^1`$ the Fermi distribution function. Here we have assumed the conduction and the impurity electrons are of different nature and take the pairing interaction $`g_{\mathrm{𝐢𝐣}}=0`$ between them, and otherwise $`g_{\mathrm{𝐢𝐣}}=g`$ to be a constant. Notice that the eigenvalues and the corresponding eigenfunction amplitudes in Eq. (4) is implicitly $`\lambda _0`$-dependent. The single occupancy constraint, as given by Eq. (4), is not sufficient to determine the values of $`b_0`$ and $`\lambda _0`$. For this purpose, we need to minimize the free energy, which is derived for the first time in terms of the quasiparticle wavefunction as:

$`F`$ $`=`$ $`2k_BT{\displaystyle \underset{n}{}}\mathrm{ln}[1f(E_n)]{\displaystyle \underset{n}{}}{\displaystyle \underset{𝐢,𝐣}{}}[\mathrm{\Delta }_{\mathrm{𝐢𝐣}}u_𝐣^nv_𝐢^n+\text{c.c.}]`$ (7)
$`+{\displaystyle \underset{n}{}}{\displaystyle \underset{𝐢,𝐣}{}}\{[\stackrel{~}{t}_{\mathrm{𝐣𝐢}}+(ϵ_d+\lambda _0)\delta _{\mathrm{𝐢𝟎}}\delta _{\mathrm{𝐢𝐣}}\mu \delta _{𝐢\mathrm{𝟎},𝐣}]v_𝐣^nv_𝐢^n+\text{c.c.}\}`$
$`+{\displaystyle \underset{\mathrm{𝐢𝐣}}{}}{\displaystyle \frac{|\mathrm{\Delta }_{\mathrm{𝐢𝐣}}|^2}{g_{\mathrm{𝐢𝐣}}}}+\lambda _0(b_0^21).`$
This minimization will give two equations for $`\lambda _0`$ and $`b_0`$, one of which is just the single occupancy constraint Eq. (4). When applied to other physical problems, this expression for the free energy is very useful in determining the stability of competing pairing states.

Equation (1) always has a trivial solution with $`\lambda _0=ϵ_d`$ and $`b_0=0`$, which describes the impurity decoupled with the conduction electrons. This trivial critical point can be shown not necessarily to have the lowest free energy. We solve Eq. (1) by taking $`\lambda _0`$ as a fixed parameter, and calculate the free energy with the obtained value of $`b_0`$ and $`\mathrm{\Delta }_{\mathrm{𝐢𝐣}}`$. By varying the value of $`\lambda _0`$ and repeating the same procedure, we obtain the free energy as a function of $`\lambda _0`$. Comparing the minimum free energy with that of $`\lambda _0=ϵ_d`$ and $`b_0=0`$, the true solution for $`\lambda _0`$ and $`b_0`$ is then determined. The numerical calculation is performed on a square lattice with size $`N_L=N_x\times N_y=35\times 35`$ and averaged over 36 wavevectors in the supercell Brillouin zone. Throughout the work, we are interested in the zero temperature limit, and take the pairing interaction $`g=2t`$, the bare impurity level $`ϵ_d=2t`$. With respect to the real materials, this value of $`g`$ leads to an exaggerated amplitude of the energy gap. Our choice of the parameters is motivated by the desire to have a short coherence length. The calculations with smaller values of energy gap have shown similar results and no qualitative changes appear. In Fig. 1, we plot $`b_0`$ and the free energy per site as a function of $`\lambda _0`$ for various values of the coupling strength $`𝒱`$. Also shown is the value of the

free energy $`F_0`$ obtained with $`\lambda _0=ϵ_d`$ and $`b_0=0`$, which is independent of $`𝒱`$. We find a critical value $`V_c0.6t`$ for both $`\mu =0.5`$ and $`0`$, below which the minimum free energy is larger than $`F_0`$. In this regime, the impurity is decoupled with conduction band, and the electron on the impurity site is described by the spinon only. When $`𝒱>𝒱_c`$, this free local spin state is unstable against the coupling to the conduction electrons. We find: For the coupling strength $`𝒱=1.0t`$, $`0.8t`$, and $`0.6t`$, respectively, $`(\lambda _0,b_0)=(2.3t,0.56)`$, $`(2.1t,0.50)`$, and $`(2t,0)`$ for $`\mu =0.5t`$, which breaks the band PH symmetry; while $`(\lambda _0,b_0)=(2.4t,0.59)`$, $`(2.1t,0.42)`$, and $`(2.0t,0)`$ for $`\mu =0`$, which conserves the band PH symmetry. Our numerical result is consistent with the conclusion by Withoff and Fradkin  based on the Kondo exchange model that if the DOS of the conduction band vanishes linearly at the Fermi energy, the Kondo effect is suppressed for small values of the coupling constant. However, for the special case where the band PH symmetry is present, the critical coupling strength is still found to be finite. At a first glance, our result seems to differ from the earlier prediction based on the Kondo spin model by Ingersent  that for the PH symmetry case the Kondo effect is absent for all values of coupling. The discrepancy with Ref.  may come from the facts that: (i) in the present model, the local PH symmetry is still broken due to the asymmetric Anderson impurity and/or the single-occupancy constraint; (ii) the existence of the Kondo effect is sensitive to the local character of conduction electrons near the impurity; and the magnetic property of the impurity and the quasiparticle properties of superconductors interplay with each other, that is, the presence of the impurity enhances the LDOS of the $`d`$-wave superconductor (as shown below), which in turn screens the impurity spin.

Once the parameter values of $`\lambda _0`$ and $`b_0`$, and the OP are obtained, the local density of states is then calculated according to

$$\rho _𝐢=2\underset{n}{}[|u_𝐢^n|^2\delta (E_nE)+|v_i^n|^2\delta (E_n+E)],$$
(8)
which is proportional to the local differential tunneling conductance as measured by the STM. Figure 2 shows the LDOS as a function of energy at the impurity site and its nearest neighbor with $`\mu =0.5t`$. Also shown in Fig. 2(b) is the bulk density of states in the absence of the impurity. In the bulk DOS, the superconducting coherent peaks are exhibited at the gap edge $`\pm \mathrm{\Delta }_g(=0.75t)`$. The peaks outside the gap comes from the Van Hove singularity. Due to the finite magnitude of the PH symmetry breaking by $`\mu 0`$, the spectrum both inside and outside the gap is asymmetric with respect to $`E=0`$. More interestingly, within the energy gap, the resonant peaks show up in the LDOS both at the impurity site and its nearest neighbors, in sharp contrast to the unitary nonmagnetic impurity case where the LDOS vanishes at the impurity site . The appearance of the double-peak structure in the LDOS directly at the impurity site also makes the result different from the weak nonmagnetic impurity case where only one peak appears below or above the Fermi energy ($`E=0`$) when the scattering is repulsive or attractive . The difference is mainly due to the fact that, in addition to a weak nonmagnetic-like impurity scattering in the Kondo resonance phase regime, there is also a disruption of the effective hopping as given by $`𝒱b_0`$ (from the point of view of the conduction electrons), immediately connecting the impurity site. The decrease of $`𝒱`$ enhances the quasiparticle scattering near the impurity by the suppression of $`𝒱b_0`$ on one hand and shifts the impurity level to $`\stackrel{~}{ϵ}_d=ϵ_d+\lambda _0`$ to the Fermi surface on

the other hand. As a consequence of this decrease: (i) At the impurity site, the Kondo resonance peaks are sharpened with their separation decreased, and finally the double peaks merge into a single peak as $`𝒱`$ approaches $`𝒱_c`$ at which the decoupling of the impurity from the conduction band takes place; (ii) At the nearest neighbor to the impurity site, the separation between double peaks is also decreased but the double-peak structure still survives when $`𝒱=𝒱_c`$, which is intrinsic to a $`d`$-wave superconductor with a short coherence when the PH symmetry is broken to some extent . When the band PH symmetry is conserved ($`\mu =0`$), as shown in Fig. 3, we see that most features of the local electronic structure are similar to the case with $`\mu 0`$ except that the spectrum is now symmetric outside the gap (but still asymmetric inside the gap due to the breaking of the local PH symmetry) and that at the nearest neighbor to the impurity only exhibits a single zero-energy peak when $`𝒱=𝒱_c`$. We have also studied the case with different values of impurity level $`ϵ_d`$ and found similar results.

Two remarks are in order: (1) Although the slave-boson mean field theory captures the essential physics in the Kondo resonance regime, where a spin-singlet state is formed, unfortunately, it is unable to identify the quantum transition from the spin-singlet state (screened impurity spin) to the spin-doublet state (unscreened impurity spin). The slave-boson mean filed solution always leads to a lower energy in the spin-singlet state than in the spin-doublet state . One must go beyond the static mean-field level and employ alternatively the non-crossing approximation to describe this transition . (2) Very recently, Polkovnikov et al.  have studied a similar problem. The difference and similarity between two works are as follows. The difference: (i) In their work, the Kondo exchange model where the impurity spin is coupled to the spin of conduction electrons. These two models should be equivalent, but the methods to solve them are different; (ii) In addition to the first nearest-neighbor hopping ($`t`$), second ($`t^{}`$) and the third ($`t^{\prime \prime }`$) nearest-neighbor hopping were also introduced in their model. Instead only the first nearest-neighbor hopping is used in the present work; (iii) In their work, the modification of the superconductivity by the Kondo impurity was neglected, while in the present work, this effect has been taken into account by our self-consistent calculation; (iv) A vanishing density of states near zero bias at the sites nearest neighboring to the impurity was obtained in their work, while in the present work, the intensity of density of states on these sites was found to be significant. The similarity: (i) Both works have shown the split peaks for the local density of states directly on the impurity site; (ii) The models in both papers are placing the magnetic moment at the impurity site. This is different from the experimental situation for high-Tc superconductors with Zn impurities where the magnetic spin would most likely be induced at sites nearest-neighboring to the Zn impurity. Therefore, these two models might be more suitable to describe the effects of Kondo impurities such as Mn, Fe, Co substituted for Cu in high-$`T_c`$ cuprates.

In summary, we have considered the problem of an Anderson impurity coupled to the conduction electrons of a $`d`$-wave superconductor. Within a slave-boson mean-field approach, we have been able to exactly diagonalize the composite Hamiltonian with the pair breaking effect taken into account. A critical coupling strength, above which the Kondo phase occurs, is found regardless of whether the PH symmetry is present. In the Kondo effect regime, the double-peak structure appears at both the impurity site and its nearest neighboring sites. Our analysis seems to indicate that the spatial pattern of the differential tunneling conductance around a Zn impurity in BSCCO  cannot be fully explained with the Kondo impurity model.

Acknowledgments: We thank G. Kotliar, S. Sachdev, and P. Schwab for useful discussions. This work was supported by the Texas Center for Superconductivity at the University of Houston and the Robert A. Welch Foundation.


text:
# Some low-temperature properties of a generalized Hubbard model with correlated hopping

## Abstract

In the present paper we study some correlation effects in a generalized Hubbard model with correlated hopping within low-temperature region using a generalized mean-field approximation. It is shown that in a series of cases the model leads to consequences deviating essentially from those of the Hubbard model. We consider the possibility of applying the result to interpret the peculiarities of physical properties of systems with narrow energy bands.

, ,

thanks: Corresponding author. E-mail: vaha@tu.edu.te.ua

Theoretical analyses, on the one hand, and available experimental data, on the other hand, point out the necessity of the Hubbard model generalization by taking into account correlated hopping . In the recent few years, similar models have been studied intensively (e.g., see Refs.  and references therein). In particular, some of these models  have been solved exactly in a wide range of parameters.

The present paper is devoted to study of correlation effects in a generalized Hubbard model with correlated hopping within low-temperature region in the case of partially filled narrow energy bands.

We start from the generalization of the Hubbard model with correlated hoppings $`T_1`$ and $`T_2`$ :

$`H=`$ $`\mu {\displaystyle \underset{i\sigma }{}}a_{i\sigma }^+a_{i\sigma }+t(n){\displaystyle \underset{ij\sigma }{}^{}}a_{i\sigma }^+a_{j\sigma }+`$
$`T_2{\displaystyle \underset{ij\sigma }{}^{}}(a_{i\sigma }^+a_{j\sigma }n_{i\overline{\sigma }}+h.c.)+U{\displaystyle \underset{i}{}}n_in_i,`$
where all the notations are usual, $`t(n)=t_0+nT_1`$, $`t_0`$ is the nearest-neighbor hopping integral, $`n`$ is the electron concentration, the prime at the sums signifies that $`ij`$.

The peculiarity of model (Some low-temperature properties of a generalized Hubbard model with correlated hopping) is taking into account the correlated hopping $`T_1`$ leading to the concentration dependence of the hopping integral $`t(n)`$ in contrast to similar models which are considered in Refs. .

In Ref.  we have obtained the single-particle Green function and the quasiparticle energy spectrum (which has the exact atomic and band limits) for arbitrary electron concentration using a variant  of the generalized mean-field approximation . For the case $`n=1`$ and $`t(n)=T_2`$ the approach recovers (see Ref. ) the exact results .

The energy gap $`\mathrm{\Delta }E`$ (which is calculated from the spectrum) is found to increase with a deviation from half-filling ($`\mathrm{\Delta }E`$ has a minimum at $`n=1`$). This explains the fact that in a metallic phase the vanadium oxides with fractional number of electrons per cation (VkO2k-1, $`k3`$) exhibit “metallicity” to lesser extent than the compounds VO2 and V2O3 (integer number of electrons per site) .

Within low-temperature region for $`U\mathrm{}`$ width of the lower Hubbard band at $`n<1`$ is $`W_l=2w\left[2/(2n)n\right]`$; and at $`n>1`$ width of the upper Hubbard band is $`W_u=2\stackrel{~}{w}(n2+2/n)`$, where $`w=z|t(n)|,\stackrel{~}{w}=z|\stackrel{~}{t}(n)|,\stackrel{~}{t}(n)=t(n)+2T_2,z`$ is the number of nearest neighbors to a site.

The concentration dependence of subband widths plotted in Fig. 1 is caused firstly, by correlation effect of narrowing of subbands, secondly, by the concentration dependence of the hopping integrals in the lower $`t(n)`$ and upper $`\stackrel{~}{t}(n)`$ Hubbard bands (this is a peculiarity of the present model). Taking into account correlated hopping leads to the essential narrowing of subband, and what is more, narrowing by this factor increases with increasing of electron concentration. Fig. 1 illustrates the non-equivalence of the cases $`n<1`$ and $`n>1`$ in the framework of the present model (in contrast to the electron-hole symmetry of the Hubbard model). One can see that at some values of the parameters $`\tau _1=T_1/|t_0|,\tau _2=T_2/|t_0|`$ correlated hopping $`W_lW_u`$, thus effective mass of current carriers within the lower Hubbard band $`m_lm_u`$ ($`m_u`$ is the effective mass of current carriers within the upper Hubbard band). So, the narrow-band materials with conductivity determinated by current carriers within the lower band can have much smaller resistivity than the compounds with conductivity determinated by current carriers within the upper band. This fact is confirmed by the experimental data : the metalooxides with less than half-filled $`3d`$-shell (Mn2O) exhibit much higher conductivity than the compounds with half or more than half-filled $`3d`$-shell (MnO, NiO).

V.H. thanks the Organizing Committee of the LT22 Conference for a grant. The authors are grateful to Prof. D. Khomskii for valuable discussions on the part of the result considered in the present work.