text:
# Quenches and crunchs: Does the system explore in aging the same part of the configuration space explored in equilibrium ?

## 1 Introduction

At the glass transition temperature $`T_g`$ (Debenedetti 1997), the characteristic relaxation time of a liquid becomes of the same order of the experimental time, thus preventing equilibrium studies at lower $`T`$. Material properties below $`T_g`$ depend on the previous history (i.e. on the preparation technique, on the cooling/compression rate and so on) as well as on the time spent in the glass-state. This time dependence, generically known as physical aging, highlights the out-of-equilibrium (OOE) condition of glasses and their extremely slow equilibration processes.

Substantial amount of work has been devoted to the understanding and to the formal description of the supercooled liquid dynamics (Götze 1999, Cummins 1999), of the physics beyond the glass transition (Angell 1995) and of physical aging (Bouchaud et al. 1998, Cugliandolo et al. 1996, Latz 2000). Numerical simulations of supercooled states, both in equilibrium and in controlled out-of-equilibrium conditions, have played an important role in the present developments (Kob 1999). Although the time scales probed by numerical simulations are very different than experimental ones ($`100ns`$ compared to $`s`$), the numerical “experiments” appear to be able to reproduce features of real materials (Utz et al. 2000, Barrat and Berthier 2001).

The description of the aging dynamics as motion in configuration space, has been very fruitful(Kob et al. 2000). It has been shown that, during aging, the system explores deeper and deeper basins of the potential energy surface. The search for deeper basins during aging resembles the exploration of deeper and deeper basins which takes place in equilibrium on cooling. This similarity has been interpreted in term of a decrease of the internal configurational temperature $`T`$ of the system under aging (Kob et al. 2000). An analysis of the curvatures of the PES basins supported the possibility that during the aging process the system visits the same type of minima as the one visited in equilibrium. This analysis supported also the possibility of a thermodynamic description of the aging system and a prediction of the internal configurational temperature in full agreement with the numerical estimates (Sciortino and Tartaglia 2001a).

The outcome of these studies, which are still far from being settled, suggest that the out-of-equilibrium glassy state, notwithstanding its out of equilibrium condition, can be uniquely determined by its kinetic temperature, its volume and the properties of the basin in which the system is trapped. These quantities allow to develop a thermodynamic description of the glass state (Davies and Jones 1953, Nieuwenhuizen 1998, Speedy 1998, Mezárd and Parisi 1999, Sciortino et al. 1999) which is currently tested in experiments (Grigera and Israeloff 1999, Knaebel et al. 2000) and further simulations(Sciortino and Tartaglia 2001b). .

From an experimental point of view, aging experiments are performed at constant pressure. Moreover, often a change in pressure (and hence in density) is used to bring the system from an equilibrium to an out-of-equilibrium state. The change in density produces a much more dramatic change in the PES than a change in temperature. This opens the possibility that the hypothesis of out-of-equilibrium dynamics as succession of quasi equilibrium states may not apply to the crunch experiments and hence that a more refined thermodynamic approach is requested for glass materials produce with a crunch route.

In this manuscript we revisit the hypothesis that during the aging process the system explores the same part of the configuration space by performing a more accurate analysis of the relation basin-depth as a function of the basin-curvature, made nowadays possible by the increased computational facilities. We also present a comparison of this relation for different quenching temperatures and for a crunch (Di Leonardo et al 2000).

## 2 The models

We consider two models, an atomic one and a molecular one.

The first microscopic model we consider is a binary (80:20) mixture of Lennard-Jones particles (BMLJ), which in the following we will call type $`A`$ and type $`B`$ particles. The interaction between two particles of type $`\alpha `$ and $`\beta `$, with $`\alpha ,\beta \{\mathrm{A},\mathrm{B}\}`$, is given by $`V_{\alpha \beta }=4ϵ_{\alpha \beta }[(\sigma _{\alpha \beta }/r)^{12}(\sigma _{\alpha \beta }/r)^6]`$. The parameters $`ϵ_{\alpha \beta }`$ and $`\sigma _{\alpha \beta }`$ are given by $`ϵ_{AA}=1.0`$, $`\sigma _{AA}=1.0`$, $`ϵ_{AB}=1.5`$, $`\sigma _{AB}=0.8`$, $`ϵ_{BB}=0.5`$, and $`\sigma _{BB}=0.88`$. The potential is truncated and shifted at $`r_{\mathrm{cut}}=2.5\sigma _{\alpha \beta }`$. $`\sigma _{\mathrm{AA}}`$ and $`ϵ_{\mathrm{AA}}`$ are chosen as the unit of length and energy, respectively (setting the Boltzmann constant $`k_\mathrm{B}=1.0`$). Time is measured in units of $`\sqrt{m\sigma _{\mathrm{AA}}^2/48ϵ_{\mathrm{AA}}}`$, where $`m`$ is the mass of the particles. 1000 particles were placed in a box of side 9.4. The studied isochore has density of 1.2. To study the quenches, we have equilibrated several independent configurations at $`T=T_i`$ in the $`NVT`$ ensemble (Nosé-Hoover thermostat). Each of the configuration has been quenched to $`T_f`$ by changing at $`t=0`$ the thermostat temperature to $`T_f`$. The thermostat constant is chosen in such a way that, within 1000 MD steps, the average kinetic energy thermalizes to $`T_f`$. To study the crunches, we have equilibrated 280 independent configurations at $`T=0.35`$ in the $`NVT`$ ensemble at density of 1.103. The coordinates of each atom has been rescaled by 0.9714 at $`t=0`$ to simulate a 10% density increase. In this way, the final density of the crunches coincides with the density of the quenches. Inherent structures of a system are calculated by conjugate gradient algorithms. The minimization procedure is iterated until the energy change is less than $`10^{15}`$. The density of states is then calculated by diagonalizing the corresponding Hessian matrix.

The second model we consider is a simple three-sites molecular model (LW), introduced by Wahnström and Lewis (Wahnström and Lewis 1993). The model is constructed by gluing in a rigid molecule three identical Lennard Jones (LJ) atoms. The shape of the molecule (an isoscele triangle) and the LJ parameters were chosen to mimic as close as possible one of the most studied glass-forming, liquid,ortho-terphenyl. The slow dynamics of this model has been recently revisited (A. Rinaldi et al. 2001) in great details.

For this model we present preliminary results for a quench, starting from $`T=480`$ K down to $`T=280`$ K at constant density $`1108.4907`$ Kg $`/`$ m3. The simulated system is composed by 343 molecules. The integration time step is $`0.01`$ ps and the aging dynamics is followed up to $`1`$ ns. Averages over more than 50 independent realizations are presented.

## 3 BMLJ model

### 3.1 Energies

