text:
# The Projector Augmented Wave Method: ab-initio molecular dynamics with full wave functions

## 1 Introduction

The main goal of electronic structure methods is to solve the Schrödinger equation for the electrons in a molecule or solid, to evaluate the resulting total energies, forces, response functions and other quantities of interest. In this paper we review the Projector augmented wave (PAW) method \[Blöchl 1994\], an electronic structure method for ab-initio molecular dynamics with full wave functions. The main goal of this paper is not to provide a particularly complete or detailed account of the methodology, but rather to lay out the underlying ideas. A more rigorous description can be found in the original paper \[Blöchl 1994\].

Density functional theory \[Hohenberg et al 1964, Kohn et al 1965\] maps a description for interacting electrons onto one of non-interacting electrons in an effective potential. The remaining one-electron Schrödinger equation still poses substantial numerical difficulties: (1) in the atomic region near the nucleus, the kinetic energy of the electrons is large, resulting in rapid oscillations of the wave function that require fine grids for an accurate numerical representation. On the other hand, the large kinetic energy makes the Schrödinger equation stiff, so that a change of the chemical environment has little effect on the shape of the wave function. Therefore, the wave function in the atomic region can be represented well already by a small basis set. (2) In the bonding region between the atoms the situation is opposite. The kinetic energy is small and the wave function is smooth. However, the wave function is flexible and responds strongly to the environment. This requires large and nearly complete basis sets.

Combining these different requirements is non-trivial and various strategies have been developed.

* Most appealing to quantum chemists has been the atomic point of view. Basis functions are chosen that resemble atomic orbitals. They exploit that the wave function in the atomic region can be described by a few basis functions, while the bonding is described by the overlapping tails of these atomic orbitals. Most techniques in this class are a compromise of a well adapted basis set with complex matrix elements on the one hand and on the other hand numerically convenient basis functions such as Gaussians, where the inadequacies are compensated by larger basis sets.
* Pseudopotentials regard an atom as a perturbation of the free electron gas. The most natural basis functions are plane waves. Plane waves are complete and well adapted to sufficiently smooth wave functions. The disadvantage of the large basis sets required is offset by the extreme simplicity to evaluate matrix elements. Finite plane wave expansions are, however, absolutely inadequate to describe the strong oscillations of the wave functions near the nucleus. In the pseudopotential approach the Pauli repulsion of the core electrons is therefore described by an effective potential that expels the valence electrons from the core region. The resulting wave functions are smooth and can be represented well by plane waves. The price to pay is that all information on the charge density and wave functions near the nucleus is lost.
* Augmented wave methods compose their basis functions out of atom-like partial waves in the atomic regions and a set of functions, called envelope functions, appropriate for the bonding in between. Space is divided accordingly into atom-centered spheres, defining the atomic region, and an interstitial region for the bonds. The partial solutions of the different regions, are matched at the interface between atomic and interstitial regions.

The projector augmented wave method is an extension of augmented wave methods and the pseudopotential approach, which combines their traditions into a unified electronic structure method.

After describing the underlying ideas of the various methods let us briefly review the history of augmented wave methods and the pseudopotential approach. We do not discuss the atomic-orbital based methods, because our focus is the PAW method and its ancestors.

The augmented wave methods have been introduced in 1937 by Slater \[Slater 1937\] and later modified by Korringa \[Korringa 1947\], Kohn and Rostokker \[Kohn et al 1954\]. They approached the electronic structure as a scattered electron problem. Consider an electron beam, represented by a plane wave, traveling through a solid. It undergoes multiple scattering at the atoms. If for some energy, the outgoing scattered waves interfere destructively, a bound state has been determined. This approach can be translated into a basis set method with energy dependent and potential dependent basis functions. In order to make the scattered wave problem tractable, a model potential had to be chosen: The so-called muffin-tin potential approximates the potential by a constant in the interstitial region and by a spherically symmetric potential in the atomic region.

The pseudopotential approach traces back to 1940 when C. Herring invented the orthogonalized plane wave method \[Herring 1940\]. Later, Phillips \[Phillips et al 1959\] and Antoncik \[Antoncik 1959\] replaced the orthogonality condition by an effective potential, that compensates the electrostatic attraction by the nucleus. In practice, the potential was modified for example by cutting off the singular potential of the nucleus at a certain value. This was done with a few parameters that have been adjusted to reproduce the measured electronic band structure of the corresponding solid.

Augmented wave and pseudopotential methods reached adulthood in the 1970s: At first O.K. Andersen \[Andersen 1975\] showed that the energy dependent basis set of Slater’s APW method can be mapped onto one with energy independent basis functions, by linearizing the partial waves for the atomic regions in energy. In the original APW approach the zeros of an energy dependent matrix had to be determined, which is problematic, if many states lie in a small energy region as for complex systems. With the new energy independent basis functions, however, the problem is reduced to the much simpler generalized eigenvalue problem, which can be solved using efficient numerical techniques. Furthermore, the introduction of well defined basis sets paved the way for full-potential calculations. In that case the muffin-tin approximation is used solely to define the basis set. The matrix elements of the Hamiltonian are evaluated with the full potential.

Hamann, Schlüter and Chiang \[Hamann et al 1979\] showed in 1979 how pseudopotentials can be constructed in such a way, that their scattering properties are identical to that of an atom to first order in energy. These first-principles pseudopotentials relieved the calculations from the restrictions of empirical parameters. Highly accurate calculations have become possible. A main disadvantage of these pseudopotentials has been the large basis set size required especially for first-row and transition metal atoms.

In 1985 R. Car and M. Parrinello published the ab-initio molecular dynamics method \[Car et al 1985\]. Simulations of the atomic motion have become possible on the basis of state-of-the-art electronic structure methods. Besides making dynamical phenomena and finite temperature effects accessible to electronic structure calculations, the ab-initio molecular dynamics method also introduced a radically new way of thinking into electronic structure methods. Diagonalization of a Hamilton matrix has been replaced by classical equations of motion for the wave function coefficients. If one applies friction, the system is quenched to the ground state. Without friction truly dynamical simulations of the atomic structure are performed. Electronic wave functions and atomic positions are treated on equal footing.

The Car-Parrinello method had been implemented first for the pseudopotential approach. There seemed to be unsurmountable barriers against combining the new technique with augmented wave methods. The main problem was related to the potential dependent basis set used sofar: the Car-Parrinello method requires a well defined and unique total energy functional of atomic positions and basis set coefficients. Therefore, it was one of the main goals of the PAW method to introduce energy and potential independent basis sets that were as accurate and numerically efficient as the previously used augmented basis sets. Other requirements have been: (1) The method should match the efficiency of the pseudopotential approach for Car-Parrinello simulations. (2) It should become an exact theory when converged and (3) its convergence should be easily controlled. We believe that these criteria have been met, which explains why the PAW method becomes increasingly wide spread today.

We would like to point out that most of these seemingly singular developments did not come out of the blue, but the ideas seemed to have evolved in the community. In the case of the PAW method, similar ideas have been developed by Vanderbilt \[Vanderbilt 1990\] in the context of ultra-soft pseudopotentials. The first dynamical simulations using a semiempirical electronic structure method have been performed by Wang and Karplus \[Wang et al 1973\] in 1973. The first ab-initio pseudopotentials have been published by Zunger \[Zunger et al 1978\] one year before Hamann, Bachelet and Schlüter \[Hamann et al 1979\].

## 2 Transformation theory

At the root of the PAW method lies a transformation, that maps the true wave functions with their complete nodal structure onto auxiliary wave functions, that are numerically convenient. We aim for smooth auxiliary wave functions, which have a rapidly convergent plane wave expansion. With such a transformation we can expand the auxiliary wave functions into a convenient basis set, and evaluate all physical properties after reconstructing the related physical (true) wave functions.

Let us denote the physical one-particle wave functions as $`|\mathrm{\Psi }_n`$ and the auxiliary wave functions as $`|\stackrel{~}{\mathrm{\Psi }}_n`$. Note that the tilde refers to the representation of smooth auxiliary wave functions. $`n`$ is the label for a one-particle state and contains a band index, a $`k`$-point and a spin index. The transformation from the auxiliary to the physical wave functions is $`𝒯`$.

$`|\mathrm{\Psi }_n=𝒯|\stackrel{~}{\mathrm{\Psi }}_n`$ (1)
We use here Dirac’s Bra and Ket notation. A wave function $`\mathrm{\Psi }_n(𝐫)`$ corresponds to a ket $`|\mathrm{\Psi }_n`$, the complex conjugate wave function $`\mathrm{\Psi }_n^{}(𝐫)`$ corresponds to a bra $`\mathrm{\Psi }_n|`$, and a scalar product $`d^3r\mathrm{\Psi }_n^{}(𝐫)\mathrm{\Psi }_m(𝐫)`$ is written as $`\mathrm{\Psi }_n|\mathrm{\Psi }_m`$. Vectors in the 3-d coordinate space are indicated by boldfaced symbols.

The electronic ground state is determined by minimizing a total energy functional $`E[\mathrm{\Psi }_n]`$ of the density functional theory. The one-particle wave functions have to be orthogonal. This constraint is implemented with the method of Lagrange multipliers. We obtain the ground state wave functions from the extremum condition for

$`F([\mathrm{\Psi }_n],\mathrm{\Lambda }_{m,n})=E[\mathrm{\Psi }_n]{\displaystyle \underset{n,m}{}}[\mathrm{\Psi }_n|\mathrm{\Psi }_m\delta _{n,m}]\mathrm{\Lambda }_{n,m}`$ (2)
with respect to the wave functions and the Lagrange multipliers $`\mathrm{\Lambda }_{n,m}`$. The extremum condition for the wave functions has the form

