text:
# Phase Transition in a Traffic Model with Passing

## I Introduction

Traffic flows on single-lane roads with no passing exhibit clustering since queues of fast cars accumulate behind slow cars. These clusters form and grow even when car density is small. The initial analysis of cluster formation was carried out in the earlier days of traffic theory, and this subject continued growing ever then. If passing is introduced, the clusters may stop growing after reaching a certain size. Indeed, previous work indicated that after a transient regime a steady state is reached. The models of Refs. assume that any car in a cluster can pass the leading car and the passing rate is independent on the location of the car within the cluster. This is certainly an oversimplification of the everyday traffic scenarios. The complementary case when only the next-to-leading cars can pass is also an idealization, yet it is closer to reality. Below we show that the latter model is also richer phenomenologically as it undergoes a dynamical phase transition.

We first comment on possible theoretical approaches. A mean-field theory is the primary candidate, and we believe that it may be very good, perhaps even exact, since clustering and passing mix positions and velocities of the cars. The Boltzmann equation approach is an appropriate mean-field scheme, and in our earlier work we indeed used it. However, the present model, where only the next-to-leading car is allowed to pass, is significantly more difficult than the model where passing was possible for all cars. Indeed, it appears impossible even to write down closed Boltzmann equations for the distribution functions like $`P(v,t)`$ and $`P_m(v,t)`$, the density of all clusters moving with velocity $`v`$, and the density of clusters of $`m`$ cars, respectively. Therefore our theoretical analysis is performed in the framework of the Maxwell approach. This scheme simplifies “collision” terms by replacing the actual collision rates, which are proportional to velocity difference of collision partners, by constants. Despite this essentially uncontrolled approximation, the Maxwell approximation is very popular in kinetic theory and it has already been used in traffic.

The important feature of our model is quenched disorder, which manifests itself in the random assignment of intrinsic velocities. Road conditions (construction zones, turns, hills, etc.) present another source of quenched randomness in real driving situations , which is ignored in our model. Quenched disorder significantly affects characteristics of many-particle systems, especially in low spatial dimensions. This general conclusion applies to the present one-dimensional traffic model as we shall show below.

## II Maxwell approximation

We now formally define the model. Free cars move with quenched intrinsic velocities randomly assigned from some distribution $`P_0(v)`$. When a car or a cluster encounters a slower one, it assumes its velocity and a larger cluster is formed. In every cluster, the next-to-leading car is allowed to pass and resume driving with its intrinsic velocity. The rate of passing is assumed to be a constant. Thus clusters move and aggregate deterministically, while passing is stochastic. The system is initialized by randomly placing single cars and assigning them uncorrelated intrinsic velocities.

Within the Maxwell approach, the joint size-velocity distribution function (the density of clusters of size $`m`$ moving with velocity $`v`$) $`P_m(v,t)`$ obeys

$`{\displaystyle \frac{P_m(v,t)}{t}}`$ $`=`$ $`\gamma (1\delta _{m,1})[P_{m+1}(v,t)P_m(v,t)]`$ (1)
$`+`$ $`\gamma \delta _{m,1}[N(v,t)+P_2(v,t)]c(t)P_m(v,t)`$ (2)
$`+`$ $`{\displaystyle _v^{\mathrm{}}}𝑑v^{}{\displaystyle \underset{i+j=m}{}}P_i(v^{},t)P_j(v,t).`$ (3)
Here $`\gamma `$ is the passing rate, so terms proportional to $`\gamma `$ account for escape, while the rest describes clustering. The escape terms are the same within Boltzmann and Maxwell approaches, and they are actually exact. The collision terms are mean-field by nature, and they are different in the Boltzmann and Maxwell approaches. For instance, in the Boltzmann case, the integral term must involve $`v^{}v`$. Eqs. (1) also contain $`c(t)`$, the total cluster density

$`c(t)={\displaystyle \underset{j1}{}}{\displaystyle _0^{\mathrm{}}}𝑑vP_j(v,t),`$ (4)
and $`N(v,t)`$, the density of clusters in which the next-to-leading car has intrinsic velocity $`v`$. This $`N(v,t)`$ causes the major trouble since it cannot be expressed through $`P_j(v,t)`$. One might try to close Eqs. (1) by introducing $`F_k(v,v^{},t)`$, the density of clusters moving with the velocity $`v^{}`$ whose $`k^{\mathrm{th}}`$ car has intrinsic velocity $`v`$. Clearly, $`N(v,t)=_0^v𝑑v^{}F_2(v,v^{},t)`$, and it appears that equations for $`F_k(v,v^{},t)`$ are closed. A more careful look, however, reveals that the governing equation for $`F_2(v,v^{},t)`$ includes three-velocity correlators.

Thus, at the first sight, the Boltzmann and Maxwell approaches appear to be equally incapable of providing closed equations for the joint size-velocity distribution function. Still, the Maxwell framework has an advantage that it does provide a closed description on the level of the cluster size distribution. Indeed, integrating Eqs. (1) over velocity and defining $`P_m(t)_0^{\mathrm{}}𝑑vP_m(v,t)`$, we find that the cluster size distribution $`P_m(t)`$ obeys

$`{\displaystyle \frac{dP_m}{dt}}=\gamma [P_{m+1}P_m]cP_m+{\displaystyle \frac{1}{2}}{\displaystyle \underset{i+j=m}{}}P_iP_j`$ (5)
for $`m2`$, and

$`{\displaystyle \frac{dP_1}{dt}}=\gamma [P_2P_1+c]cP_1.`$ (6)
Besides this formal derivation of Eqs. (5)–(6) by direct integration of Eqs. (1), it is possible to obtain these equations by enumerating all possible ways in which clusters evolve. For instance, consider Eq. (6). Collisions reduce the density of single cars, and the collision rate is clearly equal to $`c(t)`$, as it is velocity-independent in the framework of the Maxwell approach. The escape term in Eq. (6) is understood by observing that the rate of return of single cars into the system is equal to

$`\gamma \left[2P_2+{\displaystyle \underset{j3}{}}P_j\right]=\gamma \left[P_2P_1+c\right].`$
Here $`P_2(t)`$ is singled out since passing transforms it into two single cars while an escape from larger clusters produces only one freely moving car.

Eqs. (5)–(6) are closed. Mathematically similar equations were investigated previously in the context of the aggregation-fragmentation model. Therefore, we merely present essential steps of the analysis. Restricting ourselves to the steady state and introducing notations $`P_m=\gamma F_m`$, $`c_{\mathrm{}}=\gamma F`$, we recast Eqs. (5)–(6) into

$$FF_m=F_{m+1}F_m+\delta _{m,1}F+\frac{1}{2}\underset{i+j=m}{}F_iF_j.$$
(7)
These equations should be solved together with the constraints $`_{m1}P_m=c_{\mathrm{}}`$ and $`_{m1}mP_m=1`$, i.e.,

$$\underset{m1}{}F_m=F,\underset{m1}{}mF_m=\gamma ^1.$$
(8)
Note that the sum $`_{m1}mP_m(t)`$ is obviously constant due to car conservation. The constant is equal to the initial concentration $`c_0`$ as cars were initially unclustered. Here and below we always choose $`c_0=1`$.

As in Ref. , we introduce the generating function

$$(z)=\underset{m1}{}(z^m1)F_m.$$
(9)
This generating function obeys

$$\frac{1}{2}^2+\frac{1z}{z}+\frac{(1z)^2}{z}F=0,$$
(10)
with the solution

$$(z)=\frac{z1}{z}\left\{1\sqrt{12zF}\right\}.$$
(11)
The steady state solution (11) exists only when the generating function is real for all the $`0z1`$. Hence, we require that $`2F1`$. Assuming that this condition is satisfied, we expand the generating function in the powers of $`z`$ to obtain the steady state concentrations:

$$F_m=\frac{(2F)^m}{2\sqrt{\pi }}\left\{\frac{\mathrm{\Gamma }\left(m\frac{1}{2}\right)}{\mathrm{\Gamma }(m+1)}2F\frac{\mathrm{\Gamma }\left(m+\frac{1}{2}\right)}{\mathrm{\Gamma }(m+2)}\right\}.$$
(12)
This solution is still incomplete as we have not yet determined $`F`$. To find $`F`$ we use the sum rules (8). The first sum rule is manifestly obeyed, while the second sum rule yields $`mF_m=d/dz|_{z=1}=1\sqrt{12F}=\gamma ^1`$. Thus, $`F=\frac{2\gamma 1}{2\gamma ^2}`$, which translates into $`c_{\mathrm{}}=11/2\gamma `$.

The steady state solution (12) exists for sufficiently high passing rates, $`\gamma \gamma _c=1`$. For $`\gamma >1`$ and large $`m`$, the steady state size distribution simplifies to

$$P_mCm^{3/2}\left[1\left(1\gamma ^1\right)^2\right]^m,$$
(13)
with $`C=(4\pi )^{1/2}\gamma ^1(\gamma 1)^2`$. Apart from a power-law prefactor, the size distribution exhibits an exponential decay, $`P_me^{m/m^{}}`$, in the large size limit. The characteristic size diverges, $`m^{}(\gamma 1)^2`$ as the passing rate approaches the critical value $`\gamma _c=1`$. In the critical case, the size distribution has a power-law form

$$F_m=\frac{3}{4\sqrt{\pi }}\frac{\mathrm{\Gamma }\left(m\frac{1}{2}\right)}{\mathrm{\Gamma }(m+2)}m^{5/2}.$$
(14)
Let now the passing rate drops below the critical value ($`\gamma <\gamma _c`$). Since $`F`$ cannot grow beyond $`F_c=1/2`$, it stays constant. Therefore, $`F_m`$ is given by the same Eq. (14) as in the critical case, and the cluster size distribution reads $`P_m=\gamma F_m`$. This implies $`c_{\mathrm{}}=\gamma /2`$, i.e., the sum rule $`P_m=c_{\mathrm{}}`$ is valid. The second sum rule is formally violated: $`mP_m=\gamma 1`$, i.e. the cluster size distribution (14) accounts only for the fraction of all the cars present in the system. The only possible explanation is the formation of an infinite cluster that contains all the excessive cars. The second sum rule then shows that $`1\gamma `$ of all the cars in the system are in this infinite cluster.

Thus within the framework of the Maxwell approximation, our traffic model displays a phase transition which separates the disordered and jammed phases. The steady state cluster concentration has different dependence on the passing rate for these two phases:

$$c_{\mathrm{}}=\{\begin{array}{cc}11/2\gamma ,\hfill & \gamma >1\text{;}\hfill \\ \gamma /2,\hfill & \gamma <1\text{.}\hfill \end{array}$$
(15)
In the disordered phase, the size distribution decays exponentially in the large size limit. In the jammed phase, $`P_m`$ has a power law tail and in addition there is an infinite cluster which contains the following fraction of cars:

$$I=\{\begin{array}{cc}0,\hfill & \gamma >1\text{;}\hfill \\ 1\gamma ,\hfill & \gamma <1\text{.}\hfill \end{array}$$
(16)
This phase transition is similar to phase transitions in driven diffusive systems without passing and to phase transitions in aggregation-fragmentation models. Also, the mechanism of the formation of the infinite cluster has a strong formal analogy to Bose-Einstein condensation .

Turning back to the joint size-velocity distribution (1), we note that the lack of an exact expression for $`N(v)`$ in terms of $`P_m(v)`$ does not mean the lack of a mean-field relation between these quantities. Indeed, the density $`N(v)`$ of clusters in which the next-to-leading car has intrinsic velocity $`v`$, can be written as

$$N(v)=_0^v𝑑v^{}\underset{j2}{}P_j(v^{})\frac{C(v)}{_v^{}^{\mathrm{}}𝑑v^{\prime \prime }C(v^{\prime \prime })}.$$
(17)
Here $`_{j2}P_j(v^{})`$ is the density of “true” clusters (i.e., freely moving cars are excluded) moving with velocity $`v^{}`$. Then, $`C(v)=P_0(v)P(v)`$ is the density of cars with intrinsic velocity $`v`$ which are currently slowed down, i.e., they are neither single cars, nor cluster leaders. Assuming that the velocities of cars inside clusters are perfectly mixed, $`C(v)/_v^{}^{\mathrm{}}𝑑v^{\prime \prime }C(v^{\prime \prime })`$ gives the probability density that the next-to-leading car in a true $`v^{}`$-cluster has the velocity $`v`$. The product form of Eq. (17) reveals its mean-field nature, which is consistent with the spirit of our theoretical approach. One can verify that Eq. (17) agrees with the sum rule $`𝑑vN(v)=_{j2}P_j`$, thus providing a useful check of self-consistency.

Although Eqs. (1) with $`N(v)`$ given by (17) seem very complex even in the steady-state regime, several conclusions can be derived without getting their complete solution. We first simplify Eqs. (1) by introducing auxiliary functions

$$Q_m(v)=_v^{\mathrm{}}𝑑v^{}P_m(v^{}).$$
(18)
By inserting $`P_m=\frac{dQ_m}{dv}`$ into the Eqs. (1), integrating resulting equations over $`v`$, and using the boundary conditions $`Q_m(v=\mathrm{})=0`$, we find

$`\gamma \left[Q_{m+1}(v)Q_m(v)\right]`$ $``$ $`cQ_m(v)+{\displaystyle \frac{1}{2}}{\displaystyle \underset{i+j=m}{}}Q_i(v)Q_j(v)`$ (19)
$`=`$ $`\delta _{m1}\gamma q(v),`$ (20)
with

$`q(v)=Q_1(v){\displaystyle _v^{\mathrm{}}}𝑑v^{}N(v^{}).`$ (21)
Eqs. (19) are almost identical to the Eqs. (5)–(6), the velocity is just a parameter. Consequently, we anticipate qualitatively similar results, $`Q_m(v)m^{3/2}e^{m/m^{}}`$, and

$$P_m(v)m^{1/2}e^{m/m^{}}.$$
(22)
with the characteristic size $`m^{}(v,\gamma )`$ dependent on both velocity and passing rate. Our more rigorous generating function analysis, performed along the lines described above, confirms the asymptotic form (22).

## III Simulations

Now let us examine what conclusions obtained within the Maxwell approach are relevant for the original model. We first re-derive the condition for the phase transition in the complete velocity-dependent form. Let us consider a system of reference with the origin moving with the slowest car. We assume that the system is sufficiently large for the slowest car to have negligible velocity. We compare the total flux of cars clustering behind this slowest car, $`mP_mv_m`$, to the rate of escape, $`\gamma `$. Here $`v_m`$ is an average velocity of a cluster of size $`m`$. When the rate of escape becomes less than the rate of accumulation of the cars, the cluster behind the slowest car (analog of the “infinite cluster” for finite systems) grows to remove the excessive cars from the system. Hence, the phase transition point $`\stackrel{~}{\gamma }_c`$ is defined as

$$\underset{m1}{}mP_mv_m=\stackrel{~}{\gamma }_c.$$
(23)
For the Maxwell model, where $`v_m=1`$ for all $`m`$, Eq. (23) reduces to $`mP_m=\gamma _c=1`$ as obtained above. Since large clusters usually form behind slow cars, $`v_m`$ is a decreasing function of the cluster size $`m`$. In particular, $`v_m`$ is always smaller than the average car velocity $`v`$, implying $`\stackrel{~}{\gamma }_c<1`$.

For a rough estimate of $`v_m`$, consider a cluster of $`m`$ cars and assume that intrinsic velocities of the cars in the cluster are independent. The leading car has the minimal velocity, so the size-velocity distribution reads

$$P_m(v)mP_0(v)\left[_v^{\mathrm{}}𝑑v^{}P_0(v^{})\right]^{m1}P_m.$$
(24)
For concreteness, let us consider intrinsic velocity distributions which behave algebraically near the lower cutoff, $`P_0(v)v^\mu `$ as $`v0`$. Then for large clusters we get

$$P_m(v)P_m\mathrm{exp}\left(mv^{\mu +1}\right).$$
(25)
This implies that the average cluster velocity $`v_m`$ scales with $`m`$ according to $`v_mm^{1/(\mu +1)}`$, and hence $`\stackrel{~}{\gamma }_cm^{\mu /(\mu +1)}P_m`$. We conclude that the phase transition does exist in the original model, although its location is shifted towards lower passing rate compared to the Maxwell model prediction. This shift is especially significant for small $`\mu `$ ($`\mu >1`$ from the normalization requirement).

To check the relevance of other predictions of the Maxwell approach, we performed molecular dynamics simulations. We place $`N=20000`$ single cars onto the ring of length $`L=N`$, so that the average car density is equal to one. Initial positions and velocities of cars were assigned randomly. We considered linear $`P_0(v)=\frac{8}{9}v`$ ($`0<v<3/2`$), exponential $`P_0(v)=e^v`$, and $`P_0(v)=(2\pi v)^{1/2}e^{v/2}`$ velocity distributions, which correspond to $`\mu =1,0,1/2`$ for the small-$`v`$ asymptotics. All these three distributions have the average velocity equal to one.

In Fig. 1, we plot $`\mathrm{ln}[m^{3/2}P_m]`$ vs. $`m`$ for the above three velocity distributions. We take $`\gamma =1`$ which, as we concluded before, lies above the phase transition point $`\stackrel{~}{\gamma }_c`$. We expect the system to be in the disordered phase with $`P_m`$ being expressed by Eq. (13). For the exponential and $`P_0(v)=(2\pi v)^{1/2}e^{v/2}`$ intrinsic velocity distributions, there is a good agreement with the prediction of the Maxwell model (13); for the linear velocity distribution, there are some deviations for small $`m`$, but for large $`m`$ the agreement is satisfactory. The slopes of the plots decrease with $`\mu `$. Taking into account that at the point of the phase transition the slope equals to zero, this qualitatively confirms that $`\stackrel{~}{\gamma }_c`$ gets smaller when $`\mu `$ decreases.

Plots of $`P_m`$ vs. $`m`$ for intrinsic velocity distributions $`P_0(v)=e^v`$ and $`P_0(v)=(2\pi v)^{1/2}e^{v/2}`$, with passing rate $`\gamma =0.005`$ well below the phase transition point, are shown in Fig. 2a and Fig. 2b, respectively. The cluster size distribution clearly consists of two regions: almost power-law tail for smaller $`m`$ and several separate peaks for larger $`m`$. These peaks correspond to the fluctuating size of the infinite cluster, while the power-law tail describes the regular part of $`P_m`$. The apparent exponent $`\tau `$ of the power-law region $`P_mm^\tau `$ slightly varies for different passing rates and $`P_0(v)`$, though it remains confined between 3/2 and 2. It is definitely different form the value 5/2, predicted by the Maxwell model (14). The measured values of $`\tau `$ would make the total amount of cars in the system divergent, $`mP_m\mathrm{}`$, so the power-law region ends with an exponential cutoff at large $`m`$.

We now comment on the relationship of our model to earlier work. On the mean-field level, our model is similar to the models of Refs.. On the level of the process, our model reminds an asymmetric conserved-mass aggregation model (ASCMAM) where clusters undergo asymmetric diffusion, aggregation upon contact, and chipping (single-particle dissociation). Of course, our model is continuum while the ASCMAM is the lattice model. More substantial difference between the two models lies in the nature of randomness – in our model intrinsic velocities are quenched random variables, while in the ASCMAM dynamics is the only source of randomness. Nevertheless, the phenomenology of the two models appears to be quite similar. In particular, the ASCMAM undergoes a phase transition, and in the jammed phase, the cluster size distribution exibits a power law decay with the exponent close to 2. We should stress that in the jammed phase, we have not reached a scale-free critical state which must have the exponent $`\tau 2`$. Maybe quenched randomness does not allow the system to organize itself into a truly normalizable critical state. Other possible explanation relies on large fluctuations in disordered systems, i.e., our system was not large enough to ensure self-averaging.

## IV Conclusion

In this paper, we have investigated the model of traffic that involves clustering and passing of the next-to-leading car. Despite the fact that it is one of the simplest (if not the simplest) possible continuous model of one-lane traffic with passing, the model has rich kinetic behavior. Depending on the passing rate $`\gamma `$ the system organizes itself either into disordered phase where density of large clusters is exponentially suppressed, or into the jammed phase, where the cluster size distribution becomes independent on $`\gamma `$ and the infinite cluster is formed. Within the framework of Maxwell approach, which plays the role of the mean-field theory in the present context, we have shown that the model admits an analytical solution. We have argued that the Maxwell approach correctly predicts the existence of the phase transition and adequately describes the properties of the disordered phase which arises when the passing rate is high. For the jammed phase, the Maxwell approach correctly predicts that the system stores excessive cars in the infinite cluster and organizes itself into some kind of a critical state. However, the Maxwell approach cannot quantitatively describe other properties of the jammed phase. It would be interesting to design a more accurate theoretical approach which would allow to probe the characteristics of the low passing rate regime analytically. Some properties of the jammed state appear similar to the properties of the jammed state of a lattice model of Ref. which includes an asymmetric lattice diffusion, aggregation, and fragmentation. It would be interesting to gain a deeper understanding of the relationship between these models, and whether the quenched disorder is the main source of difference.

We are thankful to E. Ben-Naim and S. Redner for discussions, and to NSF and ARO for support of this work.


format: text

text:
# Self-regulated hydrodynamical process in halo stars : a possible explanation of the lithium plateau

## References

Bonifacio, P., & Molaro, P. 1997, MNRAS, 285, 847

Maeder, A., & Zahn, J.-P. 1998, A&A, 334, 1000

Mestel, L. 1953, Mon. Not. R. Astron. Soc., 113, 716

Molaro, P. 1999, preprint

Richard, O. 1999, phD Thesis, University of Toulouse

Spite, M., & Spite F. 1982, A&A, 115, 357

Spite, F., Francois, P., Nissen, P.E., & Spite, M. 1996, ApJ, 408, 262

Vauclair, S. 1999, A&A, 351, 973

Vauclair, S. 2000, this conference

Von Zeipel, H. 1924, MNRAS, 84, 665

Zahn, J.-P. 1992, A&A, 265, 115


text:
# An updated analysis of 𝜀'/𝜀 in the standard model with hadronic matrix elements from the chiral quark model

## I Introduction

The violation of $`CP`$ symmetry in the kaon system (for two recent textbooks on the subject see ref. ) is parameterized in terms of the ratios

$$\eta _{00}\frac{\pi ^0\pi ^0|_W|K_L}{\pi ^0\pi ^0|_W|K_S}\text{and}\eta _+\frac{\pi ^+\pi ^{}|_W|K_L}{\pi ^+\pi ^{}|_W|K_S}.$$
(1)
Eqs. (1) can be written as

$`\eta _{00}`$ $`=`$ $`\epsilon {\displaystyle \frac{2\epsilon ^{}}{1\omega \sqrt{2}}}\epsilon 2\epsilon ^{},`$ (2)
$`\eta _+`$ $`=`$ $`\epsilon +{\displaystyle \frac{\epsilon ^{}}{1+\omega /\sqrt{2}}}\epsilon +\epsilon ^{},`$ (3)
where $`\omega =A_2/A_0`$ is the ratio between the isospin $`I=2`$ and 0 components of the $`K\pi \pi `$ amplitudes, the anomalous smallness of which is known as the $`\mathrm{\Delta }I=1/2`$ selection rule . The complex parameters $`\epsilon `$ and $`\epsilon ^{}`$ are introduced to quantify, respectively, indirect (via $`K_L`$-$`K_S`$ mixing) and direct (in the $`K_L`$ and $`K_S`$ decays) $`CP`$ violation. They are measurable quantities, and $`\epsilon `$ has been known to be non vanishing since 1964 .

The $`\mathrm{\Delta }S=1`$ effective lagrangian $`_W`$ is given by

$$_W=\underset{i}{}C_i(\mu )Q_i(\mu ),$$
(4)
where

$$C_i(\mu )=\frac{G_F}{\sqrt{2}}V_{ud}V_{us}^{}\left[z_i(\mu )+\tau y_i(\mu )\right].$$
(5)
In (5), $`G_F`$ is the Fermi coupling, the functions $`z_i(\mu )`$ and $`y_i(\mu )`$ are the Wilson coefficients and $`V_{ij}`$ the Cabibbo-Kobayashi-Maskawa (CKM) matrix elements; $`\tau =V_{td}V_{ts}^{}/V_{ud}V_{us}^{}`$. According to the standard parameterization of the CKM matrix, in order to determine $`\epsilon ^{}/\epsilon `$ , we only need to consider the $`y_i(\mu )`$ components, which control the $`CP`$-violating part of the lagrangian. The coefficients $`y_i(\mu )`$, and $`z_i(\mu )`$ contain all the dependence of short-distance physics, and depend on the $`t,W,b,c`$ masses, the intrinsic QCD scale $`\mathrm{\Lambda }_{\mathrm{QCD}}`$, the $`\gamma _5`$-scheme used in the regularization and the renormalization scale $`\mu `$.

The $`Q_i`$ in eq. (4) are the effective four-quark operators obtained in the standard model by integrating out the vector bosons and the heavy quarks $`t,b`$ and $`c`$. A convenient and by now standard basis includes the following ten operators:

$$\begin{array}{ccc}\hfill Q_1& =& \left(\overline{s}_\alpha u_\beta \right)_{\mathrm{V}\mathrm{A}}\left(\overline{u}_\beta d_\alpha \right)_{\mathrm{V}\mathrm{A}},\hfill \\ \hfill Q_2& =& \left(\overline{s}u\right)_{\mathrm{V}\mathrm{A}}\left(\overline{u}d\right)_{\mathrm{V}\mathrm{A}},\hfill \\ \hfill Q_{3,5}& =& \left(\overline{s}d\right)_{\mathrm{V}\mathrm{A}}_q\left(\overline{q}q\right)_{\mathrm{V}\mathrm{A}},\hfill \\ \hfill Q_{4,6}& =& \left(\overline{s}_\alpha d_\beta \right)_{\mathrm{V}\mathrm{A}}_q(\overline{q}_\beta q_\alpha )_{\mathrm{V}\mathrm{A}},\hfill \\ \hfill Q_{7,9}& =& \frac{3}{2}\left(\overline{s}d\right)_{\mathrm{V}\mathrm{A}}_q\widehat{e}_q\left(\overline{q}q\right)_{\mathrm{V}\pm \mathrm{A}},\hfill \\ \hfill Q_{8,10}& =& \frac{3}{2}\left(\overline{s}_\alpha d_\beta \right)_{\mathrm{V}\mathrm{A}}_q\widehat{e}_q(\overline{q}_\beta q_\alpha )_{\mathrm{V}\pm \mathrm{A}},\hfill \end{array}$$
(6)
where $`\alpha `$, $`\beta `$ denote color indices ($`\alpha ,\beta =1,\mathrm{},N_c`$) and $`\widehat{e}_q`$ are the quark charges ($`\widehat{e}_u=2/3`$, $`\widehat{e}_d=\widehat{e}_s=1/3`$). Color indices for the color singlet operators are omitted. The labels $`(V\pm A)`$ refer to the Dirac structure $`\gamma _\mu (1\pm \gamma _5)`$.

The various operators originate from different diagrams of the fundamental theory. At the tree level, we only have the current-current operator $`Q_2`$ induced by $`W`$-exchange. Switching on QCD, the gluonic correction to tree-level $`W`$-exchange induces $`Q_1`$. Furthermore, QCD induces the gluon penguin operators $`Q_{36}`$. The penguin diagrams induce different operators because of the splitting of the color quark (vector) current into a right- and a left-handed part and the presence of color octet and singlet currents. Electroweak penguin diagrams—where the exchanged gluon is replaced by a photon or a $`Z`$-boson and box-like diagrams—induce $`Q_{7,9}`$ and also a part of $`Q_3`$. The operators $`Q_{8,10}`$ are induced by the QCD renormalization of the electroweak penguin operators $`Q_{7,9}`$.

Even though the operators in eq. (6) are not all independent, this basis is of particular interest for any numerical analysis because it has been extensively used for the the calculation of the Wilson coefficients to the next-to-leading order (NLO) order in $`\alpha _s`$ and $`\alpha _e`$ , in different renormalization schemes.

In the standard model, $`\epsilon ^{}`$ can be in principle different from zero because the $`3\times 3`$ CKM matrix $`V_{ij}`$, which appears in the weak charged currents of the quark mass eigenstate, is in general complex. On the other hand, in other models like the superweak theory , the only source of $`CP`$ violation resides in the $`K^0`$-$`\overline{K}^0`$ mixing, and $`\epsilon ^{}`$ vanishes. It is therefore of great importance for the discussion of the theoretical implications within the standard model and beyond to establish the experimental evidence and precise value of $`\epsilon ^{}`$ .

The ratio $`\epsilon ^{}/\epsilon `$ (for a review see, e.g., ) is computed as

$$\text{Re}\epsilon ^{}/\epsilon =e^{i\varphi }\frac{G_F\omega }{2\left|ϵ\right|\text{Re}A_0}\text{Im}\lambda _t\left[\mathrm{\Pi }_0\frac{1}{\omega }\mathrm{\Pi }_2\right],$$
(7)
where, the CKM combination $`\lambda _t=V_{td}V_{ts}^{}`$ and, referring to the $`\mathrm{\Delta }S=1`$ quark lagrangian of eq. (4),

$`\mathrm{\Pi }_0`$ $`=`$ $`{\displaystyle \frac{1}{\mathrm{cos}\delta _0}}{\displaystyle \underset{i}{}}y_i\text{Re}Q_i_0(1\mathrm{\Omega }_{\eta +\eta ^{}}),`$ (8)
$`\mathrm{\Pi }_2`$ $`=`$ $`{\displaystyle \frac{1}{\mathrm{cos}\delta _2}}{\displaystyle \underset{i}{}}y_i\text{Re}Q_i_2,`$ (9)
and $`Q_i_I=2\pi ,I|Q_i|K`$.

We take the phase $`\varphi =\pi /2+\delta _0\delta _2\theta _ϵ=(0\pm 4)^0`$ as vanishing , and we assume everywhere that $`CPT`$ is conserved. Therefore, $`\text{Re}\epsilon ^{}/\epsilon =\epsilon ^{}/\epsilon `$.

Notice the explicit presence of the final-state-interaction (FSI) phases $`\delta _I`$ in eqs. (8) and (9). Their presence is a consequence of writing the absolute values of the amplitudes in term of their dispersive parts. Theoretically, given that in eq. (5) $`\tau 1`$, we obtain

$$\mathrm{tan}\delta _I\frac{_iz_i\text{Im}Q_i_I}{_iz_i\text{Re}Q_i_I}.$$
(10)
Finally, $`\mathrm{\Omega }_{\eta +\eta ^{}}`$ is the isospin breaking (for $`m_um_d`$) contribution of the mixing of $`\pi `$ with $`\eta `$ and $`\eta ^{}`$.

### A Preliminary remarks

The experiments in the early 1990’s  could not establish the existence of direct $`CP`$ violation because they agreed only marginally (one of them being consistent with zero ) and did not have the required accuracy. During 1999, as we shall briefly recall in the next section, the preliminary analysis of the new run of experiments  have settled on a range of values consistent with the previous NA31 result, conclusively excluding a vanishing $`\epsilon ^{}`$. On the other hand, in order to asses precisely the value we must wait for the completion of the data analysis which will improve the present accuracy by a factor 2-3.

On the theoretical side, progress has been slow as well because of the intrinsic difficulty of a computation that spans energy scales as different as the pion and the top quark masses. Nevertheless, the estimates available before 1999  pointed to a non-vanishing and positive value, with one of them  being in the ball park of the present experimental result.

We revise our 1997 estimate  of $`\epsilon ^{}/\epsilon `$ by updating the values of the short-distance input parameters —among which the improved determination of the relevant CKM entries— and by including the Gaussian sampling of the experimental input data. We also update the values and ranges of the “long-distance” model parameters (quark and gluon condensates, and constituent quark mass), by including a larger theoretical “systematic” error ($`\pm 30\%`$) in the fit of the $`CP`$ conserving $`K\pi \pi `$ amplitudes in order to better account for the error related to the truncation of the chiral expansion.

For the sake of comparison with other approaches, we give our results in terms of the so-called $`B_i`$-parameters, in two $`\gamma _5`$ regularization schemes: ’t Hooft-Veltman (HV) and Naive Dimensional Regularization (NDR).

The combined effect of the new ranges of the input parameters and $`\text{Im}\lambda _t`$ makes the central value of $`\epsilon ^{}/\epsilon `$ slightly higher than before, while the statistical analysis of the input parameters reduces the final uncertainty. For a more conservative assessment of the error, we give the full range of uncertainty obtained by the flat span of the allowed ranges of the input parameters. The result is numerically stable as we vary the renormalization scale and scheme. We conclude by briefly reviewing other estimates of $`\epsilon ^{}`$ published in the last year.

## II The experimental status

Experimentally the ratio $`\epsilon ^{}/\epsilon `$ is extracted, by collecting $`K_L`$ and $`K_S`$ decays into pairs of $`\pi ^0`$ and $`\pi ^\pm `$, from the relation

$$\left|\frac{\eta _+}{\eta _{00}}\right|^21+6\text{Re}\frac{\epsilon ^{}}{\epsilon },$$
(11)
and the determination of the ratios $`\eta _+`$ and $`\eta _{00}`$ given in (1).

The announcement last year of the preliminary result from the KTeV collaboration (FNAL)

$$\text{Re}\epsilon /\epsilon ^{}=(2.8\pm 0.41)\times 10^3,$$
(12)
based on data collected in 1996-97, and the present result from the NA48 collaboration (CERN),

$$\text{Re}\epsilon /\epsilon ^{}=(1.4\pm 0.43)\times 10^3,$$
(13)
based on data collected in 1997  and 1998 , settle the long-standing issue of the presence of direct $`CP`$ violation in kaon decays. However, a clearcut determination of the actual value of $`\epsilon ^{}`$ at the precision of a few parts in $`10^4`$ must wait for further statistics and scrutiny of the experimental systematics.

By computing the average among the two 1992 experiments (NA31 and E731 ) and the KTeV and NA48 data we obtain

$$\text{Re}\epsilon /\epsilon ^{}=(1.9\pm 0.46)\times 10^3,$$
(14)
where the error has been inflated according to the Particle Data Group procedure ($`\sigma \sigma \times \sqrt{\chi ^2/3}`$), to be used when averaging over experimental data—in our case four sets—with substantially different central values.

The value in eq. (14) can be considered the current experimental result. Such a result will be tested within the next year by the full data analysis from KTeV and NA48 and (hopefully) the first data from KLOE at DA$`\mathrm{\Phi }`$NE (Frascati); at that time, the experimental uncertainty will be reduced to a few parts in $`10^4`$.

The most important outcome of the 1999 results is that direct $`CP`$ violation has been unambiguously observed and that the superweak scenario, in which $`\epsilon ^{}`$ = 0, can be excluded to a high degree of confidence (more than 4 $`\sigma `$’s).

## III Hadronic matrix elements

In the present analysis, we use hadronic matrix elements for all the relevant operators $`Q_{110}`$ and a parameter $`\widehat{B}_K`$ computed in the $`\chi `$QM at $`O(p^4)`$ in the chiral expansion . This approach has three model-dependent parameters which are fixed by means of a fit of the $`\mathrm{\Delta }I=1/2`$ rule. Let us here review briefly the model and how the matrix elements are computed.

The $`\chi `$QM  is an effective quark model of QCD which can be derived in the framework of the extended Nambu-Jona-Lasinio model of chiral symmetry breaking (for a review, see, e.g., ).

In the $`\chi `$QM an effective interaction between the $`u,d,s`$ quarks and the meson octet is introduced via the term

$$_{\chi \text{QM}}=M\left(\overline{q}_R\mathrm{\Sigma }q_L+\overline{q}_L\mathrm{\Sigma }^{}q_R\right),$$
(15)
which is added to an effective low-energy QCD lagrangian whose dynamical degrees of freedom are the $`u,d,s`$ quarks propagating in a soft gluon background. The matrix $`\mathrm{\Sigma }`$ in (15) is the same as that used in chiral perturbation theory and it contains the pseudo-scalar meson multiplet. The quantity $`M`$ is interpreted as the constituent quark mass in mesons (current quark masses are also included in the effective lagrangian).

In the factorization approximation, the matrix elements of the four quark operators are written in terms of better known quantities like quark currents and densities. Such matrix elements (building blocks) like the current matrix elements $`0|\overline{s}\gamma ^\mu \left(1\gamma _5\right)u|K^+(k)`$ and $`\pi ^+(p_+)|\overline{s}\gamma ^\mu \left(1\gamma _5\right)d|K^+(k)`$ and the matrix elements of densities, $`0|\overline{s}\gamma _5u|K^+(k)`$, $`\pi ^+(p_+)|\overline{s}d|K^+(k)`$, are evaluated up to $`O(p^4)`$ within the model. The model dependence in the color singlet current and density matrix elements appears (via the $`M`$ parameter) beyond the leading order in the momentum expansion, while the leading contributions agree with the well known expressions in terms of the meson decay constants and masses.

Non-factorizable contributions due to soft gluonic corrections are included by using Fierz-transformations and by calculating building block matrix elements involving the color matrix $`T^a`$:

$`0|\overline{s}\gamma ^\mu T^a\left(1\gamma _5\right)u|K^+(k),\pi ^+(p_+)|\overline{s}\gamma ^\mu T^a\left(1\gamma _5\right)d|K^+(k).`$ (16)
Such matrix elements are non-zero for emission of gluons. In contrast to the color singlet matrix elements above, they are model dependent starting with the leading order. Taking products of two such matrix elements and using the relation

$$g_s^2G_{\mu \nu }^aG_{\alpha \beta }^a=\frac{\pi ^2}{3}\frac{\alpha _s}{\pi }GG\left(\delta _{\mu \alpha }\delta _{\nu \beta }\delta _{\mu \beta }\delta _{\nu \alpha }\right)$$
(17)
makes it possible to express gluonic corrections in terms of the gluonic vacuum condensate . While the factorizable corrections are re-absorbed in the renormalization of the chiral couplings, non-factorizable contributions affect explicitly the form of the matrix elements. The model thus parameterizes all amplitudes in terms of the quantities $`M`$, $`\overline{q}q`$ , and $`\alpha _sGG/\pi `$ . Higher order gluon condensates are omitted.

The hadronic matrix elements of the operatorial basis in eq. (6) for $`K\pi \pi `$ decays, have been calculated up to $`O(p^4)`$ inclusive of chiral loops . The leading order (LO) ($`O(p^0,p^2)`$) matrix elements $`Q_i_I^{LO}`$ and the NLO ($`O(p^2,p^4)`$) corrections $`Q_i_I^{NLO}`$ for final state isospin projections $`I=0,2`$ are obtained by properly combining factorizable and non-factorizable contributions and expanding the result at the given order. The total hadronic matrix elements up to $`O(p^4)`$ can then be written as:

$$Q_i(\mu )_I=Z_\pi \sqrt{Z_K}\left[Q_i_I^{LO}+Q_i_I^{NLO}\right](\mu )+a_i^I(\mu ),$$
(18)
where $`Q_i`$ are the operators in eq. (6), and $`a_i^I(\mu )`$ are the contributions from chiral loops (which include the wave-function renormalization). The scale dependence of the $`Q_i_I^{LO,NLO}`$ comes from the perturbative running of the quark condensate and masses, the latter appearing explicitly in the NLO corrections. The quantities $`a_i^I(\mu )`$ represent the scale dependent meson-loop corrections which depend on the chiral quark model via the tree level chiral coefficients. They have been included to $`O(p^4)`$ in ref.  by applying the $`\overline{MS}`$ scheme in dimensional regularization, as for the $`\chi `$QM calculation of the tree-level chiral coefficients. The wave-function renormalizations $`Z_K`$ and $`Z_\pi `$ arise in the $`\chi `$QM from direct calculation of the $`KK`$ and $`\pi \pi `$ propagators.

The hadronic matrix elements are matched—by taking $`\mu _{SD}=\mu _{LD}`$— with the NLO Wilson coefficients at the scale $`\mathrm{\Lambda }_\chi 0.8`$ ($`m_\rho `$) as the best compromise between the range of validity of chiral perturbation and that of strong coupling expansion. The scale dependence of the amplitudes is gauged by varying $`\mu `$ between 0.8 and 1 GeV. In this range the scale dependence of $`\epsilon ^{}/\epsilon `$ remains always below 10%, thus giving a stable prediction.

### A The fit of the $`\mathrm{\Delta }I=1/2`$ rule

In order to assign the values of the model-dependent parameters $`M`$, $`\overline{q}q`$ and $`\alpha _sGG/\pi `$ , we consider the $`CP`$-conserving amplitudes in the $`\mathrm{\Delta }I=1/2`$ selection rule of $`K\pi \pi `$ decays. In practice, we compute the amplitudes

$`A_0`$ $`=`$ $`{\displaystyle \frac{G_F}{\sqrt{2}}}V_{ud}V_{us}^{}{\displaystyle \frac{1}{\mathrm{cos}\delta _0}}{\displaystyle \underset{i}{}}z_i(\mu )\text{Re}Q_i(\mu )_0,`$ (19)
$`A_2`$ $`=`$ $`{\displaystyle \frac{G_F}{\sqrt{2}}}V_{ud}V_{us}^{}{\displaystyle \frac{1}{\mathrm{cos}\delta _2}}{\displaystyle \underset{i}{}}z_i(\mu )\text{Re}Q_i(\mu )_2+\omega A_0\mathrm{\Omega }_{\eta +\eta ^{}},`$ (20)
within the $`\chi `$QM approach and vary the parameters in order to reproduce their experimental values

$$A_0(K\pi \pi )=3.3\times 10^7\text{GeV}\text{and}A_2(K\pi \pi )=1.5\times 10^8\text{GeV},$$
(21)
This procedure combines a model for low-energy QCD—which allow us to compute all hadronic matrix elements in terms of a few basic parameters—with the phenomenological determination of such parameters. In this way, some shortcomings of such a naive model (in particular, the matching between long- and short-distance components) are absorbed in the phenomenological fit. As a check, we eventually verify the stability of the computed observables against renormalization scale and scheme dependence.

The fit of the $`CP`$-conserving involves the determination of the FSI phases. The absorptive component of the hadronic matrix elements appear when chiral loops are included. In our approach the direct determination of the rescattering phases gives at $`O(p^4)`$ $`\delta _020^0`$ and $`\delta _212^0`$. Although these results show features which are in qualitative agreement with the phases extracted from pion-nucleon scattering ,

$$\delta _0=34.2^0\pm 2.2^0,\delta _2=6.9^0\pm 0.2^0,$$
(22)
the deviation from the experimental data is sizeable, especially in the $`I=0`$ component. On the other hand, at $`O(p^4)`$ the absorptive parts of the amplitudes are determined only at $`O(p^2)`$ and disagreement with the measured phases should be expected. As a matter of fact, the authors of ref.  find that at $`O(p^6)`$ the absorptive part of the hadronic matrix elements are substantially modified to give values of the rescattering phases quite close to those in eq. (22). At the same time the $`O(p^6)`$ corrections to the dispersive part of the hadronic matrix elements are very small.

This result corroborates our ansatz  of trusting the dispersive parts of the $`O(p^4)`$ matrix elements while inputting the experimental values of the rescattering phases in all parts of our analysis, which amounts to taking $`\mathrm{cos}\delta _00.8`$ and $`\mathrm{cos}\delta _21`$.

Hadronic matrix elements in the $`\chi `$QM depend on the $`\gamma _5`$-scheme utilized . Their dependence partially cancels that of the short-distance NLO Wilson coefficients. Because this compensation is only numerical, and not analytical, we take it as part of our phenomenological approach. A formal $`\gamma _5`$-scheme matching can only come from a model more complete than the $`\chi `$QM. Nevertheless, the result, as shown in Fig. 2 below, is rather convincing.

By taking

$$\mathrm{\Lambda }_{\mathrm{QCD}}^{(4)}=340\pm 40\text{MeV}$$
(23)
and fitting at the scale $`\mu =0.8`$ GeV the amplitudes in eqs. (19)–(20) to their experimental values, allowing for a $`\pm `$ 30% systematic uncertainty, we find (see Fig. 1)

$$M=195_{15}^{+25}\text{MeV},\alpha _sGG/\pi =\left(330\pm 5\text{MeV}\right)^4,\overline{q}q=\left(235\pm 25\text{MeV}\right)^3.$$
(24)
in the HV-scheme, and

$$M=195_{10}^{+15}\text{MeV},\alpha _sGG/\pi =\left(333_6^{+7}\text{MeV}\right)^4,\overline{q}q=\left(245\pm 15\text{MeV}\right)^3.$$
(25)
in the NDR-scheme.

As shown by the light (NDR) and dark (HV) curves in Fig. 2, the $`\gamma _5`$-scheme dependence is controlled by the value of $`M`$, the range of which is fixed thereby. The $`\gamma _5`$ scheme dependence of both amplitudes is minimized for $`M190200`$ MeV. The good $`\gamma _5`$-scheme stability is also shown by $`\epsilon ^{}/\epsilon `$ and $`\widehat{B}_K`$. For this reason in our previous papers  we only quoted for these observables the HV results.

The fit of the amplitude $`A_0`$ and $`A_2`$ is obtained for values of the quark and gluon condensates which are in agreement with those found in other approaches, i.e. QCD sum rules and lattice, although it is fair to say that the relation between the gluon condensate of QCD sum rules and lattice and that of the $`\chi `$QM is far from obvious. The value of the constituent quark mass $`M`$ is in good agreement with that found by fitting radiative kaon decays .

In Fig. 3 we present the anatomy of the relevant operator contributions to the $`CP`$ conserving amplitudes. It is worth noticing that, because of the NLO enhancement of the $`I=0`$ matrix elements (mainly due to the chiral loops), the gluon penguin contribution to $`A_0`$ amounts to about 20% of the amplitude.

Turning now to the $`\mathrm{\Delta }S=2`$ lagrangian,

$$_{\mathrm{\Delta }S=2}=C_{2S}(\mu )Q_{S2}(\mu ),$$
(26)
where

$$C_{2S}(\mu )=\frac{G_F^2m_W^2}{4\pi ^2}\left[\lambda _c^2\eta _1S(x_c)+\lambda _t^2\eta _2S(x_t)+2\lambda _c\lambda _t\eta _3S(x_c,x_t)\right]b(\mu )$$
(27)
where $`\lambda _j=V_{jd}V_{js}^{}`$, $`x_i=m_i^2/m_W^2`$. We denote by $`Q_{S2}`$ the $`\mathrm{\Delta }S=2`$ local four quark operator

$$Q_{S2}=(\overline{s}_L\gamma ^\mu d_L)(\overline{s}_L\gamma _\mu d_L),$$
(28)
which is the only local operator of dimension six in the standard model.

The integration of the electroweak loops leads to the Inami-Lim functions  $`S(x)`$ and $`S(x_c,x_t)`$, the exact expressions of which can be found in the reference quoted, which depend on the masses of the charm and top quarks and describe the $`\mathrm{\Delta }S=2`$ transition amplitude in the absence of strong interactions.

The short-distance QCD corrections are encoded in the coefficients $`\eta _1`$, $`\eta _2`$ and $`\eta _3`$ with a common scale- and renormalization-scheme- dependent factor $`b(\mu )`$ factorized out. They are functions of the heavy quarks masses and of the scale parameter $`\mathrm{\Lambda }_{\mathrm{QCD}}`$. These QCD corrections are available at the NLO  in the strong and electromagnetic couplings.

The scale-dependent factor of the short-distance corrections is given by

$$b(\mu )=\left[\alpha _s\left(\mu \right)\right]^{2/9}\left(1J_3\frac{\alpha _s\left(\mu \right)}{4\pi }\right),$$
(29)
where $`J_3`$ depends on the $`\gamma _5`$-scheme used in the regularization. The NDR and HV scheme yield, respectively:

$$J_3^{\mathrm{NDR}}=\frac{307}{162}\text{and}J_3^{\mathrm{HV}}=\frac{91}{162}.$$
(30)
On the long-distance side, the hadronic $`B_K`$ parameter is introduced by writing the $`\mathrm{\Delta }S=2`$ matrix element as

$$\overline{K^0}|Q_{S2}(\mu )|K^0=\frac{4}{3}f_K^2m_K^2B_K(\mu ).$$
(31)
The scale- and renormalization-scheme- independent parameter $`\widehat{B}_K`$ is then defined by means of eqs. (29)–(31) as

$$\widehat{B}_K=b(\mu )B_K(\mu ).$$
(32)
By using the input values found by fitting the $`\mathrm{\Delta }I=1/2`$ rule we obtain in both $`\gamma _5`$-schemes

$$\widehat{B}_K=1.1\pm 0.2.$$
(33)
The result (33) includes chiral corrections up to $`O(p^4)`$ and it agrees with what we found in . In the chiral limit one derives a simple expression for $`\widehat{B}_K`$ (eq. (6.3) in ref. ), which depends crucially on the value of the gluon condensate. In this limit, for central values of the parameters, we obtain $`\widehat{B}_K=0.44`$ $`(0.46)`$ in the HV (NDR) scheme. A recent calculation of $`\widehat{B}_K`$, based on QCD sum rules, finds in the chiral limit and to the NLO in the $`1/N`$ expansion $`\widehat{B}_K=0.41\pm 0.09`$ . Many calculations of $`B_K`$ have been performed on the lattice (for a recent review see ). According to ref. , the analysis of ref.  presents the most extensive study of systematic errors and gives the (quenched) value $`\widehat{B}_K=0.86\pm 0.06\pm 0.14`$, which should be taken as the present reference value, while awaiting for further progress in including dynamical quarks on the lattice.

Notice that no estimate of $`\epsilon ^{}/\epsilon `$ can be considered complete unless it also gives a value for $`\widehat{B}_K`$. The case of the $`\chi `$QM, for instance, is telling insofar as the enhancement of $`B_6`$ is partially compensated for by a large $`\widehat{B}_K`$ (and accordingly a smaller $`\text{Im}\lambda _t`$).

### B The factors $`B_i`$

The factors $`B_i`$, defined as

$$B_iQ_i^{\mathrm{model}}/Q_i^{\mathrm{VSA}}$$
(34)
have become a standard way of displaying the values of the hadronic matrix elements in order to compare them among various approaches. However they must be used with care because of their dependence on the renormalization scheme and scale, as well as on the choice of the VSA parameters.

They are given in the $`\chi `$QM in table I in the HV and NDR schemes, at $`\mu =0.8`$ GeV, for the central value of $`\mathrm{\Lambda }_{\mathrm{QCD}}^{(4)}`$. The dependence on $`\mathrm{\Lambda }_{\mathrm{QCD}}`$ enters indirectly via the fit of the $`\mathrm{\Delta }I=1/2`$ selection rule and the determination of the parameters of the model.

The uncertainty in the matrix elements of the penguin operators $`Q_{58}`$ arises from the variation of $`\overline{q}q`$. This affects mostly the $`B_{5,6}`$ parameters because of the leading linear dependence on $`\overline{q}q`$ of the $`Q_{5,6}`$ matrix elements in the $`\chi `$QM, contrasted to the quadratic dependence of the corresponding VSA matrix elements. Accordingly, $`B_{5,6}`$ scale as $`\overline{q}q^1`$, or via PCAC as $`m_q`$, and therefore are sensitive to the value chosen for these parameters. For this reason, we have reported the corresponding values of $`B_{5,6}`$ when the quark condensate in the VSA is fixed to its PCAC value.

It should however be stressed that such a dependence is not physical and is introduced by the arbitrary normalization on the VSA result. The estimate of $`\epsilon ^{}`$ is therefore almost independent of $`m_q`$, which only enters the NLO corrections and the determination of $`\widehat{B}_K`$.

The enhancement of the $`Q_{5,6}`$ matrix elements with respect to the VSA values (the conventional normalization of the VSA matrix elements corresponds to taking $`\overline{q}q`$$`(0.8\mathrm{GeV})(220\mathrm{MeV})^3`$ ) is mainly due to the NLO chiral loop contributions. Such an enhancement, due to final state interactions, has been found in $`1/N_c`$ analyses beyond LO , as well as recent dispersive studies . A large-$`N_c`$ approach, based on QCD sum rules  which reproduces the electroweak $`\pi ^+`$-$`\pi ^{}`$ mass difference and the leptonic $`\pi (\eta )`$ rare decays, disagrees in the determination of $`Q_7_{0,2}`$ at $`\mu =0.8`$ GeV, due to the sharp scale dependence found for these matrix elements. Since the operator $`Q_7`$ gives a negligible contribution to $`\epsilon ^{}/\epsilon `$, we should wait for a calculation of other matrix elements within the same framework in order to asses the extent and the impact of the disagreement with the $`\chi `$QM results.

Among the non-factorizable corrections, the gluon condensate contributions are most important for the $`CP`$-conserving $`I=2`$ amplitudes (and account for the values and uncertainties of $`B_{1,2}^{(2)}`$) but are otherwise inessential in the determination of $`\epsilon ^{}`$, for which FSI are the most relevant corrections to the LO $`1/N`$ result.

## IV Bounds on $`\text{Im}\lambda _t`$

The updated measurements for the CKM elements $`|V_{ub}/V_{cb}|`$ implies a change in the determination of the Wolfenstein parameter $`\eta `$ that enter in $`\text{Im}\lambda _t`$. This is of particular relevance because it affects proportionally the value of $`\epsilon ^{}/\epsilon `$.

The allowed values for $`\text{Im}\lambda _t\eta |V_{us}||V_{cb}|^2`$ are found by imposing the experimental constraints for $`\epsilon `$, $`|V_{ub}/V_{cb}|`$, $`\mathrm{\Delta }m_d`$ and $`\mathrm{\Delta }m_s`$ which give rise to the following equations:

$$\eta \left(1\frac{\lambda ^2}{2}\right)\left\{\left[1\rho \left(1\frac{\lambda ^2}{2}\right)\right]|V_{cb}|^2\eta _2S(x_t)+\eta _3S(x_x,x_t)\eta _1S(x_c)\right\}\frac{|V_{cb}|^2}{\lambda ^8}\widehat{B}_K=\frac{|\epsilon |}{C\lambda ^{10}},$$
(35)
where

$$C=\frac{G_F^2f_K^2m_K^2m_W^2}{3\sqrt{2}\pi ^2\mathrm{\Delta }M_{LS}}\text{and}x_i=m_i^2/m_W^2,$$
(36)
and

$`\eta ^2+\rho ^2`$ $`=`$ $`{\displaystyle \frac{1}{\lambda ^2}}{\displaystyle \frac{|V_{ub}|^2}{|V_{cb}|^2}},`$ (37)
$`\eta ^2\left(1{\displaystyle \frac{\lambda ^2}{2}}\right)^2+\left[1\rho \left(1{\displaystyle \frac{\lambda ^2}{2}}\right)\right]^2`$ $`=`$ $`{\displaystyle \frac{1}{\lambda ^2}}{\displaystyle \frac{|V_{td}|^2}{|V_{cb}|^2}},`$ (38)
where $`|V_{td}|`$ is found by means of

$$\mathrm{\Delta }m_d=\frac{G_F}{24\sqrt{2}\pi ^2}|V_{td}|^2|V_{tb}|^2m_{B_d}f_{B_d}^2B_{B_d}\eta _Bx_tF(x_t)$$
(39)
and

$$\frac{\mathrm{\Delta }m_s}{\mathrm{\Delta }m_d}=\frac{m_{B_s}f_{B_s}^2B_{B_s}|V_{ts}|^2}{m_{B_d}f_{B_d}^2B_{B_d}|V_{td}|^2}.$$
(40)
The function $`F(x)`$ in (39) in (35) can be found in .

By using the method of Parodi, Roudeau and Stocchi , who have run their program starting from the inputs listed in table II and

$$\eta _1=1.44\pm 0.18,\eta _2=0.52,\eta _3=0.45\pm 0.01,$$
(41)
which are the values we find for our inputs, it is found that

$$\text{Im}\lambda _t=(1.14\pm 0.11)\times 10^4,$$
(42)
where the error is determined by the Gaussian distribution in Fig. 4. Notice that the value thus found is roughly 10% smaller than those found in other estimates for which $`\widehat{B}_K`$ is smaller.

The effect of this updated fit is a substantial reduction in the range of $`\text{Im}\lambda _t`$ with respect to what we used in our 1997 estimate : all values smaller than $`1.0\times 10^4`$ are now excluded (as opposed as before when values as small as $`0.6\times 10^4`$ were included).

## V Estimating $`\epsilon ^{}/\epsilon `$

The value of $`\epsilon ^{}`$ computed by taking all input parameters at their central values (Table II) is shown in Fig. 5. The figure shows the contribution to $`\epsilon ^{}`$ of the various operators in two $`\gamma _5`$ renormalization schemes at $`\mu =0.8`$ GeV and $`1.0`$ GeV. The advantage of such a histogram is that, contrary to the $`B_i`$, the size of the individual contributions does not depend on some conventional normalization.

As the histogram makes it clear, the two dominant contributions come from the gluon and electroweak penguins, $`Q_6`$ and $`Q_8`$. However, the gluon penguin dominates and there is very little cancellation with the electroweak penguin operator. The dominance of the $`I=0`$ components in the $`\chi `$QM originates from the $`O(p^4)`$ chiral corrections, the detailed size of which is determined by the fit of the $`\mathrm{\Delta }I=1/2`$ rule. It is a nice feature of the approach that the renormalization scheme stability imposed on the $`CP`$ conserving amplitudes is numerically preserved in $`\epsilon ^{}/\epsilon `$. The comparison of the two figures shows also the remarkable renormalization scale stability of the central value once the perturbative running of the quark masses and the quark condensate is taken into account.

In what follows, the model-dependent parameters $`M`$, $`\alpha _sGG/\pi `$ and $`\overline{q}q`$ are uniformly varied in their given ranges (flat scanning), while the others are sampled according to their normal distributions (see Table II for the ranges used). Values of $`\epsilon ^{}/\epsilon `$ found in the HV and NDR schemes are included with equal weight.

For a given set, a distribution is obtained by collecting the values of $`\epsilon ^{}/\epsilon `$ in bins of a given range. This is shown in Fig. 6 for a particular choice of bins. The final distribution is partially skewed, with more values closer to the lower end but a longer tail toward larger values. However, because the skewness of the distribution is less than one, the mean and the standard deviation are a good estimate of the central value and the dispersion of values around it. This statistical analysis yields

$$\epsilon ^{}/\epsilon =(2.2\pm 0.8)\times 10^3.$$
(43)
A more conservative estimate of the uncertainties is obtained via the flat scanning of the input parameters, which gives

$$0.9\times 10^3<\text{Re}\epsilon ^{}/\epsilon <4.8\times 10^3.$$
(44)
In both estimates a theoretical systematic error of $`\pm 30\%`$ is included in the fit of the $`CP`$ conserving amplitudes $`A_0`$ and $`A_2`$.

The stability of the numerical outcomes is only marginally affected by shifts in the value of $`\mathrm{\Omega }_{\eta +\eta ^{}}`$ due to NLO chiral corrections  and by additional isospin breaking effects . Any effective variation of $`\mathrm{\Omega }_{\eta +\eta ^{}}`$ is anti-correlated to the value of $`\alpha _sGG/\pi `$ obtained in the fit of $`A_2`$. We have verified that this affects the calculation of $`\widehat{B}_K`$ and the consequent determination of $`\text{Im}\lambda _t`$ in such a way to compensate numerically in $`\epsilon ^{}/\epsilon `$ the change of $`\mathrm{\Omega }_{\eta +\eta ^{}}`$. Waiting for a confident assessment of the NLO isospin violating effects in the $`K\pi \pi `$ amplitudes, we have used for $`\mathrm{\Omega }_{\eta +\eta ^{}}`$ the ‘LO’ value quoted in Table II.

The weak dependence on some poorly-known parameters is a welcome outcome of the correlation among hadronic matrix elements enforced in our semi-phenomenological approach by the fit of the $`\mathrm{\Delta }I=1/2`$ rule.

We have changed the central value of $`m_c`$ from 1.4 GeV, the value we used in , to 1.3 GeV in order to make our estimate more homogeneous with others. This change affects the determination of the ranges of the model parameters, mainly increasing $`\overline{q}q`$, via the fit of the $`CP`$ conserving amplitudes. The central value of $`\epsilon ^{}/\epsilon `$ turns out to be affected below 10%.

While the $`\chi QM`$ approach to the hadronic matrix elements relevant in the computation of $`\epsilon ^{}/\epsilon `$ has many advantages over other techniques and has proved its value in the prediction of what has been then found in the experiments, it has a severe short-coming insofar as the matching scale has to be kept low, around 1 GeV and therefore the Wilson coefficients have to be run down at scales that are at the limit of the applicability of the renormalization-group equations. Moreover, the matching itself suffers of ambiguities that have not been completely solved. For these reasons we have insisted all along that the approach is semi-phenomenological and that the above shortcomings are to be absorbed in the values of the input parameters on which the fit to the $`CP`$ conserving amplitudes is based.

## VI Other estimates

Figure 7 summarizes the present status of theory versus experiment. In addition to our improved calculation (and an independent estimate similarly based on the $`\chi `$QM), we have reported five estimates of $`\epsilon ^{}/\epsilon `$ published in the last year. Trieste’s, München’s and Roma’s ranges are updates of their respective older estimates, while the other estimates are altogether new.

The estimates reported come from the following approaches:

* München’s : In the München approach (phenomenological $`1/N`$) some of the matrix elements are obtained by fitting the $`\mathrm{\Delta }I=1/2`$ rule at $`\mu =m_c=1.3`$ GeV. The relevant gluonic and electroweak penguin $`Q_6`$ and $`Q_8_2`$ remain undetermined and are taken around their leading $`1/N`$ values (which implies a scheme dependent result). In Fig. 7 the HV (left) and NDR (right) results are shown. The darker range represents the result of Gaussian treatment of the input parameters compared to flat scanning (complete range).
* Roma’s : Lattice cannot provide us at present with reliable calculations of the $`I=0`$ penguin operators relevant to $`\epsilon ^{}/\epsilon `$ as well as of the $`I=0`$ components of the hadronic matrix elements of the current-current operators (penguin contractions), which are relevant to the $`\mathrm{\Delta }I=1/2`$ rule. This is due to large renormalization uncertainties, partly related to the breaking of chiral symmetry on the lattice. In this respect, very promising is the domain-wall fermion approach which allows us to decouple the chiral symmetry from the continuum limit. On the other hand, present lattice calculations compute $`K\pi `$ matrix elements and use lowest order chiral perturbation theory to estimate the $`K\pi \pi `$ amplitude, which introduces additional (and potentially large) uncertainties. In the recent Roma re-evaluation of $`\epsilon ^{}/\epsilon `$, $`B_6`$ is taken at the VSA result varied by a 100% error. The estimate quotes the values obtained in two $`\gamma _5`$ schemes (HV and NDR). The dark (light) ranges correspond to Gaussian (flat) scan of the input parameters.
* Dortmund’s : In recent years the Dortmund group has revived and improved the approach of Bardeen, Buras and Gerard based on the $`1/N`$ expansion. Chiral loops are regularized via a cutoff and the amplitudes are arranged in a $`p^{2n}/N`$ expansion. A particular attention has been given to the matching procedure between the scale dependence of the chiral loops and that arising from the short-distance analysis. The renormalization scheme dependence remains and it is included in the final uncertainty. The $`\mathrm{\Delta }I=1/2`$ rule is reproduced, but the presence of the quadratic cutoff induces a matching scale instability (which is very large for $`B_K`$). The NLO corrections to $`Q_6`$ induce a substantial enhancement of the matrix element (right range in Fig. 7) compared to the leading order result (left). The dark range is drawn for central values of $`m_s`$, $`\mathrm{\Omega }_{\eta +\eta ^{}}`$, $`\text{Im}\lambda _t`$ and $`\mathrm{\Lambda }_{QCD}`$.
* Dubna’s : The hadronic matrix elements are computed in the ENJL framework including chiral loops up to $`O(p^6)`$ and the effects of scalar, vector and axial-vector resonances. ($`B_K`$, and therefore $`\text{Im}\lambda _t`$, is taken from ) Chiral loops are regularized via the heat-kernel method, which leaves unsolved the problem of the renormalization scheme dependence. A phenomenological fit of the $`\mathrm{\Delta }I=1/2`$ rule implies deviations up to a factor two on the calculated $`Q_6`$. The reduced (dark) range in Fig. 7 corresponds to taking the central values of the NLO chiral couplings and varying the short-distance parameters.
* Taipei’s : Generalized factorization represents an attempt to parametrize the hadronic matrix elements in the framework of factorization without a-priori assumptions . Phenomenological parameters are introduced to account for non-factorizable effects. Experimental data are used in order to extract as much information as possible on the non-factorizable parameters. This approach has been applied to the $`K\pi \pi `$ amplitudes in ref. . The effective Wilson coefficients, which include the perturbative QCD running of the quark operators, are matched to the factorized matrix elements at the scale $`\mu _F`$ which is arbitrarily chosen in the perturbative regime. A residual scale dependence remains in the penguin matrix elements via the quark mass. The analysis shows that in order to reproduce the $`\mathrm{\Delta }I=1/2`$ rule and $`\epsilon ^{}/\epsilon `$ sizable non-factorizable contributions are required both in the current-current and the penguin matrix elements. However, some assumptions on the phenomenological parameters and ad hoc subtractions of scheme-dependent terms in the Wilson coefficients make the numerical results questionable. In addition, the quoted error does not include any short-distance uncertainty.
* Trieste’s: The dark (light) ranges correspond to Gaussian (flat) scan of the input parameters. The bar on the left corresponds to the present estimate. That on the right is a new estimate , similarly based on the $`\chi `$QM hadronic matrix elements, in which however $`\epsilon ^{}/\epsilon `$ is estimated in a novel way by including the explicit computation of $`\epsilon `$ in the ratio as opposed to the usual procedure of taking its value from the experiment. This approach has the advantage of being independent of the determination of the CKM parameters $`\text{Im}\lambda _t`$ and of showing more directly the dependence on the long-distance parameter $`\widehat{B}_K`$ as determined within the model. The difference (around 10%) between the two Trieste’s estimates corresponds effectively to a larger value of $`\text{Im}\lambda _t`$, as determined from $`\epsilon `$ only, with respect to eq. (42).
* Valencia’s : The standard model estimate given by Pallante and Pich is obtained by applying the FSI correction factors obtained using a dispersive analysis à la Omnès-Mushkelishvili  to the leading (factorized) $`1/N`$ amplitudes. The detailed numerical outcome has been questioned on the basis of ambiguities related to the choice of the subtraction point at which the factorized amplitude is taken . Large corrections may also be induced by unknown local terms which are unaccounted for by the dispersive resummation of the leading chiral logs. Nevertheless, the analysis of ref.  confirms the crucial role of higher order chiral corrections for $`\epsilon ^{}/\epsilon `$, even though FSI effects alone leave the problem of reproducing the $`\mathrm{\Delta }I=1/2`$ selection rule open.
* Lund’s : The $`\mathrm{\Delta }I=1/2`$ rule and $`B_K`$ have been studied in the NJL framework and $`1/N`$ expansion by Bijnens and Prades showing an impressive scale stability when including vector and axial-vector resonances. A recent calculation of $`\epsilon ^{}/\epsilon `$ at the NLO in $`1/N`$ has been performed in ref. . The calculation is done in the chiral limit and it is eventually corrected by estimating the largest $`SU(3)`$ breaking effects. Particular attention is devoted to the matching between long- and short-distance components by use of the $`X`$-boson method . The couplings of the $`X`$-bosons are computed within the ENJL model which improves the high-energy behavior. The $`\mathrm{\Delta }I=1/2`$ rule is reproduced and the computed amplitudes show a satisfactory renormalization scale and scheme stability. A sizeable enhancement of the $`Q_6`$ matrix element is found which brings the central value of $`\epsilon ^{}/\epsilon `$ at the level of $`3\times 10^3`$.

Cut-off based approaches should also pay attention to higher-dimension operators which become relevant for matching scales below 2 GeV . The calculations based on dimensional regularization may avoid the problem if phenomenological input is used in order to encode in the hadronic matrix elements the physics at all scales.

Other attempts to reproduce the measured $`\epsilon ^{}/\epsilon `$ using the linear $`\sigma `$-model, which include the effect of a scalar resonance with $`m_\sigma 900`$ MeV, obtain the needed enhancement of $`Q_6`$ . However, the $`CP`$ conserving $`I=0`$ amplitude falls short the experimental value by a factor of two. With a lighter scalar, $`m_\sigma 600`$ MeV the $`CP`$ conserving $`I=0`$ amplitude is reproduced, but $`\epsilon ^{}/\epsilon `$ turns out more than one order of magnitude larger than currently measured .

## VII Conclusions

The present analysis updates our 1997 estimate $`\text{Re}\epsilon ^{}/\epsilon =1.7_{1.0}^{+1.4}\times 10^3`$ , which already pointed out the potential relevance of non-factorizable contributions and the importance of addressing both $`CP`$ conserving and violating data for a reliable estimate of $`\epsilon ^{}/\epsilon `$.

The increase in the central value is due to the update on the experimental inputs (mainly $`V_{ub}/V_{cb}`$). The uncertainty is reduced when using the Gaussian sampling, as opposed to the flat scan used in ref. . On the other hand the error obtained by flat scanning is larger due to the larger systematic uncertainty ($`\pm 30\%`$) used in the fit of the $`CP`$ conserving amplitudes.

Among the corrections to the leading $`1/N`$ (factorized) result, FSI play a crucial role in the determination of the gluon penguin matrix elements. Recent dispersive analyses of $`K\pi \pi `$ amplitudes show how a (partial) resummation of FSI increases substantially the size of the $`I=0`$ amplitudes, while slightly affecting the $`I=2`$ components . On the other hand, the precise size of the effect depends on boundary conditions of the factorized amplitudes which are not unambiguously known .

Finally, it is worth stressing that FSI by themselves do not account for the magnitude of the $`CP`$ conserving decay amplitudes. In our approach a combination of non-factorizable soft-gluon effects (at $`O(p^2)`$) and FSI (at the NLO) makes possible to reproduce the $`\mathrm{\Delta }I=1/2`$ selection rule. In turn, requiring the fit of the $`CP`$ conserving $`K\pi \pi `$ decays allows for the determination of the “non-perturbative” parameters of the $`\chi `$QM, which eventually leads to the detailed prediction of $`\epsilon ^{}/\epsilon `$. Confidence in our final estimate of $`\epsilon ^{}/\epsilon `$ is based on the coherent picture of kaon physics which arises from the phenomenological determination of the model parameters and the self-contained calculation of all $`\mathrm{\Delta }S=1`$ and 2 hadronic matrix elements.

###### Acknowledgements.

We thank F. Parodi and A. Stocchi for their help in the determination of $`\text{Im}\lambda _t`$ and for Fig. 4.


text:
# Sum Rule Approach to Collective Oscillations of Boson-Fermion Mixed Condensate of Alkali Atoms

## Abstract

The behavior of collective oscillations of a trapped boson-fermion mixed condensate is studied in the sum rule approach. Mixing angle of bosonic and fermionic multipole operators is introduced so that the mixing characters of the low-lying collective modes are studied as functions of the boson-fermion interaction strength. For an attractive boson-fermion interaction, the low-lying monopole mode becomes a coherent oscillation of bosons and fermions and shows a rapid decrease in the excitation energy towards the instability point of the ground state. In contrast, the low-lying quadrupole mode keeps a bosonic character over a wide range of the interaction strengths. For the dipole mode the boson-fermion in-phase oscillation remains to be the eigenmode under the external oscillator potential. For weak repulsive values of the boson-fermion interaction strengths we found that an average energy of the out-of-phase dipole mode stays lower than the in-phase oscillation. Physical origin of the behavior of the multipole modes against boson-fermion interaction strength is discussed in some detail.

Collective oscillation is one of the most prominent phenomena common to a variety of many-body systems. The realization of the Bose-Einstein condensates (BEC) for trapped Alkali atoms offers a possibility to study such phenomena of quantum systems under ideal conditions. Up to now the experimental as well as theoretical studies of collective motions of BEC have been intensively performed. Quite recently the degenerate Fermi gas of trapped 40K atoms has been realized, which motivates the study of collective motion also in Fermi gases. These developments further encourage the study of possible boson-fermion mixed condensates of trapped atoms. Now the static properties, stability conditions, and some dynamical properties of trapped boson-fermion condensates have been investigated. In the present paper, we study the behavior of collective oscillations of a boson-fermion mixed condensate at $`T=0`$ for both repulsive and attractive boson-fermion interactions. We adopt the sum rule approach that has proved to be successful in the studies of collective excitations of Bose condensates. For the mixed condensate, we introduce a mixing angle of the boson/fermion excitation operators so as to allow the in- and out-of-phase oscillations of the bose and fermi condensates.

In the sum rule approach we first calculate the energy weighted moments $`m_p=_j(E_jE_0)^p|j|F|0|^2`$ of the relevant multipole operator $`F`$, where $`|j`$’s represent the complete set of eigenstates of the Hamiltonian with energies $`E_j`$, and $`|0`$ denotes the ground state. The excitation energy is expressed as $`\mathrm{}\omega =(m_3/m_1)^{1/2}`$ which provides a useful expression for the average energy of the collective oscillation. The moments are calculated from formulae $`m_1=\frac{1}{2}0|[F^{},[H,F]]|0`$ and $`m_3=\frac{1}{2}0|[[F^{},H],[H,[H,F]]]|0`$. We consider three types of multipole operators which are defined by

$$F_\alpha ^\pm \underset{i=1}{\overset{N_b}{}}f_\alpha (\stackrel{}{r}_{bi})\pm \underset{i=1}{\overset{N_f}{}}f_\alpha (\stackrel{}{r}_{fi}),(\alpha =M,D,Q)$$
(1)
where the functions $`f_\alpha `$ are defined by $`f_M(\stackrel{}{r})=r^2`$ for monopole, $`f_D(\stackrel{}{r})=z`$ for dipole, and $`f_Q(\stackrel{}{r})=3z^2r^2`$ for quadrupole excitations. The indices $`b,f`$ denote boson/fermion, $`N_b,N_f`$ the numbers of bose/fermi particles, and $`\pm `$ correspond to the in-phase and out-of-phase oscillation of the two types of particles. We actually take a linear combination of the form

$$F_\alpha (𝐫;\theta )=F_\alpha ^+\mathrm{cos}\theta +F_\alpha ^{}\mathrm{sin}\theta (\frac{\pi }{2}<\theta \frac{\pi }{2})$$
(2)
parametrized by the mixing angle $`\theta `$. We study the value of $`\theta `$ that minimizes the calculated frequency $`\omega `$ for each $`\alpha `$.

We consider the polarized boson-fermion mixed condensate in spherically symmetric harmonic oscillator potential. The system is described by the Hamiltonian

$$H=\underset{i=1}{\overset{N_b}{}}\left\{\frac{\stackrel{}{p}_{bi}^2}{2m}+\frac{1}{2}m\omega _0^2\stackrel{}{r}_{bi}^2+\frac{1}{2}g\underset{j=1}{\overset{N_b}{}}\delta (\stackrel{}{r}_{bi}\stackrel{}{r}_{bj})\right\}+\underset{i=1}{\overset{N_f}{}}\left\{\frac{\stackrel{}{p}_{fi}^2}{2m}+\frac{1}{2}m\omega _0^2\stackrel{}{r}_{fi}^2\right\}+h\underset{i=1}{\overset{N_b}{}}\underset{j=1}{\overset{N_f}{}}\delta (\stackrel{}{r}_{bi}\stackrel{}{r}_{fj})$$
(3)
where we assume the same oscillator frequencies and masses for bosons and fermions for simplicity. The coupling constants $`g,h`$ are the boson-boson/boson-fermion interaction strengths represented by the $`s`$-wave scattering lengths $`a_{bb}`$ and $`a_{bf}`$ as $`g=4\pi \mathrm{}^2a_{bb}/m`$ and $`h=4\pi \mathrm{}^2a_{bf}/m`$. The fermion-fermion interaction has been neglected as the polarized system is considered. Following the standard calculation procedure the excitation frequencies are obtained as:

(i) Monopole

$$\frac{\omega _M(\theta )}{\omega _0}=\sqrt{2}\sqrt{1+\frac{E_{\mathrm{kin}}^++\frac{9}{4}E_{bb}+\frac{9}{4}E_{bf}+2(E_{\mathrm{kin}}^{}+\frac{9}{4}E_{bb}+\frac{3}{4}\mathrm{\Delta }^{})\mathrm{cos}\theta \mathrm{sin}\theta \mathrm{\Delta }\mathrm{sin}^2\theta }{E_{\mathrm{ho}}^++2E_{\mathrm{ho}}^{}\mathrm{cos}\theta \mathrm{sin}\theta }}$$
(4)
(ii) Quadrupole

$$\frac{\omega _Q(\theta )}{\omega _0}=\sqrt{2}\sqrt{1+\frac{E_{\mathrm{kin}}^++2E_{\mathrm{kin}}^{}\mathrm{cos}\theta \mathrm{sin}\theta \frac{2}{5}\mathrm{\Delta }\mathrm{sin}^2\theta }{E_{\mathrm{ho}}^++2E_{\mathrm{ho}}^{}\mathrm{cos}\theta \mathrm{sin}\theta }}$$
(5)
(iii) Dipole

$$\frac{\omega _D(\theta )}{\omega _0}=\sqrt{1\frac{4}{3\mathrm{}\omega _0}\frac{\mathrm{\Omega }\mathrm{sin}^2\theta }{N^++2N^{}\mathrm{cos}\theta \mathrm{sin}\theta }}$$
(6)
Here we defined $`E_{\mathrm{kin}}^\pm E_{\mathrm{kin}}^b\pm E_{\mathrm{kin}}^f,E_{\mathrm{ho}}^\pm E_{\mathrm{ho}}^b\pm E_{\mathrm{ho}}^f`$ and $`N^\pm N_b\pm N_f`$, where $`E_{\mathrm{kin}}^{\{b,f\}}`$ and $`E_{\mathrm{ho}}^{\{b,f\}}`$ are respectively the expectation values of the kinetic and oscillator potential energies for boson/fermion in the ground state. Boson-boson and boson-fermion interaction energies have been denoted by $`E_{bb}`$ and $`E_{bf}`$. The quantities $`\mathrm{\Delta },\mathrm{\Delta }^{}`$ and $`\mathrm{\Omega }`$ are given in terms of the boson/fermion densities $`n_b(r),n_f(r)`$ in the ground state by

$$\mathrm{\Delta }hd^3rr^2\frac{dn_f(r)}{dr}\frac{dn_b(r)}{dr},\mathrm{\Delta }^{}hd^3rr\left[n_f(r)\frac{dn_b(r)}{dr}\frac{dn_f(r)}{dr}n_b(r)\right],\mathrm{\Omega }h\xi ^2d^3r\frac{dn_f(r)}{dr}\frac{dn_b(r)}{dr},$$
(7)
where $`\xi =\sqrt{\mathrm{}/m\omega _0}`$. One may use the stationary condition of the ground state,

$$2E_{\mathrm{kin}}^+2E_{\mathrm{ho}}+3E_{bb}+3E_{bf}=0,2E_{\mathrm{kin}}^{}2E_{\mathrm{ho}}^{}+3E_{bb}+\mathrm{\Delta }^{}=0$$
(8)
in order to eliminate in eq.(4) the dependences on $`E_{bb},E_{bf}`$ and $`\mathrm{\Delta }^{}`$. The monopole frequency is then rewritten as

$$\frac{\omega _M(\theta )}{\omega _0}=\sqrt{5\frac{E_{\mathrm{kin}}^++2E_{\mathrm{kin}}^{}\mathrm{cos}\theta \mathrm{sin}\theta +2\mathrm{\Delta }\mathrm{sin}^2\theta }{E_{\mathrm{ho}}^++2E_{\mathrm{ho}}^{}\mathrm{cos}\theta \mathrm{sin}\theta }}.$$
(9)
We have checked that the Thomas-Fermi calculation of the ground state adopted below gives rise to a negligible difference if one evaluates either the expression (9) or (4).

We calculate the ground state energies and densities of the boson-fermion mixed system in the Thomas-Fermi approximation which is valid for $`gN_b1`$ and $`N_f1`$ except around the surface region. We take harmonic oscillator length $`\xi `$ and energy $`\mathrm{}\omega _0`$ as units, and define scaled dimensionless variables: the radial distance $`x=r/\xi `$, boson/fermion densities $`\rho _{b,f}(x)=n_{b,f}(r)\xi ^3/N_{b,f}`$, and chemical potentials $`\stackrel{~}{\mu }_{b,f}=2\mu _{b,f}/\mathrm{}\omega _0`$. We solve the coupled Thomas-Fermi equations,

$$\stackrel{~}{g}N_b\rho _b(x)+x^2+\stackrel{~}{h}N_f\rho _f(x)=\stackrel{~}{\mu _b},[6\pi ^2N_f\rho _f(x)]^{2/3}+x^2+\stackrel{~}{h}N_b\rho _b(x)=\stackrel{~}{\mu _f},$$
(10)
where $`\stackrel{~}{g}=2g/\mathrm{}\omega _0\xi ^3`$ and $`\stackrel{~}{h}=2h/\mathrm{}\omega _0\xi ^3`$.

One of the most promising candidates for the realization of the mixed condensate is the potassium isotope system. Precise values of the scattering lengths are not well known at present and different values have been reported . We take for the boson-boson interaction the parameters of the 41K-41K system in  and a trapping frequency of $`450`$Hz which gives $`\stackrel{~}{g}=0.2`$. For the boson-fermion interaction we take several values in the range $`h/g=\stackrel{~}{h}/\stackrel{~}{g}=2.373.2`$. It should be noted that the interaction strength can be controlled using Feshbach resonances.

We have performed a numerical calculation for $`N_b=N_f=10^6`$. In the ground state the fermions have a much broader distribution than bosons because of the Pauli principle. Fermions are further squeezed out of the center for a repulsive boson-fermion interaction ($`h>0`$) . They will eventually form a shell-like distribution around the surface of bosons for $`h/g1`$ and will be completely pushed away from the center ($`n_f(0)=0`$) at around $`h/g3`$. For an attractive boson-fermion interaction, on the other hand, the central densities of the bosons and fermions increase together. The system becomes unstable against collapse at around $`h/g=2.37`$ due to the strong attractive boson-fermion interaction.

Figure 1 shows the kinetic energy, the oscillator potential energy, and the interaction energy contributions to the ground state energy against the parameter $`h/g`$. The figure shows also the quantities $`\mathrm{\Delta }`$ and $`\mathrm{\Omega }`$ which represent the contributions of the boson-fermion interaction to the multipole frequencies (4)-(6) and (9). One may notice that the fermionic kinetic- and potential-energy contributions are a few times larger than the bosonic contribution in the present system. It is noted that $`\mathrm{\Delta }`$ takes large negative values at both large negative and positive regions of $`h/g`$. In the former region the bose and fermi density distributions become coherent due to the attractive interaction, and the radial integral in eq.(7) takes a large positive value. In the opposite case ($`h/g1`$), the same integral changes sign because the fermions are pushed away from the center by the repulsive boson-fermion interaction, thus giving a large negative contribution in the bosonic surface region. In the region $`0<h/g<1`$, on the other hand, the integral is slightly positive and $`\mathrm{\Delta }`$ takes a small positive value. The quantity $`\mathrm{\Omega }`$ follows the same trend as $`\mathrm{\Delta }`$, but the absolute value is much smaller than $`\mathrm{\Delta }`$, as the most important contribution to the integral comes from the surface region where $`r\xi `$.

Once the ground state parameters are obtained the frequencies (4)-(6) are minimized against $`\theta `$. Usually, the sum rule approach predicts a strength-weighted average energy of eigenstates for a given multipolarity. The calculated energy coincides with the true excitation energy if the relevant strength is concentrated in a single state. By adopting the minimization procedure we simultaneously determine the character of the low-lying collective mode and the corresponding average energy. The character of the mode is given by the value of $`\theta `$, for instance, $`\theta \pi /4`$ for the bosonic- and $`\pi /4`$ for the fermionic-modes, $`\theta 0`$ for the in-phase oscillation and $`\pi /2`$ for the out-of-phase oscillation of the two types of particles. As there are two kinds of particles involved in eq.(2), one would expect an emergence of two types of collective oscillations for each multipole. Another collective mode would have a character orthogonal to the low-lying mode as far as the phase relation of the two operators in (2) is concerned. In the present approach we calculate the frequency of the latter mode, the high-lying one, from eqs.(4)-(6), by using the operator $`F_\alpha ^+\mathrm{sin}\theta _LF_\alpha ^{}\mathrm{cos}\theta _L`$ for each $`\alpha `$, where the mixing angle $`\theta _L`$ is the one determined for the low-lying mode.

Figure 2 shows frequencies of the lower (solid lines) and the higher (dashed lines) modes for (a) monopole, (b) quadrupole and (c) dipole cases as functions of $`h/g`$. The corresponding mixing angles $`\theta `$ determined by the minimization procedure for the lower mode are plotted in Fig.3 for each multipolarity as a function of $`h/g`$. Below we discuss the behavior of the frequencies $`\omega _\alpha `$ by defining three regions of $`h/g`$: (I) $`h/g<0`$, (II) $`0<h/g1`$, (III) $`1h/g`$.

a)monopole:
For a non-interacting boson-fermion system the low-lying monopole mode is the fermionic oscillation with frequency $`\omega _M^L=2\omega _0`$, while the higher mode is the bosonic one with $`\omega _M^H=\sqrt{5}\omega _0`$ in the Thomas-Fermi approximation. Around $`h0`$ one may obtain $`\omega _M^L2\omega _0\left(1+A\stackrel{~}{h}N^{\frac{1}{6}}\right)`$, $`\omega _M^H\sqrt{5}\omega _0\left(1A^{}\stackrel{~}{h}\stackrel{~}{g}^{\frac{1}{5}}N^{\frac{3}{10}}\right)`$, where $`A=36^{\frac{1}{6}}/4\sqrt{2}\pi ^2`$ and $`A^{}=7\sqrt{3}/20\pi ^2(8\pi /15)^{\frac{1}{5}}`$. The mixing angle for the lower mode is given as $`\theta _M=\pi /4\delta _M`$ with $`\delta _M=(56^{\frac{1}{6}}/\sqrt{2}\pi ^2)\stackrel{~}{h}N^{\frac{1}{6}}`$. This behavior is seen in region (II) where the boson-fermion interaction is weakly repulsive and the lower monopole mode is of a fermionic character, see Fig.3 (solid line). Bosonic oscillation in this region is more rigid than the fermionic one because of the repulsive interaction among bosons. In region (I), the situation is quite different: The low-lying mode becomes a coherent boson-fermion oscillation as represented by the large negative value of $`\mathrm{\Delta }`$, and the excitation energy shows a sharp decrease towards the instability point $`h/g2.37`$ of the ground state, although $`\omega _M`$ does not become exactly zero within our approximation. In this region the attractive boson-fermion interaction is much more effective in the excited state than in the ground state and cancels out the increase in the kinetic energy. In region (III), too, we find that the low-lying mode becomes an in-phase oscillation. Here, the boson and the fermion densities in the ground state are somewhat separated, and the in-phase motion which keeps this separation is energetically more favorable than the out-of-phase motion as seen in the value of $`\mathrm{\Delta }`$.

b)quadrupole:
For the quadrupole excitation (Fig.2(b)), the lower (higher) energy mode is almost the pure bosonic (fermionic) oscillation over the broad range of the $`h/g`$ values studied, Fig.2 (dashed line). To the first order in $`\stackrel{~}{h}`$ the frequencies of the lower and the higher quadrupole modes are given by $`\omega _Q^L\sqrt{2}\omega _0\left(1B\stackrel{~}{h}N^{\frac{1}{6}}\right)`$, $`\omega _Q^H2\omega _0(1B^{}\stackrel{~}{h}\stackrel{~}{g}^{\frac{2}{5}}N^{\frac{7}{30}})`$, where $`B=6^{\frac{1}{6}}/4\sqrt{2}\pi ^2`$ and $`B^{}=(15/8\pi )^{\frac{2}{5}}/(146^{\frac{1}{6}}\sqrt{2}\pi ^2)`$. The corresponding mixing angle for the lower mode is given by $`\theta _Q=\pi /4+\delta _Q`$ with $`\delta _Q=2^3B^{}\stackrel{~}{h}\stackrel{~}{g}^{\frac{2}{5}}N^{\frac{7}{30}}`$. For the quadrupole mode similar mechanisms as for the monopole mode are at work concerning the dependence on $`|h/g|`$. The role of the boson-fermion interaction is however much reduced compared with the monopole case as seen by the factor 2/5 of eq.(5), which reflects that the quadrupole oscillation has five different components. Thus the quadrupole mode obtains an in-phase character only at large values of $`|h/g|`$. In the other region of $`|h/g|`$ the low-lying mode becomes a simple bosonic oscillation. This is because the fermionic mode costs a larger kinetic energy and favors $`\theta =\pi /4`$ as seen in the term $`(E_{\mathrm{kin}}^bE_{\mathrm{kin}}^f)\mathrm{sin}\theta \mathrm{cos}\theta `$ in eq.(5).