Following a quench or a crunch, the system finds itself in a region of configuration space which is not explored in equilibrium under the externally imposed temperature and volume. The motion of the system in configuration space evolves in the attempt of reaching the typical equilibrium configuration. A clear indicator of this evolution is the time evolution of the potential energy. In the PES paradigm, the potential energy of the system can be expressed as $`E=e_{IS}+E_{vib}`$, where $`e_{IS}`$ is the IS energy and $`E_{vib}`$ describes the thermal excitations about the IS (harmonic plus anharmonic vibrations). While the vibrational energy does not show a significant aging dependence, the time evolution of $`e_{IS}`$ provides a clear indication of the search for deeper and deeper potential energy minima, as shown in Fig.1.

### 3.2 Basins’ curvatures

In equilibrium, below $`T=1.`$, the system start to explore deeper and deeper basins, whose shape is depth dependent. A simple way to characterize the shape of the PES basins is to calculate the density of states in harmonic approximation, expanding the potential energy around the local minimum configuration. The resulting distribution of frequencies characterizes, in harmonic approximation, the volume in configuration space associated to the basin. The density of states allows to estimate the vibrational free energy $`f_{basin}`$ in harmonic approximation, since

$$f_{basin}(T)=\frac{k_BT}{N}\underset{i=1}{\overset{3N3}{}}log(\beta \mathrm{}\omega _i)$$
(1)
A $`T`$-independent indicator of the basin shape can be defined as $`S_{i=1}^{3N3}log(\mathrm{}\omega _i/ϵ_{AA})`$.

Studies on several models of liquids shows that along isochoric paths, the depth dependence of the average curvature of the basins is model and density dependent.

In the case of the BMLJ at the studied density, the average shape of the basins becomes wider and wider on moving to deeper and deeper basins, as shown in Fig. 2.

Fig. 2 also shows $`S(e_{IS})`$ during the aging dynamics following quenches and crunches, from different initial and bath conditions. The important observation is that during $`T`$-changes, there is an initial part of the aging process were the system explores in aging the same set of basins explored in equilibrium. For longer times, deviations start to take place and the system finds itself located in region of the configuration space which are not commonly explored in equilibrium conditions. The effect is negligible (and indeed it was not noted in previous studies) at sufficiently high $`T_f`$ values (shallow quenches), but is becomes detectable in deeper quenches. In future studies, it will become important to correlated the time at which the aging dynamics significantly separates from the equilibrium one with the different mechanism of exploration of configuration spaces which characterize also the dynamics of the system in equilibrium (i.e. saddle-dominated dynamics as compared to activated dynamics (Sciortino and Tartaglia 1997, La Nave et al 2000, Angelani et al 2000, Cavagna et al 2000)). A related interesting question is the validity of the description of the aging system as composed by two system in quasi-equilibrium at different temperatures (Sciortino and Tartaglia 2001) in conditions where the aging system explores basins which are not populated in equilibrium.

In the crunch case, in the time window accessible to numerical simulations the system never explores region of configuration space which are visited in equilibrium. The starting configuration is different from any configuration visited in equilibrium. In this respect, the ensemble of configuration of the system after an initial crunch can not be identified with an equilibrium ensemble at a different temperature, not even at infinite temperature.

Having said so, we call the reader attention on the fact that the differences observed in $`S`$ are smaller than the entire variation of $`S`$ with aging. We also recall that the time windows accessible to numerical experiment is such that the system, even after the longest waiting time simulated, has not completely forgot the initial configuration (or, in other words, the correlation functions never decay to zero in the simulated time period). In this respect, it is not that surprising that the basins explored in crunches are different from the basins explored in aging and in equilibrium. Finally we note that $`S`$, being the sum of logarithms, is driven by the low frequency values. The small frequencies are the most affected by size effects and by artificial localization of the eigenmodes. Future studies should focus on the size effect on the data reported in Fig.2.

We note on passing that some model potentials, like for example the BKS model for silica do not show, along an isochoric path any dependence of the basin shape on the potential depth(Saika-Voivod et al 2001). This class of models (and their corresponding materials) may provide simple aging and crunching dynamics compared to the BMLJ model here studied.

### 3.3 Structure of the liquid in the inherent structures

The inherent structure energy is, by definition, the integral of the pair potential over $`r`$, weighted by the radial distribution function $`g(r)`$. The changes in inherent structure energy on cooling does reflect changes in the local structure. In this section we look in detains these differences, with the aim of characterizing in a microscopic way the differences in states with the same inherent structure but different density of states.

Fig. 3 shows the AB radial distribution function in equilibrium as evaluated in the IS configurations, i.e. once the thermal distortion has been subtracted. The structure of the liquid changes in a very minor way and can be visualized only with a very fine resolution. The net effect of cooling in $`g_{AB}(r)`$ appears as a shift of less than 0.001 in the average interatomic distance.

To emphasize the $`T`$-changes we show in Fig. 4 the relation between the basin potential energy and $`g(r)`$ for a fixed $`r`$ value. The data confirm that under crunch, the local liquid order is different from the equilibrium ones.

Fig.5 shows the enlargements of $`g(r)`$ in equilibrium and contrast them with the one during aging at the same $`e_{IS}`$ value. Again, we note that the quench configuration is closer to the equilibrium one as compared with the crunch one.

## 4 LW model

In the case of the LW model for ortotherphenyl (Wahnström and Lewis 1993). , all results refer to a quench case.

### 4.1 Energies

The $`T`$-dependence of the inherent structure energy for the LW model is shown in Fig.6, together with the time evolution under aging. As for the BMLJ model, the aging dynamics is characterized by a slow progressive reduction of the IS energy.

### 4.2 Basins’ curvatures

As in the BMLJ case, the basin curvature is correlated with the basin’s depth. The $`T`$-dependence of the average frequency in equilibrium and the $`t_w`$ dependence in aging is shown in Fig.7. The frequency decreases on cooling or on aging as in the BMLJ case.

The $`e_{IS}`$ dependence of the local curvatures in equilibrium and in aging is shown in Fig.8. We note that for the case of the LW potential, in the time window explored by the numerical simulation and for the chosen $`T_f`$, the basins explored during aging coincides with the basins explored in equilibrium.

## 5 Conclusions

In this manuscript we have compared the properties of configuration space explored in equilibrium and in out of equilibrium conditions, with the aim of deepening our understanding of the physical mechanisms behind the aging process in disordered materials.

A careful analysis of the relations between curvature and depth of the potential energy basins reveals that basins which are not statistically explored in equilibrium are visited during the aging dynamics, especially during the dynamics following a crunch. In quenches, such conditions appears to hold also in the case of deep quench depth. At long times, the hypothesis that during the aging process the system explores the same part of the configuration space does not seem to hold any longer, at least on the scale of present day numerical calculations.