$`H|\mathrm{\Psi }_nf_n={\displaystyle \underset{m}{}}|\mathrm{\Psi }_m\mathrm{\Lambda }_{m,n}`$ (3)
where the $`f_n`$ are the occupation numbers and $`H=\frac{\mathrm{}^2}{2m_e}^2+v_{\mathrm{eff}}(𝐫)`$ is the effective one-particle Hamilton operator.

After a unitary transformation that diagonalizes the matrix of Lagrange multipliers $`\mathrm{\Lambda }_{m,n}`$, we obtain the Kohn-Sham equations.

$`H|\mathrm{\Psi }_n=|\mathrm{\Psi }_nϵ_n`$ (4)
The one-particle energies $`ϵ_n`$ are the eigenvalues of $`\mathrm{\Lambda }_{n,m}\frac{f_n+f_m}{2f_nf_m}`$.

Now we express the functional $`F`$ in terms of our auxiliary wave functions

$`F([𝒯\stackrel{~}{\mathrm{\Psi }}_n],\mathrm{\Lambda }_{m,n})=E[𝒯\stackrel{~}{\mathrm{\Psi }}_n]{\displaystyle \underset{n,m}{}}[\stackrel{~}{\mathrm{\Psi }}_n|𝒯^{}𝒯|\stackrel{~}{\mathrm{\Psi }}_m\delta _{n,m}]\mathrm{\Lambda }_{n,m}`$ (5)
The variational principle with respect to the auxiliary wave functions yields

$`𝒯^{}H𝒯|\stackrel{~}{\mathrm{\Psi }}_n=𝒯^{}𝒯|\stackrel{~}{\mathrm{\Psi }}_nϵ_n.`$ (6)
Again we obtain a Schrödinger-like equation, but now the Hamilton operator has a different form, $`𝒯^{}H𝒯`$, an overlap operator $`𝒯^{}𝒯`$ occurs, and the resulting auxiliary wave functions are smooth.

When we evaluate physical quantities we need to evaluate expectation values of an operator $`A`$, which can be expressed in terms of either the true or the auxiliary wave functions.

$`A`$ $`=`$ $`{\displaystyle \underset{n}{}}f_n\mathrm{\Psi }_n|A|\mathrm{\Psi }_n={\displaystyle \underset{n}{}}f_n\stackrel{~}{\mathrm{\Psi }}_n|𝒯^{}A𝒯|\stackrel{~}{\mathrm{\Psi }}_n`$ (7)
In the representation of auxiliary wave functions we need to use transformed operators $`\stackrel{~}{A}=𝒯^{}A𝒯`$. As it is, this equation only holds for the valence electrons. The core electrons are treated differently as will be shown below.

The transformation takes us conceptionally from the world of pseudopotentials to that of augmented wave methods, which deal with the full wave functions. We will see that our auxiliary wave functions, which are simply the plane wave parts of the full wave functions, translate into the wave functions of the pseudopotential approach. In the PAW method the auxiliary wave functions are used to construct the true wave functions and the total energy functional is evaluated from the latter. Thus it provides the missing link between augmented wave methods and the pseudopotential method, which can be derived as a well-defined approximation of the PAW method.

In the original paper \[Blöchl 1994\], the auxiliary wave functions have been termed pseudo wave functions and the true wave functions have been termed all-electron wave functions, in order to make the connection more evident. I avoid this notation here, because it resulted in confusion in cases, where the correspondence is not clear cut.

## 3 Transformation operator

Sofar, we have described how we can determine the auxiliary wave functions of the ground state and how to obtain physical information from them. What is missing, is a definition of the transformation operator $`𝒯`$.

The operator $`𝒯`$ has to modify the smooth auxiliary wave function in each atomic region, so that the resulting wave function has the correct nodal structure. Therefore, it makes sense to write the transformation as identity plus a sum of atomic contributions $`𝒮_R`$

$`𝒯=1+{\displaystyle \underset{R}{}}𝒮_R.`$ (8)
For every atom, $`𝒮_R`$ adds the difference between the true and the auxiliary wave function. The index $`R`$ is a label for an atomic site.

The local terms $`𝒮_R`$ are defined in terms of solutions $`|\varphi _i`$ of the Schrödinger equation for the isolated atoms. This set of partial waves $`|\varphi _i`$ will serve as a basis set so that, near the nucleus, all relevant valence wave functions can be expressed as superposition of the partial waves with yet unknown coefficients.

$`\mathrm{\Psi }(𝐫)={\displaystyle \underset{iR}{}}\varphi _i(𝐫)c_i\mathrm{for}|𝐫𝐑_R|<r_{c,R}`$ (9)
The index $`i`$ refers to a site index $`R`$, the angular momentum indices $`(\mathrm{},m)`$ and an additional index that differentiates partial waves with same angular momentum quantum numbers on the same site. With $`iR`$ we indicate those partial waves that belong to site $`R`$. $`𝐑_R`$ is the position of the nucleus of site $`R`$.

Note that the partial waves are not necessarily bound states and are therefore not normalizable, unless we truncate them beyond a certain radius $`r_{c,R}`$. The PAW method is formulated such that the final results do not depend on the location where the partial waves are truncated, as long as this is not done too close to the nucleus.

Since the core wave functions do not spread out into the neighboring atoms, we will treat them differently. Currently we use the frozen-core approximation so that density and energy of the core electrons are identical to those of the corresponding isolated atoms. The transformation $`𝒯`$ shall produce only wave functions orthogonal to the core electrons, while the core electrons are treated separately. Therefore, the set of atomic partial waves $`|\varphi _i`$ includes only valence states that are orthogonal to the core wave functions of the atom.

For each of the partial waves we choose an auxiliary partial wave $`|\stackrel{~}{\varphi }_i`$. The identity

$`|\varphi _i`$ $`=`$ $`(1+S_R)|\stackrel{~}{\varphi }_i\mathrm{for}iR`$
$`S_R|\stackrel{~}{\varphi }_i`$ $`=`$ $`|\varphi _i|\stackrel{~}{\varphi }_i`$ (10)
defines the local contribution $`𝒮_R`$ to the transformation operator. Since $`1+𝒮_R`$ shall change the wave function only locally, we require that the partial waves $`|\varphi _i`$ and their auxiliary counter parts $`|\stackrel{~}{\varphi }_i`$ are pairwise identical beyond a certain radius $`r_c`$.

$`\varphi _i(r)=\stackrel{~}{\varphi }_i(r)\mathrm{for}iR\mathrm{and}|𝐫𝐑_R|>r_{c,R}`$ (11)
In order to be able to apply the transformation operator to an arbitrary auxiliary wave function, we need to be able to expand the auxiliary wave function locally into the auxiliary partial waves.

$`\stackrel{~}{\mathrm{\Psi }}(𝐫)`$ $`=`$ $`{\displaystyle \underset{iR}{}}\stackrel{~}{\varphi }_i(𝐫)\stackrel{~}{p}_i|\stackrel{~}{\mathrm{\Psi }}\mathrm{for}|𝐫𝐑_R|<r_{c,R}`$ (12)
which defines the projector functions $`|\stackrel{~}{p}_i`$. The projector functions probe the local character of the auxiliary wave function in the atomic region. Examples of projector functions are shown in Fig. 1. From Eq. 12 we can derive $`_i|\stackrel{~}{\varphi }_i\stackrel{~}{p}_i|=1`$, which is valid within $`r_c`$. It can be shown by insertion, that the identity Eq. 12 holds for any auxiliary wave function $`|\stackrel{~}{\mathrm{\Psi }}`$ that can be expanded locally into auxiliary partial waves $`|\stackrel{~}{\varphi }_i`$, if

$$\stackrel{~}{p}_i|\stackrel{~}{\varphi }_j=\delta _{i,j}\mathrm{for}i,jR$$
(13)
Note that neither the projector functions nor the partial waves need to be orthogonal among themselves.

By combining Eq. 10 and Eq. 12, we can apply $`S_R`$ to any auxiliary wave function.

$`S_R|\stackrel{~}{\mathrm{\Psi }}`$ $`=`$ $`{\displaystyle \underset{iR}{}}S_R|\stackrel{~}{\varphi }_i\stackrel{~}{p}_i|\stackrel{~}{\mathrm{\Psi }}={\displaystyle \underset{iR}{}}\left(|\varphi _i|\stackrel{~}{\varphi }_i\right)\stackrel{~}{p}_i|\stackrel{~}{\mathrm{\Psi }}`$ (14)
Hence the transformation operator is

$`𝒯=1+{\displaystyle \underset{i}{}}\left(|\varphi _i|\stackrel{~}{\varphi }_i\right)\stackrel{~}{p}_i|`$ (15)
where the sum runs over all partial waves of all atoms. The true wave function can be expressed as

$`|\mathrm{\Psi }=|\stackrel{~}{\mathrm{\Psi }}+{\displaystyle \underset{i}{}}\left(|\varphi _i|\stackrel{~}{\varphi }_i\right)\stackrel{~}{p}_i|\stackrel{~}{\mathrm{\Psi }}=|\stackrel{~}{\mathrm{\Psi }}_n+{\displaystyle \underset{R}{}}\left(|\mathrm{\Psi }_R^1|\stackrel{~}{\mathrm{\Psi }}_R^1\right)`$ (16)
with