c)dipole:
General arguments show that for a harmonic oscillator external potential a uniform shift of the ground state density generates an eigenstate of the system, corresponding to the boson-fermion in-phase dipole oscillation with frequency $`\omega _0`$. This is evident in Fig.2c) and also in eq.(6) at $`\theta =0`$. In the regions (I) and (III) the out-of-phase oscillation is unfavorable due to the same reason as for the monopole mode: It loses the energy gained in the ground state boson-fermion configuration. For a weakly repulsive $`h`$ an interesting possibility arises: In the region (II) the out-of-phase mode of the boson-fermion oscillation lies lower than the in-phase mode. Let us first note that at $`h=0`$ the out-of-phase oscillation frequency becomes degenerate as the in-phase one because the bosonic and the fermionic dipole modes are independent. One may note that at $`h/g1`$ again the degeneracy occurs. Here the potential term for the fermion becomes almost linear to the fermion density itself (see, eq.(10)). This suggests that the fermion density is determined almost entirely by the chemical potentials and becomes nearly constant as far as the boson density is finite. A uniform dipole shift of fermions thus produces almost no effect on bosons, and results in the degeneracy of the frequency. Between $`h/g=0`$ and 1, the boson-fermion repulsion is weaker for the out-of-phase oscillation than the in-phase one (and hence the ground state) as reflected in the sign of $`\mathrm{\Omega }`$.

