text:
# Proper motions and velocity asymmetries in the RW Aur jetBased on observations collected at the Canada-France-Hawaii Telescope (CFHT), operated by the National Research Council of Canada, the Centre National de la Recherche Scientifique, and University of Hawaii.

## 1 Introduction

Jets from young stars are believed to be powered by accretion and accelerated and collimated by magnetic forces (see e.g. Königl & Pudritz ppiv (2000) for a review). However, a clear identification of the magnetic launch configuration is still needed to fully understand the role of jets in the physics of accreting protoplanetary disks. The optically revealed innermost jet regions of T Tauri stars (TTS) hold essential clues to this question, and are the subject of growing observational efforts at sub-arcsecond resolution (e.g. Lavalley-Fouquet et al. lavalley (2000); Dougados et al. paper1 (2000) (hereafter Paper I); Bacciotti, Mundt, Ray et al. bacc (2000)). HH 229, the bipolar jet from the actively accreting star RW Aur A, is a particularly interesting target. It was first identified in long-slit spectroscopy by Hirth et al. (hirth (1994)), who noted an asymmetry in radial velocity between the blueshifted and redshifted lobes of a factor $``$1.8. Further study with HST/STIS (Woitas et al. woitas (2002), hereafter Paper II) traces this asymmetry down to 0.1′′ $``$ 14 AU from the source, suggesting that it might be an intrinsic property of the jet driving engine. In addition, the RW Aur jet contains many emission knots, traced as far as 50′′-100′′ (Mundt & Eislöffel mundt (1998)) down to only 0.1′′-10′′ from the star (Paper I, Paper II), which may be used to test models of knot formation in stellar jets.

In this Letter, we present a ground-based spectro-imaging study of the RW Aur jet that brings several new results, including improved knot proper motions, jet inclination, and knot dynamical times, as well as indirect evidence for internal velocity and excitation gradients. Our observations are summarized in Sect. 2 and our results are discussed in Sect. 3.

## 2 Observations

Observations were carried out at the CFHT on December 18-20, 1998 with the integral field spectrograph OASIS. We used a spatial sampling of 0.16′′/lens and achieved an angular resolution of 0.4′′ FWHM after adaptative optics correction. Two spectral ranges were observed (6213-6534 Å and 6496-6824 Å) with a resolving power $``$ 90 km s-1. A second dataset with 0.11′′/lens was obtained over a smaller field of view. Following standard OASIS reduction (cf. Lavalley-Fouquet et al. lavalley (2000)), a specific procedure was applied to subtract at each position the strong continuum+emission spectrum of RW Aur A, scaled by the PSF, and to construct datacubes of residual forbidden line profiles. More details are given in Dougados et al. (2003, in preparation). Comparison between the two lens samplings shows that continuum-subtracted emission maps at 0.16′′/lens are reliable beyond 0.35′′ from the star, while the 0.11′′/lens data are more reliable at $`d<0.35^{\prime \prime }`$. Radial velocities are computed in the stellar reference frame (heliocentric velocity +16 km s-1; Petrov, Gahm, Gameiro et al. petrov (2001)). We adopt hereafter the notation \[OI\]$``$\[OI\]$`\lambda `$6300, \[SII\]$``$\[SII\]$`\lambda `$6731, and \[NII\]$``$\[NII\]$`\lambda `$6583.

## 3 Results

### 3.1 Proper motions of redshifted knots

Figure 1 compares our 1998 \[S II\] image of the redshifted jet with images obtained in Dec. 2000 (Paper II) and in Dec. 1997 (Paper I)111In order to replicate the correct separation between RW Aur A and B, the Dec. 1997 image published in Paper I had to be corrected for a small 5% underestimate in the FOCAM pixel size. A smooth systematic shift is seen between the 3 epochs, indicative of proper motions. Table 1 lists the knot positions and the resulting proper motions in arcsec per year. As the OASIS spatial resolution is not sufficient to detect possible offsets between different lines, we adopt for the nominal 1998 position of each knot in Table 1 an average of independent estimates using the strongest two lines (\[OI\] and \[SII\]) and the two lens samplings (0.11′′ and 0.16′′). Knots are numbered from the outside in, starting at the most distant one detected by Dougados et al. (Paper I). For the outermost knot seen with OASIS (knot R4), the 3 epochs are consistent with a constant proper motion of 0.16′′ yr-1. The proper motion for the following knot (R5) is more uncertain, as it appears double in the 2000 STIS image. The range of possible values is 0.17′′-0.27′′ yr-1. For knot R6, we find again 0.16′′ yr-1, based on the 1998 and 2000 epochs only. We note that an innermost emission peak is detected in the 0.11′′/lens dataset at $``$ 0.28′′. It may correspond to the peak at 0.25′′ in the STIS image, and would then be stationary.

### 3.2 Proper motions of blueshifted knots

Figure 2 compares \[S II\] images of the inner blueshifted jet in 1998 (this paper) and 2000 (Paper II). We lack a 1997 epoch, as deconvolution residuals in the PUEO image dominate over the faint blue jet emission out to $`2^{\prime \prime }`$ (Paper I). We find a proper motion of 0.26′′ yr-1 for the main knot (B3). A new knot (B4) has appeared at 0.6′′ in the 2000 STIS image. With a similar proper motion to B3, it would have been within 0.05′′ of the star at the date of the OASIS observations and thus unresolvable. Finally, our 0.11′′/lens dataset shows an innermost peak at $`0.19^{\prime \prime }`$, which may correspond to the 0.2′′ peak seen in the STIS image and would then be stationary.

### 3.3 Revised jet inclination

Assuming that knot proper motions trace actual fluid motions and not wave pattern speeds, the angle of inclination $`i`$ of the jet with respect to the line of sight is given by $`\mathrm{tan}(i)=(\mu /^{\prime \prime }\mathrm{yr}^1)\times (666.7/V_{\mathrm{rad}})`$, where $`V_{\mathrm{rad}}`$ is the velocity projected on the line of sight, and a distance of 140 pc to Taurus-Auriga is assumed. We adopt as $`V_{\mathrm{rad}}`$ for each knot the centroid of a gaussian fit to the \[S II\] line profile. Values from this paper and from Paper II (in italics) are given in Col. 5 of Table 1. Discrepancies in $`V_{\mathrm{rad}}`$ between the two datasets are discussed in more detail in Section 3.5. The resulting values of $`i`$ are given in Col. 7, with an uncertainty interval that includes both the uncertainty in $`V_{\mathrm{rad}}`$ and that in $`\mu `$ (see notes to Table 1). Inclination values are remarkably consistent between knots, within the estimated errors. In particular, we derive consistent inclinations for the blue and red lobes, i.e. the proper motion in the blue jet is higher than in the red jet by the same factor ($`1.6`$) as their average centroid radial velocities. This fact argues strongly for a real velocity difference between the blueshifted and redshifted jets, rather than an apparent difference due solely to slow bowshock wings dominating the emission around the redshifted jet (proper motions should then be similar in the two lobes despite differing radial velocities). A weighted average over both lobes yields222The larger inclination value $`i=53`$° found by Woitas et al. (Paper II) for knot R4 stems from a proper motion estimate 20% too high, inferred from the published PUEO image of Paper I before correction for pixel scale (see footnote 1) $`i=46`$°$`\pm `$3°.

### 3.4 Knot spacing, dynamical times, and origin

The RW Aur observations allow to measure knot spacings and dynamical timescales very close to the source and thus bring interesting constraints on the origin of knots in stellar jets. To increase the distance and time span in our analysis, we include in Table 1 more distant knots imaged with PUEO in 1997, which lie beyond the OASIS field of view and were not considered by Woitas et al. (Paper II).

Table 1 shows that the knot spacing first increases with distance from the star, from $`0.7^{\prime \prime }`$-$`0.8^{\prime \prime }`$ at $`d1.5^{\prime \prime }`$ (knots R5-R6, B3-B4) to $`1.8^{\prime \prime }`$-$`2.4^{\prime \prime }`$ at $`d3^{\prime \prime }`$ (knots R4-R5, B2-B3), and then levels off at $`\mathrm{\Delta }x3^{\prime \prime }4^{\prime \prime }`$ for $`d4^{\prime \prime }`$. Adopting as the jet characteristic radius $`r_\mathrm{j}`$ half the FWHM measured at the same distances in Paper II, the inner knot spacing is then $`\mathrm{\Delta }x810r_\mathrm{j}`$. For comparison, the axisymmetric pinch instabilities in MHD jets simulated by Lery et al. (2000) produce jet knot spacings of $`\mathrm{\Delta }x/r_\mathrm{j}3`$, while Kelvin-Helmholtz instabilities predict $`\mathrm{\Delta }x/r_\mathrm{j}V_\mathrm{j}/c_\mathrm{s}>10`$. An even more important observational constraint is that the first moving knot appears at a distance $`d_{\mathrm{min}}<\mathrm{\Delta }x`$, while the growth scale of an instability is expected to be several times its wavelength. Hence, neither kind of instability seems adequate to explain the innermost knot spacings, although a range of MHD jet configurations should be explored before they are fully ruled out.