$`|\mathrm{\Psi }_R^1`$ $`=`$ $`{\displaystyle \underset{iR}{}}|\varphi _i\stackrel{~}{p}_i|\stackrel{~}{\mathrm{\Psi }}`$ (17)
$`|\stackrel{~}{\mathrm{\Psi }}_R^1`$ $`=`$ $`{\displaystyle \underset{iR}{}}|\stackrel{~}{\varphi }_i\stackrel{~}{p}_i|\stackrel{~}{\mathrm{\Psi }}`$ (18)
In Fig. 2 the decomposition of Eq. 16 is shown for the example of the bonding p-$`\sigma `$ state of the Cl<sub>2</sub> molecule.

To understand the expression for the true wave function, Eq. 16, let us concentrate on different regions in space. (1) Far from the atoms, the partial waves are, according to Eq. 11, pairwise identical so that the auxiliary wave function is identical to the true wave function $`\mathrm{\Psi }(𝐫)=\stackrel{~}{\mathrm{\Psi }}(𝐫)`$. (2) Close to an atom, however, the true wave function $`\mathrm{\Psi }(𝐫)=\mathrm{\Psi }_R^1(𝐫)`$ is built up from partial waves that contain the proper nodal structure, because the auxiliary wave function and its partial wave expansion are equal according to Eq. 12.

In practice the partial wave expansions are truncated. Therefore, the identity of Eq. 12 does not hold strictly. As a result the plane waves also contribute to the true wave function inside the atomic region. This has the advantage that the missing terms in a truncated partial wave expansion are partly accounted for by plane waves, which explains the rapid convergence of the partial wave expansions.

Frequently, the question comes up, whether the transformation Eq. 15 of the auxiliary wave functions indeed provides the true wave function. The transformation should be considered merely as a change of representation analogous to a coordinate transform. If the total energy functional is transformed consistently, its minimum will yield an auxiliary wave function that produces a correct wave function $`|\mathrm{\Psi }`$.

## 4 Expectation values

Expectation values can be obtained either from the reconstructed true wave functions or directly from the auxiliary wave functions

$`A`$ $`=`$ $`{\displaystyle \underset{n}{}}f_n\mathrm{\Psi }_n|A|\mathrm{\Psi }_n+{\displaystyle \underset{n=1}{\overset{N_c}{}}}\varphi _n^c|A|\varphi _n^c`$ (19)
$`=`$ $`{\displaystyle \underset{n}{}}f_n\stackrel{~}{\mathrm{\Psi }}_n|𝒯^{}A𝒯|\stackrel{~}{\mathrm{\Psi }}_n+{\displaystyle \underset{n=1}{\overset{N_c}{}}}\varphi _n^c|A|\varphi _n^c`$
where $`f_n`$ are the occupations of the valence states and $`N_c`$ is the number of core states. The first sum runs over the valence states, and second over the core states $`|\varphi _n^c`$.

Now we can decompose the matrix elements into their individual contributions according to Eq. 16.

$`\mathrm{\Psi }|A|\mathrm{\Psi }`$ $`=`$ $`\stackrel{~}{\mathrm{\Psi }}+{\displaystyle \underset{R}{}}(\mathrm{\Psi }_R^1\stackrel{~}{\mathrm{\Psi }}_R^1)|A|\stackrel{~}{\mathrm{\Psi }}+{\displaystyle \underset{R^{}}{}}(\mathrm{\Psi }_R^{}^1\stackrel{~}{\mathrm{\Psi }}_R^{}^1)`$ (20)
$`=`$ $`\underset{\text{part 1}}{\underset{}{\stackrel{~}{\mathrm{\Psi }}|A|\stackrel{~}{\mathrm{\Psi }}+{\displaystyle \underset{R}{}}\left(\mathrm{\Psi }_R^1|A|\mathrm{\Psi }_R^1\stackrel{~}{\mathrm{\Psi }}_R^1|A|\stackrel{~}{\mathrm{\Psi }}_R^1\right)}}`$
$`+`$ $`\underset{\text{part 2}}{\underset{}{{\displaystyle \underset{R}{}}\left(\mathrm{\Psi }_R^1\stackrel{~}{\mathrm{\Psi }}_R^1|A|\stackrel{~}{\mathrm{\Psi }}\stackrel{~}{\mathrm{\Psi }}_R^1+\stackrel{~}{\mathrm{\Psi }}\stackrel{~}{\mathrm{\Psi }}_R^1|A|\mathrm{\Psi }_R^1\stackrel{~}{\mathrm{\Psi }}_R^1\right)}}`$
$`+`$ $`\underset{\text{part 3}}{\underset{}{{\displaystyle \underset{RR^{}}{}}\mathrm{\Psi }_R^1\stackrel{~}{\mathrm{\Psi }}_R^1|A|\mathrm{\Psi }_R^{}^1\stackrel{~}{\mathrm{\Psi }}_R^{}^1}}`$
Only the first part of Eq. 20, is evaluated explicitly, while the second and third parts of Eq. 20 are neglected, because they vanish for sufficiently local operators as long as the partial wave expansion is converged: The function $`\mathrm{\Psi }_R^1\stackrel{~}{\mathrm{\Psi }}_R^1`$ vanishes per construction beyond some augmentation region, because the partial waves are pairwise identical beyond that region. The function $`\stackrel{~}{\mathrm{\Psi }}\stackrel{~}{\mathrm{\Psi }}_R^1`$ vanishes inside the augmentation region, if the partial wave expansion is sufficiently converged. In no region of space both functions $`\mathrm{\Psi }_R^1\stackrel{~}{\mathrm{\Psi }}_R^1`$ and $`\stackrel{~}{\mathrm{\Psi }}\stackrel{~}{\mathrm{\Psi }}_R^1`$ are simultaneously nonzero. Similarly the functions $`\mathrm{\Psi }_R^1\stackrel{~}{\mathrm{\Psi }}_R^1`$ from different sites are never non-zero in the same region in space. Hence, the second and third parts of Eq. 20 vanish for operators such as the kinetic energy $`\frac{\mathrm{}^2}{2m_e}^2`$ and the real space projection operator $`|rr|`$, which produces the electron density. For truly nonlocal operators the second and third parts of Eq. 20 would have to be considered explicitly.

The expression, Eq. 20, for the expectation value can therefore be written as

$`A`$ $`=`$ $`{\displaystyle \underset{n}{}}f_n\left(\stackrel{~}{\mathrm{\Psi }}_n|A|\stackrel{~}{\mathrm{\Psi }}_n+\mathrm{\Psi }_n^1|A|\mathrm{\Psi }_n^1\stackrel{~}{\mathrm{\Psi }}_n^1|A|\stackrel{~}{\mathrm{\Psi }}_n^1\right)+{\displaystyle \underset{n=1}{\overset{N_c}{}}}\varphi _n^c|A|\varphi _n^c`$ (21)
$`=`$ $`{\displaystyle \underset{n}{}}f_n\stackrel{~}{\mathrm{\Psi }}_n|A|\stackrel{~}{\mathrm{\Psi }}_n+{\displaystyle \underset{n=1}{\overset{N_c}{}}}\stackrel{~}{\varphi }_n^c|A|\stackrel{~}{\varphi }_n^c`$
$`+`$ $`{\displaystyle \underset{R}{}}\left({\displaystyle \underset{i,jR}{}}D_{i,j}\varphi _j|A|\varphi _i+{\displaystyle \underset{nR}{\overset{N_{c,R}}{}}}\varphi _n^c|A|\varphi _n^c\right)`$
$``$ $`{\displaystyle \underset{R}{}}\left({\displaystyle \underset{i,jR}{}}D_{i,j}\stackrel{~}{\varphi }_j|A|\stackrel{~}{\varphi }_i+{\displaystyle \underset{nR}{\overset{N_{c,R}}{}}}\stackrel{~}{\varphi }_n^c|A|\stackrel{~}{\varphi }_n^c\right)`$
where $`D_{i,j}`$ is the one-center density matrix defined as

$`D_{i,j}={\displaystyle \underset{n}{}}f_n\stackrel{~}{\mathrm{\Psi }}_n|\stackrel{~}{p}_j\stackrel{~}{p}_i|\stackrel{~}{\mathrm{\Psi }}_n={\displaystyle \underset{n}{}}\stackrel{~}{p}_i|\stackrel{~}{\mathrm{\Psi }}_nf_n\stackrel{~}{\mathrm{\Psi }}_n|\stackrel{~}{p}_j`$ (22)
The auxiliary core states, $`|\stackrel{~}{\varphi }_n^c`$ allow to incorporate the tails of the core wave function into the plane wave part, and therefore assure, that the integrations of partial wave contributions cancel strictly beyond $`r_c`$. They are identical to the true core states in the tails, but are a smooth continuation inside the atomic sphere. It is not required that the auxiliary wave functions are normalized.

For example, the electron density is given by

$`n(𝐫)`$ $`=`$ $`\stackrel{~}{n}(𝐫)+{\displaystyle \underset{R}{}}\left(n_R^1(𝐫)\stackrel{~}{n}_R^1(𝐫)\right)`$ (23)
$`\stackrel{~}{n}(𝐫)`$ $`=`$ $`{\displaystyle \underset{n}{}}f_n\stackrel{~}{\mathrm{\Psi }}_n^{}(𝐫)\stackrel{~}{\mathrm{\Psi }}_n(𝐫)+\stackrel{~}{n}_c`$
$`n_R^1(𝐫)`$ $`=`$ $`{\displaystyle \underset{i,jR}{}}D_{i,j}\varphi _j^{}(𝐫)\varphi _i(𝐫)+n_{c,R}`$
$`\stackrel{~}{n}_R^1(𝐫)`$ $`=`$ $`{\displaystyle \underset{i,jR}{}}D_{i,j}\stackrel{~}{\varphi }_j^{}(𝐫)\stackrel{~}{\varphi }_i(𝐫)+\stackrel{~}{n}_{c,R}`$ (24)
where $`n_{c,R}`$ is the core density of the corresponding atom and $`\stackrel{~}{n}_{c,R}`$ is the auxiliary core density that is identical to $`n_{c,R}`$ outside the atomic region and a smooth continuation inside.