The differences in the basins appear to be located in the region of very small frequencies and may not be clearly seen if the average frequency is used as indicator of the basin curvature. We have shown here that the average of the logarithm of the frequency (a quantity which weigh more the very low frequency spectrum) is indeed a better indicator. Such quantity is important since it quantifies the basin’s vibrational free energy. On the other hand, the very low frequency part of the spectrum is very sensitive to size effects. It could be that spurious localized modes are stabilized by the limited size of the simulated system. This call for a size-dependence detailed analysis of the relation between curvature and depth in equilibrium and in aging.

Finally, we recall that the hypothesis that the part of configurations space explored in aging and in equilibrium are similar is an important element to be able to predict the value of the internal temperature of the system and the associated response of the aging system to an external perturbation. Indeed, in the case where such detailed comparison was performed, the basin shape and curvature in aging and equilibrium do coincide. Unfortunately, no evaluation of the internal temperature has been performed yet for $`t_w`$ values where the possibility of predicting the value of $`T_{int}`$ should fail. It would be very interesting to perform such accurate study (i.e. the comparison between predictions and numerical calculations) in the near future, in particular for crunches, where the breaking of the assumption is apparent even at short waiting times.

## Acknowledgements

We are supported by INFM PRA-HOP, INFM Parallel Computing Initiative and MURST-COFIN2000.


format: text

text:
# I Introduction

## I Introduction

The hadronic $`B`$ decays $`\overline{B}^0D^0\eta `$, $`\overline{B}^0D^0\eta `$, $`\overline{B}^0D^0\eta ^{}`$, $`\overline{B}^0D^0\eta ^{}`$ are of considerable interest both theoretically and experimentally. Measuring the branching ratios of these processes may shed light on the nature of the $`\eta `$$`\eta ^{}`$ mixing. On the experimental side, only limits on these decays are known at present from CLEO cleo . A theoretical estimate of $`\overline{B}^0D^{()0}\eta `$, $`\overline{B}^0D^{()0}\eta ^{}`$ is important also because these processes are a background for the semi–inclusive processes $`\overline{B}^0X_s\eta `$, $`\overline{B}^0X_s\eta ^{}`$. For a recent experimental analysis, see babarnoi .

The processes $`\overline{B}^0D^{()0}\eta `$, $`\overline{B}^0D^{()0}\eta ^{}`$ proceed via the internal spectator diagram shown in Fig. 1. They are colour suppressed since the colour of the quarks from the virtual $`W`$ has to match the colour of the produced $`c`$ quark from $`bc`$ and the spectator antiquark coming from the $`B`$ meson. Here we consider the quark component $`\eta _q`$ of the $`\eta `$ particle in the so called quark-flavour basis feldman :

$$\left(\begin{array}{c}\eta \\ \eta ^{}\end{array}\right)=\left(\begin{array}{cc}\mathrm{cos}\varphi & \mathrm{sin}\varphi \\ \mathrm{sin}\varphi & \mathrm{cos}\varphi \end{array}\right)\left(\begin{array}{c}\eta _q\\ \eta _s\end{array}\right)$$
(1)
where the Fock state description of $`\eta _{q,s}`$ is:

$`|\eta _q`$ $`=`$ $`\mathrm{\Psi }_q{\displaystyle \frac{|u\overline{u}+d\overline{d}}{\sqrt{2}}}+\mathrm{}`$ (2)
$`|\eta _s`$ $`=`$ $`\mathrm{\Psi }_s|s\overline{s}+\mathrm{}`$ (3)
The ellipses indicate higher Fock states, like glue states, while the wave functions $`\mathrm{\Psi }_{q,s}`$ are the light-cone wave functions of the parton states. Our assumption is that the $`\eta _q`$ component is the one to be considered in the factorization diagram resulting from Fig. 1 (shrinking to a point the $`W`$ propagator). We will exploit the relation between $`\eta `$, $`\eta ^{}`$ and $`\eta _q`$ to compute the amplitude for the process $`\overline{B}^0D^{()0}\eta ^{}`$.

We recall the form factor parameterization:

$`\eta (q_\eta )|V^\mu (q)|\overline{B}^0(p)`$ $`=`$ $`\left[(p+q_\eta )^\mu +{\displaystyle \frac{m_\eta ^2m_B^2}{q^2}}q^\mu \right]F_1(q^2)`$ (4)
$``$ $`\left[{\displaystyle \frac{m_\eta ^2m_B^2}{q^2}}q^\mu \right]F_0(q^2),`$
with $`F_1(0)=F_0(0)`$. Moreover we need the matrix elements:

$`\mathrm{VAC}|V^\mu |D(p)`$ $`=`$ $`if_Dp^\mu `$ (5)
$`\mathrm{VAC}|V^\mu |D^{}(ϵ,p)`$ $`=`$ $`ϵ^\mu m_D^{}f_D^{},`$ (6)
and we will use the values neub  $`f_D=200`$ MeV, $`f_D^{}=230`$ MeV.

Adopting the Bauer-Stech-Wirbel approach to factorization bsw  we can consider the diagram in Fig. 1 describing the $`\eta `$ coupling through its $`d\overline{d}`$ component. The amplitude we need to compute is:

$`G_{\overline{B}^0D^0\eta }`$ $`=`$ $`\eta D^0|H^{\mathrm{BSW}}|\overline{B}^0=\mathrm{cos}\varphi \eta _qD^0|H^{\mathrm{BSW}}|\overline{B}^0`$ (7)
$`=`$ $`{\displaystyle \frac{G_F}{\sqrt{2}}}V_{cb}V_{ud}^{}a_2f_DF_0^{(B\eta _q)}(m_D^2)(m_B^2m_{\eta _q}^2){\displaystyle \frac{\mathrm{cos}\varphi }{\sqrt{2}}},`$
where $`\varphi `$ is determined in feldman  to be $`\varphi =39.3^{}`$ and $`a_2`$ is the phenomenological coefficient of the Bauer-Stech-Wirbel effective Hamiltonian $`a_2=0.29`$. It is significantly larger than the leading–order result $`a_2^{\mathrm{L}O}0.12`$ corresponding to naive factorization. Such a large correction is not surprising. The QCD factorization formula of buch  can not be used to compute rigorously $`a_2`$ in these decays, however a rough estimate in the limit $`m_cm_b`$ gives $`a_20.25`$ Neubert:2000ag . The factor $`1/\sqrt{2}`$ in (7) accounts for the fact that only the $`d\overline{d}`$ component enters in our diagram. The strange quark component is absent.

The amplitude $`G_{\overline{B}^0D^0\eta ^{}}`$ is then connected to (7) through:

$$G_{\overline{B}^0D^0\eta ^{}}=\mathrm{tan}\varphi G_{\overline{B}^0D^0\eta }.$$
(8)
In an analogous way we find:

$`G_{\overline{B}^0D^0\eta }`$ $`=`$ $`\eta D^0|H^{\mathrm{BSW}}|\overline{B}^0=\mathrm{cos}\varphi \eta _qD^0|H^{\mathrm{BSW}}|\overline{B}^0`$ (9)
$`=`$ $`{\displaystyle \frac{G_F}{\sqrt{2}}}V_{cb}V_{ud}^{}a_2f_D^{}m_D^{}F_1^{(B\eta _q)}(m_D^{}^2)ϵ(p+q_{\eta _q}){\displaystyle \frac{\mathrm{cos}\varphi }{\sqrt{2}}}.`$
The term $`ϵ(p+q_\eta )`$ is easily evaluated summing over the polarisations of $`D^{}`$ once the square modulus of this amplitude is considered.

## II NS model

In order to get an idea about the amount of model dependence, we will consider two different models for the form factors. We start by considering a model by Neubert and Stech (NS) neub  which is a parameterization based on simple assumptions, in reasonable agreement with the available data and with most theoretical predictions. In the case of heavy-to-light transitions, the form factors no longer obey the symmetry relations valid in the heavy-to-heavy decays. One can account for this by introducing, for each form factor, a function $`\xi _i(w)`$ replacing the Isgur-Wise function:

$`F_1(q^2)`$ $`=`$ $`{\displaystyle \frac{m_B+m_M}{2\sqrt{m_Bm_M}}}\xi _1(w),`$ (10)
$`F_0(q^2)`$ $`=`$ $`{\displaystyle \frac{2\sqrt{m_Bm_M}}{m_B+m_M}}{\displaystyle \frac{w+1}{2}}\xi _0(w),`$ (11)
where

$$w=v_Bv_M=\frac{m_B^2+m_M^2q^2}{2m_Bm_M}.$$
(12)
The maximum value of $`w=w_{\mathrm{max}}`$ is obtained at $`q^2=0`$. For an estimate of the functions $`\xi _i(w)`$ a simple pole model is used:

$`\xi _1(w)`$ $`=`$ $`\sqrt{{\displaystyle \frac{2}{w+1}}}{\displaystyle \frac{1}{1+r}}{\displaystyle \frac{w_{\mathrm{max}}w_1}{ww_1}},`$ (13)
$`\xi _0(w)`$ $`=`$ $`\sqrt{{\displaystyle \frac{2}{w+1}}}{\displaystyle \frac{1}{1+r{\displaystyle \frac{w1}{w_{\mathrm{max}}1}}}}`$ (14)
where

$$r=\frac{(m_Bm_V)^2}{4m_Bm_V}\left(1+\frac{4m_Bm_V}{M_0^2(m_Bm_V)^2}\right),$$
(15)
where $`V`$ is the vector meson in the same doublet as the scalar meson involved in the decay described by the $`F_0`$ form factor (for the $`\eta `$ this corresponds to $`\omega `$). Moreover

$$w_1=\frac{m_B^2+m_M^2M_1^2}{2m_Bm_M},$$
(16)
where $`M_i`$ is the mass of the nearest resonance with the appropriate spin-parity quantum numbers (for $`M_0=M_B=5.754`$ GeV it is the $`0^+`$ pole and for $`M_1=M_B=5.325`$ GeV it is the $`1^{}`$ pole).

These form factors satisfy the relation

$$F_1(0)=F_0(0),$$
(17)
and scale as $`\xi _i(w)w^{3/2}`$ for large $`w`$, in accordance with the scaling rules obtained in Ali:1994vd .

The $`B\eta `$ form factors we obtain in this model are given in Table 1. The results for the branching ratios are given in Table 2. We also compute the corresponding $`B\pi `$ decays in order to allow a comparison with experimental data and theoretical predictions that is independent of the assumptions concerning the $`\eta `$ mixing.

## III CQM model

The computation of the form factor $`F_0`$ and $`F_1`$ can be carried out with the aid of a Constituent Quark Model (CQM) cqm . This model has shown to be particularly suitable for the study of heavy meson decays. Since its Lagrangian describes Feynman rules for the vertices (heavy meson)-(heavy quark)-(light quark) rev , transition amplitudes are computable via simple constituent quark loop diagrams in which mesons appear as external legs. We determine $`F_0^{B\eta _8}`$ considering the $`\eta _8`$ as a Goldstone boson. Therefore one has to take into account the relations between $`\eta _8`$ and $`\eta _q`$ for the process under study. To begin we observe that:

$$q^\mu \eta _8|V^\mu |\overline{B}^0=(m_B^2m_\eta ^2)F_0^{(B\eta _8)},$$
(18)
where $`q^\mu `$ is the momentum carried by the current $`V^\mu `$. For our purposes $`\eta `$ is essentially the Goldstone octet component $`\eta _8`$:

$$\eta |V^\mu |\overline{B}^0\eta _8|V^\mu |\overline{B}^0.$$
(19)
On the other hand the amplitude related to the process in Fig.1 is:

$$\eta |V^\mu |\overline{B}^0=\mathrm{cos}\varphi \eta _q|V^\mu |\overline{B}^0,$$
(20)
if we also take into account a small mixing angle, $`\theta _8`$ in the notations of feldman , we write:

$$\eta _8|V^\mu |\overline{B}^0=\frac{\mathrm{cos}\varphi }{\mathrm{cos}\theta _8}\eta _q|V^\mu |\overline{B}^0.$$
(21)
Again we stress that here we are neglecting the singlet component and we consider that the $`B\eta `$ matrix element which enters in the factorization diagram is the one containing only the $`d\overline{d}`$ components. In other terms we can write:

$$q^\mu \eta _8|V^\mu |\overline{B}^0=\frac{\mathrm{cos}\varphi }{\mathrm{cos}\theta _8}(m_B^2m_{\eta _q}^2)F_0^{(B\eta _q)},$$
(22)
in such a way that:

$$F_0^{(B\eta _q)}=\frac{\mathrm{cos}\theta _8}{\mathrm{cos}\varphi }\frac{(m_B^2m_\eta ^2)}{(m_B^2m_{\eta _q}^2)}F_0^{(B\eta _8)},$$
(23)
Notice that $`\theta _8`$ is determined in feldman  to be $`\theta _8=21.0^{}`$. For the $`F_1`$ form factor we obtain:

$$F_1^{(B\eta _q)}(m_D^{}^2)=\frac{ϵ(p+q_{\eta _8})}{ϵ(p+q_{\eta _q})}\frac{\mathrm{cos}\theta _8}{\mathrm{cos}\varphi }F_1^{(B\eta _8)}(m_D^{}^2).$$
(24)
We compute the $`\eta `$ particle via its $`\eta _8`$ component with a Goldstone like coupling in the meson-quark loop computation. The interaction Lagrangian generating this coupling has the form:

$$\overline{\psi }_a\gamma 𝐀_{ab}\gamma _5\psi _b,$$
(25)
where the indices $`a,b`$ are the $`u,d,s`$ quark indices, $`\psi `$ being a triplet of flavour-$`SU_3`$. The structure of the $`𝐀`$ matrix is:

$$𝐀_\mu =\frac{i}{2}(\xi _\mu \xi ^{}\xi ^{}_\mu \xi ),$$
(26)
and $`\xi =e^{i\pi /f}`$. The $`\pi `$ matrix has the well known form:

$$\pi =\left(\begin{array}{ccc}\frac{\pi ^0}{\sqrt{2}}+\frac{\eta _8}{\sqrt{6}}& \pi ^+& K^+\\ \pi ^{}& \frac{\pi ^0}{\sqrt{2}}+\frac{\eta _8}{\sqrt{6}}& K^0\\ K^{}& \overline{K^0}& \frac{2}{\sqrt{6}}\eta _8\end{array}\right)$$
(27)
When the $`d\overline{d}`$ component is taken into account the Feynman rule for the vertex is $`q^\mu \gamma _\mu \gamma _5/(\sqrt{6}f_\pi )`$. We isolate three contributions to the $`F_{0,1}`$ form factors: a non-derivative contribution, a direct and a polar contribution noi ; see Fig. 2. The direct and the non-derivative contributions are represented in the CQM model by loop diagrams in which the internal lines are heavy quark and light quark propagators, while the external legs contain the incoming meson and the vector current. In the non-derivative diagram the $`\eta `$ external leg is attached at the same vertex (heavy quark)-(light quark)-(current) due to the structure of the CQM interaction Lagrangian. In the direct diagram the $`\eta `$ is instead attached to the light quark internal line (see noi ). The polar diagram allows for an intermediate polar state separating the current insertion from the loop diagram. Each of these diagrams contributes to the form factors $`F_{0,1}`$ in a calculable way. Summing up these contributions and imposing the condition $`F_1(0)=F_0(0)`$, which eliminates the spurious singularity in the form factor decomposition, one obtains:

$$F_i(q^2)=\alpha _j(q^2)(F_i^{\mathrm{nd}}+F_i^{\mathrm{dir}}(q^2)+F_i^{\mathrm{pol}}(q^2)),$$
(28)
where $`i=0,1`$. The functions $`\alpha _j(q^2)`$, $`j=0,1`$ are needed to impose the condition $`F_1(0)=F_0(0)`$. This is because the sum in parenthesis reproduces the form factor on the lhs of the previous equation only at $`q_{\mathrm{max}}^2`$, i.e., $`\alpha (q_{\mathrm{max}}^2)1`$. Following noi  we choose:

$$\alpha _j(q^2)=1+\alpha _j\frac{q^2m_B^2m_\eta ^2}{2m_B\mathrm{\Lambda }_\chi }=1\alpha _jE_\eta /\mathrm{\Lambda }_\chi ,$$
(29)
where $`\mathrm{\Lambda }_\chi =1`$ GeV. In such a way we enforce that $`\alpha (q_{\mathrm{max}}^2)1`$. $`E_\eta `$ is the $`\eta `$ energy in the $`B`$ rest frame. We borrow $`\alpha _0`$ from noi  where we had $`\alpha _0=0.27`$. The form factor used in the computation of the branching ratios are of course those given in (28). The non-derivative (nd) contributions to the form factors are given by:

$$F_0^{\mathrm{nd}}=\frac{f_B}{\sqrt{6}f_\pi }F_1^{\mathrm{nd}}=\frac{f_B}{2\sqrt{6}f_\pi }.$$
(30)
The direct contributions are:

$`F_1^{\mathrm{Dir}}(q^2)`$ $`=`$ $`{\displaystyle \frac{S_1+S_2}{2}}`$ (31)
$`F_0^{\mathrm{Dir}}(q^2)`$ $`=`$ $`{\displaystyle \frac{S_1S_2+A(q^2)(S_1+S_2)}{2A(q^2)}},`$ (32)
where:

$$A(q^2)=\frac{(m_B^2m_\eta ^2)}{q^2},$$

and:

$`S_1`$ $`=`$ $`{\displaystyle \frac{2}{\sqrt{6}f_\pi }}\sqrt{Z_Hm_H}\left({\displaystyle \frac{mm_\eta ^2}{m_H}}Z+\left({\displaystyle \frac{2m}{m_H}}vq_\eta {\displaystyle \frac{m_\eta ^2}{m_H}}\right)\mathrm{\Omega }_1{\displaystyle \frac{2vq_\eta }{m_H}}\mathrm{\Omega }_4{\displaystyle \frac{2m_\eta }{m_H}}\mathrm{\Omega }_6\right)`$
$`S_2`$ $`=`$ $`{\displaystyle \frac{2}{\sqrt{6}f_\pi }}\sqrt{Z_Hm_H}(m^2Z2m\mathrm{\Omega }_1m_\eta \mathrm{\Omega }_2+2\mathrm{\Omega }_3+\mathrm{\Omega }_4\mathrm{\Omega }_5).`$
In these expressions:

$$vq_\eta =\frac{m_B^2+m_\eta ^2q^2}{2m_H}.$$
(33)
The symbols $`Z,Z_H,\mathrm{\Omega }_i`$ indicate integral expression listed in rev . These are functions of $`\mathrm{\Delta }_H`$, the difference between the heavy meson and the heavy quark constituent mass $`m_Hm_Q`$, $`\mathrm{\Delta }=\mathrm{\Delta }_Hvq_\eta `$ and $`vq_\eta `$, $`v`$ being the four velocity of the $`B`$ meson. $`q^\mu `$ is the momentum carried by the current insertion in the CQM loop diagram, $`m`$ is the constituent light quark mass, $`m=300`$ MeV. The suffix $`H`$ is referred to the $`H=(0^{},1^{})`$ multiplet predicted by heavy quark effective theory.

The polar contributions are:

$`F_1^{\mathrm{pol}}(q^2)`$ $`=`$ $`{\displaystyle \frac{\widehat{F}g}{\sqrt{6}f_\pi \sqrt{m_B}}}{\displaystyle \frac{1}{1q^2/m_B^{}^2}}`$ (34)
$`F_0^{\mathrm{pol}}(q^2)`$ $`=`$ $`{\displaystyle \frac{1}{m_B^2m_\eta ^2}}\left({\displaystyle \frac{hm_\eta \sqrt{m_B}\widehat{F}^+}{\sqrt{6}f_\pi }}\right){\displaystyle \frac{1}{1q^2/m_B^{}^2}},`$ (35)
where $`B^{}`$ is the $`0^+`$ state of the multiplet $`S=(0^+,1^+)`$ while $`B^{}`$ is the $`1^{}`$ state of the $`H`$ multiplet. $`h`$ is the coupling constant $`HS\pi `$ and in CQM one finds $`h=0.76\pm 0.13`$ while $`g=0.46\pm 0.04`$, from $`HH\pi `$. We use the central values for all these constants and the value $`\mathrm{\Delta }_H=0.4`$ GeV. Variations of $`\mathrm{\Delta }_H`$ in the range $`0.30.5`$ GeV induce few $`\%`$ variations in the form factors (see Table 2). $`\widehat{F}`$ and $`\widehat{F}^+`$ are defined as the decay constants of the $`H`$ and $`S`$ heavy mesons respectively.

