text:
# Morphology of ledge patterns during step flow growth of metal surfaces vicinal to fcc(001)

## Abstract

The morphological development of step edge patterns in the presence of meandering instability during step flow growth is studied by simulations and numerical integration of a continuum model. It is demonstrated that the kink Ehrlich-Schwoebel barrier responsible for the instability leads to an invariant shape of the step profiles. The step morphologies change with increasing coverage from a somewhat triangular shape to a more flat, invariant steady state form. The average pattern shape extracted from the simulations is shown to be in good agreement with that obtained from numerical integration of the continuum theory.

today

Epitaxial growth on vicinal surfaces is known to give rise to interesting growth instabilities under suitable conditions, e.g. to step bunching, mound formation and meandering of the step edges . The meandering instability emerges when the interlayer mass transport from the upper side of the step is reduced due to the Ehrlich-Schwoebel barrier  enhancing growth of protrusions at the step edges. This is now known as the Bales-Zangwill instability (BZI)  which tends to destabilize the ledge morphology due to terrace diffusion and asymmetric interlayer crossing. There is no diffusion along the ledges in BZI. However, recently it was found that in the case of $`1+1`$ dimensional growth there is an analogous phenomenon due to the kink Ehrlich-Schwoebel barrier for going around a kink site at the step edge. The corresponding kink Ehrlich-Schwoebel effect (KESE) leads to growth of unstable structures at the step edges with a dynamically selected wavelength . The ledge instabilities were originally found and reported experimentally on the Cu(1,1,17) vicinal surface  but attributed to the BZI scenario. More recent STM experiments on the Cu(1,1,17) surface proposed that the formation of the regular patterns is due to the KESE . Since then theoretical studies of the meandering instability have shown that the KESE indeed supersedes the BZI in the formation of the periodic patterns  and eventually leads to an in-phase motion of the step edge structures .

Instability and wavelength selection of the step edge patterns due to the KESE have been studied within the framework of a continuum step model , the solid-on-solid (SOS) lattice model , and semi-realistic Monte Carlo (MC) simulations . In particular, in the recent MC work  it was shown that on vicinal Cu(1,1,$`m`$) surfaces the observed instability is due to the KESE and the competing BZI is of no importance in the length and time scales considered. The role of dimer nucleation in determining the selected wavelength was confirmed, in good agreement with the theoretical scaling relation  and more recent SOS simulations . In the MC simulations, there was evidence of phase locking of the ledge structures at the largest coverages studied, but this was not quantitatively confirmed.

Regarding the step morphologies, a triangular shape has been predicted to occur for a strong KESE and a rounded, more flat shape for a weak KESE . However, the MC simulations of Ref.  indicate that in the case of a strong KESE there is in fact an interesting shape transition from narrow, somewhat triangular shapes in the initial stage of growth to more rounded patterns in the large coverage regime. Moreover, the MC simulations of Ref.  and SOS model results of Ref.  are in disagreement with asymptotic evolution of the step profiles as predicted by the continuum theories . In this work we study the ledge morphologies of growing steps on the vicinal Cu(1,1,$`m`$) surfaces in detail. In particular, we study the onset of the in-phase growth and the phase-locking of the step profiles in the presence of a strong KESE. Our results show that the ledge morphologies assume an invariant shape due to an interplay between various mass transport currents and phase-locking of the steps. We show how the ledge profiles from the MC simulations can be reproduced by explicitly including the relevant mass transport currents on the surface. On continuum level this indicates a delicate balance between the various currents that determines and stabilizes the invariant ledge shapes.

The model system used here is as in Ref. , based on MC simulations of a lattice gas model with energetics from the effective medium theory (for more details, see Refs. ). Our MC method is efficient enough to simulate growth of Cu up to ten monolayers (ML) under realistic temperature and flux conditions. The temperature range explored here was $`T=240310`$ K and the flux $`F=3\times 10^31.0`$ ML/s. Thus the ratio between

Fig. 1: Snapshots of typical ledge profiles with step orientations in the close packed $`\left[110\right]`$ direction at $`T=300`$ K with $`F=8\times 10^2`$ ML/s for coverages $`\theta =0.4,2.0`$, and $`10.0`$, in figures (a)–(c) (lateral and vertical scales are 1000 and 70 lattice spacings, respectively). In (d) the shape transition is shown. The profiles have been obtained by averaging over the meander periods at $`\theta =0.4`$ ML (circles) and at $`\theta =2.0`$ ML (squares). The horizontal direction is scaled with the wavelength $`\lambda =120a`$ and in the vertical direction with the roughness $`w=1.4a`$ and $`w=5.6a`$ for coverages $`0.4`$ ML and $`2.0`$ ML, respectively.

the terrace diffusion and the flux $`D/F6\times 10^59\times 10^7`$ in units of the lattice constant $`a=0.361`$nm, corresponding to a typical molecular beam epitaxy regime . The energetics of the model also specifies the important length scales controlling step flow growth. These are $`\mathrm{}_c`$, the length scale for dimer nucleation at the step edge , and the kink Schwoebel length  $`\mathrm{}_s=\mathrm{exp}[(E_sE_d)/k_BT]1`$ which is related to the energy barriers $`E_s=0.52`$ eV and $`E_d=0.26`$ eV for jumps around a kink site and along a straight edge, respectively. For the close packed $`\left[110\right]`$ ledges, $`\mathrm{}_s10^4`$ and $`\mathrm{}_c10^2`$ around room temperature corresponding to strong KESE . In Ref.  it was shown that the wavelength of the step edge patterns is given by $`\mathrm{}_c=(12D_s/FL)^\alpha `$, $`D_s`$ being the adatom diffusion constant along the straight edge, with a scaling exponent $`\alpha 0.23`$, and an effective barrier of $`E_{\mathrm{eff}}=75\pm 10`$ meV. Both are in good agreement with the exact values which give $`\alpha =1/4`$  and $`E_{\mathrm{eff}}=E_d/4=65`$ meV, respectively. Our previous study  was done on a Cu(1,1,17) surface but we have checked the results also with smaller terrace widths.

Simulation results for the step edge profiles on Cu(1,1,17) are shown in Figs. 1(a)–1(c) after deposition of $`\theta =0.4,2.0`$, and $`10.0`$ ML, respectively, at $`T=300`$ K with $`F=6\times 10^2`$ ML/s. In the beginning of growth (Fig. 1(a)) the shape of the patterns is somewhat triangular as predicted for a relatively strong KESE . The

Fig. 2: (a) Lateral moments $`M_n`$ of the step meander periods ($`n=2,4,6,8,10`$ from top to bottom) as a function of coverage in a semi-logarithmic scale. The moments approach constant values already before $`2.0`$ ML. (b) Width $`w`$ of the profiles as a function of the coverage. No saturation is observed up to $`10`$ ML. The slope of the solid line corresponds to $`\beta =1/3`$.

meandering structures are not yet completely in the same phase indicating that the diffusion field has not yet coupled the subsequent step edge trains, a typical feature for KESE dominated meandering . However, a selection of the relatively well-defined wavelength for all ledges is apparent already at this stage of growth . At larger coverages the meandering of steps begins gradually to phase-lock, seen in Fig. 1(b), and in-phase growth and phase-locking seem complete at largest studied coverage of $`10`$ ML shown in Fig 1(c). However, now the average shape of the patterns is clearly different from that at low coverages. The shape of the average patterns is more rounded, as predicted for a weak KESE. In Fig. 1(d) we show this change by comparing average ledge profiles after $`0.4`$ and $`2.0`$ ML, respectively .

From Fig. 1 it is clear that there is no coarsening of the structures when the coverage is large enough. The steady state pattern shape seems to be governed by geometric constraints which is a sign of asymmetry of the growth rates between bottom and top parts of the steps. This asymmetry is a general feature in many models of step growth with or without coarsening . Moreover, a quantitative inspection of the patterns at larger coverages suggest that the profiles have an invariant shape. This can be seen by examining the $`n`$th lateral moments $`M_n(\theta )=\zeta _i(x,\theta )x^n_{i,x}`$ of the meander periods $`\zeta _i(x,\theta )`$, shown in Fig. 2(a). The scaled even moments approach

Fig. 3: The average shape of the step patterns at $`T=300`$ K and $`\theta =8.0`$ ML with $`F=6\times 10^2`$ ML/s (squares). Changing the temperature or the flux as described in the text does not have any effects within the error bars. The solid line is the stationary profile obtained by integration of Eq. (6). Good agreement between the average and the integrated profiles is evident. The inset displays the relative differences of the profiles, $`\mathrm{\Delta }=(\zeta _{all}\zeta _i)/\zeta _{all}`$, where $`\zeta _{all}`$ is the profile with all currents included, and $`i=k,e,SB`$ denotes the solution with only a small contribution for the KESE current from the Gibbs-Thomson and symmetry breaking currents, respectively (from top to bottom in the inset).

their steady state values already at $`\theta 2`$ ML. The patterns are still changing, however, which can be seen from the roughness of the step $`w(\theta )=\sqrt{\zeta (x,\theta )^2_x}`$, where $`\zeta (x,\theta )`$ is the step profile. It does not show any sign of saturation up to the largest coverage in the simulations. Instead the roughness follows $`w(\theta )\theta ^\beta `$, with $`\beta 0.3`$ as shown in Fig. 2(b). It is interesting to note that although the roughness does not saturate the shape of the periodic structures attains an invariant form.

Our simulation results show that the profile shape is rather insensitive to deposition and temperature conditions. This suggests that the invariant shape is not dependent on the relative magnitudes of the various diffusion processes but rather is a result of geometric constraints due to crowding and in-phase evolution of the step edges. In order to justify this assumption we compare the MC profiles with continuum profiles which are obtained as stationary solutions to the dynamic equation $`_t\zeta =_xJ_{tot},`$ where $`J_{tot}`$ is the total mass current at the step edge. The most important partial currents which we take into account in the total current here, when expressed in terms of the variable $`m(x)=(_x\zeta )/\sqrt{1+(_x\zeta )^2}`$ and appropriately scaled, are the mass current due to the destabilizing strong KESE

$$J_k=\frac{m(\sqrt{1m^2}\left|m\right|)\sqrt{1m^2}}{(\left|m\right|+L_c^1\sqrt{1m^2})^2},$$
(1)
Fig. 4: The mass currents Eqs. (1)-(4) are shown using the integrated profile as an input. Note the difference between the vertical scales for (c) and (d).

the stabilizing current due to the Gibbs-Thomson effect and edge diffusion

$$J_e=\frac{2D_S\stackrel{~}{\mathrm{\Gamma }}}{F}\left(\sqrt{1m^2}+\frac{D_L}{D_SL}\right)(_{xx}m)\sqrt{1m^2},$$
(2)
the current into the step edge from the deposition flux

$$J_d=Lm\sqrt{1m^2},$$
(3)
and the front-back symmetry breaking current

$`J_{SB}`$ $`=`$ $`{\displaystyle \frac{D_S\stackrel{~}{\mathrm{\Gamma }}L}{F}}(_xm)(_{xx}m)\sqrt{1m^2}`$ (5)
$`+{\displaystyle \frac{L^2}{3}}m(_xm)(3m^2).`$
In these expressions $`L`$ is the terrace width, $`D_S`$ is the macroscopic diffusion constant on the terrace, $`D_L`$ is the macroscopic diffusion constant along the step edge, and $`\stackrel{~}{\mathrm{\Gamma }}`$ is the step stiffness (see Refs.  for the definitions and experimental values of the parameters, respectively). All length scales are given in the units of the lattice constant. By requiring the condition of stationarity

$$J_{tot}J_k+J_e+J_d+J_{SB}=0,$$
(6)
we obtain a second order differential equation for $`m(x)`$. The stationary profiles are obtained by solving Eq. (6) numerically for given initial conditions $`m(\pm 1)=\pm m_0`$.