Before we continue, let us discuss a special point: The matrix element of a general operator with the auxiliary wave functions may be slowly converging with the plane wave expansion, because the operator $`A`$ may not be well behaved. An example for such an operator is the singular electrostatic potential of a nucleus. This problem can be alleviated by adding an intelligent zero: If an operator $`B`$ is purely localized within an atomic region, we can use the identity between the auxiliary wave function and its own partial wave expansion

$`0`$ $`=`$ $`\stackrel{~}{\mathrm{\Psi }}_n|B|\stackrel{~}{\mathrm{\Psi }}_n\stackrel{~}{\mathrm{\Psi }}_n^1|B|\stackrel{~}{\mathrm{\Psi }}_n^1`$ (25)
Now we choose an operator $`B`$ so that it cancels the problematic behavior of the operator $`A`$, but is localized in a single atomic region. By adding $`B`$ to the plane wave part and the matrix elements with its one-center expansions, the plane wave convergence can be improved without affecting the converged result.

## 5 Total Energy

Like wave functions and expectation values also the total energy can be divided into three parts.

$`E([\stackrel{~}{\mathrm{\Psi }}_n],R_i)`$ $`=`$ $`\stackrel{~}{E}+{\displaystyle \underset{R}{}}\left(E_R^1\stackrel{~}{E}_R^1\right)`$ (26)
The plane-wave part $`\stackrel{~}{E}`$ involves only smooth functions and is evaluated on equi-spaced grids in real and reciprocal space. This part is computationally most demanding, and is similar to the expressions in the pseudopotential approach.

$`\stackrel{~}{E}`$ $`=`$ $`{\displaystyle \underset{n}{}}\stackrel{~}{\mathrm{\Psi }}_n|{\displaystyle \frac{\mathrm{}^2}{2m_e}}^2|\stackrel{~}{\mathrm{\Psi }}_n`$ (27)
$`+`$ $`{\displaystyle \frac{e^2}{8\pi ϵ_0}}{\displaystyle d^3rd^3r^{}\frac{[\stackrel{~}{n}(𝐫)+\stackrel{~}{Z}(𝐫)][\stackrel{~}{n}(𝐫^{})+\stackrel{~}{Z}(𝐫^{})]}{|𝐫𝐫^{}|}}`$
$`+`$ $`{\displaystyle d^3r\stackrel{~}{n}(𝐫)ϵ_{xc}(𝐫,[\stackrel{~}{n}])}+{\displaystyle d^3r\overline{v}(𝐫)\stackrel{~}{n}(𝐫)},`$
where $`\stackrel{~}{Z}(𝐫)`$ is an angular dependent core-like density that will be described in detail below. The remaining parts can be evaluated on radial grids in a spherical harmonics expansion. The nodal structure of the wave functions can be properly described on a logarithmic radial grid that becomes very fine near nucleus,

$`E_R^1`$ $`=`$ $`{\displaystyle \underset{i,jR}{}}D_{i,j}\varphi _j|{\displaystyle \frac{\mathrm{}^2}{2m_e}}^2|\varphi _i+{\displaystyle \underset{nR}{\overset{N_{c,R}}{}}}\varphi _n^c|{\displaystyle \frac{\mathrm{}^2}{2m_e}}^2|\varphi _n^c`$ (28)
$`+`$ $`{\displaystyle \frac{e^2}{8\pi ϵ_0}}{\displaystyle d^3rd^3r^{}\frac{[n^1(𝐫)+Z(𝐫)][n^1(𝐫^{})+Z(𝐫^{})]}{|𝐫𝐫^{}|}}`$
$`+`$ $`{\displaystyle d^3rn^1(𝐫)ϵ_{xc}(𝐫,[n^1])}`$
$`\stackrel{~}{E}_R^1`$ $`=`$ $`{\displaystyle \underset{i,jR}{}}D_{i,j}\stackrel{~}{\varphi }_j|{\displaystyle \frac{\mathrm{}^2}{2m_e}}^2|\stackrel{~}{\varphi }_i`$ (29)
$`+`$ $`{\displaystyle \frac{e^2}{8\pi ϵ_0}}{\displaystyle d^3rd^3r^{}\frac{[\stackrel{~}{n}^1(𝐫)+\stackrel{~}{Z}(𝐫)][\stackrel{~}{n}^1(𝐫^{})+\stackrel{~}{Z}(𝐫^{})]}{|𝐫𝐫^{}|}}`$
$`+`$ $`{\displaystyle d^3r\stackrel{~}{n}^1(𝐫)ϵ_{xc}(𝐫,[\stackrel{~}{n}^1])}+{\displaystyle d^3r\overline{v}(𝐫)\stackrel{~}{n}^1(𝐫)}`$
The nuclear charge density $`eZ(𝐫)`$ is defined as a sum of $`\delta `$-functions on the nuclear sites, $`Z(𝐫)=_R𝒵_R\delta (𝐫𝐑)`$, with the atomic numbers $`𝒵_R`$. Note that the self energy of a point charge is infinite and must be subtracted out.

The compensation density $`\stackrel{~}{Z}(𝐫)=_R\stackrel{~}{Z}_R(𝐫)`$ is given as a sum of angular momentum dependent Gauss functions, which have an analytical Fourier transform. A similar term occurs also in the pseudopotential approach. In contrast to the norm-conserving pseudopotential approach however, the compensation charge is non-spherical and it is constantly adapting to the instantaneous environment. It is constructed such that the augmentation charge densities

$`n_R^1(𝐫)+Z_R(𝐫)\stackrel{~}{n}_R^1(𝐫)\stackrel{~}{Z}_R(𝐫)`$ (30)
have vanishing electrostatic multi-pole moments for each atomic site. As a result the sum of all one-center contributions from one atom does not produce an electrostatic potential outside their own atomic region. This is the reason that the electrostatic interaction of the one-center parts between different sites vanish.

The compensation charge density as given here is still localized within the atomic regions, but a technique similar to an Ewald summation allows to replace it by a very extended charge density. Thus we can achieve, that all functions in $`\stackrel{~}{E}`$ converge as fast as the auxiliary density itself.

The potential $`\overline{v}`$, which occurs in Eqs. 27 and 29 enters the total energy in the form of a zero described in Eq. 25

$`{\displaystyle \underset{n}{}}f_n\stackrel{~}{\mathrm{\Psi }}_n|\left(\overline{v}{\displaystyle \underset{i,j}{}}|\stackrel{~}{p}_i\stackrel{~}{\varphi }_i|\overline{v}|\stackrel{~}{\varphi }_j\stackrel{~}{p}_j|\right)|\stackrel{~}{\mathrm{\Psi }}_n`$ (31)
The main reason for introducing this potential is that the self-consistent potential resulting from the plane wave part is not necessarily optimally smooth. The potential $`\overline{v}`$ allows to influence the plane wave convergence beneficially, without changing the converged result. $`\overline{v}`$ must be localized within the augmentation region, where equation 12 holds.

## 6 Approximations

Once the total energy functional provided in the previous section has been defined, everything else follows: Forces are partial derivatives with respect to atomic positions. The potential is the derivative of the potential energy with respect to the density, and the Hamiltonian follows from derivatives with respect to wave functions. The fictitious Lagrangian approach of Car and Parrinello \[Car et al 1985\] does not allow any freedom in the way these derivatives are obtained. Anything else than analytic derivatives will violate energy conservation in a dynamical simulation. Since the expressions are straightforward, even though rather involved, we will not discuss them here.

All approximations are incorporated already in the total energy functional of the PAW method. What are those approximations?

* Firstly we use the frozen core approximation. In principle this approximation can be overcome.
* The plane wave expansion for the auxiliary wave functions must be complete. The plane wave expansion is controlled easily by increasing the plane wave cutoff defined as $`E_{PW}=\frac{1}{2}\mathrm{}^2G_{max}^2`$. Typically we use a plane wave cutoff of 30 Ry.
* The partial wave expansions must be converged. Typically we use one or two partial waves per angular momentum $`(\mathrm{},m)`$ and site. It should be noted that the partial wave expansion is not variational, because the partial wave expansion changes the total energy functional and not only the basis set.

We do not discuss here numerical approximations such as the choice of the radial grid, since those are easily controlled.

We mentioned earlier that the pseudopotential approach can be derived as a well defined approximation from the PAW method: The augmentation part $`\mathrm{\Delta }E=E^1\stackrel{~}{E}^1`$ is a functional of the one-center density matrix $`D_{i,j}`$ defined in Eq. 22. The pseudopotential approach can be recovered if we truncate a Taylor expansion of $`\mathrm{\Delta }E`$ about the atomic density matrix after the linear term. The term linear to $`D_{i,j}`$ is the energy related to the nonlocal pseudopotential.