Another possible origin for the knots, which can easily yield $`d_{\mathrm{min}}<\mathrm{\Delta }x`$, is time variability in the jet ejection speed (cf. Equ. (3) in Raga et al. raga1 (1990)). The last column of Table 1 gives the dynamical ages of each knot assuming constant proper motion (with $`t=0`$ arbitrarily set at our OASIS observation date of Dec.19, 1998). For knots lacking proper motion measurements, we adopt the typical $`\mu `$ in the same lobe (in brackets in Col. 6 of Table 1). We find that the timescale between knots initially increases with distance, from 3 to 10 yr, and then appears to level off at $``$ 20 yr beyond 4′′. Furthermore, the inner knots do not appear to be “synchronized” between the two jet sides, while there is tentative evidence of synchronization when the 20 yr period appears (knots R3 and B1). We conclude that there would be two distinct processes of time variability in RW Aur: one that is irregular and asymmetric — possibly random (random time variability of the jet speed produces asymptotically an increase in knot spacing with distance; Raga raga2 (1992)) and another more regular and symmetrical with a 20 yr period — possibly related to a stellar magnetic cycle. We suggest that asynchronous variability might occur, e.g., through local magnetic reconnexion events that affect only one jet side.

The apparently stationary knots at $`0.2^{\prime \prime }`$-$`0.3^{\prime \prime }`$ = 30-45 AU from the star need to have a different origin, possibly recollimation shocks (e.g. Gomez de Castro & Pudritz gomez (1993)) or temperature gradients in a progressively heated wind (e.g. Garcia et al. garcia (2001)), but predicted distances are presently too model-dependent to yield useful constraints.

### 3.5 Velocity gradients and velocity asymmetry

Figure 3 plots the \[S II\] centroid radial velocity as a function of distance from the star for the OASIS and STIS datasets (the latter is shifted to correct for proper motions). In the blue jet, the various datasets agree to within 10 km $`\mathrm{s}^1`$, comparable to expected instrumental and calibration errors. The situation is different in the redshifted jet: the 0.11′′/lens and 0.16′′/lens datasets still agree with each other within $``$ 10 km $`\mathrm{s}^1`$, but they lie typically 20 km $`\mathrm{s}^1`$ below the STIS values out to the position of knot R5. Further out, however, at the position of knot R4, all data converge to a very similar $`V_{\mathrm{rad}}`$.

These discrepancies can only be understood if the faster gas seen on the redshifted jet axis by STIS (0.1′′ wide slit) suffers from beam-smearing in the OASIS data. In other words, sharp transverse velocity gradients must occur on scales $`0.4^{\prime \prime }`$ within the redshifted jet, out to distances $`<3^{\prime \prime }`$ from the star333such gradients seem confirmed by preliminary analysis of transverse STIS slits (Woitas et al., work in progress). We note that \[N II\] centroids from OASIS (shown as crosses in Fig. 3) agree better with the \[S II\] STIS centroids, despite an angular resolution identical to the \[S II\] OASIS data. This fact further indicates that the faster axial gas contributes a greater fraction of the total flux in \[N II\] than in \[S II\], i.e. that it has higher excitation. Interestingly, we find no evidence for such velocity and excitation gradients in the blueshifted jet. If confirmed, this difference might be related to the blue/red velocity asymmetry between the two jet lobes; e.g., bowshock wings or turbulent boundary-layers could be somehow enhanced in the slower lobe (thus also explaining its higher brightness). Further theoretical and observational work is needed to explore this possibility.

###### Acknowledgements.

LLM is grateful to Observatoire de Paris and to Ministère de l’Education Nationale for a postdoctoral fellowship.


format: text

text:
# Scaling of impact fragmentation near the critical point

## I Introduction

From asteroids to nuclei, fragmentation phenomena can be seen everywhere. Many scientists and engineers have been fascinated by the fragmentation process, and much effort has been made to understand fragmentation Beysens1 . In particular, impact fragmentation of brittle solids has been investigated by simulation and experimentally Hayakawa1 ; Ishii1 ; Oddershede1 ; Meibom1 ; Kadono1 . The results have shown that cumulative fragment mass (or size) distributions exhibit a power-law dependence, in which the exponent depends on the dimensionality of the fractured objects. These results do not depend on the details of fragmentation, such as the material the fractured object is made from or the fragmentation method. Consequently, Oddershede et al. concluded that brittle fragmentation is a self-organized critical phenomenon Oddershede1 ; Bak1 .

Cumulative distributions obey power-law dependence when the imparted energy is sufficiently large, but brittle solids do not break when the imparted energy is very small. Therefore, one might assume that there is a transition in the fragmentation phenomenon. Recently, Kun and Herrmann examined critical behavior in fragmentation transition and concluded that fragmentation transition belongs to the percolation transition universality class Kun1 . In addition, Campi has shown that nuclei fragmentation has the same statistical properties as percolation Campi1 . By contrast, Åström et al. performed other fragmentation model simulations and obtained a scaling law whose critical exponents differ from those of percolation Astrom1 .

While there are some simulation results, there have been no experiments on critical fragmentation. For a proper analysis, the existing simulation results must be compared to experimental data. Here, we experimentally investigate brittle fragmentation as a kind of critical phenomenon, and we report experimental results for glass tube fragmentation near the critical point.

## II Experimental procedure

Figure 1 is a schematic illustration of the experimental apparatus. We used glass tubes ($`50`$ mm outside diameter, $`2`$ mm thick, and $`150`$ mm, $`100`$ mm, or $`50`$ mm long) as the fractured objects. They have 2-D geometry, and it is easy to impart the impact vertically because they can stand by themselves. Kun and Herrmann Kun1  and Åström et al. Astrom1  performed numerical simulation of 2-D samples.

A glass tube was put in a plastic bag, which was in turn placed between a hard stainless steel stage (dimensions $`400`$ mm $`\times `$ $`500`$ mm $`\times `$ $`50`$ mm) and a stainless steel plate ($`10`$ mm thick). The stage was fixed by placing heavy weights on it. A cylindrical brass weight with a flat bottom ($`3.77`$ kg) was dropped on the stainless steel plate vertically, and a planar failure wave propagated to the glass tube from the impact circle (cross section of the glass tube). Consequently, the glass tube was cleaved by the impact. This is a type of 2-D fragmentation. Four poles ($`10`$ mm diameter) were used as guide poles for the falling brass weight. Stoppers for the stainless steel plate were placed at each guide pole to prevent a secondary impact. While the samples were not annealed after the cutting processing, they had almost the flat end. We used a level in order to check the verticality between the weight and the stainless steel plate at each fragmentation. The variable in this experiment was the height that the weight was dropped from.

After fragmentation, we collected the fragments in the plastic bag and measured the mass of each fragment with an electronic balance. The data were analyzed, and cumulative fragment mass distributions were obtained. We fractured a total of $`60`$ samples ($`34:150`$ mm, $`15:100`$ mm, and $`11:50`$ mm length). We defined the falling height $`h`$ as a length from the top of stainless steel plate to the bottom of the falling weight. In our experiments, it was set in the range $`61`$ mm $`h`$ $`407`$ mm. This range corresponds to $`0.050ϵ0.251`$ (Nm/g) in terms of the imparted energy (released potential energy of the dropped weight) per unit sample mass, $`ϵ`$.

## III Results

When the height was too low to break the glass tube, no visible cracks were seen. Once a visible crack was produced, fragmentation occurred suddenly. We defined the point at which a sample began to cleave as the fragmentation critical point. We could not obtain samples that had only visible macro-cracks, but did not fragment. We show some photos of fragments in Fig. 2. We observed that with a flat impact, low impact energy fragmentation events produced vertical main cracks. Therefore, vertical cleaving produces a few, large fragments (Fig. 2(a)). In the mediate range, some large fragments remained (Fig. 2(b)). For large impact energy events, the largest fragment was smaller, and it was impossible to distinguish the main cracks from the many small fragments by observation (Fig. 2(c)).

### III.1 Cumulative distribution

First, we plotted the cumulative mass distribution of each fragmentation. Figure 3 shows some typical cumulative fragment mass distributions, $`N(m)`$, for $`150`$-mm-long samples. $`N(m)`$ is defined as follows:

$$N(m)=_m^{\mathrm{}}n(m^{})𝑑m^{},$$
(1)
where $`m`$ and $`n(m)`$ are the fragment mass and the number of fragments of mass $`m`$, respectively. The distribution curves for samples of different sizes were similar to the example in Fig. 3. One can see the clear power-law dependence for fully fragmented samples (far from the fragmentation transition point). The distribution curves near the critical point do not show clear power-law dependence, and they have a rather flat part. This tendency is quite different from that of percolation universality. In percolation universality, the clearest power-law dependence should appear beneath the critical point. In the fully fragmented data curves, we can observe a region that satisfies $`N(m)m^{(\tau 1)}`$, i.e., $`n(m)m^\tau `$. The obtained value of the characteristic exponent $`\tau `$ is nearly $`1.5`$ ($`<2`$). This result conflicts with the scaling ansatz of percolation Stauffer1 . Therefore, we believe fragmentation criticality to differ from that of percolation. Conversely, our result $`\tau 10.5`$ concurs with Åström’s result Astrom1 . In addition, the value of $`\tau `$ obtained in this experiment is consistent with Hayakawa’s scaling for 2-D fragmentation with a planar failure wave Hayakawa1 .

### III.2 The weighted mean mass scaling

Since fragmentation events near the fragmentation critical point have less fragments, they cannot be analyzed statistically using only $`N(m)`$. Therefore, we used the $`k`$-th order moment of the fragment mass distribution, $`M_k`$, defined as