Note that the form factors calculated in this section already include the factor taking into account that the process we are considering only involves the $`d\overline{d}`$ component of $`\eta `$ through the numerical coefficient of the Feynman rule. Therefore when combining (7) and (23) one has to be careful to remove a factor $`1/\sqrt{2}`$ from (7):

$$G_{\overline{B}^0D^0\eta }=\frac{G_F}{\sqrt{2}}V_{cb}V_{ud}^{}a_2f_DF_0^{(B\eta _8)}(m_D^2)(m_B^2m_\eta ^2)\mathrm{cos}\theta _8$$
(36)
in terms of the $`\eta _8`$ state. From (24) and (9) we get in a similar way:

$$G_{\overline{B}^0D^0\eta }=\frac{G_F}{\sqrt{2}}V_{cb}V_{ud}^{}a_2f_D^{}m_D^{}F_1^{(B\eta _8)}(m_D^{}^2)ϵ(p+q_\eta )\mathrm{cos}\theta _8.$$
(37)
Numerical results are given in Table 1 and 2. The two models give similar results for the form factors and branching ratios. The amount of model dependence due to the form factors is difficult to quantify. By comparing the two models one can conclude that they affect the branching ratios at the level of $`\pm 0.1\times 10^4`$. The $`a_2`$ coefficient is considered fixed at the value 0.29 in the calculation. As it enters as $`a_2^2`$ in the branching ratio, by changing its value to 0.25 may add to the error on the branching ratio up to $`\pm 0.2\times 10^4`$. However if the quark picture used in the calculation is correct such a large variation is not expected as $`\overline{B}^0D^0\eta `$ is in this case on the same footing as $`\overline{B}^0D^0\pi `$.

## IV Conclusion

In order to test the idea of relating the $`\eta `$ and $`\eta ^{}`$ using the quark–flavour basis in the decays considered in this work, experimental data will have to measure the ratio of Br’s $`\eta ^{}/\eta `$. If different from what estimated, it probably means that the quark mechanism (see Fig. 1) is not sufficient to explain the decays into $`\eta ^{}`$. This would suggest to look at some other mechanism, like the gluon anomaly explored in atw  where the large observed production of $`\eta ^{}`$ in $`BX_s\eta ^{}`$ decays has been studied. The mechanism suggested there is based on the subprocess $`bsg^{}s\eta ^{}g`$ where the virtual gluon emerging from the standard model penguin couples to $`\eta ^{}`$ via a gluon anomaly vertex $`g^{}g\eta ^{}`$. The possibility that the other gluon g is emitted by the light quark inside a $`B`$ meson has been examined in ahmadi . The $`BK\eta ^{}`$ decay has been recently examined in the context of perturbative QCD in sanda . In the case of $`BD\eta ^{}`$ studied here it is more difficult to imagine some gluonic mechanism of the kind described before. Anyway we cannot exclude this possibility since the $`\eta ^{}`$ can have a large glue component ahmadi ; see also norge . Under the assumptions considered in this work $`\overline{B}^0D^0\eta `$ is similar to $`\overline{B}^0D^0\pi `$. The calculation is related to the one for $`\overline{B}^0D^0\eta ^{}`$ using the quark-flavour basis. We have checked that two different ways of computing the form factors give similar results.

## Note added in proofs

Very recently there was a phenomenological determination of the $`a_2`$ coefficient Cheng:2001sc  based on new data from Belle newdata . This seems to be in contrast with all previous determinations based on the factorization approach buch  addressing the problem of final–state interaction in the $`BD\pi `$ channel. Such a problem could arise also in this analysis where the $`BD\eta `$ is considered as the starting point, and the factorization value for $`a_2`$ is assumed. We will consider this aspect in future work.

## Acknowledgements

We would like to thank R. Gatto and A. Hicheur for useful discussions. ADP acknowledges support from EU-TMR programme, contract CT98-0169. Institut de Physique Nucléaire de Lyon (IPN Lyon) is UMR 5822.


text:
# Photoinduced Changes of Reflectivity in Single Crystals of 𝑌𝐵𝑎₂𝐶𝑢₃𝑂_6.5 (Ortho II)

## Abstract

We report measurements of the photoinduced change in reflectivity of an untwinned single crystal of $`YBa_2Cu_3O_{6.5}`$ in the ortho II structure. The decay rate of the transient change in reflectivity is found to decrease rapidly with decreasing temperature and, below $`T_c`$, with decreasing laser intensity. We interpret the decay as a process of thermalization of antinodal quasiparticles, whose rate is determined by an inelastic scattering rate of quasiparticle pairs.

Much of current researchon cuprate superconductivity focuses on the low-energy electronic structure. Perhaps the most pressing question is whether the pseudogap is a precursor of the superconducting gap, or has an entirely different origin. Another key question is whether quasiparticles live in a two-dimensional environment, or one with pronounced one-dimensional properties on short length and time scales.

The relaxation dynamics of non-equilibrium quasiparticles can be extremely sensitive to such local electronic properties . A promising technique to study quasiparticle relaxation is optical pump and probe spectroscopy. In this type of experiment a pulse of light stimulates a non-equilibrium state with excess energy and quasiparticle density. A weaker pulse, which is delayed relative to the pump, probes the nonequilibrium state by detecting changes in the absorptance or reflectance. Dynamics can be measured with sub-picosecond time resolution while anisotropy can be explored by varying the polarization of the pump and probe beams.

In this paper we report pump and probe measurements performed on an untwinned single crystal of $`YBa_2Cu_3O_{6.5}`$ (Ortho II). In this material the charge reservoir layer is an ordered array of alternating filled and empty copper oxygen chains. Together with $`YBa_2Cu_4O_8`$ it is one of two underdoped materials in which doping does not introduce disorder. We find, in contrast with an earlier report  on underdoped films, that the decay rate of the photoinduced change in reflectivity depends strongly on the excitation density and temperature. Because it vanishes in the limit that both tend to zero, we interpret the decay rate as an inelastic scattering rate of quasiparticle pairs. We present evidence that the decay corresponds to a process of thermalization, rather than recombination, of the nonequilibrium quasiparticle population. The dynamics of thermalization change suddenly upon crossing from the pseudogap regime to the superconducting state.