$`\mathrm{\Delta }E(D_{i,j})`$ $`=`$ $`\mathrm{\Delta }E(D_{i,j}^{at})+{\displaystyle \underset{i,j}{}}(D_{i,j}D_{i,j}^{at}){\displaystyle \frac{\mathrm{\Delta }E}{D_{i,j}}}+O(D_{i,j}D_{i,j}^{at})^2`$ (32)
$`=`$ $`E_{self}+{\displaystyle \underset{n}{}}f_n\stackrel{~}{\mathrm{\Psi }}_n|v_{nl}|\stackrel{~}{\mathrm{\Psi }}_n+O(D_{i,j}D_{i,j}^{at})^2`$
Thus we can look at the PAW method also as a pseudopotential method with a pseudopotential that adapts to the instantaneous electronic environment, because the explicit nonlinear dependence of the total energy on the one-center density matrix is properly taken into account.

What are the main advantages of the PAW method compared to the pseudopotential approach?

Firstly all errors can be systematically controlled so that there are no transferability errors. As shown by Watson \[Watson et al 1998\] and Kresse \[Kresse et al 1999\], most pseudopotentials fail for high spin atoms such as Cr. While it is probably true that pseudopotentials can be constructed that cope even with this situation, a failure can not be known beforehand, so that some empiricism remains in practice: A pseudopotential constructed from an isolated atom is not guaranteed to be accurate for a molecule. In contrast, the converged results of the PAW method do not depend on a reference system such as an isolated atom, because it uses the full density and potential.

The PAW method provides access to the full charge and spin density, which is relevant for hyperfine parameters. Hyperfine parameters are sensitive probes of the electron density near the nucleus. In many situations they are the only information available that allows to deduce atomic structure and chemical environment of an atom. There are reconstruction techniques for the pseudopotential approach, which however, are poor man’s versions \[Van de Walle et al 1993\] of the PAW method.

The plane wave convergence is more rapid than in norm-conserving pseudopotentials and should in principle be equivalent to that of ultra-soft pseudopotentials \[Vanderbilt 1990\]. Compared to the ultra-soft pseudopotentials, however, the PAW method has the advantage that the total energy expression is less complex and therefore is expected to be more efficient.

The construction of pseudopotentials requires to determine a number of parameters. As they influence the results, their choice is critical. Also the PAW methods provides some flexibility in the choice of auxiliary partial waves. However, this choice does not influence the converged results.

## 7 Recent Developments

Since the first implementation of the PAW method in the CP-PAW code, a number of groups have adopted the PAW method. The second implementation was done by the group of Holzwarth \[Holzwarth et al 1997\]. The resulting PWPAW code is freely available \[Tackett et al 2001\]. This code is also used as a basis for the PAW implementation in the AbInit project \[AbInit\]. An independent PAW code has been developed by Valiev and Weare \[Valiev et al 1999\]. Recently the PAW method has been implemented into the VASP code \[Kresse et al 1999\]. The PAW method has also been implemented by W. Kromen into the ESTCoMPP code of Blügel and Schröder \[Blügel et al 2001\].

Another branch of methods uses the reconstruction of the PAW method, without taking into account the full wave functions in the self-consistency. Following chemist notation this approach could be termed “post-pseudopotential PAW”. This development began with the evaluation for hyperfine parameters from a pseudopotential calculation using the PAW reconstruction operator \[Van de Walle et al 1993\] and is now used in the pseudopotenital approach to calculate properties that require the correct wave functions

The implementation by Kresse and Joubert \[Kresse et al 1999\] has been particularly useful as they had an implementation of PAW in the same code as the ultra-soft pseudopotentials, so that they could critically compare the two approaches with each other and LAPW calculations. Their conclusion is that both methods compare well in most cases, but they found that magnetic energies are seriously – by a factor two – in error in the pseudopotential approach, while the results of the PAW method were in line with other all-electron calculations using the linear augmented plane wave method. As a short note, Kresse and Joubert incorrectly claim that their implementation is superior as it includes a term that is analogous to the non-linear core correction of pseudopotentials \[Louie et al 1982\]: this term however is already included in the original version in the form of the pseudized core density.

Several extensions of the PAW have been done in the recent years: For applications in chemistry truly isolated systems are often of great interest. As any plane-wave based method introduces periodic images, the electrostatic interaction between these images can cause serious errors. The problem has been solved by mapping the charge density onto a point charge model, so that the electrostatic interaction could be subtracted out in a self-consistent manner \[Blöchl 1995\]. In order to include the influence of the environment, the latter was simulated by simpler force fields using the molecular-mechanics-quantum-mechanics (QM-MM) approach \[Woo et al 2000\].

In order to overcome the limitations of the density functional theory several extensions have been performed. Bengone \[Bengone et al 2000\] implemented the LDA+U approach \[Anisimov et al 1991\] into the CP-PAW code. Soon after this, Arnaud \[Arnaud et al 2000\] accomplished the implementation of the GW approximation into the CP-PAW code. The VASP-version of PAW \[Hobbs et al 2000\] and the CP-PAW code have now been extended to include a noncollinear description of the magnetic moments. In a non-collinear description the Schrödinger equation is replaced by the Pauli equation with two-component spinor wave functions

The PAW method has proven useful to evaluate electric field gradients \[Petrilli et al 1998\] and magnetic hyperfine parameters with high accuracy \[Blöchl et al 2000\]. Invaluable will be the prediction of NMR chemical shifts using the GIPAW method of Pickard and Mauri \[Pickard et al 2001\], which is based on their earlier work \[Mauri et al 1996\]. While the GIPAW is implemented in a post-pseudopotential manner, the extension to a self-consistent PAW calculation should be straightforward. An post-pseudopotential approach has also been used to evaluate core level spectra \[Jayawardane et al 2001\] and momentum matrix elements \[Kageshima et al 1997\].

## 8 Acknowledgment

We are grateful for carefully reading the manuscript to Dr. J. Noffke, Dr. R. Schmidt and P. Walther and to Prof. K. Schwarz for his continuous support. This work has benefited from the collaborations within the ESF Programme on ’Electronic Structure Calculations for Elucidating the Complex Atomistic Behavior of Solids and Surfaces’.


format: text

text:
# Electroweak Symmetry Breaking by Strong Dynamics and the Collider Phenomenology1footnote 11footnote 1This is the report of the Strong Electroweak Symmetry Breaking Working Group at 2001 Snowmass summer studies.

## I Introduction

Particle physics is on the verge of major discovery. General arguments indicate that new physics in the electroweak symmetry breaking sector must show up below the scale of 1 TeV. The experiments at the Tevatron and next generation high energy colliders such as the LHC and a TeV $`e^+e^{}`$ linear collider will fully explore the new physics at the electroweak scale.

In theories of dynamical electroweak symmetry breaking, the electroweak interactions are broken to electromagnetism by the vacuum expectation value of a fermion bilinear. These theories may thereby avoid the introduction of fundamental scalar particles, of which we have no examples in nature thus far. Prominent examples include Technicolor, topcolor, and related models. If the new dynamical scale is somewhat higher than 1 TeV, then the low energy effects or the early signature at collider experiments may be anomalous gauge boson interactions, enhanced $`WW`$ scattering signals, or contact 4-fermion interactions. In this report, we first briefly introduce the dynamical electroweak symmetry breaking models and parameterization of the anomalous couplings. We then summarize the collider sensitivities to probe the new dynamics at future $`e^+e^{}`$ linear colliders in Sec. II, and at hadron colliders in Sec. III.

### I.1 Technicolor

The earliest modelsWeinberg (1979); Susskind (1979) of dynamical electroweak symmetry breakingLane ; Simmons (2001) include a new non-Abelian gauge theory (“Technicolor”) and additional massless fermions (“technifermions”) which feel this new force. The global chiral symmetry of the fermions is spontaneously broken by the formation of a technifermion condensate, just as the approximate chiral $`SU(2)\times SU(2)`$ symmetry in QCD is broken down to $`SU(2)`$ isospin by the formation of a quark condensate. If the quantum numbers of the technifermions are chosen correctly (e.g. by choosing technifermions in the fundamental representation of an SU$`(N)`$ Technicolor gauge group, with the left-handed technifermions being weak doublets and the right-handed ones weak singlets) this condensate can break the electroweak interactions down to electromagnetism.

The breaking of the global chiral symmetries implies the existence of Goldstone bosons, the “technipions” ($`\pi _T`$). Through the Higgs mechanism, three of the Goldstone bosons become the longitudinal components of the $`W`$ and $`Z`$, and the weak gauge bosons acquire a mass proportional to the technipion decay constant (the analog of $`f_\pi `$ in QCD). The quantum numbers and masses of any remaining technipions are model dependent. There may be technipions which are colored (octets and triplets) as well as those carrying electroweak quantum numbers, and some color-singlet technipions are too lightEichten and Lane (1981); Farhi and Susskind (1981) unless additional sources of chiral-symmetry breaking are introduced. The next lightest Technicolor resonances are expected to be the analogs of the vector mesons in QCD. The technivector mesons can also have color and electroweak quantum numbers and, for a theory with a small number of technifermions, are expected to have a mass in the TeV rangeDimopoulos et al. (1981).

While Technicolor chiral symmetry breaking can give mass to the $`W`$ and $`Z`$ particles, additional interactions must be introduced to produce the masses of the standard model fermions. The most thoroughly studied mechanism for this invokes “extended Technicolor” (ETC) gauge interactionsEichten and Lane (1981); Dimopoulos and Susskind (1979). In ETC, Technicolor, color and flavor are embedded into a larger gauge group which is broken to Technicolor and color at an energy scale of 100s to 1000s of TeV. The massive gauge bosons associated with this breaking mediate transitions between quarks/leptons and technifermions, giving rise to the couplings necessary to produce fermion masses. The ETC gauge bosons also mediate transitions among technifermions themselves, leading to interactions which can explicitly break unwanted chiral symmetries and raise the masses of any light technipions. The ETC interactions connecting technifermions to quarks/leptons also mediate technipion decays to ordinary fermion pairs. Since these interactions are responsible for fermion masses, one generally expects technipions to decay to the heaviest fermions kinematically allowed (though this need not hold in all models).