$$M_k=\underset{m}{}m^kn(m).$$
(2)
Since we are interested in fragmentation criticality, the height of the falling weight was controlled near the critical point. Despite our control of the height (imparted energy), the $`M_k`$ data fluctuate too much due to individual differences in the glass tube samples, e.g., the initial density of micro-cracks and residual stress distributions. When we considered the imparted energy as a control parameter, we could not find clear critical behavior quantitatively. For example, we show a log-log plot of $`\frac{M_2}{M_1}`$ vs. $`ϵ`$ in Fig. 4(a). Where, $`\frac{M_2}{M_1}`$ is the weighted mean fragment mass. It is hard to say the data in Fig. 4(a) are particularly convincing about showing linear behavior. We could not consider that power-law fitting was appropriate with confidence. Alternatively, we used multiplicity $`\mu `$ analysis, which Campi introduced to the analysis of nuclei fragmentation Campi1 . Multiplicity is defined as

$$\mu =m_{\mathrm{min}}\frac{M_0}{M_1},$$
(3)
where $`M_0`$, $`M_1`$, and $`m_{\mathrm{min}}`$ correspond to the number of fragments, the total mass of the fragments, and the smallest cut-off fragment mass, respectively. We fix $`m_{\mathrm{min}}=0.01\text{g}`$ for all the analyses. This value usually corresponds to the smallest 2-dimensionality for the samples. According to this definition, the fragmentation critical point corresponds to $`\mu =0`$. While the multiplicity is a resultant parameter, rather than a control parameter, it can be considered the control parameter that incorporates individual differences in the samples. Furthermore, because the multiplicity is easy to measure and calculate, it has been applied to other fragmentation systems Campi1 . We use $`\mu `$ as a pseudo-control-parameter in this paper. Certainly, $`\mu `$ and the imparted energy are positively correlated, roughly. The larger the impact energy imparted, the larger $`\mu `$ becomes.

In Fig. 4(b), we plot $`\frac{M_2}{M_1}`$ vs. $`\mu `$ of all the fragmentation event data. It indicates the tendency toward divergence near the critical point. The criticality appears independent of sample size. All the data match the same power-law line. Hence, we assume that the relationship between the weighted mean fragment mass $`\frac{M_2}{M_1}`$ and the multiplicity, $`\mu `$, holds as,

$$\frac{M_2}{M_1}\mu ^\sigma .$$
(4)
The obtained value of $`\sigma `$ is $`0.84\pm 0.05`$. This value is nontrivial because $`\sigma 1`$.

### III.3 Scaling function

We discuss the possibility of scaling $`N(m)`$ curves using $`\sigma `$. For large $`\mu `$ data, $`N(m)`$ curves have the same power-law, in which the power is $`\tau 10.5`$. This value is independent of the size of the samples. This similarity allows us to collapse the $`N(m)`$ curves of different-sized samples into a unified scaling function. We expect the following scaling for distribution curves:

$$P(m)f(m/\mu ^\sigma )m^{(\tau 1)}g(m/\mu ^\sigma ).$$
(5)
Where $`P(m)`$ is the probability distribution of fragments whose masses are larger than $`m`$. We show the scaling result in Fig. 5 with the approximation $`P(m)=\frac{N(m)}{M_0}`$. The curves in Fig. 5 include the results for samples of three different sizes (11 curves of large $`\mu `$ data). Figure 5 shows good collapse of $`N(m)`$ curves into the scaling function $`f`$ and $`g`$. The function $`f`$ has a clear power-law regime, in which the power is $`0.5`$, and a unified cut-off scale. The former corresponds to the scaling region and the latter corresponds to exponential decay due to the finite size effect. Accordingly, the scaling function $`g`$ has a constant part and a rapidly decaying part. The flat part in $`g`$ exists only for the fully fragmented data, because clear power-law form $`N(m)m^{(\tau 1)}`$ can be satisfied with well fractured samples. The form of scaling function might change with $`\mu `$, however, we consider they will not depend on the sample size. We could not confirm the clear collapse for small $`\mu `$. More detailed experiments are necessary to verify the scaling functions in wide $`\mu `$ range.

Bak et al. reported a similar scaling function for earthquake statistics Bak2 . While they discussed the time correlation of earthquakes, their scaling function is analogous to ours. They concluded that the unified scaling function for earthquakes depends only on the number of earthquakes occurring within the area and period considered. In Eq. 5, our scaling function depends on the variable $`x=\frac{m}{\mu ^\sigma }`$, which is the normalized mass of the fragments. Moreover, the scaling Eq. 5 resembles the scaling function Åström et al. obtained Astrom1 .

### III.4 Multi-scaling and simple biased cascade model

Next, we discuss the behavior of the higher-order moments of the fragment mass distributions. We can easily expand Eq. 4 and introduce the higher-order weighted mean fragment mass scaling as,

$$\frac{M_{k+1}}{M_k}\mu ^{\sigma _k}.$$
(6)
The open circles in Fig. 6 show the $`\sigma _k`$ values obtained from fitting the all experimental data like Fig. 4(b). While the value $`\sigma _{k=0}`$ is the trivial value $`1`$, $`\sigma _k\mathrm{}`$ approaches a nontrivial value of approximately $`0.6`$.

In order to understand this behavior, we consider a simple biased cascade model. Let us set a unit mass fragment initially. First, this unit mass is divided into two pieces of mass $`a`$ and $`1a`$ $`(0a1)`$. The same $`a`$-biased partition occurs for each fragment in the next step. This cascade continues until the imparted energy dissipates. From the above definition and Eq. 6, the following equation holds at each step:

$$\frac{a^{k+1}+(1a)^{k+1}}{a^k+(1a)^k}=2^{\sigma _k}.$$
(7)
Substituting $`a=\frac{1}{2}`$ into Eq. 7, the trivial value $`\sigma _k=1`$ is obtained for all $`k`$. This case corresponds to mono-scaling critical fragmentation. In Fig. 6, however, $`\sigma _k`$ varies with $`k`$. This indicates a multi-scaling property of critical impact fragmentation. The solid curve in Fig. 6 is the value calculated using Eq. 7 with $`a=\frac{1}{3}`$. The curve and experimental data roughly agree within the error bars. Although this model is very simple, it explains the trend of the higher-order behavior. We conclude that fragmentation near the critical point has a multi-scaling nature. However, this only alters the problem. The origin of the symmetry breaking given with $`a`$ and the relation between $`a`$ (or $`\sigma _k`$) and $`\tau `$ remain unsolved. The cascade model does not seem to yield the appropriate value of $`\tau `$ ($`1.5`$) directly. Since this model is too simple, it needs to be improved using more detailed analyses and experiments.

Similar multiplicative model was proposed for the energy cascade of eddy by Meneveau and Sreenivasan Meneveau1 . According to their result, the cascade of division into $`0.3`$ and $`0.7`$ can explain the energy cascade of turbulent flow. However, the origin of the value $`0.3`$ has remained unsolved. The value $`0.3`$ slightly agrees the value $`\frac{1}{3}`$ we obtained.

## IV Discussion

Kun and Herrmann have analyzed the largest fragment mass as an order parameter Kun1 . In our experiments, the largest fragment mass had a very wide distribution, as seen in Fig. 5, while the weighted mean fragment mass behaved more calmly. Figure 3 suggests that the power-law exponent might change with $`\mu `$. Ching et al. reported the impact energy dependence of the power-law exponent Ching1 . Their interest focused on the very large imparted energy region and not on the neighborhood of the fragmentation critical point. In addition, if we fit the power-law form to all $`N(m)`$ curves, the results involve large uncertainty, and do not show a clear trend. Ishii and Matsuhita reported that the distribution function for the low impact energy region has a log-normal form Ishii1 . We have not investigated the distribution of the concrete function form. The relationship between the two energy regions and the study of the concrete function form remain open to study.

Our experiments were focused on the 2-D fragmentation with the flat impact. Kun and Herrmann investigated the point impact 2-D fragmentation Kun1 . This might be a reason of the discrepancy between our results and their ones. Hayakawa suggested that the universality of the fragmentation depends not only on the dimensionality of fractured object, but also on the propagation dimensionality of the failure wave Hayakawa1 . Åström et al. performed numerical simulations of 2-D fragmentation with the flat expansion Astrom1 . Their results resemble ours in many points as described above. We believe that our results are valid for 2-D fragmentation by the flat impact to one side. Other dimensional objects or other failure wave propagation manner may produce other values of scaling exponents such as $`\tau `$ and $`\sigma _k`$. Moreover, scaling-law among these exponents may exist. As written in Sec. I, the exponent $`\tau `$ is independent of the material of fractured objects in fully fragmented state. We guess the other critical exponents are also independent of that, however, it should be examined. These problems are open questions.

In general, the low impact energy region is difficult to study. Small differences of initial condition are enlarged by nonlinearity of fragmentation dynamics. Figure 4(a) presents this difficulty well. Even in Fig. 4(b), the fluctuating data can be seen. Some unrealized experimental parameters might not be controled. We think that the coincidence between the experimental result and the model shown in Fig. 6 is rough due to this reason.