We measured the decay of the change in reflectivity at a photon energy of 1.5 eV due to excitation with photons of the same energy. The optical pulses, produced by a mode-locked Ti:Sapphire laser, have duration of 100 fs, repetition rate 12 ns, and center wavelength 800 nm. Both the pump and probe beams were focused onto the sample with a 20 cm focal length lens, yielding a spot size of 75 $`\mu `$ diameter. A photoelastic modulator varied the pump intensity at 100 kHz and a vibrating mirror oscillating at 40 Hz varied the time delay between pump and probe. This double modulation provides sensitivity to the fractional change in reflectivity $`\mathrm{\Delta }R/R`$ of $`10^7`$.

Fig. 1a shows $`\mathrm{\Delta }R/R`$ due to absorption of a pump pulse of energy 0.07 nJ vs. time ($`t`$), at several temperatures. These reflectivity transients were measured with pump and probe polarization each aligned perpendicular to the chain axis but the decay rate is the same when measured with polarization along the chains. The amplitude of $`\mathrm{\Delta }R/R`$ increases and the decay rate decreases as the temperature is lowered. Fig. 1b shows the temperature dependence of the maximum amplitude $`(\mathrm{\Delta }R/R)_{max}`$. The reflectivity change is smaller when measured with pump and probe parallel to the chains, particularly above $`T_c`$.

The $`T`$ dependence shown in Fig. 1 is consistent with earlier reports of $`\mathrm{\Delta }R/R`$ in underdoped (twinned) thin films excited at 1.5 eV and probed at either near-infrared  or far-infrared  frequencies. These data have revealed a strong correlation of $`\mathrm{\Delta }R/R`$ with the strength of spin and current fluctuations at the characteristic energy scale for superconductivity in the cuprates. As emphasized in Ref., the increase of $`\mathrm{\Delta }R/R`$ with decreasing $`T`$ tracks the growth of a spin resonance peak seen by neutron scattering  and a gap-like feature seen in optical conductivity . Furthermore, $`\mathrm{\Delta }R/R`$ shows the same behavior with varying carrier concentration as do these spectroscopic features, namely abrupt onset at $`T_c`$ for near optimal doping  and a slow onset above $`T_c`$ for underdoped materials . Such considerations suggest that 1.5 eV photons ultimately generate low-energy excitations, despite the fact that they initially create an electron-hole pair whose energy is much larger than the maximum gap energy $`\mathrm{\Delta }_0`$.

It is generally agreed that the generation of low energy quasiparticles occurs via a rapid cascade process, in which the parent electron-hole pair creates a large number of excitations at much lower energy. What has remained less clear is the mechanism by which this nonequilibrium population of quasiparticles causes a change in reflectivity in the visible region of the spectrum. However, recent measurements  of the photoinduced change in far-IR reflectivity point to a mechanism for changes in reflectivity at 1.5 eV. Because it is crucial to understanding the decay of $`\mathrm{\Delta }R/R`$ with time, we describe this mechanism below.

The far-IR measurements reveal that photoexcitation shifts conductivity spectral weight from low energies (most likely the condensate contribution at zero energy) to $``$100 meV. This change corresponds to filling in the gap feature in the optical conductivity. The Kramers-Kronig relations dictate that this redistribution of spectral weight changes the real part of the dielectric function at much higher frequencies $`(\omega )`$ according to:

$$\mathrm{\Delta }ϵ_1(\omega )=\frac{8}{\omega ^4}𝑑\omega ^{}\mathrm{\Delta }\sigma _1(\omega ^{})\omega ^2$$
(1)
where $`\mathrm{\Delta }\sigma _1`$ is the change in the real part of the optical conductivity. Eq. 1 shows that $`\mathrm{\Delta }ϵ_1`$ is proportional to the second moment of the change in the conductivity, or $`\mathrm{\Sigma }^{(2)}`$. Because the spectral weight shifts from zero to $``$100 meV, we can relate the second moment to the spectral weight shift, $`i.e.`$, $`\mathrm{\Sigma }^{(2)}(100meV)^2_{0^+}𝑑\omega \mathrm{\Delta }\sigma _1(\omega )`$. Finally, using literature  values $`ϵ_1=1.1`$ and $`ϵ_2=2.7`$ at 1.5 eV, we can establish an approximate calibration between $`\mathrm{\Delta }R/R`$ and $`_{0^+}𝑑\omega \mathrm{\Delta }\sigma _1(\omega )`$.

The data shown in Fig. 2 provide a test of this mechanism for reflectivity change at 1.5 eV. We plot $`(\mathrm{\Delta }R/R)_{max}`$ at 15 K vs. excitation density expressed in units of photons per planar Cu. The right-hand scale shows the spectral weight shift that, according to the preceding analysis, generates the change in reflectivity shown on the left-hand scale. The spectral weight is normalized to the condensate spectral weight assuming a penetration depth of 1600 Å. Notice that $`\mathrm{\Delta }R/R`$ departs from linearity at $`1.5\times 10^4`$, which corresponds to a shift of 30$`\%`$ of the condensate spectral weight. This reasonable value for the onset of nonlinearity lends credibility to the calibration. It is important to note that a measurable reflectivity change occurs only if spectral weight shifts to frequencies of the order of $`\mathrm{\Delta }_0`$. If instead the spectral weight were shifted to the quasiparticle Drude peak, $`\mathrm{\Delta }R/R`$ would be smaller by the factor $`(1/\tau \mathrm{\Delta }_0)^2`$, which is $`10^6`$ for these exceptionally clean ortho II crystals .

Given this interpretation of the change in reflectivity, we now turn to the central result of this paper, namely the rate of decay of $`\mathrm{\Delta }R/R`$ measured as a function of $`T`$ and pump intensity $`I`$. We have seen in Fig. 1 that the decay of $`\mathrm{\Delta }R/R`$ becomes slower as $`T`$ decreases. Fig. 3a, a plot of the reflectivity change at 15 K created by pulses of varying intensity, shows that the decay rate also slows as the excitation density decreases.

The decays shown in Fig. 3a are non-exponential in that the instantaneous rate of decay, $`\gamma (t)dln(\mathrm{\Delta }R)/dt`$, decreases as $`t`$ increases. To characterize the intensity dependence of the decay rate in a relatively simple way, we focus on the initial decay rate, or $`\gamma (0)`$. We estimate $`\gamma (0)`$ from the slope of a linear fit to $`\mathrm{\Delta }R/R`$ between 100$`\%`$ and 75$`\%`$ of $`(\mathrm{\Delta }R/R)_{max}`$, normalized by the average value of $`\mathrm{\Delta }R(t)/R`$ in that interval. Given our time resolution, this procedure provides a good estimate of $`\gamma (0)`$ for rates slower than $`1ps^1`$. Fig. 3b shows $`\gamma (0)`$ as a function of $`(\mathrm{\Delta }R/R)_{max}`$, obtained by analysis of the curves similar to the ones in Fig. 3a. $`\gamma (0)`$ decreases linearly with $`(\mathrm{\Delta }R/R)_{max}`$ and extrapolates to a nonzero intercept as the excitation density approaches zero.