In the present paper, we studied collective oscillations of trapped boson-fermion mixed condensates using sum rule approach. We introduced a mixing angle of bosonic and fermionic multipole operators so as to study if the in- or out-of-phase motion of those particles is favored as a function of the boson-fermion interaction strength. For the monopole and quadrupole cases, the coupling of the bose- and fermi-type oscillation is not large for moderate values of the coupling strength $`h`$. At large values of $`|h/g|`$, the low-lying modes become an in-phase oscillation of bosons and fermions. This is especially so for the monopole oscillation at attractive boson-fermion interaction: The excitation energy of this mode almost vanishes at the instability point of the ground state. In the case of the dipole mode, in contrast, the in-phase oscillation remains an exact eigenmode with a fixed energy for harmonic oscillator potentials, while the average energy of the out-of-phase oscillation is strongly dependent on the boson-fermion interaction. We found that at weak repulsive values of the interaction the out-of-phase motion stays lower than the in-phase oscillation. In this paper we calculated also the frequencies of the high-lying modes for each multipole, by adopting the operators orthogonal to the low-lying modes. These modes, too, are collective in character and, in the present framework, their average frequencies showed rather strong dependences on the boson-fermion coupling strength. Deeper insight into the collective modes studied in this paper will require a detailed investigation of the solutions of, e.g., the self-consistent RPA type equations that allow an arbitrary radial dependence of the excitation operators. Studies in this direction are now in progress.