In conclusion, we investigated the 2-D impact fragmentation of brittle solids experimentally. The measured characteristic exponents imply that the impact fragmentation of brittle solids does not belong to the percolation transition universality class. Instead, we found that the weighted mean fragment mass $`\frac{M_2}{M_1}`$ and the multiplicity parameter $`\mu `$ are related to the nontrivial scaling exponent $`\sigma `$. The cumulative fragment mass distributions of different-sized samples could be collapsed into the unified scaling function in the large $`\mu `$ regime. We calculated the generalized weighted mean fragment mass scaling $`\frac{M_{k+1}}{M_k}\mu ^{\sigma _k}`$. The behavior of the multi-scaling exponent $`\sigma _k`$ was modeled using a simple biased cascade partition model. The critical behavior of other fragmentation systems should be studied in order to understand the universality of fragmentation phenomena in detail. Experiments with other materials, such as ceramics, or with samples of different dimensions may prove interesting.

We thank Professor H. Sakaguchi and Professor K. Nomura for their useful comments.


text:
# Theory of Analogous Force on Number Sets

## I Introduction

There are suggestive parallels between various aspects of number theory and physics phenomena. The connection often involves affinity with functions arising in statistical mechanics Ste92 . For example, a connection between given number sets and the thermodynamic formalism has been established for the multifractal spectrum via a Legendre transform Lee88  and the statistical mechanics applications for prime numbers via the Riemann partition function Wol99 . Bose condensation has been indirectly linked to the zeta function in Jul90 . Primes as such have also attracted recent applications in acoustics and dynamical systems Dah01 . Furthermore, the existence of an imaginary quantum mechanical potential for the prime numbers has been idealized in Mus97  after some energy considerations. In a similar context, additive number theory in terms of scale invariance has been used to study the origin and form of the uneven distribution of integers of many random data sets Hil98 ; Tol00 ; Adh68 ; Math ; Got02 ; Sny02 .

In this paper we introduce an alternative approach which considers sets of $`x`$-natural numbers to play the role of particles of known type in an isolated elastic system and assume that their discrete probability distribution can be split as the product of some functions. By applying thermodynamics techniques, we develop an analogous statistical mechanics version of the average energy and energy states of arbitrary number sets.

From the Helmoltz free energy in a reversible isothermal process, a relation for an analogous tensile force $`F>0`$ is derived by treating $`x`$ as variations in length in the direction of $`F`$ and considering the number system as a elastic specimen (body). In particular, our theory allows to give an alternative physics answer to the uneven distribution of integers in terms of a screened inverse square force which ideally may act on these numbers. It also allows to establish second order polynomial recursions to predict sequences of integers following generalized Fibonacci numbers. We also analyze a possible application of the distinct forces to primes.

Hereafter, the word x-number is used to refer to quantities which are integers and positive.

## II Underlying Theory

### II.1 A priori probability postulate

Let us define the discrete probability distribution of finding a certain positive natural number $`x`$ as

$$p_x=g_1(x)g_2(x),$$
(1)
with $`g_1(x)>0`$ and $`0<g_2(x)1`$ and satisfying the normalization $`_xp_x=1`$. Such discrete probability distribution functions (and also discrete mixtures of them), used to describe complex phenomena, are given by the product of two different functions. The real function $`g_2(x)`$ here represents the set of digital frequencies which remain fixed under base changes given by the real function $`g_1(x)`$.

This factorization may be seen somewhat arbitrary, however it is based on composed functions which include the Binomial(A,B), Furry or Geometrical, inversed Logarithmic, Pascal or Polya and Poisson distributions. Furthermore, a probability decomposition into (two) functions is no new. This type of factorization appears formally, e.g., in the analysis of stochastic processes on graphs. According to the well-known Hammersley-Clifford Theorem, a probability function $`p(x)`$ is formed by a normalized product of positive functions on cliques (i.e., fully connected subgraph) of the graph $`𝒢`$ when the stochastic process $`x`$ is Markov relative to $`𝒢`$.

To simplify notation in all later sums and products the variable $`x`$ represents an upper bound (or cutoff).

### II.2 Thermodynamic Analogy

Using the above probability postulate the information theoretical entropy $`S_I=_xp_xlnp_x`$ then becomes

$$S_I=\underset{x}{}g_1(x)g_2(x)\{lng_1(x)+lng_2(x)\}.$$
(2)
In order to make a thermodynamic analogy, let us consider next the Helmholtz free energy $`A=ETS`$ assuming an analogous temperature $`T`$ and volume $`V`$ to be independent variables, with $`E`$ the analogous internal energy of number sets. Because of the close relationship of Eq.(2) to the thermodynamic entropy $`S`$, we can readily identify

$`A/kT`$ $``$ $`{\displaystyle \underset{x}{}}p_xlng_1(x),`$ (3)
$`E/kT`$ $``$ $`{\displaystyle \underset{x}{}}p_xlng_2(x).`$ (4)
with $`k`$ the Boltzmann constant. Reduced units are used throughout. The choice of the latter quantities may be seen quite arbitrary, but as mentioned in the conclusion these can be justified in analogy with standard thermodynamics results.

In this case, for the correct application of conventional thermodynamics the normalization of $`p_x`$ is used to define both the Boltzmann and the Shannon entropies.

For ’non-interacting’ number sets $`x>1`$, we can get an approximated statistical mechanics version of the average energy as the sum of individual contributions $`E\overline{E}=_x\overline{n}_xϵ_x`$. The total number of positive natural number $`x`$’s is given by $`N=_x\overline{n}_x`$ and $`ϵ_x`$ represents their analogous energy states. Thus, from Eqs.(1) and (4) we further identify

$`\overline{n}_x`$ $``$ $`lng_2(x)0,`$ (5)
$`ϵ_x/kT`$ $``$ $`g_1(x)g_2(x)=p_x0.`$ (6)
The number of states available to a molecular system is given by the partition function $`\mathrm{\Omega }`$ which satisfies the relation $`S=kln\mathrm{\Omega }(N,V,E)`$ in the microcanonical ensemble and relates the canonical partition function by $`Q(N,V,T)=\mathrm{\Omega }(N,V,E)e^{E/kT}`$. Therefore, from Eqs.(2), (4) and (6) we also have in analogy

$$Q_x=\underset{x}{}[g_1(x)]^{p_x}=\underset{x}{}[g_1(x)]^{ϵ_x/kT}.$$
(7)
As expected, the form of our canonical partition function depends on the number probability distribution of Eq.(1).

### II.3 Analogous Force

To derive a relation for an analogous external tension which we may ideally apply on number sets in a gedanken experiment, let us treat $`x`$ as a variation in length in the direction of a tensile force $`f_x`$ and the number system as an elastic body. As a first approximation, the required relation is then obtainable directly from the first and second laws of thermodynamics in a reversible isothermal process Tre75 . For a change at constant $`T`$, it follows that $`dA=dW=fd\mathrm{}`$, where $`W`$ is the work done on the system by its surroundings in a small elongation $`d\mathrm{}`$. Hence, by partial differentiation of Eq.(3),

$$f_x=\left(\frac{A}{x}\right)_T=\frac{ϵ_x}{g_1(x)}\left(\frac{g_1(x)}{x}\right)_T+lng_1(x)\left(\frac{ϵ_x}{x}\right)_T.$$
(8)
Since the free energy depends on the upper cutoff $`x`$ then the derivative was taken in terms of that variable, ıi.e. the derivative of the last term in the sums. If $`f_x<0`$, this corresponds to an analogous compressive force on the number system. We shall show next that this new tensile force may indeed be useful to characterize number sets.

## III Probability Distributions and Related Functions

### III.1 First Digits Patterns: $`g_1(x)c_0`$ , $`g_2(x)f[lnx]`$

In terms of the first digits of unbiased random samples from random distributions, the so-called Benford’s law follows the logarithmic distribution

$$p_x=log_{10}(1+\frac{1}{x})c_0ln\gamma (x),$$
(9)
with $`1/c_0ln\mathrm{\hspace{0.33em}10}`$, $`\gamma (x)=1+\frac{1}{x}`$ and $`x=1,2,\mathrm{}`$. This empirical law deals with the distribution of first digits in numbers in which significant positive digits are not uniformly distributed, but smaller digits and smaller combinations of significant digits are favored.

The sum of probabilities in this case is

$$\underset{x}{}p_x\beta _x=c_0\underset{x}{}ln\gamma (x)=c_0ln\underset{x}{}(\frac{x+1}{x})=c_0ln(1+x).$$
(10)
This result implies that $`p_x`$ is not a proper probability since it imposes some maximum value Hil98 .

From Eq.(7) the evaluation of the analogous partition function is straight forward, namely

$$Q_x=\underset{x}{}\gamma (x)^{c_0lnc_0}=(1+x)^{c_0lnc_0}.$$
(11)
According to Eq.(8), our analytical expression for the analogous force acting on this number set becomes

$`{\displaystyle \frac{f_x}{c_0lnc_0}}`$ $`=`$ $`{\displaystyle \frac{1}{x+x^2}}.`$ (12)
It is from this base-invariant behavior of a distinct $`f_x`$ that we can make an attempt to understand the significant digit distribution of Eq.(9) in terms of an screened inverse square force.

### III.2 Relevance to Primes: $`g_1(x)e^1`$ , $`g_2(x)x^\sigma `$

If we generalize the above to a power law behavior for the sequence of primes only greater than one with an exponent $`\sigma >1`$ , i.e., $`p_{x_p}=1/(ex_p^\sigma )`$, then the class of analogous forces on prime numbers becomes

$$f_{x_p}=\frac{\sigma }{x_p^{\sigma +1}},$$
(13)
-a fingerprint for primes.