The stationary solution is found using the value $`m_00.97`$ as the boundary condition in order to match the end points with the slopes of the patterns obtained from the MC simulations . The other parameter values of the integration are based on known energetics of Cu, yielding $`\mathrm{}_c=7001600`$, $`\mathrm{}_s=2\times 10^42\times 10^5`$, $`D_L/(D_SL)=250700`$, and $`D_S\stackrel{~}{\mathrm{\Gamma }}/F=0.54000`$ in the range $`T=240300`$K and $`F=3\times 10^310^1`$ ML/s. For the step stiffness we used the expression $`\stackrel{~}{\mathrm{\Gamma }}=\mathrm{exp}[E_k/k_BT]/2`$, where $`E_k=0.13`$ eV is the kink energy . In all cases we set $`L=10`$ for the terrace width. The resulting profiles are shown in Fig. 3 with various values of the parameters. The shape is rather independent of the details of the currents in agreement with simulations. In Fig. 3 the average shapes obtained from the simulations are plotted with a few different flux rates. In the inset we show how the resulting profile deviates from the complete one when each of the mass currents is forced to be small.

In Fig. 4 we show the mass currents using the integrated profile as an input. It is now seen that for the invariant profile there is a delicate compensation of the currents, the Gibbs-Thomson current compensated by the sum of the KESE, the deposition, and the symmetry breaking currents. This compensation happens for the specific shape of the profile, and cannot take place e.g. in the case of a triangular shaped profile as obtained in the initial stages of growth. In determination of the stationary profile shape the front-back symmetry breaking and the geometric constraints contained implicitly in the initial conditions are crucial.

In summary, the MC  and the SOS  simulations have proven that the KESE is the dominant mechanism behind the meandering instability and that it leads to the selection of the dominant wavelength determined by dimer nucleation at step edges. In this work we have shown that the KESE also induces an invariant shape of the step profiles during in-phase growth. This occurs even though the overall roughness of the step structures $`w(\theta )`$ shows no signs of saturation. The value of the corresponding scaling exponent $`\beta 0.3`$ is consistent with the case of an isolated step . The SOS model gives for the strong KESE an exponent $`\beta 0.57`$ , while for a collection of steps in the phase-locking regime $`\beta =1/2`$ . This puzzling behavior of dynamical scaling is apparently related to the strict in-phase growth and consequent formation of the invariant shape of the profile. The fact that the shape remains invariant although the roughness does not show any sign of saturation indicates a subtle coupling of the step edge currents with the stationary morphology. By numerically integrating the continuum equation we have shown how the interplay between various surface currents determines the invariant step shapes.

Acknowledgements: We wish to thank J. Kallunki and Z. Chvoj for helpful discussions and the Center of Scientific Computing, Ltd. for computing time. TSR and CG thank colleagues at Fyslab, HUT for their wonderful hospitality during their stay in Helsinki. This work has been supported by the Academy of Finland, in part through its Center of Excellence program. We also acknowledge partial support from the National Science Foundation, USA under Grant EEC-0085604 (including International Supplement).


format: text

text:
# V344 Lyr: an unusual large-amplitude SU UMa-type dwarf nova with a short supercycle

## 1 Introduction

Dwarf novae are a class of cataclysmic variables (CVs), which are close binary systems consisting of a white dwarf and a red dwarf secondary transferring matter via the Roche lobe overflow. Thermal instability in the resultant accretion disk is now widely believed to cause dwarf nova-type outbursts (see Osaki (1996) for a review). In dwarf novae whose binary mass-ratios ($`q`$=$`M_2`$/$`M_1`$) are small enough, a different kind of instability – tidal instability – occurs (Whitehurst, 1988), which is now widely believed to be a cause of superoutburst and superhumps in SU UMa-type dwarf novae (for a recent review of SU UMa-type stars and their observational properties, see Warner (1995)).

Among SU UMa-type dwarf novae, some objects show large-amplitude outbursts. Although there is evidence that outburst amplitudes of dwarf novae are known to comprise a continuum from small to large amplitudes (see a discussion in Patterson et al. (1996)), these large-amplitude objects are sometimes symbolically known as TOADs (Tremendous Outburst Amplitude Dwarf Novae: Howell et al. (1995)). At the extreme end, there exist so-called WZ Sge-type dwarf novae (originally proposed by Bailey (1979); see also Downes & Margon (1981) and O’Donoghue et al. (1991)). WZ Sge-type dwarf novae (and related systems) are known to show peculiar outburst characteristics among SU UMa-type dwarf novae: 1) long interval between outbursts (WZ Sge itself has a recurrence time of 23 to 33 yr, which is the longest among all known dwarf novae), 2) lack or very low frequency of normal outbursts, which are more abundant than superoutbursts in usual SU UMa-type dwarf novae (no normal outburst has ever been observed in WZ Sge itself) and 3) extremely long (up to $``$100 d) and bright superoutbursts.

V344 Lyr is a dwarf nova discovered by Hoffmeister (1966). Kato (1993) detected superhumps with a period of 0.09145(2) d, which confirms that V344 Lyr is a member of SU UMa-type dwarf novae with long orbital periods. The reported range of variability (13.8 – $`[`$20 $`V`$) makes V344 Lyr a candidate for rare large-amplitude dwarf novae with long orbital periods. Since many of large-amplitude SU UMa-type dwarf novae are known to be most infrequently outbursting systems, the proposed supercycle of 240 d (Kato, 1993) would make V344 Lyr an exceptional object. In order to clarify the situation, we examined our observations and the observations between 1994 July and 2001 June, made by the VSNET Collaboration 111http://www.kusastro.kyoto-u.ac.jp/vsnet/.

## 2 Observations and results

### 2.1 Outburst cycle length and statistics

The observations by the authors and those from VSNET were made visually, using $`V`$-magnitude calibrated comparison stars. The typical error of visual estimates was 0.2 mag, which will not affect the following discussion. The total number of observations was 1089. The number of positive observations was 209, corresponding to the outburst duty cycle of 19%. Although this value may have suffered from some degree of selection bias, the large duty cycle is already unusual for a system with a large outburst amplitude. We selected outbursts from these data, which are summarized in Table 1 and Figure 1.

Most outbursts were unambiguously classified as either normal outbursts or superoutbursts based on their duration and peak magnitudes, but the faintness of the object and sometimes unfavourable observing conditions made some of classifications slightly ambiguous. Such outbursts are flagged as “?”.

### 2.2 Outburst amplitude

Since the outburst amplitude is sensitive to outburst and quiescent magnitudes, we have tried to recalibrate these values using the modern $`V`$ scale. Since recent outburst observations have been made using modern $`V`$-magnitude comparison stars, we have adopted $`V`$=14.0 as the maximum magnitude. \[This value is in agreement with the $`V`$=14.2 obtained by Kato (1993). Even a better agreement would be achieved when considering that the CCD observation by Kato (1993) started a few days after the visual maximum.\]

Regarding the quiescent magnitude, Hoffmeister (1966) stated that “the object is invisible on the Palomar Sky Survey plates”, from which an upper limit of 20 has been adopted in the literature. We have confirmed this finding by a direct inspection of paper reproductions of the POSS I plates. The object is also missing from the USNO A1.0 and A2.0 catalogs, which reach a limiting magnitude of 19.5. We have also found that the object is faintly seen on the digitized POSS I plates, available at the USNOFS Image and Catalogue Archive, and yielded a magnitude of $`R`$=19.5$`\pm `$0.3, by comparison with the modern comparison stars. Since short-period dwarf novae have colors close to $`VR`$=0 even in quiescence, this value is considered to be a good approximation of the quiescent $`V`$ magnitude.

We have also inspected the available images, on which V344 Lyr was probably recorded at or close to its minimum (the object was recorded in outburst on POSS II plates). The images taken at Ouda Station, Kyoto University (Ohtani et al., 1992), on 1991 February 23 have yielded $`I_\mathrm{c}`$=17.8$`\pm `$0.2. An unfiltered CCD image on 1996 July 17 shows no hint of the object down to 19.0 mag. The Guide Star Catalog plates, as given in Downes & Shara (1993), show the object at 18–19 mag. The above results are summarized in Table 2. These results may suggest a considerable variation of the quiescent magnitude, but it would not be surprising that the object may not have completely reached quiescence in some observations, when one takes the large outburst duty cycle of 19 % into account. From these observations, we adopted 5.5$`\pm `$0.3 mag for the outburst amplitude of V344 Lyr.

## 3 Discussion

### 3.1 Length of supercycle

The most prominent sequence of superoutbursts can be found between JD 2450882 and 2451094 (Figure 2). The detection of three successive superoutbursts establishes the supercycle of $``$110 d, which is half of the estimate by Kato (1993). The other superoutbursts and candidate superoutbursts are well expressed by this period, but the supercycle does not seem to be extremely stable: the supercycle seems to be longer ($``$120–130 d) in 1994–1995. Using the best sampled interval (JD 2450312 – 2451628), a regression of the observed epochs of superoutbursts yielded a supercycle of 109.6 d. This supercycle length is one of the shortest among known SU UMa-type dwarf novae (cf. Nogami et al. (1997)), except ER UMa stars (Kato & Kunjaya (1995); Robertson et al. (1995); Misselt & Shafter (1995); Nogami et al. (1995)) which are unusual SU UMa-type dwarf novae with an exceptionally large outburst frequency (supercycle lengths are shorter than 60 d; normal outbursts occur every 3–5 d) and a small outburst amplitude (for an observational review and theoretical explanations of ER UMa stars, see Kato et al. (1999) and Osaki (1995), respectively).

### 3.2 Outburst amplitude versus supercycle length

As shown in Section 3.1, the supercycle of V344 Lyr is 109.6 d, which is exceptionally short for an SU UMa-type dwarf nova with a large outburst amplitude. Table 3 lists the properties of SU UMa-type dwarf novae with well-established short supercycle lengths ($`T_\mathrm{s}<300`$ d)222 The limit has been chosen for two reasons: (1) $`T_\mathrm{s}`$ tends to be stable in frequently outbursting systems (i.e. short $`T_\mathrm{s}`$ systems), the extreme cases being ER UMa stars. Systems with superoutbursts less than once per year often do not have a fixed $`T_\mathrm{s}`$. (2) Due to unavoidable seasonal observational gaps, long $`T_\mathrm{s}`$ systems have uncertainties in unambiguously identifying $`T_\mathrm{s}`$. . The source data are from Nogami et al. (1997), with recent revisions and additions mentioned in the individual notes. Excluding SX LMi (unusual low-amplitude system as discussed by Nogami et al. (1997)), Z Cha (high-inclination eclipsing system) at lower right, and V344 Lyr, we found the following good correlation between $`T_\mathrm{s}`$ and outburst amplitude ($`A`$) for systems $`T_\mathrm{s}>80`$ d (i.e. normal SU UMa stars).

$$A=2.04\mathrm{log}(T_\mathrm{s})0.4.$$
(1)
Most remarkable is that the outburst amplitudes of SU UMa-type dwarf novae are confined to a narrow range (3.1 to 4.0 mag) for systems between $`T_\mathrm{s}`$=84 and $`T_\mathrm{s}`$=174, except V344 Lyr. This distribution has a mean amplitude of 3.5 mag and a standard deviation of 0.3 mag. The observed amplitude of V344 Lyr (5.5$`\pm `$0.3 mag) is 7$`\pm `$1$`\sigma `$ larger than the average. Even adopting a conservative upper limit of the minimum magnitude of 18.5 (from GSC plate scan), the value is still 3$`\sigma `$ above the average.

Warner (1987) presented, with an assumption of an optically thick disk, a formulation of the dependence of the absolute $`V`$ magnitude on the system inclination. A disk seen at nearly edge-on (i=85) is 3.5 mag brighter than the pole-on (i=0) view. Although this inclination effect explains the existence of low-amplitude systems among eclipsing systems (e.g. Z Cha), it would be difficult to explain the existence of high-amplitude systems, since a pole-on disk is expected to be only 0.4 mag brighter than systems with moderate inclinations (Warner, 1987). Such a large outburst amplitude as in V344 Lyr would thus be difficult to reproduce simply by the geometrical effect. Since superoutbursts of V344 Lyr have quite typical properties of usual SU UMa-type dwarf novae (Kato, 1993), it is unlikely the superoutbursts of V344 Lyr are intrinsically brighter, which implies that the quiescent accretion disk of V344 Lyr is fainter than the average.

The standard disk instability theory expects that the quiescent luminosity is a strong function of the mass-transfer rate and the quiescent viscosity (cf. Warner (1987); Smak (1989); Cannizzo (1998)). Since the same theory expects that the supercycle length is inversely proportional to the mass-transfer rate (Ichikawa & Osaki, 1994), the mass-transfer rate in V344 Lyr is expected to be higher than the average, which would work to reduce the outburst amplitude. The modification of the quiescent viscosity could be a viable idea, but is against the observed high frequency of normal outbursts. V344 Lyr thus shows a considerable departure from known SU UMa-type systems, and may require an additional mechanism to effectively reduce the quiescent luminosity or to increase the outburst frequency. Further detailed observations of V344 Lyr, especially orbital parameters and long-term quiescent photometry are therefore necessary.