text:
# Cyclotron lines in X-ray pulsars as a probe of relativistic plasmas in superstrong magnetic fields

## Introduction

Accreting magnetized neutron stars are an ideal cosmic laboratory for high energy relativistic physics. Cyclotron resonance features are the signature of the presence of a superstrong magnetic field, following the first discovery in Her X-1 (Trümper et al. trump ). These features are due to the discrete Landau energy levels for motion of free electrons perpendicular to the field in presence of a locally uniform superstrong magnetic field. A slight deviation from a pure harmonic relationship in the spacing of the different levels is expected due to relativistic effects ($`\frac{\omega _n}{m_e}=((1+2n\frac{B}{B_{\mathrm{crit}}}\mathrm{sin}^2\theta )^{\frac{1}{2}}1)/\mathrm{sin}^2\theta `$, e.g. Araya and Harding araya ). Therefore the detection of these features in the emitted X–ray spectra is in principle a direct measure of the field intensity.
As the number of sensitive measurements in the hard X–ray interval (above $``$ 10 keV) is continuously growing, a sample is available to search for possible correlations between the observed parameters. A detailed modeling is difficult and a parametrized shape of the continuum still is not available from theoretical models, but substantial advances in our understanding of the radiation transport in strongly magnetized atmospheres were done in the last decade (e.g. Alexander et al. alex1 , Alexander and Mészáros alex2 , Araya and Harding araya , Isenberg et al. isenba ; isenbb , Nelson et al. nelson ). Some of these new results focused on the properties of the cyclotron resonance features observed in the spectra of accreting X–ray pulsars. In this report we discuss the current status of the measurements of cyclotron lines, with emphasis on the possible correlations between observable parameters.