Using Eq.(7), the analogous partition function approaches

$$\frac{Q_{x_p}}{e^{1/e}}\underset{x_p>1}{}\frac{1}{1x_p^\sigma }=\underset{x}{}\frac{1}{x^\sigma }.$$
(14)
Another possible relation of our statistical thermodynamic approach to prime numbers follows by considering the average system energy. In the limit of Eq.(12), this takes the form

$$\frac{E}{c_0}\underset{x}{}\frac{lnx}{x}.$$
(15)
In principle this could be asymptotic to $`_x\frac{1}{\pi (x)}`$, where the counting function $`\pi (x)`$ denotes the number of prime numbers $`x_p`$ that do not exceed $`x`$ Ten00 .

## IV Discussion

### IV.1 Benford’s law, Analogous Force and Thermodynamics

As an example, in Fig.1 we plot the probability distribution and our distinct force on (up to 300) numbers from Fibonacci series and compare results with the present analogous statistical thermodynamics approach. The theoretical calculation of the analogous force is obtained using Eq.(12) (and its large-$`x`$ behavior). The numerical calculation of $`f_x`$ for Fibonacci numbers is also rather simple. This is derived using Eq.(8) and from the fact that $`c_x`$ is constant. By the two-point central formula we approximate

$$f_x/lnc_0=\left(\frac{p_x}{x}\right)_T\frac{1}{2}(p_{x+1}p_{x1}).$$
(16)
This means that the probability distribution $`p_x`$ is related to an analogous local potential field at each number (i.e., particle) as in the well known diffusion limited aggregation model Wit81 . Assuming a sort of superposition principle, the total analogous force $`F_xf_x`$ can also be readily evaluated. From Eq.(12) we obtain $`F/c_0=lnc_0/\gamma (x)`$.

Since $`c_00.4343`$ in logarithm base $`10`$, we then have an analogous system entropy $`S>0`$, the analogous Helmholtz free energy $`A<0`$, the analogous partition function $`Q_x>0`$, an analogous expansive force on the integers $`f_x>0`$ and the sum of probabilities $`\beta _x>0`$. Interestingly, from Eq.(10) we also have that this sum only becomes normalized at $`x=9`$ (i.e., a decadelong). On the other hand, the increase in the analogous entropy as $`x`$ increases follows the natural tendency of thermodynamic systems to move toward disordered states having more microscopic states characterized by our analogous macroscopic variables.

Besides, the particular $`\gamma (x)`$ under consideration allows to derive a measure for an analogous inverse temperature of the integers system. From Eqs.(2), (4) and (9), we thus approximate

$$\frac{\mathrm{\Delta }(S/k)}{\mathrm{\Delta }(E/kT)}=1+\frac{lnc_0}{ln(ln\gamma (x))},$$
(17)
where $`\mathrm{\Delta }`$ represents functions difference between consecutive number fluctuations.

For $`x>>1`$, Eq.(12) implies that $`p_xc_0x^1`$. Such a power law behavior with exponent one is the solution of the functional form $`p_x^{}=p_{ax}=A_ap_x`$ which may correspond to an analogous scale invariance Math .

The results shown in the figure allows to argue that our simple physics-based approach is not merely speculative. Indeed it yields to a different and more complete understanding of Benford’s law than that of previous attempts Adh68 ; Math ; Got02 ; Sny02 . It includes scale invariance, frequencies behavior under base changes and explanation for a logarithmic law distribution besides the analogous thermodynamic characterization of given number sets such as time series or primes.

### IV.2 Screening

For $`x>>1`$, the force $`\frac{f_x}{c_0lnc_0}`$ can be approximated as $`\frac{e^{\frac{1}{x}}}{x^2}`$, so we could consider the exponential term as being screening. The presence of $`exp(1/x)`$ screening in this analogous inverse-squared force implies that the second leading digits (i.e., $`10x99`$) are much more uniformly distributed than the first digits, and the third than the second, and so on according to observations Hil98 . This is reflected by the numerical differences $`f_1f_9=0.17708`$, $`f_{10}f_{99}=0.00326`$, $`f_{100}f_{999}=0.00004`$. The main message is that the asymmetry in favor of smaller significant digits or combinations of significant digits can be understood in terms of analogous screened inverse square force acting on the number sets.

Indeed screening could be a way to unveil the fact that significant digits are not independent, but that instead knowledge of one digit affects the likelihood of another. Furthermore, since screening is a consequence of having the base invariant probability distribution of Eq.(1), then the choice of digital frequencies bounded by $`1<1+\frac{1}{x}e`$ leads to support a logarithmic law distribution.

### IV.3 Recursions: A Possible Scenario for Number Prediction

From a physical point of view, the equilibrium of forces within defect-free plate under external shear can be approximated by the discrete Laplace equationTak85 ; Can91 ; Can92

$$\underset{k=1}{\overset{4}{}}G_k\mathrm{\Delta }u_k=0$$
(18)
where $`\mathrm{\Delta }u_k=u_k(i)u(i)`$ denotes relative displacements of the $`i`$th lattice point with respect to its nearest neighbors and $`G_k`$ is the corresponding modulus of rigidity of the bond connecting the $`i`$th lattice point to its neighbor.

From Eq.(18) it follows then that for a 1D system

$$u(i)=\widehat{G}_1u_1(i)+\widehat{G}_2u_2(i),$$
(19)
where $`\widehat{G}_i=G_i/_iG_i(i=1,2)`$.

It is not unreasonable to assume that the moduli of rigidity of our analogous elastic system are proportional to each other, at least prior to any equivalent breakdown condition which may prevail. Therefore if we consider $`u_1(i)(G_2/G_1)u_2(i1)`$, then we have that all displacements $`u_2`$, $`u_1`$ and $`u`$ follow generalized fibonacci sequences of the form Gol77

$$𝒰(i+2)=𝒰(i)+𝒰(i+1).$$
(20)
Within our thermodynamics framework, the sequences of $`x`$-integers played the role of analogous ’non-interacting’ particles. Besides, their magnitude may be seen as representing the above effective (increasing) displacements, say expected values $`x=[𝒰]`$, at different times $`i`$. The analogous force may be in turn seen as an applied shear acting on our 1D isolated elastic system of finite size $`x`$, which is being pulled out (or in).

Therefore following these ideas from fracture phenomena, a possible scenario to deduce second order polynomial recursions leading to predict sequences of numbers becomes possible.

## V Concluding Remarks

In this paper we have introduced a quite abstract statistical mechanical system whose states are represented by natural numbers. By considering some explicit probability distributions for the natural numbers, we claim that they lead to a better understanding of probabilistic laws associated with number theory.

Within this kind of thermodynamic formalism we have investigated probability distributions which lead to Benford’s law. Sequences of numbers have been treated as the size measure of finite sets and, in analogy with statistical thermodynamics, they have been assumed to play the role of particles of known type in an isolated elastic system.

The relation to physics is rather in the thermodynamical formalism and terminology for number theoretical problems. The assignment of statistical quantities here proposed are not arbitrary. Our analogies for the considered functions, e.g., Helmoltz free energy, etc, follows the similar thermodynamical approach as proposed within Multifractal physics. The expression of the Shannon entropy in Eq.(2) split into two terms with the identification of one of the terms with an analogous free energy and the other with an analogous internal energy can be justified (although not rigorously) in terms of the form of our canonical partition function which depends on the number probability distribution of Eq.(1) as may be expected. Furthermore as derived from such identification, we can argue that the relation in Eq.(17) tends to unity in the large $`x`$-limit in analogy with the standard thermodynamic definition $`T=\left(\frac{E}{S}\right)_{V,N}`$.

Though this scheme may be seen as rather abstract, we have shown that it leads to very interesting connections with statistical thermodynamics and allows to interpret physically the observed probability distribution of first digits in terms of a stretching force (tension) which ideally may pull integers number sets along analogous lengths.

The introduction of an analogous force within our formalism has been done on a formal level. The new insights gained by such a formal introduction and its utility include, in particular, the study of the dynamics of given number sets by establishing second order polynomial recursions to predict sequences of generalized Fibonacci numbers. This gives an answer to the interesting theoretical question of the appearance of the Benford’s law in Fibonacci numbers.

Following recent interest in Benford’s law Hil98 , we believe the new force may find application to detect fraud in financial data Nig96 , by looking for anomalies in their analogous macroscopic thermodynamic variables, in the use of logarithmic computers to speed calculations Knu81 , by analyzing the screening behavior of the analogous inverse square forces on floating points, or in the characterization of the hopcount in Internet Van00 .

## Acknowledgement

The author would like to thank the ARPL Group of the Abdus Salam ICTP, Trieste, for hospitality and support.


text:
# 1 Introduction

## 1 Introduction

It is a matter of common belief in the high-energy physics community that there must be physics beyond the Standard Model (SM) around the TeV scale. This belief is founded in the fact that some of the parameters of the SM can be made stable against quantum corrections only by invoking new physics at the TeV scale. With this in mind, an enormous international effort is being poured into the construction of the 14-TeV LHC at CERN, which is expected to run by 2008 and come up with signals for new physics. However, as the LHC is a hadronic machine, it is unlikely to yield signals which are unequivocal as to the underlying theory. It is therefore widely deemed necessary to build a high-energy $`e^+e^{}`$ collider to confirm the nature and detailed properties of this new physics.