In addition to quark masses, ETC interactions must also give rise to quark mixing. One expects, therefore, that there are ETC interactions coupling quarks of the same charge from different generations. A stringent limit on these flavor-changing neutral current interactions comes from $`K^0`$$`\overline{K}^0`$ mixingEichten and Lane (1981). These force the scale of ETC breaking and the corresponding ETC gauge boson masses to be in the 100-1000 TeV range (at least insofar as ETC interactions of first two generations are concerned). To obtain quark and technipion masses that are large enough then requires an enhancement of the technifermion condensate over that expected naively by scaling from QCD. Such an enhancement can occur if the Technicolor gauge coupling runs very slowly, or “walks”Holdom (1981). Many technifermions typically are needed to make the TC coupling walk, implying that the Technicolor scale and, in particular, the technivector mesons may be much lighter than 1 TeVLane ; Eichten and Lane (1989). It should also be noted that there is no reliable calculation of electroweak parameters in a walking Technicolor theory, and the values of precisely measured electroweak quantitiesLangacker and Erler  cannot directly be used to constrain the models.

In existing colliders, technivector mesons are dominantly produced when an off-shell standard model gauge-boson “resonates” into a technivector meson with the same quantum numbersEichten et al. (1984). The technivector mesons may then decay, in analogy with $`\rho \pi \pi `$, to pairs of technipions. However, in walking Technicolor the technipion masses may be increased to the point that the decay of a technirho to pairs of technipions is kinematically forbiddenEichten and Lane (1989). In this case the decay to a technipion and a longitudinally polarized weak boson (an “eaten” Goldstone boson) may be preferred, and the technivector meson would be very narrow. Alternatively, the technivector may also decay, in analogy with the decay $`\rho \pi \gamma `$, to a technipion plus a photon, gluon, or transversely polarized weak gauge boson. Finally, in analogy with the decay $`\rho e^+e^{}`$, the technivector meson may resonate back to an off-shell gluon or electroweak gauge boson, leading to a decay into a pair of leptons, quarks, or gluons.

### I.2 Top Condensate and Related Models

The top quark is much heavier than other fermions and must be more strongly coupled to the symmetry-breaking sector. It is natural to consider whether some or all of electroweak-symmetry breaking is due to a condensate of top quarksLane ; Miransky et al. (1989). Top-quark condensation alone, without additional fermions, seems to produce a top-quark mass largerBardeen et al. (1990) than observed experimentally, and is therefore not favored. Topcolor assisted TechnicolorHill (1995) combines Technicolor and top-condensation. In addition to Technicolor, which provides the bulk of electroweak symmetry breaking, top-condensation and the top quark mass arise predominantly from “topcolor,” a new QCD-like interaction which couples strongly to the third generation of quarks. An additional, strong, U(1) interaction (giving rise to a topcolor $`Z^{}`$) precludes the formation of a $`b`$-quark condensate.

The top-quark seesaw model of electroweak symmetry breakingDobrescu and Hill (1998) is a variant of the original top-condensate idea which reconciles top-condensation with a lighter top-quark mass. Such a model can easily be consistent with precision electroweak tests, either because the spectrum includes a light composite HiggsChivukula et al. (1999); Chivukula and Hoelbling (2001) or because additional interactions allow for a heavier HiggsCollins et al. ; He et al. (2001). Such theories may arise naturally from gauge fields propagating in compact extra spatial dimensionsDobrescu (1999).

A variant of topcolor-assisted Technicolor is flavor-universal, in which the topcolor SU(3) gauge bosons, called colorons, couple equally to all quarksChivukula et al. (1996); Popovic and Simmons (1998). Flavor-universal versions of the seesaw modelSimmons (1989); Burdman and Evans (1999) incorporating a gauged flavor symmetry are also possible. In these models all left-handed quarks (and possibly leptons as well) participate in electroweak symmetry-breaking condensates with separate (one for each flavor) right-handed weak singlets, and the different fermion masses arise by adjusting the parameters which control the mixing of each fermion with the corresponding condensate. A prediction of these flavor-universal models, is the existence of new heavy gauge bosons, coupling to color or flavor, at relatively low mass scales. A mass limit of between 0.8 and 3.5 TeV is setBertram and Simmons (1998) depending on the coloron-gluon mixing angle. Precision electroweak measurements constrainBurdman et al. (2000) the masses of these new gauge bosons to be greater than 1–3 TeV in a variety of models, for strong couplings.

### I.3 Enhanced gauge-boson couplings and fermion compositeness

If the new strong dynamics scale is somewhat higher than that accessible to the next generation of colliders, the expected signature would be enhanced gauge-boson self-interactions conventionally parameterized by the “anomalous couplings”Hagiwara et al. (1987); Dawson and Valencia (1991); Dobado et al. (1991); Appelquist and Wu (1993), and the fermion contact interactions the so-called “fermion compositeness”Eichten et al. (1983) at a scale $`\mathrm{\Lambda }`$.

Although the current LEP and Tevatron experiments have put stringent bounds on the anomalous gauge-boson self-interactions, the anticipated size of those couplings due to new strong dynamics may be of order $`v^2/\mathrm{\Lambda }^21/16\pi ^2<10^3`$, smaller than the current bounds. Experiments at future colliders will reach sensitivity to this level. In particular, high energy scattering of longitudinal gauge-bosons $`W_L,Z_L`$ as the electroweak Goldstone bosons should be the most direct probe to the electroweak symmetry breaking sector. General arguments such as unitarity Lee et al. (1977); Chanowitz and Gaillard (1985) indicate that new physics associated with the electroweak symmetry breaking must show up in some form at the scale of TeV, which can be accessible most likely only at higher energy colliders of next generation. Regarding the fermion compositeness, higher sensitivity will be reached at higher energies due to the energy-dependent nature of the dimension 6-operators Eichten et al. (1983). In the next two sections, we will summarize the studies of the above physics scenarios at future colliders.

## II Strong Dynamics at $`e^+e^{}`$ Linear Colliders

An $`e^+e^{}`$ linear collider with $`\sqrt{s}=0.51.5`$ TeV and a luminosity of $`5001000`$ pb-1 can be a very effective probe of strong electroweak symmetry breaking. Production mechanisms and backgrounds are limited to electroweak processes, so that signal and background cross sections can be calculated exactly. The initial state is well defined not only in terms of four-momentum, but also in terms of electron (and possibly positron) helicity. Also, complete final state helicity analyses are possible, due to the fact that most if not all of the final state kinematic variables can be reconstructed.

In this section we review the $`e^+e^{}`$ collider phenomenology of strong $`W^+W^{}`$ interactions which appear when there is no light Higgs particle with large couplings to vector gauge bosons. Detection of directly produced narrow-width spinless particles such as technipions Casalbuoni et al. (1999) and top-pions Yue et al. (2000) is straightforward up to the kinematic limit, and will not be discussed further.

### II.1 $`e^+e^{}\nu \overline{\nu }W^+W^{}`$, $`\nu \overline{\nu }ZZ,W^+W^{}Z,ZZZ,\nu \overline{\nu }t\overline{t}`$

The first step in studying the reaction $`e^+e^{}\nu \overline{\nu }W^+W^{}`$ is to separate the scattering of a pair of longitudinally polarized $`W`$’s, denoted by $`W_\mathrm{L}W_\mathrm{L}`$, from transversely polarized $`W`$’s and background such as $`e^+e^{}e^+e^{}W^+W^{}`$ and $`e^{}\overline{\nu }W^+Z`$. Studies have shown that simple cutsBarger et al. (1995) can be used to achieve this separation in $`e^+e^{}\nu \overline{\nu }W^+W^{}`$, $`\nu \overline{\nu }ZZ`$ at $`\sqrt{s}=1000`$ GeV, and that the signalsBoos et al. (1998, 2000) are comparable to those obtained at the LHCBarger et al. (1990); Bagger et al. (1994); Chanowitz and Kilgore (1994); Bagger et al. (1995). Furthermore, by analyzing the gauge boson production and decay angles it is possible to use these reactions to measure chiral Lagrangian parameters with an accuracy greater than that which can be achieved at the LHC Chierici et al. (2001).

The chiral Lagrangian parameters associated with quartic gauge boson couplings can also be measured with the triple gauge boson production processes $`e^+e^{}W^+W^{}Z`$ and $`e^+e^{}ZZZ`$ Barger and Han (1988); Barger et al. (1989); Han et al. (1998). These measurements complement the $`W^+W^{}`$ fusion measurements, and they will play a crucial role in multi-parameter chiral Lagrangian analyses.

The reaction $`e^+e^{}\nu \overline{\nu }t\overline{t}`$ provides unique access to $`W^+W^{}t\overline{t}`$ since this process is overwhelmed by the background $`ggt\overline{t}`$ at the LHC. Techniques similar to those employed to isolate $`W_\mathrm{L}W_\mathrm{L}W^+W^{},ZZ`$ can be used to measure the enhancement in $`W_\mathrm{L}W_\mathrm{L}t\overline{t}`$ productionBarklow (1996); Larios et al. (1998); Ruiz Morales and Peskin (1999); Han et al. (2001). Even in the absence of a resonance it will be possible to establish a clear signal. The ratio $`S/\sqrt{B}`$ is expected to be 12 for a linear collider with $`\sqrt{s}=1`$ TeV, 1000 fb-1 and $`80\%/0\%`$ electron/positron beam polarization, increasing to 22 for the same luminosity and beam polarization at $`\sqrt{s}=1.5`$ TeV.