## The data

The BeppoSAX satellite has observed all the bright persistent and three bright transient (recurrent) accreting X–ray pulsars. Apart from the case of X Persei (Di Salvo et al. robba ), a source with a luminosity substantially lower than the other sources in the sample, the spectra observed by BeppoSAX can be empirically described using the classical power–law–plus–cutoff spectral function by White et al. white . The sensitive broad band BeppoSAX observations also allowed the detailed characterization of low energy components below a few keV (like in Her X–1, Dal Fiume et al. herx1 , and in 4U1626–67, Orlandini et al. 1626 ) and the detection of absorption features in the hard X–ray range of the spectra, interpreted as cyclotron resonance features.

A summary of the properties of the broad band spectra and of the cyclotron lines as measured with BeppoSAX is reported in Dal Fiume et al. ddf . From these measurements we obtained evidence of a correlation between the centroid energy of the feature and its width. This correlation is presented and discussed elsewhere (Dal Fiume et al. ddf ; ddf98 ).

### Transparency in the line

A straightforward parameter to be obtained from observations is the transparency in the line, defined as the ratio between the transmitted observed flux and the integrated flux from the continuum without the absorption feature. This ratio likely depends on the harmonic number of the feature we are observing (e.g. Wang et al. wang ) and on the physical parameters of the specific accretion column. From an observational point of view, this ratio is strongly affected by the modelization of the “continuum” shape, that is by the spectral shape used to describe the differential broad band photon number spectrum. In Figure 1 we report the observed transparencies obtained dividing the observed by the expected flux, both integrated in a $`\pm 2\sigma `$ interval around the line centroid (here $`\sigma `$ is the Gaussian width of the measured cyclotron feature). To further emphasize the uncertainty in this estimate, we added a 10% error to the data. The purely statistical uncertainties are substantially smaller.