A high-energy $`e^+e^{}`$ collider operating in the range 500 GeV to a few TeV must be of the linac type to avoid large synchrotron radiation losses. Such a machine would require a very high luminosity in order to be competitive with other accelerators. This is, in fact, an essential feature of the design of proposed linear colliders. The relatively clean environment of an $`e^+e^{}`$ collider allows the complete reconstruction of individual events, making it possible to carry out precision tests of new physics. Moreover, one has the facility of polarizing the beams, which can play a significant role in background elimination. It is also possible to run such a machine in the $`e^{}e^{}`$, $`e\gamma `$ and $`\gamma \gamma `$ modes. In this paper we consider the possibility of using a high energy $`e^+e^{}`$ collider to detect signatures of low-scale quantum gravity in a model with large extra dimensions. In particular, we focus on Bhabha scattering with radiated gravitons.

One of the most exciting theoretical developments of recent years has been the idea that there could be one or more extra spatial dimensions and the observable Universe could be confined to a four-dimensional hypersurface in a higher-dimensional ‘bulk’ spacetime. Such ideas, which can be motivated by superstring theories, give rise to elegant solutions to the well-known gauge hierarchy problem of high energy physics, which is just the instability against quantum corrections mentioned above. What is even more interesting, perhaps, is the suggestion that there could be observable signals of quantum gravity at current and future accelerator experiments, and this possibility has spawned a vast and increasing body of work over the past five years. This relatively new set of ideas, commonly dubbed ‘Brane World Phenomenology’, bases itself on two main principles: the concept of hidden compact dimensions and the string-theoretic idea of $`D_p`$-branes. The simplest brane-world scenario is the so-called Arkani-Hamed—Dimopoulos—Dvali (ADD) model, in which there are $`d`$ extra spatial dimensions, compactified on a $`d`$-torus of radius $`R_c`$ each way. Together with the four canonical Minkowski dimensions, this constitutes the ‘bulk’ spacetime. In this scenario the radius $`R_c`$ of the extra dimensions can be as large as a quarter of a millimeter. However, the SM fields are confined to a four-dimensional slice of spacetime, with thickness not more than $`10^{17}`$ cm, which is called the ‘brane’. If the ADD model is embedded in a string-theoretic framework, the ‘brane’ is, in fact a $`D_3`$-brane, i.e. a 3+1 dimensional hypersurface on which the ends of open strings are confined. However, it is not absolutely essential to embed the model in a string theory, and the word ‘brane’ or ‘wall’ is then used simply to denote the hypersurface (or thin slice) where the SM fields are confined. A crucial feature of this model is that gravity, which is a property of spacetime itself, is free to propagate in the bulk. As a result

* Planck’s constant in the bulk $`M_S`$, identified with the ‘string scale’, and fixed by the string tension $`\alpha ^{}`$, is related to Planck’s constant on the brane $`M_P`$ ($`1.2\times 10^{19}`$ GeV) by

$$(M_S)^{2+d}=(4\pi )^{d/2}\mathrm{\Gamma }(d/2)M_P^2(R_c)^d$$
(1)
which means that for $`R_c0.2`$ mm, it is possible to have $`M_S`$ as low as a TeV for $`d2`$. (The normalization of reference  has been adopted here). This solves the gauge hierarchy problem simply by bringing down the scale of new physics (i.e. strong gravity in this case) to about a TeV and thereby providing a natural cut-off to the SM, since the string scale $`M_S`$ now controls graviton-induced processes on the brane.
* There are a huge number of massive Kaluza-Klein excitations of the (bulk) graviton field on the brane, with masses $`m_n=n/R_c`$, and these collectively produce gravitational excitations of electroweak strength, which may be observable at current experiments and those planned in the near future.

It is only fair to mention that a major drawback of the ADD model is that it creates a new hierarchy between the ‘string scale’ $`M_s1`$ TeV and the size of the extra dimensions $`R_c^11`$ $`\mu `$eV. In fact, the huge size of the extra dimensions (compared to the Planck length) is not stable under quantum corrections, which tend to shrink it down until $`M_sR_c^1M_P10^{19}`$ GeV, at which stage the original hierarchy problem is reinstated. Nevertheless, there are several variants of the ADD model which address this problem in various ways, and some of these ideas may not be far from the truth. From a phenomenological point of view, it is, therefore, reasonable to postpone addressing the stability issue, and proceed to study the minimal ADD model and its consequences for experiment.

The experimental consequences of ADD gravity have been mainly studied in the context of of ($`a`$) real graviton emission and ($`b`$) virtual graviton exchange. In the former case the final state gravitons in the ADD model are ‘invisible’, escaping the detector because of their feeble individual interactions ($`M_P^1`$) with matter. The final state, involving missing energy due to gravitons, will be built up by making an incoherent sum over the tower of graviton modes. In the case of virtual gravitons, the final state is built up by making a coherent sum. In either case, it may be shown that, after summing, the Planck mass $`M_P`$ cancels out of the cross-section, leaving an interaction of near-electroweak strength.

At an $`e^+e^{}`$ collider, the most-frequently discussed signal for ADD gravitons is the process $`e^+e^{}\gamma G_n`$, where $`G_n`$ is the $`n`$th Kaluza-Klein excitation of the graviton field. This process leads to single photon events with missing energy and momentum, and is expected to be among the earliest signals of ADD gravity at an $`e^+e^{}`$ collider. In addition, the process $`e^+e^{}\mu ^+\mu ^{}G_n`$ has also been studied. In this work, we have studied the process

$$e^+e^{}e^+e^{}G_n$$

which is simply Bhabha scattering with a radiated graviton. The final state would contain an electron-positron pair with substantial missing energy. The reasons for taking up a study of this process are as follows.

* At an $`e^+e^{}`$ collider, final states with an $`e^+e^{}`$ pair will be one of the first things to be analyzed, since this is needed for beam calibration. Among such states, it should be a simple matter to select events with large missing energy as well.
* Unlike $`e^+e^{}\mu ^+\mu ^{}G_n`$, this process does not suffer from large $`s`$-channel suppression, especially if one considers final state electrons (positrons) in relatively forward (backward) directions.
* Taken in conjunction with the process $`e^+e^{}\gamma G_n`$, this process could help to distinguish signals of the ADD scenario from other physics models beyond the SM. This is elaborated in the subsequent discussions.
* In some ways, as we shall see, this process has signatures more distinct from other new physics effects than $`e^+e^{}\gamma G_n`$.

In the following section we describe the salient features of our calculation of the signal process. Backgrounds and some strategies for their elimination are discussed in section 3. Section 4 contains the results of our numerical analysis. In section 5 we discuss how this process can be used to pinpoint the model if signatures of the suggested type are indeed observed. We summarize and conclude in the final section.

## 2 The Process $`e^+e^{}e^+e^{}G_n`$

As gravity couples to the energy-momentum tensor, each component of the tower of ADD gravitons couples to all SM fields, as well as to each SM interaction vertex. The relevant Feynman rules can be found, for example, in Refs.. Using these, it can be shown that the process $`e^+e^{}e^+e^{}G_n`$ is driven, at tree-level, by 28 Feynman diagrams. Some representative diagrams are shown in Fig. 1, which are obtained by ‘dressing’ an $`s`$-channel process $`e^+e^{}e^+e^{}`$ with graviton radiation. Fig. 1($`a`$) represents a graviton emission from one of the external electron legs; there will be 4 such diagrams with the graviton emitted from each leg in turn. Fig. 1($`b`$) represents a graviton emission from the gauge boson propagator. Fig. 1($`c`$) represents a graviton emission from one vertex; there will be another such with graviton emission from the other vertex. Thus Fig. 1 encompasses 14 diagrams, 7 each with $`Z`$ and photon exchanges. In addition, there will be another set of 14 diagrams with the $`\gamma ,Z`$ exchange in the $`t`$-channel, making 28 diagrams in all. This latter set is absent in the case of $`e^+e^{}\mu ^+\mu ^{}G_n`$.

Figure 1. Representative Feynman diagrams for the process $`e^+e^{}e^+e^{}G_n`$.

Considering the process

$$e^{}(k_1,\lambda _1)e^+(k_2,\lambda _2)G_n(p_1)e^+(p_2)e^{}(p_3)$$

with a single Kaluza-Klein mode $`G_n(p_1)`$ in the final state, we obtain the unpolarized cross-section

$$\sigma (m_n)=\frac{1}{2s}\frac{d^3\stackrel{}{p_1}}{(2\pi )^32E_1}\frac{d^3\stackrel{}{p_2}}{(2\pi )^32E_2}\frac{d^3\stackrel{}{p_3}}{(2\pi )^32E_3}\delta ^4(k_1+k_2p_1p_2p_3)\frac{1}{4}\underset{\lambda _1,\lambda _2}{}|M_n(\lambda _1,\lambda _2)|^2$$
(2)
where the final state helicities are summed over. This cross-section is a function of the mass $`m_n`$ of the graviton excitation.

The squared matrix element $`|M_n(\lambda _1,\lambda _2)|^2`$ for the process can be written as the sum

$$|M_n(\lambda _1,\lambda _2)|^2=\left|\underset{i=1}{\overset{28}{}}M_n^{(i)}(\lambda _1,\lambda _2)\right|^2$$
(3)
where $`|M_n^{(i)}(\lambda _1,\lambda _2)|^2`$ represents the squared helicity amplitude arising from the $`i`$th diagram. In our calculation, we make use of the helicity amplitude technique to write down amplitudes for all these 28 diagrams and evaluate the individual terms in the above sum using the subroutine HELAS.