### Notes on individual objects in Table 3:

– Large change of $`T_\mathrm{s}`$ has been inferred (Fried et al., 1999). The value given in the table corresponds to the most stable outbursting period (Kato et al. (1996); Kato et al. in preparation).

– Newly determined $`T_\mathrm{s}`$ from Kato (2001a).

– Newly discovered ER UMa-type dwarf nova. The data are from Ishioka et al. (2001).

$`T_\mathrm{s}`$ from Kato et al. (2000). The maximum magnitude of $`V`$=12.6 (Mason et al., 1982) was probably in error; no such bright outburst has been observed in plate searches (see also Richter (1989); Kato et al. (1998)). The maximum magnitude is taken from the recent VSNET observations.

$`T_\mathrm{s}`$ and maximum magnitude have been revised using VSNET observations.

$`T_\mathrm{s}`$ listed in Kato & Nogami (1997) and Nogami et al. (1997) was erroneous. The new $`T_\mathrm{s}`$ is determined from VSNET observations, which show very regular supercycles, resembling those of YZ Cnc. The maximum and minimum magnitudes are also from VSNET observations, based on the newly calibrated $`V`$-magnitude sequence.

– Typical $`T_\mathrm{s}`$ is given. The system sometimes shows reduced outburst and supuroutburst frequencies (Rosenzweig et al., 2000).

– Typical $`T_\mathrm{s}`$, as measured from VSNET observations between 1996 and 1998, is given. The system is known to sometimes show a reduced frequency of superoutbursts.

– Revised maximum and minimum magnitudes from VSNET observations.

– Revised $`T_\mathrm{s}`$ based on the analysis of the VSOLJ data. Revised superhump period from Kato (2001b). The maximum and minimum magnitudes are from VSNET observations and citetmis96sequence, respectively.

– Mean $`T_\mathrm{s}`$ from Kato (2001c), which reported a variation of supercycles between 169 and 229 d. Kato (2001c) reported rather unusual outburst characteristics, having a low ratio of (normal outbursts)/(superoutbursts). The minimum magnitude given is our new calibration of the DSS 1, while Liu et al. (1999) reported a rough estimate of $`V`$=20.8 from their spectroscopic observation. The system requires further detailed study in quiescence.

– Revised $`T_\mathrm{s}`$ and maximum magnitude from VSNET observations.

– Revised $`T_\mathrm{s}`$ and maximum magnitude from VSNET observations. The minimum magnitude is taken from Misselt (1996).

– Mean $`T_\mathrm{s}`$ from Kato (2001d), which showed a variation of supercycles between 250 and 312 d.

The following objects in Nogami et al. (1997) have been excluded from the list:

Nogami & Kato (1997) listed a possible supercycle of 140 d, while recent reports to VSNET show that the intervals are typically longer and outbursts occur rather irregularly (see also Kato et al. (2000)).

Nogami et al. (1997) listed a possible supercycle of 290 d. Recent study by Nogami et al. (2001) indicates that supercycles vary more irregularly.

## 4 Summary

We observed the large-amplitude SU UMa-type dwarf nova V344 Lyr. Combined with reports to VSNET, we have succeeded in determining its mean supercycle length as 109.6 d. This value is one of the smallest among known SU UMa-type dwarf novae (except unusual ER UMa-type dwarf novae). The observed outburst amplitude (5.5$`\pm `$0.3 mag) in V344 Lyr is found to be exceptionally large for a system with such a short supercycle. Since a short supercycle strongly suggests a relatively high mass-transfer rate, the extreme outburst parameters of V344 Lyr would require a mechanism to effectively reduce the quiescent luminosity or to increase the outburst frequency.

## Acknowledgments

The authors are grateful to VSNET members, especially to L. Szentasko, T. Vanmunster, J. Pietz, L. T. Jensen, E. Broens, M. Reszelski for providing vital observations. The authors are grateful to the staff of the VSOLJ (Variable Star Observers League) for their on-line database. We are also grateful to an anonymous referee, whose comments significantly improved the paper. This research has made use of the USNOFS Image and Catalogue Archive operated by the United States Naval Observatory, Flagstaff Station (http://www.nofs.navy.mil/data/fchpix/).


text:
# TWO BODY 𝐵 DECAYS, FACTORIZATION AND Λ_{𝑄𝐶𝐷}/𝑚_𝑏 CORRECTIONSTalk given by G. Martinelli

## 1 Introduction

The advent of the $`B`$ factories, Babar babar  and Belle belle , opens new perspectives for very precise measurements of non-leptonic $`B`$-decays ferroni  and calls for a significant improvement of the theoretical predictions. In this respect, important progress has been recently achieved by systematic studies of factorization in $`B`$ decays li ; beneke . These studies confirmed the physical idea previous  that factorization holds in heavy hadron decays, $`m_Q\mathrm{\Lambda }_{QCD}`$, for the leading terms of the $`\mathrm{\Lambda }_{QCD}/m_Q`$ expansion. They leave open, however, the central question of whether i) the leading terms predict with sufficient accuracy the relevant $`B`$-meson decay rates or ii) the power-suppressed corrections, which cannot be evaluated without some model assumption, are phenomenologically important. This problem has been addressed in a series of papers charming lorolast . In particular, the main conclusion of ref. noilast  is that non-perturbative $`\mathrm{\Lambda }_{QCD}/m_b`$ corrections from the leading operators of the effective weak Hamiltonian, conventionally called Charming (or GIM) penguins, are very important in cases where the factorized amplitudes are either colour or Cabibbo suppressed. One of the consequences of this observation is that factorization at the leading order in $`\mathrm{\Lambda }_{QCD}/m_b`$ is unable to reproduce the observed $`B\pi \pi `$ and $`BK\pi `$ $`BR`$s even taking into account the uncertainties of the input parameters.

The results of ref. noilast  seem to be at variance with the conclusions of ref. lorolast  where it is stated that, in the QCD factorization approach, “an acceptable fit to the branching fractions is obtained even if we impose $`\gamma <90^o`$ as implied by the standard constraints on the unitarity triangle” 111 Indeed the analysis of ref. Ciuchini:2000de  gives $`\gamma =(54.8\pm 6.2)^o`$ and ref. hocker  quotes $`34^o<\gamma <82^o`$..

Given the contradictory conclusions of the different studies, it is very important to clarify the situation by comparing input parameters, methods of analysis and results. In order to help the debate on this subject, we discuss in the following:

* The theoretical framework of factorization and the calculation of the amplitudes at the leading order of the $`\mathrm{\Lambda }_{QCD}/m_b`$ expansion;
* A comparison of the leading terms, including their uncertainties, with the measured $`B\pi \pi `$ and $`BK\pi `$ branching ratios;
* Results for the $`B\pi \pi `$ and $`BK\pi `$ $`BR`$s including charming (and GIM) penguins;
* A comparison of the results obtained with charming and GIM penguins noilast  with the model of the $`\mathrm{\Lambda }_{QCD}/m_b`$ corrections adopted in lorolast .

## 2 Factorization

In this section we recall the basic ingredients necessary to compute the relevant amplitudes either whithin factorization or including the corrections which arise at higher order in the $`\mathrm{\Lambda }_{QCD}/m_b`$ expansion.

The physical amplitudes for $`BK\pi `$ and $`B\pi \pi `$ can be conveniently written in terms of RG-invariant parameters built using the Wick contractions of the effective Hamiltonian BS . In the heavy quark limit, following the approach of ref. beneke , it is possible to compute these RG invariant parameters using factorization. The formalism of ref. beneke  has been developed so that it is also possible to include the perturbative corrections to order $`\alpha _s`$, i.e. at the next-to-leading order in perturbation theory 222 An alternative framework, provided by the approach of ref. li , will not be discussed here. For a recent discussion see also lastchris ..

For the sake of discussion, it is instructive to start from the explicit expression of the $`B_dK^+\pi ^{}`$ amplitude. In terms of the parameters defined in BS , this amplitude reads

$`𝒜(B_dK^+\pi ^{})=`$ $``$ $`V_{us}V_{ub}^{}\left(E_1(s,u,u;B_d,K^+,\pi ^{})P_1^{\mathrm{GIM}}(s,u;B_d,K^+,\pi ^{})\right)`$ (1)
$`+`$ $`V_{ts}V_{tb}^{}P_1(s,u;B_d,K^+,\pi ^{}).`$
We have

$`E_1(s,u,u;B_d,K^+,\pi ^{})`$ $`=`$ $`a_1^u(K\pi )Q_1^u_{\mathrm{fact}}+a_2^u(K\pi )Q_2^u_{\mathrm{fact}}+\stackrel{~}{E}_1`$
$`P_1(s,u;B_d,K^+,\pi ^{})`$ $`=`$ $`{\displaystyle \underset{i=3}{\overset{10}{}}}a_i^c(K\pi )Q_i_{\mathrm{fact}}+\stackrel{~}{P}_1`$
$`P_1^{\mathrm{GIM}}(s,u;B_d,K^+,\pi ^{})`$ $`=`$ $`{\displaystyle \underset{i=3}{\overset{10}{}}}(a_i^c(K\pi )a_i^u(K\pi ))Q_i_{\mathrm{fact}}+\stackrel{~}{P}_1^{\mathrm{GIM}},`$ (2)
where $`Q_i_{\mathrm{fact}}`$ denotes the factorized matrix elements and $`a_i^f`$ the parameters introduced in beneke . Eqs. (1) and (2) are exact. Unfortunately, similar equations require the knowledge of several non-perturbative parameters, which at present cannot be extracted from the data. To be predictive, we will then use our physical intuition to reduce their number. In eqs. (1) and (2):

1. The CKM matrix elements can either be taken from other experimental measurements or from a fit to the non-leptonic $`BR`$s, assuming that factorization is accurate enough;
2. The coefficients $`a_i^f(M_1M_2)`$ (e.g. $`a_1^u(K\pi )`$) are renormalization group invariant, as the corresponding factorized matrix elements, and have been computed perturbatively at the NLO in refs. beneke ; lorolast ;
3. The coefficients $`a_6^f`$ and $`a_8^f`$, which have also been computed to one-loop order, are instead scheme dependent. Their scheme-dependence is cancelled by the hadronic matrix elements of the penguin operators $`Q_6`$ and $`Q_8`$ respectively. Assuming factorization for the chirally enhanced contributions beneke , the latter can be expressed in terms of the ratio

$$r_\chi ^K(\mu )=\frac{2m_K^2}{\overline{m}_b(\mu )(\overline{m}_s(\mu )+\overline{m}_q(\mu ))},$$
(3)
which is formally of $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$ but numerically important (an analogous parameter $`r_\chi ^\pi (\mu )`$ can be defined for $`\pi \pi `$ decays). We will discuss the rôle of these terms below.
4. The leading amplitudes $`Q_i_{\mathrm{fact}}`$ are computed in terms of decay constants and semi-leptonic form factors. The form factors can either be taken from theoretical calculations latticeff ; qcdsrff  or fitted from the experimental $`BR`$s (the possibility of extracting them from the corresponding $`B\pi `$ semileptonic form factor at small momentum transfer is at present rather remote).
5. The tilded parameters, namely $`\stackrel{~}{P}_1`$, $`\stackrel{~}{P}_1^{\mathrm{GIM}}`$ and $`\stackrel{~}{E}_1`$, are genuine, non-perturbative $`\mathrm{\Lambda }_{QCD}/m_b`$ corrections which cannot be computed at present. If we neglect Zweig-suppressed contributions, by $`SU(2)`$ symmetry one can show that all the Cabibbo-enhanced $`\mathrm{\Lambda }_{QCD}/m_b`$ corrections to $`BK\pi `$ decays can be reabsorbed in $`\stackrel{~}{P}_1`$. Several corrections are contained in $`\stackrel{~}{P}_1`$: this parameter includes not only the charming penguin contributions, but also annihilation and penguin contractions of penguin operators. It does not include leading emission amplitudes of penguin operators ($`Q_3`$$`Q_{10}`$) which have been explicitly evaluated using factorization. Had we included these terms, this contribution would exactly correspond to the parameter $`P_1`$ of ref. BS . The parameter $`\stackrel{~}{P}_1`$ ($`P_1`$) encodes automatically not only the effect of the annihilation diagrams considered in keum , but all the other contributions of $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$ with the same quantum numbers of the charming penguins. In this respect it is the most general parameterization of all the perturbative and non-perturbative contributions of the operators $`Q_5`$ and $`Q_6`$ ($`Q_3`$ and $`Q_4`$), including the worrying higher-twist infrared divergent contribution to annihilation discussed in ref. lorolast ; beneke2 . The parameter $`\stackrel{~}{P}_1`$ has the same quantum numbers and physical effects as the original charming penguins proposed in charming , although it has a more general meaning.
6. If one also includes $`B\pi \pi `$ decays we have several other parameters, for example $`P_1^{\mathrm{GIM}}`$ and $`P_3`$, in the formalism of ref. BS . A closer look to $`P_3`$ shows that this term is due either to Zweig suppressed annihilation diagrams (called CPA and DPA in ref. charming ) or to annihilation diagrams which are colour suppressed with respect to those entering $`\stackrel{~}{P}_1`$. For this reason in ref. noilast  $`P_3`$ was taken to be zero. $`\stackrel{~}{P}_1`$ is equal to the corresponding parameter in $`K\pi `$ decays if $`SU(3)`$ symmetry is assumed. In our analysis we have used the same value of $`\stackrel{~}{P}_1`$ for all $`K\pi `$ and $`\pi \pi `$ channels. $`P_1^{\mathrm{GIM}}`$ will be discussed later on.