### II.2 $`e^+e^{}W^+W^{}`$

Strong gauge boson interactions induce anomalous triple gauge couplings (TGC’s) at tree-levelHagiwara et al. (1987); Appelquist and Wu (1993); Dawson and Valencia (1991); Dobado et al. (1991):

$`\kappa _\gamma `$ $`=`$ $`1+{\displaystyle \frac{e^2}{32\pi ^2s_w^2}}\left(L_{9L}+L_{9R}\right)`$
$`\kappa _Z`$ $`=`$ $`1+{\displaystyle \frac{e^2}{32\pi ^2s_w^2}}\left(L_{9L}{\displaystyle \frac{s_w^2}{c_w^2}}L_{9R}\right)`$
$`g_1^Z`$ $`=`$ $`1+{\displaystyle \frac{e^2}{32\pi ^2s_w^2}}{\displaystyle \frac{L_{9L}}{c_w^2}}.`$
where $`\kappa _\gamma `$, $`\kappa _Z`$, and $`g_1^Z`$ are TGC’s, $`s_w^2=\mathrm{sin}^2\theta _w`$, $`c_w^2=\mathrm{cos}^2\theta _w`$, and $`L_{9L}`$ and $`L_{9R}`$ are chiral Lagrangian parametersBagger et al. (1993). Assuming QCD values for $`L_{9L}`$ and $`L_{9R}`$, $`\kappa _\gamma `$ is shifted by $`\mathrm{\Delta }\kappa _\gamma 3\times 10^3`$.

Table 1 contains the estimates of the TGC precision that can be obtained at $`\sqrt{s}=500`$ and 1000 GeV for the CP-conserving couplings $`g_1^V`$, $`\kappa _V`$, and $`\lambda _V`$ Abe et al. (2001). These estimates are derived from one-parameter fits in which all other TGC parameters are kept fixed at their tree-level SM values. The $`4\times 10^4`$ precision for the TGCs $`\kappa _\gamma `$ and $`\kappa _Z`$ at $`\sqrt{s}=500`$ GeV can be interpreted as a precision of $`0.26`$ for the chiral Lagrangian parameters $`L_{9L}`$ and $`L_{9R}`$. Assuming naive dimensional analysisManohar and Georgi (1984) such a measurement would provide a $`8\sigma `$ ($`5\sigma `$) signal for $`L_{9L}`$ and $`L_{9R}`$ if the strong symmetry breaking energy scale were 3 TeV (4 TeV).

When $`W^+W^{}`$ scattering becomes strong the amplitude for $`e^+e^{}W_\mathrm{L}W_\mathrm{L}`$ develops a complex form factor $`F_T`$ in analogy with the pion form factor in $`e^+e^{}\pi ^+\pi ^{}`$Peskin (1984); Iddir et al. (1990). To evaluate the size of this effect the following expression for $`F_T`$ can be used:

$$F_T=\mathrm{exp}\left[\frac{1}{\pi }_0^{\mathrm{}}𝑑s^{}\delta (s^{},M_\rho ,\mathrm{\Gamma }_\rho )\{\frac{1}{s^{}siϵ}\frac{1}{s^{}}\}\right]$$

where

$$\delta (s,M_\rho ,\mathrm{\Gamma }_\rho )=\frac{1}{96\pi }\frac{s}{v^2}+\frac{3\pi }{8}\left[\mathrm{tanh}(\frac{sM_\rho ^2}{M_\rho \mathrm{\Gamma }_\rho })+1\right].$$

Here $`M_\rho ,\mathrm{\Gamma }_\rho `$ are the mass and width respectively of a vector resonance in $`W_\mathrm{L}W_\mathrm{L}`$ scattering. The term

$$\delta (s)=\frac{1}{96\pi }\frac{s}{v^2}$$

is the Low Energy Theorem (LET) amplitude for $`W_\mathrm{L}W_\mathrm{L}`$ scattering at energies below a resonance. Below the resonance, the real part of $`F_T`$ is proportional to $`L_{9L}+L_{9R}`$ and can therefore be interpreted as a TGC. The imaginary part, however, is a distinct new effect.

The expected $`95\%`$ confidence level limits for $`F_T`$ for $`\sqrt{s}=500`$ GeV and a luminosity of 500 $`fb^1`$ are shown in Figure 1, along with the predicted values of $`F_T`$ for various masses $`M_\rho `$ of a vector resonance in $`W_\mathrm{L}W_\mathrm{L}`$ scattering. The signal significances obtained by combining the results for $`e^+e^{}\nu \overline{\nu }W^+W^{}`$, $`\nu \overline{\nu }ZZ`$Barger et al. (1995); Boos et al. (1998) with the $`F_T`$ analysis of $`W^+W^{}`$ Barklow (2000) are displayed in Fig. 2 along with the results expected from the LHCATLAS Detector and Physics Performance Technical Design Report. Vol. 2 (1999). At all values of the center-of-mass energy a linear collider provides a larger direct strong symmetry breaking signal than the LHC for vector resonance masses of 1200, 1600 and 2500 GeV. Only when the vector resonance disappears altogether (the LET case in the lower right-hand plot in Fig. 2 ) does the direct strong symmetry breaking signal from the $`\sqrt{s}=500`$ GeV linear collider drop below the LHC signal. At higher $`e^+e^{}`$ center-of-mass energies the linear collider signal exceeds the LHC signal.

## III Strong Dynamics at Hadron Colliders

Hadron colliders offer exciting possibilities for searches for new particles and other signs of new strong dynamics and compositeness. High luminosity $`pp`$ and $`p\overline{p}`$ machines should copiously produce proposed strongly-coupled resonances including technihadrons and excited quarks. They also probe contact interactions and vector boson scattering at extremely high energy scales. In this section we describe the expected physics reach of hadron colliders that exist (the Tevatron), are under construction (the LHC) and are being designed (the VLHC).

### III.1 The Tevatron

The Tevatron at Fermilab has taken approximately 100 pb-1 of $`p\overline{p}`$ collision data at $`\sqrt{s}=1.8`$ TeV (Run I). In March 2001 Run II began, with an increased energy ($`\sqrt{s}=1.96`$ TeV) and a planned integrated luminosity of 2 fb-1 (Run IIa), followed by extended high luminosity running for a total in excess of 15 fb-1 per experiment. In Tables 25 announced results from Run I are tabulated along with extrapolations to RunIIa and a possible 30 fb-1 complete RunII.

### III.2 The LHC

Despite the challenge at hadron colliders in the search for new strong dynamics at the TeV scale, much theoretical work has been performed at the LHCBarger et al. (1990); Bagger et al. (1994); Chanowitz and Kilgore (1994); Bagger et al. (1995). Many studies of strong EWSB at ATLAS and CMS have been performed and summarized in several placesATLAS (1999); ATLAS and CMS (2001); p1p (1999); Butterworth et al. . An expected “low luminosity” period will collect 30 fb-1 of data at $`\sqrt{s}=14`$ TeV over the first three years of operation, and will be followed by a similar “high luminosity” period collecting up to 300 fb-1. High luminosity running (up to $`10^{34}`$ cm-2s-1) presents many experimental challenges with an average of 20 collisions per beam crossing, degrading tracking and electron identification capabilities particularly in the forward region.

As an example of a Technicolor resonance search, ATLAS have considered the production of 500 GeV technirho in a multiscale Technicolor model and its signal in the channel $`\rho _T^\pm WZl^\pm \nu l^+l^{}`$ ATLAS (1999). This study assumes the 30 fb-1 of low luminosity data and hence the full lepton ID and tracking capabilities of the detector. The expected signal significance is strongly dependent on the input model parameters: a narrow resonance ($`\mathrm{\Gamma }_{\rho _T}=1.1`$ GeV) which is not allowed to decay to $`\pi _T\pi _T`$ ($`m_{\pi _T}>m_{\rho _T}/2`$) could have $`S/\sqrt{B}80`$; but for $`\mathrm{\Gamma }_{\rho _T}=110`$ GeV and $`m_{\pi _T}=110`$ GeV this would drop to an indiscernable $`S/\sqrt{B}0.3`$.

The masses of observable resonances at the LHC are expected to be 5-10$`\times `$ those at the Tevatron. A $`Z^{}`$ with couplings similar to those of the Standard Model $`Z`$ should be observable up to $`m_Z^{}5`$ TeV and direct observation of excited quarks of $`m_q^{}6`$ TeV is possibleHarris (1996c). The reach for compositeness scales is similarly enhanced, with 300 fb-1 of dijet data being sensitive to $`\mathrm{\Lambda }40`$ TeV.

A further possibility at the LHC is that as $`\sqrt{\widehat{s}}`$ begins to exceed 1 TeV, strong interaction effects in $`WW`$ scattering could become detectable. If jets can be reliably tagged in the forward region at high luminosities, a signal should be observable with the full 300 fb-1.

### III.3 The Super-LHC