Recalling that the observed final state, namely, $`e^+e^{}`$ plus missing energy, consists of an incoherent sum over the tower of graviton modes, we get

$`\sigma (e^+e^{}e^+e^{}+\overline{)}E)`$ $`=`$ $`{\displaystyle \underset{n}{}}\sigma (m_n)`$ (4)
$``$ $`{\displaystyle _0^\sqrt{s}}𝑑m\rho (m)\sigma (m)`$
approximating the discrete (but closely-spaced) tower of states by a continuum. The density of states is given by

$$\rho (m)=\frac{2R_c^dm^{d1}}{(4\pi )^{d/2}\mathrm{\Gamma }(d/2)},$$
(5)
and the integration is cut off at the kinematic limit $`\sqrt{s}`$.

As the calculation, even using HELAS, is long and cumbersome, some checks on the numerical results are called for. The most useful check is provided by the Ward identities arising from general coordinate invariance, which constitute an essential feature of any theory involving gravity. We can write the amplitude for the emission of any graviton in the form

$$M_n(\lambda _1,\lambda _2)=T^{\mu \nu }(\lambda _1,\lambda _2)ϵ_{\mu \nu }^{(n)}(p_1)$$
(6)
where $`ϵ_{\mu \nu }^{(n)}(p_1)`$ is the polarization tensor for the $`n`$th (massive) graviton mode. The tensor $`T^{\mu \nu }(\lambda _1,\lambda _2)`$ is the same for every mode, including the massless mode $`ϵ_{\mu \nu }^{(0)}(p_1)`$, which is the usual graviton of four-dimensional Einstein gravity. This must now satisfy the Ward identities

$$p_1^\mu T_{\mu \nu }(\lambda _1,\lambda _2)=p_1^\nu T_{\mu \nu }(\lambda _1,\lambda _2)=0,$$
(7)
where we note that

$`T^{\mu \nu }(\lambda _1,\lambda _2)={\displaystyle \underset{i=1}{\overset{28}{}}}T_i^{\mu \nu }(\lambda _1,\lambda _2)`$ (8)
with $`i`$ indicating the $`i`$th diagram, as above. The consistency check therefore requires a perfect cancellation, for each choice of $`\lambda _1`$ and $`\lambda _2`$, between all such terms222In fact the $`s`$ and $`t`$-channel diagrams form two independent gauge-invariant sets, so that the actual cancellation takes place between 14 diagrams at a time., which is highly sensitive to errors in signs and factors. We have also checked that if only the $`s`$-channel results are taken our numerical results are in close agreement with those of Ref., and this, of course, implies that the $`t`$-channel contributions can be easily checked using crossing symmetry.

In our numerical analysis of this problem, we have also considered diagrams mediated by gravitons in place of $`\gamma ,Z`$. Such diagrams, in fact, are of the same order in perturbation theory after summation over the graviton propagators, but it turns out that the contributions are very small — typically less than 1% for experimentally-allowed values of $`M_S`$ — due to the smallness of the effective graviton coupling after summation.

It is also worth mentioning that contributions due to graviscalars in the ADD model can also be neglected. It can be shown that diagrams with graviscalars coupling to external fermion legs and vertices undergo large cancellations and the final contribution is proportional to the electron mass. This leaves only the diagrams with graviscalars emitted from the massive $`Z`$-propagator, which is again suppressed at high energies by a factor of $`M_Z^2/s`$ compared to the gravitensor contributions, the latter being jacked up by derivative couplings.

## 3 Background Elimination

The actual signal for the process considered in the previous section consists of a hard electron and a hard positron, with substantial missing energy. The SM backgrounds arise from all processes of the form $`e^+e^{}e^+e^{}\nu \overline{\nu }`$, where the neutrinos can be of any flavor. At the lowest order, there are 96 Feynman diagrams which give rise to this final state. Some of the principal sub-processes are listed below.

$`e^+e^{}`$ $``$ $`e^\pm W^{}\nu _e^{()}e^+e^{}\nu _e\overline{\nu }_e`$
$`e^+e^{}`$ $``$ $`e^+e^{}Ze^+e^{}(\nu \overline{\nu })`$
$`e^+e^{}`$ $``$ $`\nu _e\overline{\nu }_eZ\nu _e\overline{\nu }_e(e^+e^{})`$
$`e^+e^{}`$ $``$ $`W^+W^{}(e^+\nu _e)(e^{}\overline{\nu }_e)`$
$`e^+e^{}`$ $``$ $`ZZ(e^+e^{})(\nu \overline{\nu })`$
$`e^+e^{}`$ $``$ $`\gamma ^{}Z(e^+e^{})(\nu \overline{\nu })`$ (9)
Of course, for every process with real $`W/Z`$, there will also be many diagrams with off-shell particles. Taken together, all these diagrams constitute a formidable background to the suggested signal. A judiciously chosen set of kinematic cuts, however, enable us to reduce these backgrounds enormously. For example, diagrams with the neutrino pair (missing energy) arising from a real $`Z`$ boson can be easily removed by putting a cut on the $`\nu \overline{\nu }`$ (missing) invariant mass. A similar cut on the $`e^+e^{}`$ pairs would seem obvious but it is not advisable, since it would remove a large part of the signal arising from radiative return to the $`Z`$-pole through graviton emission. This, however, is not a serious problem, as the background contribution due to $`e^+e^{}\nu _e\overline{\nu }_eZ\nu _e\overline{\nu }_e(e^+e^{})`$ is not very large. On the other hand, the background due to the process $`e^+e^{}e^\pm W^{}\nu _e^{()}`$ is more problematic, as most cuts tending to reduce this also affect the signal adversely.

In principle, an additional source of background can be the ‘two-photon’ process $`e^+e^{}e^+e^{}e^+e^{}`$ where one electron-positron pair escapes detection by being emitted close to the beam pipe. Here one has to remember that such forward electrons, being rather hard (with energy on the order of $`100`$ GeV or above), should be detectable in the end caps of the electromagnetic calorimeter, although their energy resolution will be rather poor. Thus it should be possible to impose a veto on hard electrons up to within a few degrees of the beam pipe, whereby the above background can be virtually eliminated. Events with the forward electrons even closer to the beam pipe are removed via a cut on missing $`p_T`$.

Taking all these considerations into account, we focus on an $`e^+e^{}`$ collider operating at a center-of-mass energy of 500 GeV (1 TeV). The initial set of kinematic cuts used in our analysis are as follows:

* The final state electron (positron) should be at least $`10^0`$ away from the beam pipe. This tames the collinear singularities arising from $`t`$-channel photon exchange. At the same time, it also ensures that any background effects from beamstrahlung are mostly eliminated.
* The electron (positron) should have a transverse momentum $`p_T^e>10`$ GeV.
* We demand a missing transverse momentum $`p_T^{\mathrm{miss}}>15(25)`$ GeV for $`\sqrt{s}=`$ 500 GeV (1 TeV). As mentioned above, this also helps in reducing the two-photon background.
* The tracks due to electron and positron must be well-separated, with $`\mathrm{\Delta }R>0.2`$, where $`R=\sqrt{\mathrm{\Delta }\eta ^2+\mathrm{\Delta }\varphi ^2}`$, in terms of the pseudorapidity $`\eta `$ and the azimuthal angle $`\varphi `$.
* The opening angle between the electron and positron tracks is required to be limited by $`5^0<\theta _{e^+e^{}}<175^0`$. This ensures elimination of possible cosmic ray backgrounds, and also ensures sufficient missing energy.
* The missing invariant mass $`M_{\mathrm{miss}}`$ should satisfy the cut $`|M_{\mathrm{miss}}M_Z|>10`$ GeV. This important cut eliminates many SM backgrounds in which two final state neutrinos arise from a $`Z`$-decay.

Incorporating all these cuts into a Monte Carlo event generator, we obtain the distribution in missing invariant mass shown in Figure 2. It is clear from the figure that the missing invariant mass of the signal is much harder than that of the background. This is easy to understand if we realize that the missing mass is a measure of the energy of the emitted graviton. Now it is well-known that the higher the energy, the more strongly the graviton is coupled to matter and consequently the higher the production cross-section. Moreover, a higher energy leads to a higher density of states, and this further enhances the cross-section. In Figure 2, we have also shown the value of the signal when only $`s`$-channel diagrams are considered, which would be the case if we were to study processes like $`e^+e^{}\mu ^+\mu ^{}G_n`$ or $`e^+e^{}\tau ^+\tau ^{}G_n`$. Comparison with the signal of interest makes it obvious that there is a substantial contribution to the Bhabha scattering signal from the $`t`$-channel, which makes it more viable as a signal of ADD gravity than the other two. In fact, the emission of a high-energy graviton (corresponding to large missing invariant mass) reduces the energy of the underlying Bhabha scattering process, causing a large enhancement due to the $`Z`$-pole. This is the reason why the signal for large $`M_{\mathrm{miss}}`$ is completely dominated by the $`s`$-channel contribution. At the same time, when the graviton is relatively soft, the $`s`$-channel diagrams are strongly suppressed, but the $`t`$-channel contributions are not; as a result the signal is dominated by the $`t`$-channel contribution for low $`M_{\mathrm{miss}}`$. We can therefore obtain a clear separation of signal from background by imposing a cut