We are now ready to discuss factorization for the leading terms of the $`\mathrm{\Lambda }_{QCD}/m_b`$ expansion. Factorization is the theory of non leptonic decays which is obtained in the limit $`m_b\mathrm{}`$. Thus it consists in neglecting all terms of $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$ ($`\stackrel{~}{E}_1=0`$, $`\stackrel{~}{P}_1=0`$, $`r_\chi ^{K,\pi }=0`$, etc.). At lowest order in perturbation theory, called also naïve factorization, the $`a_i^f`$ are simple combinations of the Wilson coefficients and do not depend on the hadron wave functions. The inconvenience of naïve factorization is that physical amplitudes still have a marked dependence on the renormalization scale because, contrary to the Wilson coefficients, the factorized matrix elements are scale independent.

The scale dependence is reduced by working at $`𝒪(\alpha _s)`$, both for the Wilson coefficients and the matrix elements. In refs. beneke  it has been shown that, at this order, all the dangerous infrared divergences can be reabsorbed in the definition of the hadronic wave functions. For the leading terms in the $`\mathrm{\Lambda }_{QCD}/m_b`$ expansion there are strong arguments to support the idea that this will remain true at all orders in $`\alpha _s`$, see also ref. allorders . Thus, in the limit $`m_b\mathrm{}`$, it is likely that factorization is preserved by strong interactions. At $`𝒪(\alpha _s)`$ or higher, the coefficients $`a_i^f`$ depend on the specific detail of the hadron wave-fuctions. For this reason, the uncertainties relative to the wave functions, as the residual renormalization scale dependence, must be taken into account in the evaluation of the uncertainties for the theoretical predictions. The approximation in which we neglect all the $`\mathrm{\Lambda }_{QCD}/m_b`$ corrections, but include the perturbative corrections to the leading contribution, is called QCD factorization or simply factorization. Factorization implies an important consequence: predictions of non-leptonic decay rates are model independent to the extent that the few relevant hadronic parameters, namely the kaon and pion decay constants, $`f_{K,\pi }`$, the semileptonic form factors, $`f_{K,\pi }(0)`$, and the hadronic wave functions are known.

In cases like $`BK\pi `$ decays, where the factorized amplitudes are Cabibbo suppressed, the corrections of $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$, which unfortunately are model dependent, become important. At lowest order in $`\alpha _s`$, the chirally enhanced terms proportional to $`r_\chi ^K`$ ($`r_\chi ^\pi `$) are computable by assuming that factorization can be applied beyond the leading order. A substantial difficulty arises, however, at $`𝒪(\alpha _s\mathrm{\Lambda }_{QCD}/m_b)`$. Although the chirally enhanced corrections from $`Q_{6,8}`$ are infrared finite, other contributions of the same order from different operators are infrared divergent, signaling that they belong to the class of the non-perturbative contributions which appear beyond factorization. These cannot be predicted using the same hadronic quantities of the factorized amplitudes. For this reason, any phenomenological analysis which aims at including in a coeherent way the terms of $`𝒪(\alpha _s\mathrm{\Lambda }_{QCD}/m_b)`$ is forced to introduce extra model-dependent non-perturbative parameters besides $`f_{K,\pi }`$, $`f_{K,\pi }(0)`$ and the hadronic wave functions. This implies that, at $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$, model dependence is unavoidable (even in the subsector of the chirally enhanced contributions) and it is present in both the analyses of ref. noilast  and ref. lorolast , which we will compare below.

Different model-dependent assumptions were made in the two approaches:

1. In ref. noilast , $`\stackrel{~}{E}_1`$ and $`\stackrel{~}{P}_3`$ were neglected and $`SU(2)`$ symmetry was assumed for $`\stackrel{~}{P}_1`$ ($`SU(3)`$ when it was used for $`B\pi \pi `$ decays). The same approximations were made for $`\stackrel{~}{P}_1^{\mathrm{GIM}}`$. The complex parameter $`\stackrel{~}{P}_1`$ ($`\stackrel{~}{P}_1^{\mathrm{GIM}}`$) was then fitted to reproduce the experimental $`BR`$s.
2. In ref. lorolast  the effects of the chirally enhanced $`\mathrm{\Lambda }_{QCD}/m_b`$ corrections were either computed perturbatively or encoded in two complex phenomenological parameters called $`X_H`$ and $`X_A`$. An uncertainty of 100 $`\%`$ to the “default” values (e.g. $`X_H=2.4`$) was assigned to these parameters in order to determine allowed bands for the predicted $`BR`$s. The bands include all other sources of uncertainties.

We now compare the leading amplitudes to the experimental results in order to test factorization. Before using predictions based on factorization to test the Standard Model and look for signals of new physics, it is crucial to check how large are the errors induced by our ignorance of the $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$ corrections which we are unable to compute. Our position, indeed, is that we have more confidence in the SM rather than in factorization. In order to test factorization, we ought to use all the information that we have from other measurements. Thus, for example, whereas the size of the error on $`|V_{cb}|`$ and $`|V_{ub}|`$ can be debated, there is no question that these experimental inputs must be included in any analysis that aims at testing (or using) factorization. We also stress that the value of $`|V_{ub}|`$ is not expected to be affected by the presence of new physics beyond the SM.

The CP parameter $`\gamma `$ does in general change if there is physics beyond the SM. It remains an interesting exercise, however, to verify whether, by taking the value of $`\gamma `$ from the Unitarity Triangle Analysis (UTA) in the SM, the predicted $`BR`$s are in agreement with the data. If, by using $`\gamma `$ from the UTA, one is unable to reproduce the experimental $`B\pi \pi `$ and $`BK\pi `$ $`BR`$s, this implies that “either there is new physics or $`\mathrm{\Lambda }_{QCD}/m_b`$ corrections are important” noilast .

In our analysis we have used the likelihood method which has been described in all details in ref. Ciuchini:2000de . Without entering in the “ideological” controversy about frequentistic and bayesian methods, we only note here that in Ciuchini:2000de  it has been shown that, at 95 % C.L., the Bayesian analysis give the same results as the frequentistic Babar Scanning method (and its variations) when the same inputs are used. Thus we will present our contour plots, corresponding to fig. 17 of ref. lorolast , both with factorization plus chirally enhanced contributions and with the non-factorizable charming (and GIM) penguin corrections. Besides this, we will also give tables with the relevant $`BR`$s, both in the factorization approximation with chirally enhanced terms and with the charming (and GIM) penguin corrections included.

We end this section with some remarks. In our approach, we have first checked that, within factorization and the SM, it is impossible to fit the experimental $`BR`$s. The $`\mathrm{\Lambda }_{QCD}/m_b`$ terms, that we are then forced to include in order to reproduce the experimental results, are non-perturbative quantities, infrared divergent in perturbation theory, on which we do not have any knowledge a priori. For this reason we decided to fit them on the data. The experimental numbers are nicely reproduced and the corrections to factorization are well consistent with the expected size (i.e. $`\stackrel{~}{P}_1`$ is of $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$ with respect to the leading contributions).

In ref. lorolast  the subleading power corrections are varied in predefined intervals and the change in the predicted $`BR`$s is interpreted as uncertainty on the factorized amplitudes. In our opinion the uncertanties on factorization are only those coming from the CKM matrix elements or the form factors, etc. The $`\mathrm{\Lambda }_{QCD}/m_b`$ terms instead are really contributions beyond factorization: if they are necessary to reproduce the data then it is not possible to make model independent predictions. This would remain true even if we knew without any uncertainty the hadronic parameters entering at the lowest order of the $`\mathrm{\Lambda }_{QCD}/m_b`$ expansion 333 Note that some $`BR`$s are dominated by the $`\mathrm{\Lambda }_{QCD}/m_b`$ corrections..

Thus we are in the Bermuda triangle: i) without the $`𝒪(r_\chi )`$ and $`𝒪(\alpha _sr_\chi )`$ terms, that is within QCD factorization, we cannot reproduce the data; ii) the inclusion of the computable subset of $`𝒪(\alpha _sr_\chi )`$ terms only is inconsistent since there is no reason to exclude the other non-perturbative non-computable contributions of the same order. In any case, we will show that, by using all the available experimental information, also this case is very difficult to reconcile with the data; iii) the complete set of $`𝒪(\alpha _sr_\chi )`$ corrections leads us beyond factorization and the results are model dependent. Indeed there is no proof that the one-loop finite chirally enhanced terms remain infrared finite at higher orders in $`\alpha _s`$. Moreover, if the corrections of $`𝒪(\alpha _sr_\chi \alpha _s\mathrm{\Lambda }_{QCD}/m_b)`$ are phenomenologicaly important, it is difficult to understand why, at the same level of numerical accuracy, other non-chirally enhanced non-perturbative $`\mathrm{\Lambda }_{QCD}/m_b`$ terms should not also be taken into account.

The sad conclusion is that the very nice theory of QCD factorization developed in beneke  is insufficient to fit the data because power corrections, which are model dependent, are important in $`BK\pi `$ and $`B\pi \pi `$ decays. Finally, model dependence does not implies that we are unable to make any prediction. If the assumptions made in our approach are reasonable 444 For example that we may neglect $`\stackrel{~}{E}_1`$., by fitting $`\stackrel{~}{P}_1`$ to the data and with the increasing experimental precision we may hope to extract also the value of $`\gamma `$ or to constrain $`\mathrm{sin}2\alpha `$.

## 3 Results

In this section we present our analysis and a detailed comparison with the results of ref. lorolast . In order to obtain our results we used the likelihood method as in noilast . The input parameters, given in table 1, are also the same but for a few differences:

* We added the $`𝒪(\alpha _s)`$ corrections to the coefficients of the penguin and electropenguin operators computed in lorolast , which appeared after the completion of noilast . In this respect the criticisms of ref. beuc  do not apply to the present analysis. We will see below that, even including these new ingredients, the main physics conclusions of ref. noilast  are confirmed.
* As for the matrix elements of $`Q_{6,8}`$, we include the “computable” factorizable chirally enhanced terms in the definition of $`\stackrel{~}{P}_1`$, which in this way contains all the possible Cabibbo enhanced, model-dependent corrections of $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$. In practice this corresponds to fit $`\stackrel{~}{P}_1`$ as in ref. noilast  with $`r_\chi ^{K,\pi }=0`$. There is no substantial difference between the old choice and the new one, since the contribution of the chirally enhanced terms and of $`\stackrel{~}{P}_1`$ have exactly the same quantum numbers.
* Although this is a minor source of uncertainty we also allow a variation of the renormalization scale and of the parameters of the hadronic wave functions.

We have also verified that by using the inputs of ref. lorolast  we obtain essentially the same results, and hence arrive to the same physics conclusions.

In this section we present:

1. A brief discussion of the results obtained in QCD factorization, namely with all the terms of $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$ set to zero.
2. The results obtained in our approach by fitting $`\stackrel{~}{P}_1`$ and using $`\gamma `$ as determined from the Unitary Triangle Analysis in ref. Ciuchini:2000de

$$\gamma =(54.8\pm 6.2)^o.$$
(4)
On the basis of this study we predict the $`B`$$`\overline{B}`$ asymmetries of the $`BR`$s, such as

$$𝒜(B_d\pi ^+\pi ^{})=\frac{BR(\overline{B}_d^0\pi ^+\pi ^{})BR(B_d^0\pi ^+\pi ^{})}{BR(\overline{B}_d^0\pi ^+\pi ^{})+BR(B_d^0\pi ^+\pi ^{})};$$
(5)
3. The results obtained by letting $`\gamma `$ free and a comparison of our results with those of ref. lorolast .

### Results in QCD factorization

In this subsection we compare the model independent results obtained with QCD factorization, namely including only the terms which survive when $`m_b\mathrm{}`$, with the experimental data. In this case there is still a residual model dependence due to our ignorance of the semileptonic form factors $`f_{K,\pi }(0)`$ and, at order $`\alpha _s`$, to the ignorance of the hadron wave functions. We vary the semileptonic form factors with flat p.d.f. in the intervals given in table 1 and the parameters of the hadron wave functions in the intervals given in table 2 of ref. lorolast . We also vary the renormalization scale between $`m_b/2`$ and $`2m_b`$, see table 1 of lorolast . The data show a generalized disagreement with the QCD factorization predictions. In particular the allowed region in the $`\rho `$$`\eta `$ plane and the value of $`\gamma `$ do not have any overlap with the corresponding ones from the unitarity triangle analysis Ciuchini:2000de . This remains true even if we double the uncertainty on $`|V_{ub}|`$. We conclude that QCD factorization cannot be reconciled with data.

### Factorization with Charming and GIM penguins

We now discuss the effects of charming and GIM penguins, parameterized by $`\stackrel{~}{P}_1`$ and $`\stackrel{~}{P}_1^{\mathrm{GIM}}`$. $`\stackrel{~}{P}_1`$ is a complex amplitude that we fit on the $`BK\pi `$ $`BR`$s. In order to have a reference scale for its size, we introduce a suitable “Bag” parameter, $`\widehat{B}_1`$, by writing

$$\stackrel{~}{P}_1=\frac{G_F}{\sqrt{2}}f_\pi f_\pi (0)g_1\widehat{B}_1,$$
(6)
where $`G_F`$ is the Fermi constant, $`g_1`$ is a Clebsh-Gordan parameter depending on the final $`K\pi `$ ($`\pi \pi `$) channel and $`B_1=|B_1|\mathrm{exp}(i\varphi )`$. Note that $`\widehat{B}_1`$ differs from the parameter defined in ref. noilast  because it now includes all the chirally enhanced $`\mathrm{\Lambda }_{QCD}/m_b`$ corrections, part of which were previously explicitly calculated using factorization. In a similar way we introduce $`\widehat{B}_1^{\mathrm{GIM}}`$.

We fit all the $`BR`$s given in table 1 with GIM and charming penguins included and taking the value of $`\gamma `$ determined from the UTA, see eq. (4). We find

$`|\widehat{B}_1|`$ $`=`$ $`0.13\pm 0.02,\varphi =(188\pm 82)^o,`$
$`|\widehat{B}_1^{\mathrm{GIM}}|`$ $`=`$ $`0.17\pm 0.08,\varphi ^{\mathrm{GIM}}=(181\pm 59)^o,`$ (7)
where the notation is self-explaining. Note that the size of the charming and GIM penguin effects is of the expected magnitude. Note that $`|\widehat{B}_1^{\mathrm{GIM}}|`$ is very poorly determined. The reason is that only $`BR(B_d^0\pi ^+\pi ^{})`$ in the fit is sensitive to GIM penguins. Thus, in practice, we are trying to fit two parameters, namely $`|\widehat{B}_1^{\mathrm{GIM}}|`$ and $`\varphi ^{\mathrm{GIM}}`$, to a single $`BR`$.

The results for the $`BR`$s can be found in table 2 with the label “Charming+GIM”. They show that the extra charming and GIM parameters radically improve the agreement for the measured $`BK\pi `$ and $`B\pi \pi `$ $`BR`$s. We do not claim, however, to be able to predict $`BR(B_d\pi ^+\pi ^{})`$ since many effects of the same order besides charming and GIM contributions, which in this case are not Cabibbo enhanced, were ignored: our results instead show that accurate predictions for $`B_d\pi \pi `$ decays can only be obtained by controlling quantitatively the $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$ corrections, which is presently far beyond the theoretical reach. To give a complete information, and for comparison with ref. lorolast  we also fit the data by letting $`\gamma `$ free. In this case we obtain $`\gamma =(89\pm 42)^o`$. At present, the precision of the data and the number of free parameters does not allow a useful determination of $`\gamma `$.

The large absolute values of $`\varphi `$, and the sizable effects that penguins have on the $`BR`$s, stimulated us to consider whether we could find observable particle-antiparticle asymmetries as the one defined in eq. (5). We find large effects in $`BR(B^+K^+\pi ^0)`$, $`BR(B_dK^+\pi ^{})`$ and $`BR(B_d\pi ^+\pi ^{})`$. As discussed before, for $`BR(B_d\pi ^+\pi ^{})`$ our predictions suffer from very large uncertainties due to contributions which cannot be fixed theoretically. For this reason, the values of the asymmetry reported in table 3 are only an indication that a large asymmetry could be observed also in this channel. There is a sign ambiguity in $`𝒜\mathrm{sin}\gamma \mathrm{sin}\varphi `$. This ambiguity can be solved only by an experimental measurement or, but this is extremely remote, by a theoretical calculation of the relevant amplitudes. For each channel, we give the absolute value of the asymmetry in table 3. Note that within factorization all asymmetries would be unobservably small, since the strong phase is a perturbative effect of $`𝒪(\alpha _s)`$ beneke . The possibility of observing large asymmetries in these decays opens new perspectives. These points will be the subject of a future study.

### Comparison with BBNS

In this section we make a critical comparison with the latest analysis of BBNS lorolast . We recall that we added all the perturbative corrections and allowed the variations of the non-perturbative parameters which were implemented by BBNS.

We start by discussing the rôle of chirally enhanced terms that do not suffer from infrared divergences at $`𝒪(\alpha _s)`$. This is an istructive case since, if not for other infrared-divergent contributions of the same order, the inclusion of only these terms, although not justified, would still allow model independent predictions, in the sense discussed before. Thus the question is whether these chirally enhanced terms alone, part of which have the same effect as $`\stackrel{~}{P}_1`$, can describe the data. Thus we have repeated the analysis in the UTA case with only factorization and chirally enhanced infrared-finite terms with the results for the $`BR`$s given in table 2 with the label “chirally enhanced”. The combined probability that all the predicted values for the $`K\pi `$ channels are within two $`\sigma `$s from the experimental numbers is 6%. If we relax the constraint on $`\gamma `$ we obtain $`\gamma =(127\pm 20)^o`$ with a probability of 0.3%that $`\gamma <80^o`$ and of 2% that $`\gamma <90^o`$. The corresponding $`\rho `$$`\eta `$ contour plot is given in fig. 1 (upper–left). We conclude that this model (model in the sense of including only a chosen subset of the chirally enhanced terms) is strongly disfavoured by the data. The reader may be surprised of the difference between fig. 1 and fig. 17 of ref. lorolast  (where the non-perturbative parameters $`X_{A,H}`$ were however included). The difference is explained by the fact that we used the experimental measurements of $`|V_{ub}|`$. If, following ref. lorolast ; beuc , we let both $`\rho `$ and $`\eta `$ free, the $`\rho `$$`\eta `$ contour plot changes and becomes that shown in fig. 1 (upper-right), which is very similar to fig. 17 of ref. lorolast . For completeness, we also give in fig. 1 the p.d.f. of $`\gamma `$ in the two cases, together with the p.d.f. obtained with the UTA analysis of ref. Ciuchini:2000de . When the experimental information on $`|V_{ub}|`$ is used the two p.d.f. have no overlap (lower-left). In the other case, it is still possible to find values of $`\gamma `$ compatible with the UTA analysis (lower-right), at the price of a rather low value of $`|V_{ub}|`$. The situation is well illustrated by fig. 2, where the contour plot in the $`\gamma `$$`|V_{ub}|`$ plane of the joint p.d.f. from non-leptonic decays is compared with the UTA range for $`\gamma `$ and the allowed interval for the measured $`|V_{ub}|`$, at 1$`\sigma `$. This figure demonstrates that, given the strong correlation between $`\gamma `$ and $`|V_{ub}|`$, it is crucial to take into account the experimental knowledge of $`|V_{ub}|`$. An important theoretical remark is in order at this point. The coefficient $`a_6^f`$ is enhanced by the $`𝒪(\alpha _s)`$ corrections and may play the same rôle of charming penguins. It is very scaring that its actual value is strongly affected by the contribution of the chromomagnetic operator computed at tree level lorolast . It is very hard to believe that this contribution which, besides all general considerations, is also non-local, will remain infrared-finite in higher orders and can be really evaluated in this way.

We now discuss the effect of $`X_{A,H}`$ on the final results. The addition of these parameters, in the range proposed in ref. lorolast , leaves the situation substantially unaltered as shown in fig. 3. This holds true also for the $`BR`$s, which are given in table 2 with the label BBNS. The reason is that on the one hand $`X_A`$ and $`X_H`$ do not have the same quantum numbers, and hence effects, of charming penguins, on the other the range chosen a priori, on the basis on one-loop perturbation theory, is not large enough to improve the agreement with the measured $`BR`$s (it essentially increases the uncertainty on the predictions). Of course by choosing a low value of the strange quark mass and of $`|V_{ub}|`$ (without using the experimental information coming from its measurements), a large value of $`f_K(0)`$ and a small value for $`f_\pi (0)`$ etc. it is still possible to find some point in parameter space where the $`\chi ^2`$ is good. That this can be used to fit $`\gamma `$ and test the SM is hard to believe though.

The situation would be different if $`X_{A,H}`$ are let free to vary and fitted to the data. We have done this exercise and found that the preferred value of $`\rho _A`$ is much larger than the values allowed in the interval chosen in lorolast , whereas $`\rho _H`$ is not determined by the fit. We conclude that since $`X_A`$ and $`X_H`$ are infrared divergent quantities, the value of which cannot be predicted, and since without the inclusion of non-perturbative contributions of $`𝒪(\mathrm{\Lambda }_{QCD}/m_b)`$ is not possible to reproduce the experimental data, we are bound to use model dependent assumptions in the analysis of non-leptonic $`BK\pi `$ and $`B\pi \pi `$ decays.

## 4 Conclusion

We have analyzed the predictions of QCD factorization for $`B\pi \pi `$ and $`BK\pi `$ decays. Even taking into account the uncertainties of the input parameters, we find that QCD factorization is unable to reproduce the observed $`BR`$s. The introduction of charming and GIM penguins charming  allows to reconcile the theoretical predictions with the data. Instead of varying the non-perturbative phenomenological parameters in preassigned ranges, we prefer to try to fit them on the data. With the present theoretical and experimental accuracy, we find that it is still not possible to determine the CP violation angle $`\gamma `$. The situation is expected to improve in the near future with more accurate experimental measurements of the relevant $`BR`$s.

Contrary to factorization, we predict large asymmetries for several of the particle–antiparticle $`BR`$s, in particular $`BR(B^+K^+\pi ^0)`$, $`BR(B_dK^+\pi ^{})`$ and, possibly, $`BR(B_d\pi ^+\pi ^{})`$. This opens new perspectives for the study of CP violation in $`B`$ systems.

## Acknowledgments

We thank C. Sachrajda for useful discussions on this work. G.M. and L.S. thanks the TH division at CERN where part of this work has been done.


text:
# The Chiral Critical Point in 3-Flavour QCD The work has been supported by the TMR network ERBFMRX-CT-970122 and by the DFG under grant FOR 339/1-2.

## 1 INTRODUCTION

The phase transition in 3-flavour QCD is first order in the chiral limit ($`m0`$). This first order transition will continue to persist for small but non-zero values of $`m`$ up to a critical value, $`\overline{m}`$, of the quark mass. It is expected that QCD at this chiral critical point belongs to the universality class of the 3d Ising model . We will present here evidence for this scenario .