There has been some discussion of upgrading the LHC in luminosity and energy after the 300 fb-1 run is complete. A possible (though unlikely) doubling of the energy has been considered along with a tenfold increase in instantaneous luminosity. Since the LHC detectors were not designed for these conditions only jet and muon information is likely to be usefull. Such an upgrade could double the reach for a $`Z^{}`$ ($`m_Z^{}10`$ TeV) and compositeness ($`\mathrm{\Lambda }80`$ TeV), and significantly increase the sensitivity for excited quarks ($`m_q^{}9`$ TeV) and the scale of $`WW`$ scattering available ($`\sqrt{\widehat{s}}1.5`$ TeV, assuming that forward jet tagging is still possible). Unfortunately, most of these gains come from the energy increase, which is less plausible than a simple luminosity upgrade.

### III.4 The VLHC

A staged 40-175 TeV $`p\overline{p}`$ collider operating at luminosities comparable to the LHC (1-2$`\times 10^{34}`$ cm-2s-1) has been proposedGroup (2001). Studies of such a machine’s physics reach are in progress (see also the E4 Working Group report), but the direct reach for excited quark resonances is expected to be $`m_q^{}25`$ TeV for 10 fb-1 at $`\sqrt{s}=100`$ TeVHarris (1996c), and $`WW`$ scattering could be probed at the scale of $`23`$ TeV.

New signatures could become detectable at such high center-of-mass energies. For example, in topcolor models, direct $`\chi `$ pair production and subsequent decays $`\chi htt\overline{t}t`$ could occurHe et al. (2001), with a $`6t`$ final state. Such a heavy state may only be copiously produced. The cross section for this process with $`m_\chi =1`$ TeV would be $`10`$ pb, as shown in Fig. 3.

In interactions with $`\sqrt{\widehat{s}}\mathrm{\Lambda }_{TC}`$, it is possible (in analogy with QCD) that asymptotically free techniquarks could be produced that subsequently hadronize into technijets consisting of weak vector bosons and technihadrons. A technijet would manifest itself as an extremely massive but significantly boosted (and hence not necessarily wide) jet in a VLHC detector. The production rate for such a process can be significant: For $`m_{Q_T}=400`$ GeV with $`\sqrt{s}=100`$ TeV the dijet differential cross section for technijets exceeds that for $`t\overline{t}`$ for dijet masses $`>900`$ GeV. Exploration of technijets could provide the ultimate determination of the TC dynamics. As shown in Fig. 4, a representative techni-quark may decay subsequently into multiple jets and the separation between any two jets may be small enough so that experimental signature would be a very massive (but not too fat) jet.

###### Acknowledgements.

We would like to thank the participants in the Strong Electroweak Symmetry Breaking Working group for the contribution during the Snowmass workshop, on which the materials summarized in this report are based. T.L.B. was supported by Department of Energy contract DE-AC03-76SF00515. R.S.C. was supported in part by a DOE grant DE-FG02-91ER40676. T.H. was supported in part by a DOE grant No. DE-FG02-95ER40896 and in part by the Wisconsin Alumni Research Foundation.


text:
# Carrier relaxation dynamics in heavy fermion compounds

## I Introduction

Femtosecond time-resolved optical spectroscopy is an excellent experimental alternative to conventional spectroscopic methods that probe the low energy electronic structure in strongly correlated electron systems Kabanov . In particular, it has been shown that carrier relaxation dynamics are very sensitive to changes in the low energy density of states (e.g. associated with the formation of a low energy gap or pseudogap) providing new insights into the low energy electronic structure in these materials. In this report we present the first studies of carrier relaxation dynamics in heavy fermion (HF) systems by means of femtosecond time-resolved optical spectroscopy. Our results show that the carrier relaxation dynamics, below the Kondo temperature (TK), are extremely sensitive to the low energy density of states (DOS) near the Fermi level to which localized f-moments contribute. Specifically, we have performed measurements of the photoinduced reflectivity $`\mathrm{\Delta }`$R/R dynamics as a function of temperature and excitation intensity on the series of HF compounds YbXCu4 (X=Ag, Cd, In) Graf  in comparison to their non-magnetic counterparts LuXCu4.

## II Experimental

The experiments were performed on freshly polished flux-grown single crystals. We used a standard pump-probe set-up with a mode-locked Ti:Sapphire laser (30 fs pulses at 800 nm and 80 MHz repetition rate) for both pump and probe pulse trains. The photoinduced (PI) change in reflectivity $`\mathrm{\Delta }`$R/R was measured using photodiode and lock-in detection. The pump laser fluence was kept below 10 $`\mu `$J/cm2 and the pump/probe intensity ratio was $``$30. The steady-state heating effects were accounted for as described in Mihailovic  yielding an uncertainty in temperature of $`\pm `$3 K.

## III Results and discussion

In Fig. 1a we show the PI reflectivity traces as a function of temperature on the HF system YbCdCu4 (TK $``$100K, low temperature Sommerfeld coefficient $`\gamma `$ $``$200 mJ/mol K2). The relaxation time is virtually temperature independent above $``$140 K and shows a quasi-divergence (hundred-fold increase) as T = 0 K is approached.

Similar behavior (see Fig. 1b) is observed also in YbAgCu4 (TK $``$100K, $`\gamma `$ $``$200 mJ/mol K2) Graf . Importantly, YbAgCu4 is a much better metal than YbCdCu4 revealing that the observed dynamics are not merely due to changes in the carrier concentration. Furthermore, the relaxation time of the non-HF compound LuAgCu4 behaves similar to regular metals Groenevald , showing a gradual decrease upon cooling. This, coupled with the data on the HF compounds in Fig. 1b, indicates that the upturn in the relaxation time at $``$TK and the subsequent divergence as T $``$ 0 K derives from the low energy electronic structure of these HF compounds.

An important and unresolved issue in HF physics Degiorgi  is whether the electronic ground state is a many-body Kondo state that forms only below TK, or is more simply described in terms of the thermal population of a narrow gap system with the gap arising from hybridization between the conduction electrons and localized *f* electrons (note that the gap is present at all temperatures above and below TK in this case). As we have shown, the ultrafast dynamics are clearly sensitive to the low energy electronic state and such experiments may provide insight into the true nature of the HF ground state. However, a complete analysis (currently in progress) requires many-body calculations beyond the scope of this paper. In the following, we present a heuristic model, which provides some physical insight into the low temperature electron dynamics in HF compounds.

The two-temperature model is commonly employed to interpret the ultrafast dynamics in metals Groenevald  with the sequence of relaxation events after photoexcitation as follows: the photoexcited carriers first rapidly thermalize via electron-electron scattering ($`\tau _{ee}`$100 fs), giving rise to a rapid change in reflectivity on a sub-picosecond timescale. This is followed by energy transfer to the lattice with a characteristic electron-phonon relaxation time $`\tau _{eph}`$ on the order of a picosecond Groenevald . At low excitation fluences, as used here, the increase in the electronic temperature after photoexcitation is small compared to the lattice temperature, T. In this limit the T-dependence of $`\tau _{eph}`$ in the two-temperature model is given by $`\tau _{eph}=\gamma (T)T/g(T)`$, where $`g(T)`$ is the electron-phonon coupling function Groenevald . Since $`g(T)`$ is constant above the Debye temperature $`\mathrm{\Theta }_D`$, and proportional to T4 at T$`<`$$`<`$ $`\mathrm{\Theta }_D`$, $`\tau _{eph}`$ is expected to decrease linearly upon cooling from room temperature, showing a minimum with a subsequent increase going as T-3 at low temperatures. However, this low-T upturn in $`\tau _R`$ is not observed in LuAgCu4 (Fig 1b) - as it was not observed in conventional Fermi liquid metals like Au and Ag Groenevald . The lack of the low-T quasi-divergence in $`\tau _R`$ in conventional metals was shown to result due to the fact that $`\tau _{ee}`$ and $`\tau _{eph}`$ are comparable below $``$300-600 K Groenevald  in the experimental configuration of low laser powers. A non-thermal electron model Groenevald  that takes into account a ”finite” $`\tau _{ee}`$ can account for a) the increased room temperature $`\tau _{eph}`$ and b) the low temperature plateau in the observed temperature dependence of $`\tau _{eph}`$ in Au and Ag Groenevald  as well as in LuAgCu4.

However, in the HF compounds we have studied $`\gamma `$ is, below TK, two orders of magnitude larger than in conventional metals. This is expected to give rise to a hundred-fold increase in $`\tau _{eph}`$ in these systems and the two-temperature model (without considering the initial non-thermal electron distribution) is expected to be valid. This would give rise to $`\tau _R`$ $`T^3`$ at low-T as experimentally observed (Fig. 1b). Furthermore, since $`g(T)`$ depends only on $`\mathrm{\Theta }_D`$ and a prefactor that can be estimated from the high-T value of $`\tau _R`$ Groenevald  (at 300K the value of $`\tau _R`$ is close to the estimated $`\tau _{eph}`$), one can extract the Sommerfeld constant via $`\gamma (T)=g(T)\tau _R/T`$. We find a remarkable agreement between the low temperature $`\gamma (T)`$ extracted in this way from $`\tau _R(T)`$ and the theoretical calculation Rajan  for the case of single impurity model with impurity spin 7/2 - Fig 1c, giving further support to this analysis.

## IV Conclusions

In conclusion, we have utilized ultrafast optical techniques to study the dynamics of quasiparticles in HF compounds YbXCu4. We have observed a power-law divergence in the relaxation time at low temperatures. The dramatic hundred-fold increase in the relaxation time at low temperatures in YbXCu4 (and the lack of this quasi-divergence in the non-HF LuXCu4 analogs) clearly demonstrates that our measurements probe the heavy fermion ground state.

* phone: 1-505-665-8839, fax 1-505-665-7652, e-mail: jdemsar@lanl.gov