The dependence of $`\gamma (0)`$ on excitation density is consistent with two-particle kinetics, for example a process in which a pair of quasiparticles scatter off of each other. For high excitation densities the rate of such scattering is proportional to the nonequilibrium density. As the excitation density decreases, becoming smaller than the thermal population, $`\gamma (0)`$ approaches the equilibrium scattering rate at that temperature .

Fig. 4 presents the intensity dependence of the decay rate for a range of temperatures, in the form of a logarithmic plot of $`\gamma (0)`$ vs. $`T`$. Points with the same (open) symbol show $`\gamma (0)`$ measured at the same excitation density, but at different temperatures. The data sets coincide at high temperatures but diverge suddenly at $`T_c`$, showing that the dependence of $`\gamma (0)`$ on $`I`$ starts abruptly upon crossing from the pseudogap to the superconducting state. The solid symbols show the decay rate in the limit of zero excitation density, $`\gamma _{th}`$, as determined from the intercept of plots such as in Fig. 3b. The same quantity is plotted on linear axes in the inset. With decreasing temperature $`\gamma _{th}`$ begins to decrease in the normal state and is proportional to $`T^3`$ for temperatures less than $``$40 K.

Several aspects of the decay rate of $`\mathrm{\Delta }R/R`$ suggest a strong connection to quasiparticle dynamics, as well as structure in the low-energy density of states. For example, the decay rate appears sensitive to the opening of the pseudogap, decreasing near $`T^{}`$ in a manner which is remarkably similar to the nuclear spin relaxation rate in underdoped YBCO . The onset of superconductivity also affects the decay dynamics — the decay rate depends on excitation density only below $`T_c`$. Finally the $`T^3`$ dependence of $`\gamma _{th}`$ is consistent with the prediction of a theoretical analysis  of quasiparticle scattering and recombination in the presence of antiferromagnetic correlations. Although the rates were calculated within the random phase approximation, the asymptotic $`T^3`$ behavior is a universal feature of quasiparticle scattering in a d-wave superconductor .

The simplest interpretation of the connection between $`\gamma (0)`$ and quasiparticle dynamics is that $`\mathrm{\Delta }R(t)/R`$ reflects the decay of the nonequilibrium quasiparticle density. According to this view, $`\gamma (0)`$ is the rate at which quasiparticle pairs scatter into the condensate and $`\gamma _{th}`$ is its value in thermal equilibrium. Although such an interpretation appears to be consistent with the prediction of Ref. 16, this description of the decay of $`\mathrm{\Delta }R/R`$ is not internally consistent. The inconsistency arises because quasiparticle interactions alone cannot reduce either the excess particle or energy density. For quasiparticles to scatter into the condensate they must shed their excess energy, which for antiferromagnetic interactions appears as a spin fluctuation. However, in the low-energy regime ($`\omega <<J`$) the spin fluctuation is essentially the same entity as a quasiparticle-hole pair. Therefore such a ’recombination’ process actually scatters one pair into another, with no net change in the number or energy density. This is analogous to a problem in the interpretation of transport experiments, where quasiparticle scattering alone cannot reduce the total momentum .

The ’true’ recombination rate, in the sense of recovery of the of equilibrium quasiparticle density, would likely differ markedly from the reflectivity decay rates shown in Fig. 4. Relaxation of the quasiparticle density requires that excess energy flow irreversibly from the electrons to the lattice . We expect that process to be extremely slow on the time scale of quasiparticle scattering and to exhibit single rather than two-particle kinetics. This expectation is based on the physics of the ’phonon bottleneck’ in s-wave superconductors, an effect which may be even more restrictive in d-wave superconductors.

The bottleneck arises because the process in which a quasiparticle pair emits a phonon and scatters into the condensate is reversible — the emitted phonon can immediately regenerate the quasiparticle pair. The quasiparticle density relaxes at the rate at which the typical phonon energy decreases below the threshold for pair creation. Not only is this rate slower than the quasiparticle scattering rate, it is also approximately independent of the excess quasiparticle density.

The relaxation of the quasiparticle density may be even slower in the cuprates than in conventional superconductors. In s-wave superconductors the phonons decouple from the electrons when their energy becomes less than $`2\mathrm{\Delta }`$. In the d-wave superconductors, where there is no threshold for pair creation, the typical phonon energy must degrade to $`k_BT`$ in order for the electronic system to reach equilibrium. The extremely weak coupling between nodal quasiparticles and phonons compounds this effect. The decay of a quasiparticle pair into an acoustic phonon is kinematically forbidden because the quasiparticle velocity is larger than the sound velocity .

If recombination is extremely slow, why does $`\mathrm{\Delta }R/R`$ decay on a picosecond time scale? The remaining possibility is that the decay corresponds to thermalization, rather than recombination. If $`\mathrm{\Delta }R/R`$ is not a direct measure of the quasiparticle density, but proportional instead to the $`\mathrm{\Sigma }^{(2)}`$ as we have argued previously, it can decrease during thermalization, despite a constant, or even increasing, density of quasiparticles. Microwave measurements show that $`\mathrm{\Sigma }^{(2)}`$ for nodal particles is small, because their spectral weight is confined to a very narrow ($``$10 $`\mu `$eV) Drude peak . The observed $``$100 meV spectral weight shifts , and consequently large values of $`\mathrm{\Sigma }^{(2)}`$, are likely due to the antinodal quasiparticles. These considerations suggest that the decay of $`\mathrm{\Delta }R/R`$ measured at 1.5 eV tracks the conversion of quasiparticles from a nonthermal distribution, which has substantial $`\mathrm{\Sigma }^{(2)}`$, to a quasithermal distribution in which $`\mathrm{\Sigma }^{(2)}`$ is a factor $`10^6`$ smaller. Measurements in which the nonequilibrium state is probed in the frequency range of the Drude conductivity  show much longer decay times (typically milliseconds), consistent with an experiment which selectively probes the nodal quasiparticle density.

In summary, we present measurements of the photoinduced change in reflectivity in an underdoped crystal of YBCO, and discuss a new perspective on this effect. We interpret the decay rate of $`\mathrm{\Delta }R/R`$ as the inelastic scattering rate of quasiparticle pairs, which sets the time scale for thermalization of the nonequilibrium population. The thermalization rate decreases with the development of the pseudogap, becomes excitation density dependent in the superconducting state, and varies as $`T^3`$ in the low excitation limit. Coordinated measurements of photoinduced changes in optical response over a broad spectral range, from terahertz to visible frequencies, are needed to test this perspective.

We acknowledge valuable discussions with D.J. Scalapino. This work was supported under NSF-9870258, DOE-DE-AC03-76SF00098, Canadian Institute for Advanced Research and NSERC.