The lattice formulation of QCD with three degenerate quark masses depends on the bare quark mass $`m`$ and the gauge coupling $`\beta 6/g^2`$. The thermodynamics is described in terms of the partition function $`Z(\beta ,m)=𝒟Ue^{S(\{U\},\beta ,m)}`$. The Euclidean action, $`S`$, is given in terms of the pure gauge ($`S_G`$) and the fermion action ($`S_F`$)

$$S(\{U\},\beta ,m)=\beta S_G(\{U\})S_F(\{U\},m),$$
(1)
where $`S_F=0.75\text{Tr}\mathrm{ln}M`$ is expressed in terms of the fermion matrix $`M(\{U\})`$. We will mainly use standard staggered fermions and the standard Wilson gauge actions.

The universal properties at the chiral critical point are controlled by an effective Hamiltonian,

$$_{eff}(\tau ,\xi )=\tau +\xi ,$$
(2)
with a global $`Z(2)`$ symmetry. It controls the critical behaviour at $`(\beta _c(\overline{m}),\overline{m})`$ and can be expressed in terms of two operators $`,`$, i.e. the energy-like and ordering-field like operators that couple to two relevant scaling fields $`\tau `$ and $`\xi `$. As the symmetry of $`_{eff}`$ is not shared in any obvious way by the QCD Lagrangian one may expect that in the vicinity of the chiral critical point the operators and couplings appearing in the QCD Lagrangian are linear combinations of $``$, $``$ and $`\tau `$, $`\xi `$ respectively. One may use a linear ansatz

$`\tau `$ $`=`$ $`\beta \beta _c+A(m\overline{m}),`$
$`\xi `$ $`=`$ $`m\overline{m}+B(\beta \beta _c),`$ (3)
as well as

$$=S_G+r\overline{\psi }\psi ,=\overline{\psi }\psi +sS_G.$$
(4)
Here $`\overline{\psi }\psi 0.75\text{Tr}M^1`$ denotes the chiral condensate expressed in terms of bosonic observables. As $`S_G`$, $`\overline{\psi }\psi `$ or related observables like the Polyakov loop $`L`$ are mixtures of $``$ and $``$, the corresponding susceptibilities will all receive contributions from fluctuations of $``$ as well as $``$. Asymptotically therefore all of them will show identical finite size scaling behaviour which will be dominated by the largest universal exponent $`\gamma /\nu `$. Unfortunately this is quite similar for the 3d-$`Z(2)`$ and $`O(N)`$ spin models, so that a finite size scaling analysis of susceptibilities is not a good indicator for the universality class. The situation is different for cumulants constructed from linear combinations of $`\overline{\psi }\psi `$ and $`S_G`$,

$$B_4(x)=\frac{\left(\delta M(x)\right)^4}{\left(\delta M(x)\right)^2^2},M(x)=\overline{\psi }\psi +xS_G.$$
(5)
Here $`\delta X`$ denotes $`XX`$. For all values of $`x`$ different from $`1/r`$ the cumulants behave asymptotically like the Binder cumulant for the order parameter; cumulants calculated on different size lattices for different quark masses will intersect at some value of the quark mass. In the infinite volume limit these intersection points will converge to a universal value which is characteristic for the universality class, in fact, it is quite different for the classes of 3d-$`Z(2)`$ and $`O(N)`$ symmetric spin models; e.g. $`B_4=1.604`$ for $`Z(2)`$ , 1.242(2) for $`O(2)`$  and 1.092(3) for $`O(4)`$ . The cumulants $`B_4(x)`$ thus seem to be appropriate observables to locate the chiral critical point as well as to determine its universality class without knowing in detail the correct scaling fields.

## 2 LOCATING THE CRITICAL POINT

We make use of the finite size scaling properties of Binder cumulants $`B_4(x)`$ evaluated at $`\beta _{pc}(m)`$. From previous studies one knows that the endpoint in 3-flavour QCD is located close to $`m=0.035`$ . We have performed calculations on lattices of size $`N_\sigma ^3\times 4`$, with $`N_\sigma =8`$, 12 and 16. We have used four values of the quark mass in the interval $`m[0.03,0.04]`$ and for each of these masses we calculated thermodynamic observables for 3 to 4 different values of the gauge coupling $`\beta `$. In general we collected for each pair of couplings $`(13)10^4`$ configurations generated with the hybrid-R algorithm111We used trajectories of length $`\tau =0.675`$ generated with a discrete step size $`\delta \tau =0.015`$.. Interpolations between results from different $`\beta `$-values have been performed using the Ferrenberg-Swendsen multi-histogram technique. The cumulant of the chiral condensate, $`B_4(0)`$, is shown in Fig. 1. We note that the cumulants calculated on different size lattices intersect at a quark mass close to $`m=0.035`$. The value of $`B_4(0)`$ at the intersection point is compatible with the universal value of the Binder cumulant for the 3d Ising model. We have fitted the cumulants on a given lattice size using a linear ansatz in the quark mass. From this we find for the critical mass and the cumulant,

$$\overline{m}=0.033(1),B_4(0)=1.639(24).$$
(6)
This result strongly suggests that the chiral critical point indeed belongs to the universality class of the 3d Ising model.

As discussed above the intersection point of Binder cumulants constructed from $`M(x)`$ according to Eq. 5 will be independent of the choice of $`x`$ in the infinite volume limit. The optimal choice $`x=s`$, will minimize the finite volume effects, which are quadratic in $`\mathrm{\Delta }=xs`$. The determination of the $`x`$-value for which the intersection point of the Binder cumulant is closest to the $`Z(2)`$ value provides a finite volume estimate for the mixing parameter $`s`$. One finds $`s_{\mathrm{min}}=0.430(23)`$.

## 3 THE MIXING PARAMETERS

The parameters $`r`$ and $`s`$ are fixed by demanding that $``$ should obey basic properties of the order parameter for spontaneous symmetry breaking at the critical point. We obtain two conditions for the mixing parameters $`r`$, $`s`$ in terms of the parameter $`B`$ and expectation values of $`S_G`$, $`\overline{\psi }\psi `$,

$`r`$ $`=`$ $`B,`$
$`s`$ $`=`$ $`{\displaystyle \frac{\delta \overline{\psi }\psi \delta S_GB(\delta \overline{\psi }\psi )^2}{(\delta S_G)^2B\delta \overline{\psi }\psi \delta S_G}}.`$ (7)
The parameter $`B`$ is needed to define lines of constant $`\xi `$. The line of first order phase transitions defines the zero external field line ($`\xi =0`$) of the effective Hamiltonian. We thus can extract $`B`$ from the quark mass dependence of the pseudo-critical couplings. Using Eq. 3 one obtains

$$B^1=\left(\frac{\mathrm{d}\beta _c(m)}{\mathrm{d}m}\right)|_{m=\overline{m}}.$$
(8)
Knowing $`B`$ we also know $`r`$ and can construct $`s`$ using Eq. 7. A first estimate may be given using our data on the largest lattice ($`16^3\times 4`$). The slope of $`\beta _{pc}(m)`$ defines $`r^1`$. We estimate $`r=0.51(2)`$ from a straight line fit. The mixing parameter $`r`$ is large and definitely non-zero. The importance of choosing the correct mixing parameter $`r`$ becomes apparent from an analysis of joint probability distributions for $`\delta `$ and $`\delta `$. These are little affected by changes of $`s`$ but they strongly depend on $`r`$. In order to improve our determination of $`r`$, we demand that fluctuations in $``$ vanish for any fixed value of $`\delta `$. This maximizes the $`Z(2)`$ symmetry of the contour plot shown in Fig.2.

One obtains in this way $`r=0.550(7)`$. A determination of the parameter $`s`$ from Eq. 7 then yields $`s=0.41(51)`$. Although this value has large errors it is consistent with the result found from the intersection point of the Binder cumulant. As a best estimate we therefore obtain

$$r=\mathrm{\hspace{0.33em}0.550}\pm 0.007,s=\mathrm{\hspace{0.33em}0.430}\pm 0.023.$$
(9)
## 4 THE PHYSICAL SCALE

In order to determine a physical scale for the chiral critical point we have calculated hadron masses on a $`16^4`$ lattice at $`(\beta _c,\overline{m})`$. For the pseudo-scalar meson mass we find $`m_{\mathrm{ps}}/T_c=1.853(1)`$. Using estimates for the critical temperature  we thus estimate for the pseudo-scalar meson mass at the chiral critical point $`m_{\mathrm{ps}}290`$ MeV. The entire analysis of the critical point discussed so far has been performed with unimproved gauge and fermion actions on rather coarse lattices. Improved actions are not expected to modify the results on the universal properties of the critical point. They may, however, well influence the quantitative determination of the critical point. We therefore have investigated the critical point also in calculations with improved gauge and staggered fermion actions, the p4-action . We have performed calculations on a lattice of size $`12^3\times 4`$ with a bare quark mass $`m=0.005`$ and on a lattice of size $`16^3\times 4`$ with $`m=0.01`$. Again a Ferrenberg-Swendsen reweighting is used to determine the pseudo-critical couplings and Binder cumulants at these couplings. The smaller quark mass leads to a Binder cumulant below the $`Z(2)`$ value ($`B_4(0)=1.31(12)`$) and for the larger quark mass, the Binder cumulant is significantly above the $`Z(2)`$ value ($`B_4(0)=2.14(10)`$). We thus have obtained an upper and lower limit for the critical point, suggesting a critical bare quark mass of $`\overline{m}=0.0075(25)`$. Extrapolating the meson masses calculated in  to this value of the bare quark mass, we estimate for the pseudo-scalar meson mass at the critical point $`m_{\mathrm{ps}}192(25)`$ MeV. The critical mass thus is considerably smaller than estimated from calculations with unimproved gauge and fermion actions.


text:
# Color Evaporation and Elastic Υ Photoproduction at DESY HERA

## Acknowledgments

This work was partially financed by CNPq, CAPES and FAPERGS, BRAZIL.


text:
# Resubmitted to Physical Review Letters January 2002 Statically Transformed Autoregressive Process and Surrogate Data Test for Nonlinearity

## Abstract

The key feature for the successful implementation of the surrogate data test for nonlinearity on a scalar time series is the generation of surrogate data that represent exactly the null hypothesis (statically transformed normal stochastic process), i.e. they possess the sample autocorrelation and amplitude distribution of the given data. A new conceptual approach and algorithm for the generation of surrogate data is proposed, called statically transformed autoregressive process (STAP). It identifies a normal autoregressive process and a monotonic static transform, so that the transformed realisations of this process fulfill exactly both conditions and do not suffer from bias in autocorrelation as the surrogate data generated by other algorithms. The appropriateness of STAP is demonstrated with simulated and real world data.

The surrogate data test for nonlinearity has been widely used in real applications in order to establish statistically the existence of nonlinear dynamics and justify the use of nonlinear tools in the analysis of time series Kantz and Schreiber (1997); Diks (2000). The most general null hypothesis $`\text{H}_0`$ of this test so far is that the examined time series $`𝐱=[x_1,\mathrm{},x_n]^{}`$ is a realisation of a normal (linear) process $`𝐬=[s_1,\mathrm{},s_n]^{}`$ undergoing a possibly nonlinear static transform $`h`$, i.e. $`x_i=h(s_i)`$, $`i=1,\mathrm{},n`$. To test the null hypothesis, an estimate $`q`$ from a nonlinear method applied to the original data, $`q_0`$, is compared to the estimates, $`q_1,\mathrm{},q_M`$ on $`M`$ surrogate time series representing the null hypothesis Schreiber and Schmitz (2000); Kugiumtzis (2002). A properly designed surrogate time series $`𝐳`$ should possess the same autocorrelation as the original data, $`r_z(\tau )=r_x(\tau )`$ for a range of lags $`\tau `$, and the same amplitude distribution, $`F_z(z_i)=F_x(x_i)`$ ($`F_x(x_i)`$ is the cumulative density function (cdf) of $`x_i`$), and be otherwise random.