This measured transparency is related to a simple physical parameter, the opacity to photons with energy approximately equal to the cyclotron resonance energy. However the emerging integral flux and the shape of the line itself are non trivially related to the radiation transport in this energy interval, a rather difficult problem to be solved.
From the phenomenological point of view, one can observe that the measured transparencies cluster around 0.5–0.6, with the notable exception of Cen X–3.

### Magnetic field intensity and spectral hardness

The influence of the magnetic field intensity on the broad band spectral shape is debated. Early attempts to estimate a possible dependence of electron temperature, and therefore of broad band spectral shape, on the magnetic field intensity were done by Harding et al hard84 . Actually they conclude that “the equilibrium atmospheres have temperatures and optical depths that are very sensitive to the strength of the surface magnetic field”. If this is the case and if the broad band spectral hardness is related, as one could naively assume, to the temperature of the atmosphere, some correlation between this hardness and the cyclotron line energy should appear in data.

This seems to be the case shown in Figure 2. Here we report the ratio between photon fluxes in two “hard” bands (the flux in 20–100 keV divided by the flux in 7–15 keV) versus the cyclotron line centroid. The ratio between the two fluxes is affected by the choice of the continuum, as in Figure 1. We therefore also in this case added 10% error bars that indicate this uncertainty. The statistical errors are substantially smaller.
The number of sources in this plot is still very limited and therefore one cannot exclude that this apparent correlation is merely due to the limited size of the sample. Nevertheless the apparent correlation is in the right direction, i.e. harder spectra are observed for higher field intensities.
We parenthetically add that no cyclotron resonance feature was observed in the pulse–phase averaged spectra of the two hardest sources of this class observed with BeppoSAX (GX1+4 Israel et al. israel  and GS1843+00 Piraino et al. piraino ). If this correlation proves to be correct, this may suggest that cyclotron resonance features in these two sources should be searched at the upper limit of the BeppoSAX energy band or beyond.

### Conclusions

In conclusion, even if no complete, parametric theoretical approach to model the observed spectra of accreting X–ray pulsars is still available, some quantitative measures of parameters of hot plasmas in superstrong magnetic fields are possible.
Modeling the transparency in the cyclotron resonance feature is a complex problem. Further information will be extracted from maps of this transparency as a function of pulse phase.
The correlation between spectral hardness and field intensity is in agreement with theoretical models. This correlation, if confirmed, can be used as a rough estimate of the magnetic field intensity from the measured spectral hardness.

###### Acknowledgements.

This research is supported by Agenzia Spaziale Italiana (ASI) and Consiglio Nazionale delle Ricerche (CNR) of Italy. BeppoSAX is a joint program of ASI and of the Netherlands Agency for Aerospace Programs (NIVR).