$$M_{\mathrm{miss}}>350(450)\mathrm{GeV}\mathrm{for}\sqrt{s}=500\mathrm{GeV}(1\mathrm{TeV})$$

These choices are optimal, after taking into account the fact that the signal falls rapidly with increasing $`M_s`$ and $`d`$.

An important strategy for background reduction is the use of the beam polarization facility at a linear collider. If the electron (positron) beam has a right (left) polarization efficiency $`𝒫_e`$ ($`𝒫_p`$), the cross-section formula corresponding to equation 2 is obtained by the replacement

$`{\displaystyle \underset{\lambda _1,\lambda _2}{}}|M_n(\lambda _1,\lambda _2)|^2`$ $``$ $`(1+𝒫_e)(1𝒫_p)|M_n(+,+)|^2+(1+𝒫_e)(1+𝒫_p)|M_n(+,)|^2`$
$`+(1𝒫_e)(1𝒫_p)|M_n(,+)|^2+(1𝒫_e)(1+𝒫_p)|M_n(,)|^2`$
Typical values for the polarization efficiencies used in this analysis are $`𝒫_e=0.8`$ and $`𝒫_p=0.6`$, which tend to favor the first two terms on the right side of the above expression. The use of polarized beams leads to a drastic reduction in the SM background, chiefly because it causes suppression of the $`W`$-induced diagrams. The numerical effects are presented in the next section.

## 4 Results and Discussions

After numerical evaluation of signal and background, we find it convenient to present our results, using the significance $`S/\sqrt{B}`$, where $`S(B)=\sigma _{S(B)}`$, for given values of $`d`$ and $`M_S`$, Here $`\sigma _{S(B)}`$ denotes the cross-section for the signal (background), while $``$ is the integrated luminosity. Our numerical results are presented assuming $`=500`$ fb-1; however, it is a simple matter to scale the significance values for different values of luminosity.

In Figure 3($`a`$) we show curves showing the variation of the significance with the string scale $`M_S`$, for $`d=2,3,4,5`$ and $`6`$ respectively, for $`\sqrt{s}=500`$ GeV. Figure 3($`b`$) shows a similar plot for $`\sqrt{s}=1`$ TeV. In both cases, we have considered unpolarized electron and positron beams.

As has been already mentioned, the use of beam polarization improves changes of detecting the signal quite dramatically. This is illustrated in Figure 4, which is similar to Figure 3, except that the polarization efficiencies have been taken to be $`𝒫_e=0.8`$ and $`𝒫_p=0.6`$ for the electron and positron respectively.

From the above graphs, it is straightforward to find out the maximum value of $`M_s`$ that can be probed at the linear collider for any given value of $`d`$, both for polarized and unpolarized beams. In Table 1, we present the values of $`M_s`$ that one can probe at 99% confidence level.

$$\begin{array}{cccccccc}& & & & & & & \\ \sqrt{s}& \hfill d=& 2& 3& 4& 5& 6& \\ & & & & & & & \\ & & & & & & & \\ 0.5& \hfill \mathrm{polarized}& 4.53& 2.79& 2.01& 1.60& 1.35& \\ & \hfill \mathrm{unpolarized}& 3.16& 2.09& 1.59& 1.30& 1.12& \\ & & & & & & & \\ & & & & & & & \\ 1.0& \hfill \mathrm{polarized}& 7.40& 4.63& 3.41& 2.76& 2.36& \\ & \hfill \mathrm{unpolarized}& 5.36& 3.57& 2.74& 2.28& 2.36& \end{array}$$

Table 1. 99% C.L. discovery limits on the string scale $`M_S`$ for an ADD scenario using (graviton) radiative Bhabha scattering, for different numbers of extra dimensions. An integrated luminosity of $`500fb^1`$ has been assumed. All numbers are in TeV.

If one compares these with the corresponding reach of the process $`e^+e^{}\gamma +\overline{)}E`$, one will notice that search limits are approximately of the same order. However, the comparison should not be made too literally, because of several reasons. First, the polarization efficiencies in our case, while matching with those of , are different from those of , while the center-of-mass energy at which reference  has calculated the effects is slightly different from ours. In addition, one needs to match the event selection criteria more carefully for a full-fledged comparison. Finally, it should be borne in mind that in our analysis we have adopted a normalization of the string scale $`M_s`$, which, though widely used, is not necessarily uniform in the literature. However, in spite of such non-uniformities, Figures 3 and 4, together with Table 1, can perhaps be taken as faithful indications of the fact that the predictions on radiative Bhabha scattering are comparable to those on the photon-graviton channel, so far as the limits of the probe on the string scale in ADD models are concerned.

## 5 Pinning Down the Model

A few comments are in order on how to distinguish graviton signals involving missing energy from similar ones arising from other kinds of new physics. For example, let us consider the well-known process $`e^+e^{}\gamma G_n`$, leading to a single-photon-plus-missing-energy signal. Such signals can be obtained in several other models. A model with an extra $`Z^{}`$ boson could lead to a process $`e^+e^{}\gamma Z^{}`$ with the $`Z^{}`$ decaying to neutrinos. Similarly, in supersymmetric models, a process like $`e^+e^{}\gamma \stackrel{~}{\chi }_1^0\stackrel{~}{\chi }_1^0`$ (where $`\stackrel{~}{\chi }_1^0`$ is the lightest neutralino) would also lead to a single photon signal. Although the $`Z^{}`$ signal can be differentiated by looking for a resonant peak in the missing invariant mass, the supersymmetric signal is more difficult to disentangle.

The situation for the process considered in this paper is somewhat more encouraging. It is true that both types of new physics considered in the preceding paragraph can produce the same signal. Processes of the form $`e^+e^{}e^+e^{}Z^{}`$, with the $`Z^{}`$ decaying to neutrinos are one possibility. In supersymmetry, too, it is possible to pair-produce either sleptons or charginos, followed by decays to electrons (positrons) and invisible particles. The $`Z^{}`$ signal can again be differentiated by a peak in the missing invariant mass. The supersymmetric signal, in this case, always arises from production of a pair of real sparticles, which should emerge in opposite hemispheres and, if light enough, will be highly boosted. As a result, the observed electron and positron should lie, most of the time, in opposite hemispheres. For the graviton signal, however, there is no such correlation; in fact the distribution in the opening angle between $`e^+`$ and $`e^{}`$ turns out to be practically flat. This distinction clearly does not work very well when the produced sparticles are heavy.

A more elegant way of distinguishing these signals from those other forms of new physics is to compare the cross-sections arising from both the processes $`e^+e^{}\gamma G_n`$ and $`e^+e^{}e^+e^{}G_n`$, both of which are determined by the two parameters $`d`$ and $`M_S`$, and hence, will have some correlation. In Figure 5, we have plotted the cross-sections $`\sigma (\gamma G_n)`$ versus $`\sigma (e^+e^{}G_n)`$ for different values of $`d`$. In view of the discussion in the previous section, we would like to note that here the rates for $`e^+e^{}\gamma G_n`$ have been computed with the same normalization of the string scale as that used for the Bhabha scattering process. Each curve corresponds to variation of $`M_S`$ over the range allowed by current experimental constraints. Once experimental data are available from a high energy $`e^+e^{}`$ machine, we should be able to pinpoint a small region (ideally a point) on the graph(s) in Figure 5. The position of this would immediately direct us to the value of $`d`$; comparison of the cross-sections would now yield a measurement of the value of $`M_S`$. The identification of the model parameters should be unambiguous since the curves for different values of $`d`$ do not intersect except at the origin. In any case, the area in the vicinity of the origin is limited by the search limits indicated in Table 1 since it corresponds to very high values of $`M_S`$. The shaded region corresponds to $`M_S=\sqrt{s}`$ beyond which the theoretical calculations are unreliable. It should, of course be remembered that further information on the fundamental parameters of theory can be extracted from various kinematic distributions which are sensitive to these parameters.

## 6 Summary and Conclusions

In the paper we have considered the process $`e^+e^{}e^+e^{}G_n`$, which is essentially Bhabha scattering with radiated gravitons in the ADD model. After identifying suitable event selection criteria, we find that this process can act as an effective probe of large extra dimensions at a high-energy $`e^+e^{}`$ collider, especially with a center-of-mass energy of order 1 TeV and with polarized beams. The string scale $`M_S`$ that can be probed in this channel is found to be comparable to that which is accessible through the alternative process $`e^+e^{}\gamma G_n`$. It is also shown that a study of the (graviton) radiative Bhabha scattering process may provide some further handle on the essential characteristics of ADD-like theories. And finally, taking a cue from the wisdom that, in coming to any conclusion on new physics possibilities, it is always advantageous to have more than one type of data, we have demonstrated how our predictions can be combined with those on the photon-graviton channel to obtain rather trustworthy revelations on models with large extra dimensions.

Acknowledgments

The authors acknowledge useful discussions with D.Choudhury, R.M.Godbole, U.Mahanta, and S.K.Rai. This work was initiated as part of the activity of the Indian Linear Collider Working Group (Project No. SP/S2/K-01/2000-II of the Department of Science and Technology, Government of India). SD thanks the SERC, Department of Science and Technology, Government of India for partial support. The work of BM was partially supported by the Board of Research in Nuclear Sciences (BRNS), Government of India. SR thanks the Harish-Chandra Research Institute for hospitality while this paper was being written.