It has been reported that erroneous results are likely to occur mainly due to the insufficiency of the algorithms to generate surrogate data that preserve the original linear correlations Kugiumtzis (2001). The prominent algorithm of amplitude adjusted Fourier transform (AAFT) Theiler et al. (1992), used in most applications so far, is built based on the assumption of monotonicity of $`h`$. When $`h`$ is not monotonic, the AAFT algorithm is found to favour the rejection of $`\text{H}_0`$ due to the mismatch of the original linear correlations Kugiumtzis (1999). The iterated AAFT (IAAFT) algorithm improves the match of the autocorrelation of AAFT Schreiber and Schmitz (1996), but with about the same accuracy for all the surrogates, so that the small variance, combined with the small bias, may be another cause for false rejections in some cases Kugiumtzis (1999). Another algorithm making use of simulation annealing seems to perform similarly to IAAFT Schreiber (1998). Recently, a correction of AAFT (CAAFT) that results in unbiased match of the linear correlations was proposed in Kugiumtzis (2000).

In this paper, we develop further the somehow profound rationale of the CAAFT algorithm and formulate a new conceptual approach for the generation of surrogate data consistent with $`\text{H}_0`$, which we then solve analytically. The new algorithm, called STAP, generates the surrogate data as realisations of a suitable statically transformed autoregressive process (STAP), i.e. the process under $`\text{H}_0`$ is designed as a static transform of a suitable normal process.

The main idea behind the STAP algorithm is that for any stationary process $`𝐗=[X_1,X_2,\mathrm{}]^{}`$, for which we have finite scalar measurements $`𝐱`$, there is a scalar linear stochastic process $`𝐙`$ with the same autocorrelation $`\rho _X`$ and marginal cdf $`\mathrm{\Phi }_X`$ as for the observed process ($`\rho _Z(\tau )=\rho _X(\tau )`$ and $`\mathrm{\Phi }_Z(z_i)=\mathrm{\Phi }_X(x_i)`$), i.e. $`𝐙`$ is a scalar “linear copy” of the observed $`𝐗`$. The objective is to derive $`𝐙`$ through a static monotonic transform $`g`$ on a scalar normal process $`𝐔`$ with a proper autocorrelation $`\rho _U`$, i.e. $`Z_i=g(U_i)`$. Thus $`g`$ and $`\rho _U`$ have to be properly selected, so that $`𝐙`$ has the desired properties. In practice, the surrogate data set is a finite realisation $`𝐳=[z_1,\mathrm{},z_n]^{}`$ of the process $`𝐙`$ and $`g`$ and $`\rho _U`$ are estimated based solely on $`𝐱`$. Note that $`g`$ and $`𝐔`$ are in general different from $`h`$ and $`𝐒`$, respectively, of $`\text{H}_0`$ (they are the same if $`h`$ is monotonic Kugiumtzis (2000)). Thus with this approach $`\text{H}_0`$ can be formulated more generally that the time series is generated by a linear stochastic process.

Let $`\mathrm{\Phi }_0`$ be the marginal cdf of a standard normal process $`𝐔`$. A suitable choice for $`g`$, so that $`\mathrm{\Phi }_Z(z_i)=\mathrm{\Phi }_X(x_i)`$, is defined as Papoulis (1991)

$$Z_i=g(U_i)=\mathrm{\Phi }_X^1(\mathrm{\Phi }_0(U_i)),$$
(1)
where $`g`$ is monotonic by construction. Assuming that $`\mathrm{\Phi }_X(x_i)`$ is continuous and strictly increasing and that $`1<\rho _U<1`$, which are both true for all practical purposes, there is a function $`\varphi `$ depending on $`g`$, such that $`\rho _Z=\varphi (\rho _U)`$ for any lag $`\tau `$ Johnson and Kotz (1972); Der Kiureghian and Liu (1986). If $`g`$ has an analytic form, then it may be possible to find an analytic expression for $`\varphi `$ as well. In that case, given that $`\rho _X`$ is known and by setting $`\rho _Z\rho _X`$, one can invert $`\varphi `$ to find $`\rho _U=\varphi ^1(\rho _X)`$, if $`\varphi ^1`$ exists.

In general, the function $`g`$, as defined in eq.(1), does not have an analytic form because $`\mathrm{\Phi }_X`$ is not known analytically, but it can be approximated by an analytic function, e.g. a polynomial $`p_m`$ of degree $`m`$,

$$Z_i=g(U_i)p_m(U_i)=a_0+\underset{j=1}{\overset{m}{}}a_jU_i^j.$$
(2)
Low degree polynomials have been used to approximate such transforms Elphinstone (1983); Hawkins (1994). Then using the definition for the autocorrelation, the approximate expression for $`\varphi `$ reads

$$\rho _X=\varphi (\rho _U)=\frac{_{s=1}^m_{t=1}^ma_sa_t(\mu _{s,t}\mu _s\mu _t)}{_{s=1}^m_{t=1}^ma_sa_t(\mu _{s+t}\mu _s\mu _t)},$$
(3)
where an arbitrary lag $`\tau `$ is implied as argument for the autocorrelations, $`\mu _s`$ is the $`s`$th central moment of $`U_i`$ being $`\mu _{2k+1}=0`$, $`\mu _{2k}=13\mathrm{}(2k1)`$, $`k0`$, and $`\mu _{s,t}`$ is the $`s`$th-$`t`$th central joint moment of the bivariate standard normal distribution of $`(U_i,U_{i\tau })`$, defined as follows Hutchinson and Lai (1990)

$`\mu _{2k,2l}`$ $`=`$ $`{\displaystyle \frac{(2k)!(2l)!}{2^{k+l}}}{\displaystyle \underset{j=0}{\overset{\nu }{}}}{\displaystyle \frac{(2\rho _U)^{2j}}{(kj)!(lj)!(2j)!}}`$
$`\mu _{2k+1,2l+1}`$ $`=`$ $`{\displaystyle \frac{(2k+1)!(2l+1)!}{2^{k+l+1}}}{\displaystyle \underset{j=0}{\overset{\nu }{}}}{\displaystyle \frac{(2\rho _U)^{2j+1}}{(kj)!(lj)!(2j+1)!}}`$
$`\mu _{i,j}`$ $`=`$ $`0\text{if}k+l=\text{odd},`$
where $`\nu =\mathrm{min}(k,l)`$. By substituting the expression for the moments in eq.(3), the expression for $`\varphi `$ can be brought to a polynomial form of the same order $`m`$

$$\rho _X=\varphi (\rho _U)=\underset{j=1}{\overset{m}{}}c_j\rho _U^j,$$
(4)
where the vector of coefficients $`𝐜=[c_1,\mathrm{},c_m]^{}`$ is expressed only in terms of $`𝐚=[a_1,\mathrm{},a_m]^{}`$ (the expressions are rather involved and therefore not presented here). Simpler expressions can be derived using the Tchebycheff-Hermite polynomials Bhatt and Dave (1964). Thus an analytic expression for $`\rho _U`$ is possible if eq.(4) can be solved with respect to $`\rho _U`$. From our simulations, we conjecture that if $`g`$ is monotonic then $`\varphi `$ is also monotonic in $`[1,1]`$. Then $`\varphi ^1`$ exists and a unique solution for $`\rho _U`$ can be found from eq.(4). The proper standard normal process $`𝐔`$ is completely defined by $`\rho _U`$ and applying the transform $`g`$ of eq.(1) to the components of $`𝐔`$, the “linear copy” $`𝐙`$ of the given process $`𝐗`$ is constructed.

Note that the solution for $`\rho _U`$ is given analytically from the polynomial approximation of $`g`$ and it requires only the knowledge of the coefficients $`𝐚`$ of the polynomial and the autocorrelation $`\rho _X`$.

In practice, we operate with a single time series $`𝐱`$ rather than a process $`𝐗`$ and with the sample estimates $`F_x`$ and $`r_x`$ for $`\mathrm{\Phi }_X`$ and $`\rho _X`$, respectively. The steps of the algorithm are as follows:

1. Estimate the vector of coefficients $`𝐚=[a_1,\mathrm{},a_m]^{}`$ of the polynomial $`p_m`$ from the graph of $`x_i=F_x^1(F_0(w_i))`$, i.e. the graph of $`𝐱`$ vs $`𝐰`$ after their ranks are matched, where $`𝐰=[w_1,\mathrm{},w_n]^{}`$ is standard white normal noise.
2. Compute $`𝐜=[c_1,\mathrm{},c_m]^{}`$ for the given $`𝐚`$ from eqs.(34).
3. Find $`r_u`$ from eq.(4) for the given $`𝐜`$ and $`r_x`$ using the sample estimates $`r`$ instead of $`\rho `$. The common practice is that the solution exists and it is unique. If this is not the case, repeat the steps 1–3 for a new $`𝐰`$ until a unique solution is obtained.
4. Generate a realisation $`𝐮`$ of a standard normal process with autocorrelation $`r_u`$. We choose to do this simply by means of an autoregressive model of some order $`p`$, AR($`p`$). The parameters $`𝐛=[b_0,b_1,\mathrm{},b_p]^{}`$ of AR($`p`$) are found from $`r_u`$ using the normal equations solved effectively by the Levinson algorithm Brockwell and Davis (1991). The AR($`p`$) model is run to generate $`𝐮`$

$$u_{i+1}=b_0+\underset{j=1}{\overset{p}{}}b_ju_{i(j1)}+e_i,e_i\text{N}(0,1).$$
5. Transform $`𝐮`$ to $`𝐳`$ by reordering $`𝐱`$ to match the rank order of $`𝐮`$, i.e. $`z_i=F_x^1(F_0(u_i))`$.

Note that $`𝐮`$ possesses the sample normal marginal cdf $`F_0`$ and the proper $`r_u`$, so that $`𝐳`$ possesses $`F_z=F_x`$, $`r_z=r_x`$, and is otherwise random, as desired. In practice however, the equality $`r_z=r_x`$ is not exact and $`r_z`$ may vary substantially around $`r_x`$. Two possible reasons for this are the insufficient approximation of $`g`$ in step 1 and the inevitable variation of the sample autocorrelation of the generated $`𝐮`$ in step 4, which decreases with the increase of data size. The former is due to the limited power of polynomials in approximating monotonic functions and this shortcoming causes also occasional repetitions of the first steps of the algorithm as stated in step 3 111Other model classes may give better approximations of $`g`$ and should be investigated, but polynomials were used here to derive a simple analytic form for the transform of the autocorrelations.. The latter constitutes an inherent property of the so-called “typical realisation” approach (i.e. a model is used to generate the surrogate data) and cannot be controlled. However, less variation in the autocorrelation is achieved when the AR($`p`$) model is optimised making the following steps, in the same way as for the CAAFT algorithm Kugiumtzis (2000):

1. Apply the algorithm presented above $`K`$ times and get $`𝐳^1,\mathrm{},𝐳^K`$ surrogate time series.
2. Compute $`r_{z^1},\mathrm{},r_{z^K}`$ and find the one, $`r_{z^l}`$ closest to $`r_x`$ 222We found that a robust way to achieve this is to compute the error $`_{i=1}^\tau (r_x(i)r_{z^j}(i))^2`$ for $`j=1,\mathrm{},K`$ and $`\tau =1,\mathrm{},\tau _{\text{max}}`$ and then select the trial $`l`$ that gives the minimum error most times, where $`\tau _{\text{max}}`$ minimums are totally computed, each time over the $`K`$ trials..
3. Use the parameters $`𝐛`$ of the $`l`$-repetition to generate the $`M`$ surrogate data (steps 4–5 of the algorithm above).

The $`K`$ repetitions above as well as the occasional repetitions of steps 1–3 of the first part of the algorithm may slow down the algorithm if the time series is long, but they have no impact on the principal function of the algorithm. Simply, some realisations of white noise $`𝐰`$ are discarded in the search of the parameters $`𝐛`$ of the most suitable AR model that generates the surrogate data (through the $`g`$ transform).

The free parameters of the STAP algorithm are the degree $`m`$ of the polynomial approximation of $`g`$, the order $`p`$ of the AR model, the number $`K`$ of repetitions for the optimisation of AR($`p`$) and the maximum lag $`\tau _{\text{max}}`$, used to compare $`r_{z^1},\mathrm{},r_{z^K}`$ to $`r_x`$. Usually, a small $`m`$ ($`m10`$) is sufficient. For $`p`$, there is no optimal range of values but it may vary with the shape of $`r_x`$, e.g. a slowly decaying $`r_x`$ may be better modelled by a larger $`p`$. In all our simulations, we set $`K=M=40`$ and $`\tau _{\text{max}}=p`$.

The proper performance and superiority of STAP over AAFT and IAAFT was confirmed from simulations on different toy models. CAAFT was found to perform very similarly to STAP. We show in Fig. 1 comparative results for AAFT, IAAFT and STAP for three representative synthetic systems: the cube of an AR(1) process ($`s_{i+1}=0.3+0.8s_i+e_i,e_i\text{N}(0,1),x_i=s_i^3`$), the square of the same AR(1) ($`x_i=s_i^2`$), both being consistent with $`\text{H}_0`$, and the $`x`$ variable of the Rössler system Rössler (1976), not consistent with $`\text{H}_0`$.

For each system, we generate 100 time series of $`2048`$ samples each and for each realisation we generate $`M=40`$ surrogate data of each type.

As discriminating statistics $`q^i`$ we choose the correlation coefficient (CC) of the fit with the series of Volterra polynomials of degree 2 and order $`v=5`$. The polynomials for the first $`i`$ terms, where $`i=1,\mathrm{},5`$, are linear and for terms $`i=6,\mathrm{},20`$ are nonlinear (see also Kugiumtzis (2000)). To quantify the discrimination we use the significance $`S^i=\frac{|q_0^iq^i|}{\sigma _q^i}`$ for each polynomial of $`i`$ terms, where $`q_0^i`$ is the statistic $`q^i`$ on the original data, $`q^i`$ and $`\sigma _q^i`$ are the mean and standard deviation of the statistic $`q^i`$ on the $`M`$ surrogate data. The null hypothesis $`\text{H}_0`$ is formally rejected at the 0.05 significance level when $`S>1.96`$, under the assumption of normality for the statistic $`q`$, which turns out to hold in general. The percentages of rejections for each of the three systems are shown in Fig. 1. Very similar results were found when deciding the rejection from the rank ordering of $`q_0^i,q_1^i,\mathrm{},q_M^i`$.

For the linear statistics, STAP gives consistently and correctly no rejections, i.e. unbiased match of the linear correlations, whereas AAFT and IAAFT give very large percentages of rejections for all but the first case where AAFT gives about $`5\%`$ rejections, as expected. For IAAFT, the rejections occur because $`\sigma _q^i`$ is very small ($`10`$ to $`20`$ times smaller than for AAFT and STAP), though the bias $`q_0^iq^i`$ is smaller than for STAP.

For the nonlinear statistics, the feature of the linear statistics persists for the two first systems (consistent with $`\text{H}_0`$) and for all three algorithms, i.e. correctly no rejections with STAP and erronous rejections with AAFT and IAAFT. We cannot explain why the polynomial fit for the IAAFT surrogates improves with the addition of nonlinear terms (see Fig. 1). For the nonlinear system, STAP properly converges with the addition of few nonlinear terms to the correct $`100\%`$ rejection level, which AAFT and IAAFT possessed already with linear statistics.

Next, we verify the three algorithms on two human EEG data sets, one recorded many hours before an epileptic seizure accounting for normal brain activity and another recorded during an epileptic seizure. The epileptic EEG seems to exhibit a pattern of oscillations indicating a more deterministic character than the non-epileptic EEG, and this constitutes a well established result in physiology. This is demonstrated also with the fit of Volterra polynomials in Fig. 2. The fit improves with the inclusion of the first nonlinear terms for the epileptic EEG but not for the non-epileptic EEG. These findings are confirmed statistically by the test with STAP surrogate data while the results of the test with AAFT and IAAFT are more or less confusing.

In particular, the $`\text{H}_0`$ on the normal EEG is erroneously rejected with AAFT because the difference in the fit between original and surrogate data is about the same for the linear and nonlinear polynomials (see Fig. 2a and c). The same test result is obtained with IAAFT for large nonlinear polynomials, whereas again there is significant difference in the linear fits between original and IAAFT surrogates (not easily discernible as both bias and variance are very small). Using STAP surrogates, the $`\text{H}_0`$ is not rejected for both linear and nonlinear fits.

For the epileptic EEG, there is again a clear difference in the linear fit between original data and AAFT surrogates and a smaller but equally significant difference between original data and IAAFT surrogates (see Fig. 2b and c). The significance $`S`$ for both AAFT and IAAFT increases with the addition of the first couple of nonlinear terms, much more for IAAFT due to the small variance of CC. However, the deviation in the linear fit does not support reliable rejection of $`\text{H}_0`$. On the other hand, using STAP surrogates, the $`\text{H}_0`$ is properly rejected only for the nonlinear statistics and with high confidence ($`S<2`$ for the linear fit and $`S5`$ for the nonlinear fit).

In general, the test with STAP surrogates tends to be more conservative, i.e. small discriminations are found less significant, as the data size decreases. For example, the test on 296 sunspot samples (for which a small leap of the polynomial fit with the addition of nonlinear terms was observed) gave rejection of $`\text{H}_0`$ for AAFT and IAAFT but not for STAP (not shown here, see also Kugiumtzis (2000)). However, this should not be considered as a drawback of the STAP algorithm, as one expects that the power of the test reduces with the decrease of data size.

A new algorithm that generates surrogate data for the test for nonlinearity has been presented, called statically transformed autoregressive process (STAP). The key feature of STAP is that it finds analytically the autocorrelation of an appropriate underlying normal process for the test. This is the main difference of the STAP algorithm from the corrected AAFT (CAAFT) algorithm, where the autocorrelation is estimated numerically. Both CAAFT and STAP algorithms do not suffer from the severe drawback of the AAFT algorithm, i.e. bias in the match of the original autocorrelation. The AAFT algorithm is essentially impractical for real applications because it favours the rejection of $`\text{H}_0`$ as a result of the bias in the autocorrelation. From the numerical simulations, it turns out that the IAAFT algortihm may also give small bias in the linear correlations, favouring also the rejection of $`\text{H}_0`$. On the other hand, the STAP algorithm performs properly and gives reliable rejections of $`\text{H}_0`$, only whenever this appears to be the case.


text:
# One-loop renormalization of heavy-light currentsTalk presented by J. Harada.

## 1 INTRODUCTION

Although lattice QCD offers a nonperturbative method of calculating weak matrix elements from first principles, in practice a perturbative renormalization is also required to extract the continuum quantities for heavy-light systems. In this talk we discuss the renormalization of heavy-light vector and axial vector currents. These currents are needed for heavy quark phenomenology, such as the calculation of the decay constants and semi-leptonic form factors of heavy-light mesons. Here we calculate explicitly the mass dependent renormalization factors of heavy-light currents at one-loop order, when the heavy quark is treated with the Fermilab formalism . Results for the Wilson action have been obtained first in Ref.  and preliminary results for clover action have been reported in previous lattice conferences . For tree-level improvement at order $`1/m_Q`$, we include so-called rotation term here. Tadpole-improved renormalization factors are also presented. We also calculate mass dependent Brodsky-Lepage-Mackenzie scale $`q^{}`$ , with and without tadpole-improvement. More details of this work will be given in Ref. .

## 2 ONE-LOOP RESULTS

The renormalization factors $`Z_{J_\mathrm{\Gamma }}`$ of heavy-light currents are simply the ratio of the lattice and continuum radiative corrections:

$`Z_{J_\mathrm{\Gamma }}={\displaystyle \frac{[Z_{2h}^{1/2}\mathrm{\Lambda }_\mathrm{\Gamma }Z_{2l}^{1/2}]^{cont}}{[Z_{2h}^{1/2}\mathrm{\Lambda }_\mathrm{\Gamma }Z_{2l}^{1/2}]^{lat}}},`$ (1)
where $`Z_{2h}`$ and $`Z_{2l}`$ are wave-function renormalization factors of the heavy and light quarks, and the vertex function $`\mathrm{\Lambda }_\mathrm{\Gamma }`$ is the sum of one-particle irreducible three-point diagrams. We calculate explicitly $`Z_A`$ and $`Z_V`$ at one-loop order of perturbation theory.

In view of the mass dependence, we write

$`e^{m_1^{[0]}a/2}Z_{J_\mathrm{\Gamma }}=1+{\displaystyle \underset{l=1}{\overset{\mathrm{}}{}}}g_0^{2l}Z_{J_\mathrm{\Gamma }}^{[l]},`$ (2)
so that the $`Z_{J_\mathrm{\Gamma }}^{[l]}`$ are only mildly mass dependent. Fig. 1 plots the full mass dependence of the renormalization factors for the axial vector current $`Z_{A_4}^{[1]}`$. These numerical results are for the SW action with and without rotation, and also for Wilson action without rotation. Our results agree with those previously obtained, for $`c_{SW}=0`$  and for $`c_{SW}=1,d_1=0`$  111The coefficient $`d_1`$ is field rotation parameter. See Ref. .. We find that in both cases our results smoothly interpolate from the static to massless limit. The resulting analytical expressions are in Ref. . Fig. 2 plots the tadpole-improved renormalization factor for $`Z_{A_4}^{[1]}`$. From this figure, we can see that tadpole-improvement significantly reduces the one-loop coefficients of renormalization factors. Results for $`Z_{A_i}`$ and $`Z_{V_{4,i}}`$ are given in Ref. .

The slope of our mass-dependent renormalization factors in the massless limit is related to the improvement coefficients $`b_J`$ and $`c_J`$ . We find

$`b_V^{[1]}`$ $`=`$ $`0.153239(14),`$ (3)
$`b_A^{[1]}`$ $`=`$ $`0.152189(14),`$ (4)
$`c_V^{[1]}`$ $`=`$ $`0.016332(7),`$ (5)
$`c_A^{[1]}`$ $`=`$ $`0.0075741(15).`$ (6)
These results agree perfectly with Ref. . We also obtain by subtracting the integrands first,

$`b_V^{[1]}b_A^{[1]}=0.0010444(16)`$ (7)
which is more accurate than the difference of the two numbers quoted above. We find our one-loop result of $`b_Vb_A`$ are far from nonperturbative calculations .

## 3 SETTING THE SCALE

The typical gluon momentum $`q^{}`$ in the $`V`$-scheme, as suggested by Brodsky, Lepage and Mackenzie (BLM), is defined by

$$\mathrm{ln}(q_{}^{}{}_{}{}^{2})\frac{d^4qf(q)\mathrm{ln}(q^2)}{d^4qf(q)},$$
(8)
where $`q`$ is the momentum of gluon, and the form $`d^4qf(q)`$ is the one-loop integral for a particular renormalization constant, for example, $`d^4qf(q)=Z_{J_\mathrm{\Gamma }}^{[1]}`$. Previously $`q^{}`$ has been calculated for the light-light current   and the static-light current . Here we calculate the mass dependent $`q^{}`$ for the heavy-light current. Results are plotted in Fig. 3. For Wilson action case, our results agree with Ref.  in the massless limit. From Fig. 3, we can see that the mass dependence of $`q^{}`$ is weak from massless limit to $`m_0a1`$, especially for clover with rotation case. The original BLM prescription of $`q^{}`$ breaks down at larger masses, because $`Z_{J_\mathrm{\Gamma }}^{[1]}`$(denominator in Eq. (8)) goes through zero at there. A prescription for $`q^{}`$ in this case is given in Ref. . We also calculate tadpole-improved $`q^{}`$ and results are plotted in Fig. 4. We can see that plaquette tadpole-improvement significantly reduces $`q^{}`$, on the other hand, the reduction is rather small for $`\kappa _c`$ tadpole-improvement. We summarize the results in the massless limit in Table. 1.

We can also obtain the BLM scale for improvement coefficients $`b_J`$ and $`c_J`$ . Then it is interesting to compare BLM perturbation theory with non-perturbative calculations of these coefficients . We will present these results for $`q^{}`$ and the mentioned comparison in another publication .

## 4 CONCLUSIONS

We have obtained one-loop results of $`Z_A`$ and $`Z_V`$ with tree-level rotation, which should be useful for lattice calculations of $`f_B`$ and of form factors for $`B\pi l\nu `$. We have also obtained the BLM scale $`q^{}`$ for arbitrary masses, which should reduce the uncertainty of one-loop calculations.

## Acknowledgments

S.H. and T.O. are supported by the Grant-in-Aid of the Japanese Ministry of Education, (Nos.11740162, No.12640279). K.-I.I. and N.Y. are supported by the JSPS Research Fellowships. Fermilab is operated by Universities Research Association Inc., under contract with the U.S. Department of Energy.