text:
# Adsorption in non interconnected pores open at one or at both ends: a reconsideration of the origin of the hysteresis phenomenon.

## Abstract

We report on an experimental study of adsorption isotherm of nitrogen onto porous silicon with non interconnected pores open at one or at both ends in order to check for the first time the old (1938) but always current idea based on Cohan’s description which suggests that the adsorption of gas should occur reversibly in the first case and irreversibly in the second one. Hysteresis loops, the shape of which is usually associated to interconnections in porous media, are observed whether the pores are open at one or at both ends in contradiction with Cohan’s model.

Adsorption isotherms in mesoporous materials exhibit two main typical features: (i) a sharp increase in the amount of gas adsorbed at pressure below the saturated vapor pressure of the bulk gas which is attributed to capillary condensation in pores of confining sizes; (ii) the desorption process fails to retrace the path of the adsorption one, so that a hysteresis loop appears revealing the irreversibility of the phenomenon.

The first explanation of the capillary condensation was proposed by Zsigmondy  who, for this purpose, refered to the macroscopic Kelvin equation. In order to explain the hysteresis phenomenon, Cohan  went further and proposed different scenari for filling and emptying of a pore which he also described with the macroscopic concept of meniscus. In a cylindrical pore open at both ends, the shape of gas/adsorbate interface is different during the adsorption and desorption processes leading to a hysteresis phenomenon. On the contrary, in a cylindrical pore closed at one end, the adsorption is reversible since the same meniscus is present during the filling and the emptying of the pore.

The irreversibility of the adsorption process in a cylindrical pore open at both ends, expected by Cohan on the basis of thermodynamic considerations, has been predicted by new numerical approaches. Thus, Grand Canonical Monte-Carlo (GCMC) simulations , analytical approaches  together with Density Functional Theory (DFT)  have also proposed that the hysteresis phenomenon is an intrinsic property of the phase transition in a single idealized pore and arises from the existence of metastable gas-like and/or liquid-like states. It is noteworthy that no calculations have been yet performed for cylindrical pores closed at one end.

Macroscopic descriptions based on Cohan’s model are currently used to describe the adsorption in non interconnected mesoporous media and in particular to determine the Pore Size Distribution (PSD) via the Barrett-Joyner-Halenda (BJH) method .

However, the experimental validation of Cohan’s model has never been checked since the synthesis of ordered mesoporous media with non interconnected pores open at one or at both ends was unsuccesful for a long time. Using porous silicon, we report in this paper, for the first time, on an experimental investigation of Cohan’s predictions. For this check, this mesoporous material exhibits two essential properties: (i) the tubular pores are non interconnected and do not exhibit narrow sections and (ii) it exists under two forms with pores closed at one end when the porous layer is supported by the substrate or open at both ends when the layer is removed from the substrate.

Porous silicon is prepared by electrochemical etching of Si single crystal in a HF solution. The morphology of the porous Si layers depends on the type and concentration of the dopant atoms . The adsorption studies presented here have been performed with porous layers prepared from highly boron doped  Si substrate. The porosity of the layers is controlled by the current density and HF concentration, and the thickness by the duration of the anodic dissolution. The pore volume is determined with high precision by weighting the sample before and after the electrochemical dissolution. To calibrate the porosity and the thickness of the layer we need a third weighing after the dissolution of the porous part in a NaOH solution. It is noteworthy that the porous layers are perfectly reproduced. The porous layers under study had a porosity of 51% and the length of the pores is 21 $`\mu `$m.
The morphology of these layers (see Fig. 1(a) and 1(b)) can be described as honeycomb-like structure with pores perpendicular to the Si substrate separated from each other by Si single crystal walls (inset Fig. 1(b)) of constant thickness ($``$ 5 nm).

The pores, with sections of polygonal shape, are of the same length and closed at one end by the Si substrate. At the end of the porous layer formation which takes about 20 minutes, an increase of the current density leads to the electropolishing regime during which the Si walls at the bottom of the pores are dissolved so that the porous layer comes off the Si substrate. It is noteworthy that such a process takes only a few seconds and thus that the PSD of the porous layer is not changed. Playing with this possibility, porous Si layer or porous Si membrane can be obtained with pores open at one or at both ends, respectively.

FIG. 1.: Bright field TEM images in plane view (a) and cross-section (b) of a porous Si layer with 51% porosity prepared from $`p^+`$ type ($`3.10^3\mathrm{\Omega }`$.cm)  Si substrate. Observation axis . Pores (white) are separated by Si walls (black). The pore density is $`\mathrm{3.2.10}^{11}/\mathrm{cm}^2`$. Inset 1a: PSD. Inset 1b: transmission electronic diffraction pattern. We have shown that there is no lateral interconnection between the pores in these layers. The method involves deposition by evaporation at ambient temperature, of an aluminum layer, 500 nm thick, on a part ($`\mathrm{A}_\mathrm{c}`$ and $`\mathrm{B}_\mathrm{c}`$) of the top of two porous samples A and B, following by a thermal annealing under neutral atmosphere of argon at $`450^o`$C for half a hour to improve the contact between the Al deposit and the top of the layers. We have shown by Rutherford BackScattering analysis (RBS) with 2 MeV $`\alpha `$ particles that, the Al deposited on the samples is well localized on the top of the pores and constitutes a cap. One of the two samples (A) was thermally oxidized in $`\mathrm{O}_2`$ enriched in $`{}_{}{}^{18}\mathrm{O}`$ (99%), at $`300^o`$C, 12 mbar, for one hour. This treatment leads to the formation of a thin silicon oxide ($``$1 nm) on the silicon walls. The $`{}_{}{}^{18}\mathrm{O}`$ contents in the porous layers under the Al cap ($`\mathrm{A}_\mathrm{c}`$ and $`\mathrm{B}_\mathrm{c}`$) and beside it ($`\mathrm{A}_\mathrm{o}`$) were determined by the nuclear reaction $`{}_{}{}^{18}\mathrm{O}`$(p,$`\alpha `$)$`{}_{}{}^{15}\mathrm{N}`$ at a proton energy near the narrow resonance at $`\mathrm{E}_\mathrm{p}`$= 629 keV ($`\mathrm{\Gamma }`$ = 2 keV). For the measurements under the Al caps the energy of the protons beam energy was increased to take into account the energy loss of the incident protons in the Al layer. The nuclear reaction spectra corresponding to $`\mathrm{A}_\mathrm{o}`$ and $`\mathrm{A}_\mathrm{c}`$ are shown in Fig. 2. The $`\alpha `$ peak from $`\mathrm{A}_\mathrm{o}`$ indicates of the presence in the porous layer of $`{}_{}{}^{18}\mathrm{O}`$ atoms coming from the thermal oxidation. Due to the presence of Al cap, the $`\alpha `$ peak from $`\mathrm{A}_\mathrm{c}`$ is expected at an energy lower by 100 keV than that from the $`\alpha `$ peak from $`\mathrm{A}_\mathrm{o}`$. The counts from $`\mathrm{A}_\mathrm{c}`$ and $`\mathrm{B}_\mathrm{c}`$ which are one order of magnitude lower than that from $`\mathrm{A}_\mathrm{o}`$ are equal, and hence corresponds only to a background coming from nuclear reactions with silicon, boron atoms, etc. This experiment reveals that this medium is composed of non interconnected single pores, the top of which is in direct contact with the gas reservoir during adsorption. To our knowledge, this is the only porous material for which the absence of interconnection has been experimentally proven.
FIG. 2.: $`\alpha `$ spectra from $`{}_{}{}^{18}\mathrm{O}`$(p,$`\alpha `$)$`{}_{}{}^{15}\mathrm{N}`$: $`\mathrm{E}_\mathrm{p}`$= 629 keV. White circles: porous Si layer after thermal oxidation in $`{}_{}{}^{18}\mathrm{O}_{2}^{}`$ ($`\mathrm{A}_\mathrm{o}`$). Black circles: part of the layer protected by Al cap ($`\mathrm{A}_\mathrm{c}`$). The arrow indicates the energy of the expected $`\alpha `$ peak from $`\mathrm{A}_\mathrm{c}`$.  Another important property of this medium is the constancy of the Si walls thickness revealed by the TEM plane view performed on a porous layer thinned down to 1000 Å, i.e. about 10 times the size of the pores. This is a strong argument against the presence of narrow sections inside the pores, the presence of which would lead, through shadow effects, to an inhomogeneous apparent thickness of the Si walls.

The PSD of the porous Si layers have been estimated from a numerical treatment of the TEM images. The negative TEM photographs are digitized by means of a CCD monochrome camera. Then, the value of the threshold to make the image binary is chosen to reproduce the porosity of the sample. The analysis of this binary image allows one to obtain the cross section area and perimeter of each pore. The PSD, which corresponds to cylindrical pores having the same surface area as the polygonal pores, is shown in inset of Fig. 1(a). The mean diameter of the pores is 13nm with a standard deviation of $`\pm `$6nm.

Adsorption measurements were made by a volumetric technique using a commercial apparatus under secondary vacuum of around $`10^6`$ torr. Nitrogen adsorption isotherms were measured up to the saturating vapor pressure of bulk nitrogen (P0) and at 77 K. Porous Si layers were introduced in a cell just after their preparation. The samples were outgassed at ambient temperature for a few hours. This outgassing does not change the nearly complete passivation of the inner surface by a monolayer of $`\mathrm{SiH}_{1,2,3}`$ groups which takes place during the electrochemical etching and which is stable for a few days after the formation .
FIG. 3.: (a) N2 adsorption isotherms at 77K for porous Si with pores open at one (black diamonds) or at both (white diamonds) ends. (b) Isotherms calculated by introducing in Eqs. (1) and (2) the experimental t($`\mathrm{P}`$) values (see text) and the PSD shown in Fig. 1.

Figure 3(a) shows N2 adsorption isotherms at 77 K corresponding to porous Si layers with pores open at one or at both ends, the PSD of which is represented in inset of Fig. 1(a). These two curves exhibit wide and asymmetrical hysteresis loops with shapes corresponding to type H2 with reference to the 1985 IUPAC classification (or type E in that of de Boer ).

The adsorbed amount of gas sharply increases between 0.55 and 0.8 P0 and then reaches a plateau region where all the pores are filled with the dense phase. The steep desorption process occurs over a pressure range about three times narrower than that of the adsorption.

On the assumption that the density of the dense phase equals that of the bulk liquid nitrogen, we found that the pore volume extracted from the value of the adsorbed amount on the plateau region is in a very good agreement with that measured by gravimetry.

Cohan’s model predicts that, for a cylindrical pore open at both ends, the condensation takes place at a pressure $`\mathrm{P}_\mathrm{c}`$ given by an equation similar to the classical Kelvin equation where the curvature of the meniscus is replaced by a cylindrical curvature :

$$\mathrm{ln}\left(\frac{\mathrm{P}_\mathrm{c}}{\mathrm{P}_0}\right)=\frac{\sigma \mathrm{V}_\mathrm{L}}{\mathrm{RT}(\mathrm{r}\mathrm{t}(\mathrm{P}_\mathrm{c}))}$$
(1)
where VL is the molar liquid volume, $`\sigma `$ is the liquid-gas surface tension at temperature T, r is the radius of the pore and t($`\mathrm{P}`$) the thickness of the adsorbed layer at P. In the case of a cylindrical pore closed at one end, surface thermodynamics considerations predict a preferential adsorption on the bottom of the pore; a hemispherical meniscus is then present provided that there is perfect wetting, so that the condensation occurs at a pressure $`\mathrm{P}_\mathrm{h}`$ given by the following modified Kelvin equation:

$$\mathrm{ln}\left(\frac{\mathrm{P}_\mathrm{h}}{\mathrm{P}_0}\right)=\frac{2\sigma \mathrm{V}_L}{\mathrm{RT}(\mathrm{r}\mathrm{t}(\mathrm{P}_\mathrm{h}))}$$
(2)
The evaporation of the wetting fluid proceeds through the same hemispherical meniscus in both cases i.e. at the same pressure $`\mathrm{P}_\mathrm{h}`$.

Thus, according to Cohan’s model : (i) the adsorption phenomenon should be reversible in pore open at one end while it should be irreversible in pore open at both ends; (ii) the adsorption branch for a pore open at both ends should be located at a pressure ($`\mathrm{P}_\mathrm{c}`$) higher than that for a pore closed at one end ($`\mathrm{P}_\mathrm{h}`$); (iii) for a given size distribution of cylindrical pores open at both ends, the adsorption branch should be steeper than that of the desorption branch (provided that the diameter exceeds 3-4 nm) and the shape of the hysteresis loop becomes rather symmetrical as the pore diameter increases, as shown by numerical calculations based on Eqs. (1) and (2) in which we have introduced a currently used t($`\mathrm{P}`$) proposed by Halsey  or the t($`\mathrm{P}`$) reproducing our experimental values prior to the capillary condensation.

To extract the latter t($`\mathrm{P}`$) values, we have prepared a porous Si layer with a PSD having a mean value as high as possible (about 50 nm) in order to shift the capillary condensation pressures to the highest possible values. The obtained isotherm exhibits a hysteresis loop located above $`0.8\mathrm{P}_0`$. We have thus determined the numerical values of t($`\mathrm{P}`$) up to $`0.8\mathrm{P}_0`$ by distributing homogeneously at each pressure the adsorbed amount among the pores, the PSD of which was extracted as explained above. In order to draw the t($`\mathrm{P}`$) above $`0.8\mathrm{P}_0`$, we have chosen to extrapolate the polynomial law which fits the experimental t($`\mathrm{P}`$).

The isotherm calculated by introducing in Eqs. (1) and (2) the t($`\mathrm{P}`$) and the PSD shown in inset Fig. 1(a) is represented Fig. 3(b).

We will successively discuss the predictions (i) and (iii) which are in large discrepancy with the experimental results (Fig. 3(a)); the prediction (ii) which is in agreement with experiment and finally, the large shift between the absolute position of the experimental isotherms and the calculated ones.

(i) The adsorption-desorption cycle observed in the case of pores closed at one end is irreversible. The interconnectivity related pore blocking effect described by Mason  to explain the origin of the hysteresis phenomenon could not occur because the porous Si layers are composed of non interconnected pores. Moreover, we reject the ‘ink bottle’ effect, previously introduced by Everett , as a possible explanation of the hysteresis loop we observe on the adsorption isotherms because of the absence of narrow sections inside the pores.

(iii) The hysteresis loops observed in the porous layer and porous membrane exhibit the same shape of type H2. We have observed, by carrying out adsorption isotherms on porous Si layers with a PSD having mean values varying from 13 nm up to 50 nm, that the position on the pressure axis of the corresponding adsorption and desorption branches are shifted towards $`\mathrm{P}_0`$ whereas the shape of the hysteresis loops remains of type H2.

Such an asymmetrical shape is not expected by the single pore models. Indeed, numerical calculations based on Cohan’s model or on DFT  predict a symmetrical shape of the hysteresis loop for such large pores in contradiction with what we observe. Furthermore, although porous Si exhibits non-interconnected pores, the hysteresis shape (type H2) we observe is the same as those systematically obtained for highly interconnected mesoporous solids such as controlled pore glass or Vycor . From a theoretical point of view, both macroscopic  and microscopic  approaches described this shape as a signature of the disorder introduced by the interconnectivity between pores. We must conclude that a hysteresis loop of type H2 is not characteristic of interconnected porous materials and we thus reject a previous conclusion made by Ball and Evans according to which a de Boer’s shape of hysteresis is a consequence of interconnected network effect.

We believe that the desorption process we observe results from another effect. At the end of the adsorption process, when all the pores are filled by the dense phase, the external surface of the porous layer is covered by a film which connects the pores. A possible interpretation of the observed delay to the emptying of pores could be the absence of formation of a concave meniscus on top of each pore. Thus, question arises whether metastable states exist which would correspond to a flat interface between the film and gas. In such an idea, during the desorption process, the film thickness could decrease homogeneously until the system jumps to the stable state corresponding to a concave meniscus on top of each pore. The steepness of the desorption branch could be explained if the jump is sufficiently delayed until a pressure at which the emptying of the major part of pores occurs. This idea needs of course to be checked by numerical calculations.

Note that the analysis of the desorption branch using BJH method yields a PSD with a width 0.5 nm about which is very far from that of the actual PSD ($`\pm `$ 6nm). It is conventionally admitted  that the desorption process in porous materials occurs at coexistence of the high and low dense phases. On such an assumption, a steep desorption branch is interpretated as the signature of the presence of pores of unique size. Our results obviously show that it is not necessarily the case.

(ii) The sign of the relative positions of the two adsorption branches ($`\mathrm{P}_\mathrm{h}<\mathrm{P}_\mathrm{c}`$) is in agreement with Cohan’s model. These experiments reveal two different scenari for the filling of pores open at one or at both ends. Preferential adsorption in bottom of closed pores is confirmed by analytical calculations performed in an oblique corner formed by the intersection of two planar surfaces . We are investigating pore closure effect on adsorption in cylindrical geometry by means of GCMC simulations and DFT calculations.

As shown in Fig. 3(b), the curve which corresponds to the modified Kelvin equation (2) is largely outside of the experimental hysteresis loops. For pore sizes as large as in our samples, this equation should account for the equilibrium phase transition. Indeed, microscopic calculations based on DFT performed by Ball and Evans  have shown that, in the case of a distribution of cylindrical pores open at both ends centered on 7.2nm, which is twice smaller than the mean size of our pores, the coexistence curve is located at a pressure which is close to that predicted by equation (2). To our knowledge, it is the first time that, for so large pore sizes, such a large discrepancy is found . This may be due to the use of a cylindrical geometry in the calculations which do not account for the actual polygonal shape of pores. GCMC simulations and calculations based on the DFT are in progress.

In summary, we have studied adsorption phenomenon in porous Si medium, a well ordered material with tubular pores without narrow sections and for which we have experimentally proven the absence of interconnection. Hysteresis loops of type H2 are observed, and that, whether the pores are open at one or at both ends. These results are in contradiction with Cohan’s model and show that the only observation of hysteresis loop of type H2 is not the signature of the presence of interconnections in porous material as it is generally admitted. The steep desorption process cannot be described by the macroscopic Kelvin equation in which we introduce the experimental values of the thickness of adsorbed amount prior to the capillary condensation. We note that the application of the widely used BJH method to extract the PSD from this branch is not systematically justified.


format: text

text:
# JLAB-THY-01-38 Excited nucleon spectrum using a non-perturbatively improved clover fermion action

## 1 INTRODUCTION

The $`N^{}`$ spectrum exhibits many features that are emblematic of QCD; that they contain three valence quarks is a manifestation of SU(3) symmetry, while the linearly rising Regge trajectories are evidence for a flux-tube-like confining force. Finally, the fine and hyperfine structure reveals details of the interquark interactions. Thus the study of the spectrum of excited nucleon resonances has always been a vital component of the experimental hadronic physics program, as typified by the Hall B experimental effort at Jefferson Laboratory.

The motivation for studying the $`N^{}`$ spectrum in many ways mirrors that for the hybrids, which also provide an arena in which to explore flux tubes. One distinction is that, in the case of nucleons, all spin-parity combinations are accessible to the quark model, and thus there are no nucleon “spin exotics”. However, both hybrids and excited nucleons share the feature that they are sensitive to the presence of excited glue.

The observed excited nucleon spectrum already poses many questions. There are “missing resonances”, predicted by the quark model but not yet observed experimentally. Most strikingly, there is the so-called Roper resonance, and the lightest negative-parity $`\mathrm{\Lambda }`$ state, the $`\mathrm{\Lambda }(1405)`$ with anomalously light masses suggesting that they may be “molecular” states, rather than true three-quark resonances. Lattice gauge theory calculations have a vital role to play in resolving these questions.

In an attempt to address these issues, there have been several recent calculations of the lowest-lying excited nucleon masses, emphasising in particular the extent to which the parity partners are accessible to lattice calculation. The first calculation employed the highly-improved $`D_{234}`$ fermion action, whilst a second calculation employed domain-wall fermions. Both calculations exhibited a clear splitting between the masses of the $`N^{1/2+}(937)`$ and $`N^{1/2}(1535)`$ states, with an indication that the lowest-lying radial excitation of the nucleon has a mass considerably larger than that of the $`N^{}(1535)`$, casting doubt on the nature of the Roper resonance.

In this talk, we present a calculation of the lowest lying negative parity nucleon mass using an $`𝒪(a)`$-improved Sheikholeslami-Wohlert (SW), or clover, fermion action; preliminary results were presented in ref. , and a comprehensive analysis appeared in ref. . By choosing the coefficient of the improvement term appropriately, all $`𝒪(a)`$ discretisation uncertainties are removed, ensuring that the continuum limit is approached with a rate proportional to $`a^2`$. For a subset of our lattices, masses of the first radial excitation of the nucleon and of the parity partner of the $`\mathrm{\Delta }`$ are presented.

The calculation of the excited nucleon spectrum places particularly heavy demands on lattice spectroscopy. The excited nucleon states are expected to be large; the size of a state is expected to double with each increase in orbital angular momentum. Thus a lattice study of the excited nucleon spectrum requires large lattice volumes, with correspondingly large computational requirements. Furthermore, the states are relatively massive, requiring a fine lattice spacing, at least in the temporal direction. These requirements could be satisfied with much greater economy using the clover fermion action than using the domain-wall or overlap formulation. Thus it is important to establish that the negative parity states are indeed accessible to calculations using the clover action. Finally, by comparing the masses obtained using the clover action with a calculation, at a single quark mass, using the Wilson fermion action, we also gain insight into the nature of the interaction responsible for the splitting in the parity doublet.

The layout of the remainder of this talk is as follows. In the next section, I shall describe the construction of hadronic operators, and in particular the projection onto positive- and negative-parity states, and summarises other computational details. Section 3 contains the results for the masses of the lowest-lying positive- and negative-parity nucleon using the SW fermion action, and for a subset of the lattices the masses obtained for the Roper resonance. The talk concludes with a discussion, a comparison with other lattice and phenomenological calculations, and prospects for future studies.

## 2 SIMULATION DETAILS

### 2.1 Baryon operators

Historically, there has been relatively little study of the negative-parity baryon states, though there had been numerous computations of the nucleon radial excitations. It is therefore useful to review the rôle that parity plays in the construction of baryon operators, and how this rôle is manifest in the case of lattice calculations. We illustrate the discussion through consideration of the usual nucleon interpolation operators:

$`N_1^{1/2+}`$ $`=`$ $`ϵ_{ijk}(u_i^TC\gamma _5d_j)u_k,`$ (1)
$`N_2^{1/2+}`$ $`=`$ $`ϵ_{ijk}(u_i^TCd_j)\gamma _5u_k,`$ (2)
$`N_3^{1/2+}`$ $`=`$ $`ϵ_{ijk}(u_i^TC\gamma _4\gamma _5d_j)u_k.`$ (3)
These operators have an overlap with particle states of both positive and negative parities. In order to construct operators of definite parity, we introduce the parity-projection operator:

$$P=(1\pm \gamma _4),$$
(4)
and construct the correlators

$$C_\pm (t)=\underset{\stackrel{}{x}}{}0|TN(\stackrel{}{x},t)(1\pm \gamma _4)\overline{N}(0)|0>,$$
(5)
where $`N`$ is a baryon interpolating operator. On a periodic or anti-periodic lattice, the two time orderings result in forward propagating states of positive (negative) parity, and backward propagating states which are anti-particles of negative (positive) parity states if the projection operator is chosen with positive (negative) sign. Thus the best delineation that may be constructed is between a forward propagating particle, and the backward propagating anti-particle of the parity partner. At large distances, where $`t1`$ and $`N_tt1`$, the correlators assume the time behaviour

$`C_+(t)`$ $``$ $`A^+e^{m^+t}+A^{}e^{m^{}(N_tt)}`$ (6)
$`C_{}(t)`$ $``$ $`A^{}e^{m^{}t}+A^+e^{m^+(N_tt)}`$ (7)
where $`m_i^+`$ and $`m_i^{}`$ are the lightest positive- and negative-parity masses respectively.

For the case of the nucleon ground state, the “diquark” part of both $`N_1`$ and $`N_3`$ couples upper (large) spinor components, while that in $`N_2`$ involves both an upper and lower spinor component and thus vanishes in the non-relativistic limit. Thus our expectation is that the operators $`N_1`$, $`N_3`$ should give a better overlap with the positive-parity ground state nucleon compared with operator $`N_2`$, and this is confirmed in lattice simulations. For some of the lattices used in the calculation, the positive-parity nucleon mass is obtained using the non-relativistic quark operators , defined by

$$\psi \psi ^{\mathrm{NR}}=\frac{1}{2}(1+\gamma _4)\psi ,\overline{\psi }^{\mathrm{NR}}=\overline{\psi }\frac{1}{2}(1+\gamma _4).$$
(8)
corresponding to the “large” components in the non-relativistic limit; this results in a factor of two reduction in computational effort at the expense of some loss of statistics.

### 2.2 Simulation Details

The calculation is performed in the quenched approximation to QCD. The standard Wilson gluon action is employed, whilst the quark propagators are computed using the Sheikholeslami-Wohlert clover fermion action, using the non-perturbative coefficient of the clover term determined in refs.  and . Thus the hadron masses are free of all $`𝒪(a^2)`$ discretisation errors.

The propagators are computed using both local sources, and spatially extended “smeared” sources, using either the fuzzed-source method (F) or Jacobi-smeared sources (J) as described in references  and  respectively. The simulation parameters are summarized in Tables 1. The errors on the fitted masses are computed using a bootstrap procedure. Different numbers of configurations are used at different quark masses at some of our $`\beta `$ values, precluding the use of correlated fits in some of the chiral extrapolations. Thus a simple uncorrelated $`\chi ^2`$ fit is performed, with the uncertainties computed from the variation in the $`\chi ^2`$. Calculations of the light hadron spectrum using these lattices were presented in references  and .

The clover term removes the leading chiral-symmetry-breaking effects from hadron masses at finite $`a`$. Since the lack of degeneracy between the positive- and negative-parity baryon states is a consequence of the spontaneous breaking of chiral symmetry, we also compare the measurements of the mass splitting at a single quark mass with those obtained using the Wilson fermion action.

## 3 RESULTS

The masses of the lowest-lying $`N^{1/2+}`$ and $`N^{1/2}`$ states are obtained from a simultaneous, four-parameter fit to the positive- and negative-parity correlators $`C_+(t)`$ and $`C_{}(t)`$ constructed from the “good” nucleon operator $`N_1`$ of eqn. 1. A clear signal for the mass of the lightest negative-parity state notably higher than that for the positive-parity state is observed, and is illustrated for one of our ensembles in Figure 1.

For the chiral extrapolation of the hadron masses, we adopt the ansatz

$`(am_X)^2=(aM_X)^2+b_2(am_\pi )^2,`$ (9)
where we use upper-case letters to denote masses obtained in the chiral limit and $`X`$ is either $`N^{1/2+}`$ or $`N^{1/2}`$. We include data at different volumes, but at the same $`\beta `$, in the chiral extrapolations, but treat the fuzzed and Jacobi-smeared data independently; they give slightly different masses, reflecting the different contributions of the excited states in the two cases. Note that the fit to $`(am_X)^2`$ is formally the same as one to $`am_X`$ at small pion mass, but in practice provides a notably better description of the data.

It is known that there are non-analytic terms contributing to the chiral extrapolation of physical quantities, and the contribution of the non-analytic terms arising from the pion-induced baryon self energies has been explored in reference . In the quenched approximation, the leading, non-analytic terms can assume a different form. From quenched chiral perturbation theory, one expects leading non-analytic terms which are linear in $`m_\pi `$ ; the coefficient of the leading term is predicted to be negative. We therefore also attempted to fit the data at $`\beta =6.4`$ to the form

$$(am_X)^2=(aM_X)^2+b_1(am_\pi )+b_2(am_\pi )^2.$$
(10)
Fitting the $`N^{1/2+}`$ state we found a positive value of $`b_1`$, though with a very large error that would still accommodate a negative value. Thus it is unclear whether eqn. 10 provides a reliable form with which to extrapolate data obtained with quark masses around that of the strange to the chiral limit, emphasising once again the need for data at the smaller quark masses at which the pion cloud emerges. The chirally extrapolated masses obtained using this procedure, with $`b_1`$ unconstrained, differ from those obtained using eqn. 9 by around 5%. In view of the difficulties discussed above, and because on some of our ensembles we have insufficient data points to perform a three-parameter fit, we quote our final results from fits to eqn. 9. An example of such a chiral extrapolation for the data at $`\beta =6.4`$ is shown in Figure 2.

Ideally, the forms eqn. 9 and 10 should be applied to the baryons masses obtained in the infinite volume limit. Indeed, the talk by Stuart Wright at this workshop emphasised the importance of correctly accounting for finite-volume effects in chiral extrapolations . Our data do not admit this procedure, but an investigation of the masses of the negative-parity state obtained on the smaller and larger of our lattices at $`\beta =6.0`$ and $`6.2`$ suggests that the masses of the negative-parity state on the larger lattices are smaller than those on the smaller lattices by around 5%.

In order to look at the discretisation uncertainties in our data, we show in Figure 3 the masses in units of $`r_0`$ against the $`a^2/r_0^2`$, where $`r_0=0.5\text{fm}`$ is the hadronic scale .

On a subset of our lattices, we have extracted the ground-state positive- and negative-parity masses using the operator $`N_2`$ of eqn. 2. The effective mass in the positive-parity channel using this operator is shown in Figure 1; it indicates an effective ground state mass considerably higher than that of the lightest negative-parity state. The signal for both the positive- and negative-parity ground states degrades appreciably with decreasing pion mass, as can be seen in Figure 4 where we show the extracted lowest-lying masses of both parities using each of the $`N_1`$ and $`N_2`$ interpolating operators. A better determination of the first radial excitation would require the measurement of the cross-correlator between $`N_1`$ and $`N_2`$, as done in reference . However, the conclusion that the operator $`N_2`$ has a negligible overlap with the ground state nucleon at these values of the pseudoscalar mass seems inescapable.

The measurement of heavier excitations is inevitably noisier because of the reduced signal-to-noise ratio. On some of our lattices we have measured the correlators for the $`I=3/2`$ $`\mathrm{\Delta }`$, using an interpolating operator

$$\mathrm{\Delta }^{3/2,1/2}=ϵ_{ijk}(u_i^TC\gamma _\mu u_j)u_k.$$
(11)
This has an overlap onto both spin-$`3/2`$ and spin-$`1/2`$ states, but these can be distinguished using a suitable spin projection, and the two parities delineated as described above. The masses of the ground state of both positive and negative parity are shown in Figure 5. Further success at determining the higher spin states has been made in a recent determination of the $`N^{3/2}`$ masses .

## 4 CONCLUSIONS

We have seen that the use of the clover-improved fermion action enables the delineation of the positive- and negative-parity excited state masses, in accord with calculations using highly improved Wilson, and domain-wall fermion actions. On a subset of our lattices, we also computed propagators using the unimproved Wilson fermion action at a single quark mass, and found a parity mass splitting consistent with that of the clover fermion action. On reflection, this is hardly surprising; under $`SU(6)`$ spin-flavour symmetry, the low-lying positive-parity baryons can be assigned to an $`l=0`$ multiplet, whilst the low-lying negative-parity baryons can be assigned to an $`l=1`$ multiplet. Thus the mass-splitting can be likened to the $`PS`$ splitting in the meson sector, which is well-known to be faithfully reproduced using the Wilson fermion action.

The striking feature common to all the lattice calculations is the apparent inversion of the ordering of mass of the lowest radial excitation of the nucleon, and the mass of the parity partner of the nucleon. All the calculations, using different gluonic and fermionic actions, suggest a splitting between the masses of the first positive-parity radial excitation and the lightest negative-parity state comparable to that between the masses of the negative-parity state and the nucleon. In contrast, the observed lightest radial excitation of the nucleon is the Roper resonance with a mass of $`1440\mathrm{M}\mathrm{e}\mathrm{V}`$, whilst the lightest negative-parity state is the heavier $`N^{}(1535)`$. We can see that lattice calculations of the spectrum are already posing interesting questions about the nature of the Roper resonance.

I have not discussed the issue of the nature of the $`\mathrm{\Lambda }(1405)`$ in this talk. In the case of the positive-parity states, lattice calculations, such as that of ref. , suggest that the quark mass behaviour of the baryon masses, in the range over which they are measured, is well described by terms corresponding to the “centres of gravity” of the quark masses within the hadron. This expectation has been confirmed in ref. , where the lightest negative-parity $`\mathrm{\Lambda }`$ resonance has a mass considerably in excess of $`1405\text{MeV}`$. It should be noted, however, that a recent large $`N_C`$ calculation has been able to accommodate a light $`\mathrm{\Lambda }(1405)`$ .

We have seen that lattice calculations are already addressing some of the salient issues in baryon spectroscopy. The calculation of the radial excitations and the higher spin states will present further challenges; the need to overcome the reduced signal-to-noise ratio afflicting heavier excitations, the construction of suitable lattice operators to isolate higher excitations, and the exploitation of improved data analysis methods. As these issues are tackled, we will be able to tackle the more ambitious projects of transition form factors and hadronic decays.

## ACKNOWLEDGEMENTS

This work was supported in part by DOE contract DE-AC05-84ER40150 under which the Southeastern Universities Research Association (SURA) operates the Thomas Jefferson National Accelerator Facility, by EPSRC grant GR/K41663, and PPARC grants GR/L29927 and GR/L56336. CMM acknowledges PPARC grant PPA/P/S/1998/00255. MG acknowledges financial support from the DFG (Schwerpunkt “Elektromagnetische Sonden”). Propagators were computed using the T3D at Edinburgh, the T3E at ZIB (Berlin) and NIC(Jülich), the APE100 at NIC (Zeuthen), and the Calico cluster at Jefferson Laboratory.

We are grateful for fruitful discussion with S. Dytman, R. Edwards, J. Goity, N. Isgur, R. Lebed, F. Lee, C. Michael, C. Morningstar, C. Schat, A. Thomas and S. Wallace.


text:
# References

LU TP 01-37

December 2001

Comment on the Pion Pole Part of the Light-by-Light Contribution to the Muon $`g2`$

Johan Bijnens
Department of Theoretical Physics 2, University of Lund

Sölvegatan 14A, S 223 62 Lund, Sweden

Elisabetta Pallante
S.I.S.S.A., Via Beirut 2-4, 34014 Trieste, Italy

Joaquim Prades
Centro Andaluz de Física de las Partículas Elementales (CAFPE) and

Departamento de Física Teórica y del Cosmos, Universidad de Granada,

Campus de Fuente Nueva, E-18002 Granada, Spain

## Abstract

We comment on the recent calculations of the pion pole part of the light-by-light contribution to the muon anomalous magnetic moment and we point out where the analysis in our previous work was mistaken.

In the recent paper , the pion pole part of the light-by-light contribution to the muon anomalous magnetic moment was reevaluated and analytical expressions for a large class of $`\pi ^0\gamma ^{}\gamma ^{}`$ form factors fulfilling OPE and large–$`N_c`$ QCD constraints were obtained. In  the leading double logarithmic dependence in the case of point-like pion-photon-photon couplings was derived analytically. This result has been confirmed in . The absolute value of the contribution was found to be in good agreement with the evaluations done previously , but opposite in sign.

We have gone back to our old programs and notes in order to find the source of the discrepancy. The pion pole contribution, as well as all the others, was calculated twice independently inside our collaboration and afterwards the results were compared. Both analyses agreed and had the correct sign analytically, but unfortunately both analyses made the same sign mistake in the overall normalization factor in putting it into the integration program. Since both computations agreed and also agreed with previously known results we saw no reason to doubt our results for the pion pole. More studies of variations in the form factor  did also not indicate that there was a major problem with the size of the form-factor. The same sign problem also affects the axial-vector contribution, but not the others we discussed. The scalar exchange contribution in our approach is linked by chiral symmetry constraints to the quark-loop low-energy contribution and had the correct sign.

We have also checked that the analytical formulas for the amplitudes in  are all correct with the exception of an overall minus sign missing in Eq. (2.7), that was an errata in the manuscript and not the origin of the sign mistake in the pion pole contribution.

The often quoted argument that there was no matching with the short-distance quark-loop contribution was tested in our work, in particular in Section 7.1 of . There, it was found that the short-distance quark loop matched perfectly well numerically when added together with another contribution that our treatment of the low-energy part required. For us, there was thus no reason to doubt that we had a reasonable matching with the expected short-distance behaviour even with the pion-pole and the quark-loop contributions having a different sign.

One should also not forget that, while the pion and the other pseudoscalar pole contributions are the dominant ones, there are other contributions at the 20 to 30 percent level that also need to be estimated. In the works  they have been evaluated and were found to cancel to a large extent. We intend to reevaluate all contributions and come back to them in a later more extended note. Our number corrected for the sign mistakes becomes

$`a_\mu ^{\mathrm{LbL}}`$ $`=`$ $`(2.1\pm 0.3)10^{10}[\text{quark-loop}]+(0.68\pm 0.2)10^{10}[\text{scalar}]`$ (1)
$`+`$ $`(8.5\pm 1.3)10^{10}[\text{pseudoscalar}]+(0.25\pm 0.1)10^{10}[\text{axial-vector}]`$
$`+`$ $`(1.9\pm 1.3)10^{10}[\pi K\text{-loop}]`$
$`=`$ $`(8.3\pm 3.2)10^{10}.`$
We have also checked that using the pointlike coupling our cut-off dependence numerically agrees with the one calculated analytically in .

A similar note  has appeared by the authors of the other full evaluation .


text:
# Flat Histogram Method of Wang-Landau and N-Fold Way

## 1 Introduction

While Monte Carlo methods in statistical thermodynamics already find broad application , the standard approach using the Metropolis algorithm  has the important disadvantage that the entropy of the simulated model system is not an output of the calculation. Furthermore in systems with a rugged landscape of the (coarse-grained) free energy the convergence of the method in practice often is problematic - the system may get “stuck” in one valley of this free energy landscape for a long time.
In order to overcome these difficulties, many interesting approaches (e.g., umbrella sampling , multicanonical Monte Carlo , expanded ensemble methods ) have been proposed in the literature, and some of these techniques emphasize the idea of directly sampling the energy density of states (e.g., ). Such approaches are clearly very promising, and in particular the method of Wang and Landau  has the additional merit that it is straightforward to implement. However, it is a difficult matter to judge the efficiency of the various methods, and estimate in beforehand the computational effort that is necessary to reach a desired level of accuracy for a chosen linear dimension $`L`$ of the (lattice) model under study. In particular, a recent comparative study of Wang and Swendsen  cautioned against a too optimistic view on these matters, and suggested that the Wang-Landau method  for medium and large lattice sizes $`L`$ is considerably less efficient than various versions of transition matrix Monte Carlo methods , which are more complicated to implement however.
In the present paper, we reconsider the Wang-Landau algorithm  and modify it by combining it with the N-fold way algorithm of Bortz et al. . It has been known that this latter algorithm is far superior to the single spin flip Metropolis algorithm  at low temperatures. However, as will be demonstrated here, our new algorithm performs significantly better near the right edge of an energy interval, since for large $`L`$ it is necessary to split the total energy range in many subintervals (which are slightly overlapping, of course, so that they can be joined unambigously), which then can be sampled in parallel on the processors of a multi-processor machine. This improvement of accuracy on the right boundaries of energy intervals is crucial for obtaining a very good overall accuracy. Thus, we are able to show that this particular criticism does no longer apply with respect to this new algorithm and a performance of the algorithm comparable in efficiency to the transition matrix Monte Carlo methods is reached.
In sec. 2, the flat histogram method of Wang and Landau  is briefly reviewed, as well as the N-fold way algorithm , and our combination of these two methods is presented. Sec. 3 then describes the implementation and results for the two-dimensional Ising model, treated also in refs. , and a comparative analysis of errors in the density of states is presented. Sec. 4 then gives, as an example of an application to a nontrivial model of current interest  results for a three-dimensional Ising model in a $`L\times L\times D`$ thin film geometry ($`L=32`$, $`D=6`$), where competing boundary fields ($`h_1=h_2`$) act on the two free $`L\times L`$ surfaces. It is shown that the present method produces results compatible with those of the conventional single spin flip algorithm  and reaches a better accuracy with less effort in computing time. Sec. 5 finally summarizes some conclusions.

## 2 The flat histogram method of Wang-Landau

### 2.1 Single spin flip

Recently, F. Wang and Landau  proposed a Monte Carlo algorithm for classical statistical models which uses a random walk in energy space in order to obtain an accurate estimate for the density of states $`g(E)`$. This method is based upon the fact that a flat energy histogram $`H(E)`$ is produced if the probability for the transition to a state of energy $`E`$ is proportional to $`1/g(E)`$.
This observation is utilized in the following way. Initially, $`g(E)`$ is set equal to one for all energies. A spin is then chosen at random and flipped with probability $`\mathrm{min}(1,g(E)/g(E^{}))`$ whereby $`E^{}`$ is the energy of the system with the chosen spin being overturned. The density of states $`g(E)`$ is not constant during the random walk, but is updated according to $`g(E)g(E)f`$ after each spin flip trial whether the spin is flipped or not. A histogram H(E) records how often a state of energy $`E`$ is visited. In the beginning of the random walk the modification factor $`f`$ can be as large as $`e2.7182818`$.111If one chooses to sample $`S(E)=\mathrm{log}g(E)`$ the modification factor becomes an increment $`S(E)S(E)+\mathrm{log}f`$ with $`\mathrm{log}f\mathrm{log}e0.4342944`$.Each time the energy histogram satisfies a certain flatness criterion, $`f`$ is reduced according to $`ff^{\frac{1}{2}}`$ and $`H(E)`$ is set to zero for all energies. The histogram is considered as flat if

$$H(E)ϵH(E)$$
(1)
for all $`E`$ where $`ϵ`$ is usually between $`0.7`$ and $`0.95`$ and $`H(E)`$ is the average histogram. The simulation is ended if $`f`$ is close enough to one, i.e., smaller than some predetermined final modification factor $`f_{final}`$. To speed up simulations it is possible to perform several random walks on adjacent energy intervals on independent processors. Disadvantages of the single spin flip version of this algorithm are the small acceptance rates for energy intervals which contain the groundstate and low-lying excited states and the relatively large errors of $`g(E)`$ at right edges of energy intervals as reported in . Since this flat histogram method produces only a relative density of states $`g(E)`$, one has to normalize $`g(E)`$ in order to get the absolute density of states $`\widehat{g}(E)`$. This can be done by using known values of the density of states, for example, the groundstate degeneracy or other constraints on $`\widehat{g}(E)`$. In case of the two-dimensional nearest neighbor Ising model

$$=J\underset{i,j}{}s_is_j\text{with}s_i=\pm 1,$$
(2)
one has

$`\widehat{g}(2JN)`$ $`=`$ $`2,`$ (3)
$`\widehat{g}(E)`$ $`=`$ $`\widehat{g}(E),`$ (4)
$`{\displaystyle \underset{E}{}}\widehat{g}(E)`$ $`=`$ $`2^N,`$ (5)
where $`N`$ is the number of spins. If one uses the groundstate to normalize $`g(E)`$ one can calculate canonical averages from the density of states for all temperatures which become exact for $`T0`$. It is thus highly desirable to estimate the density of states for low-lying energies with sufficient accuracy. Since we deal with low acceptance rates of the single spin flip approach in these energy-ranges one can reduce the simulational effort enormously by using a rejection-free update scheme. In  J. S. Wang and Swendsen presented an efficient Monte Carlo method using the N-fold way to obtain an estimate for the transition matrix from which the density of states can be determined via an optimization procedure. They have tested their algorithm with the two-dimensional Ising model and compared it to various other methods from which the density of states can be calculated. In particular, they found their method to perform better than the single spin flip Wang-Landau flat histogram method. Adopting the idea of Wang and Landau in the context of transition matrix they improved the efficency of the original algorithm only for small system sizes up to $`L=8`$, but for larger system sizes they reported that one has the problem of sticking to a Gaussian distribution for the histogram. In the next subsection we will discuss how the algorithm of Wang and Landau can be successfully combined with the rejection-free N-fold way in order to enhance its performance.

### 2.2 N-fold way

In the N-fold way  a flip occurs at each step of the algorithm and one then calculates the life-time of the resulting state. The usual single spin flip dynamics is thus preserved and observables become life-time weighted averages over the generated states. We will now describe the method in the context of the algorithm of Wang and Landau. Since the density of states $`\widehat{g}(E)`$ can be very large especially for large system sizes we consider $`S(E)=\mathrm{log}g(E)`$ during simulations. In the beginning $`S(E)`$ is set to zero for all $`E`$. Initially the system may be in the state $`\sigma `$ with energy $`EI=[E_{min},E_{max}]`$, whereby $`I`$ denotes the energy-range for which one wants to estimate $`g(E)`$. One then partitions all spins into classes according to their energetic local environment, i.e., the energy difference $`\mathrm{\Delta }E_i`$ a spin flip will cause. For a two-dimensional nearest neighbor Ising model each spin belongs to one of only $`M=10`$ classes. The total probability $`P`$ of any spin of class $`i`$ being overturned is given by

$$P(\mathrm{\Delta }E_i)=n(\sigma ,\mathrm{\Delta }E_i)p(EE+\mathrm{\Delta }E_i),i=1,\mathrm{},M,$$
(6)
with $`n(\sigma ,\mathrm{\Delta }E_i)`$ being the number of spins of state $`\sigma `$ which belong to class $`i`$ and $`p(EE+\mathrm{\Delta }E_i)`$ is given by

$$p(EE+\mathrm{\Delta }E_i)=\{\begin{array}{cc}\hfill \mathrm{min}(1,g(E)/g(E+\mathrm{\Delta }E_i))\text{if}& E+\mathrm{\Delta }E_iI\hfill \\ \hfill 0\text{if}& E+\mathrm{\Delta }E_iI.\hfill \end{array}$$
(7)
In order to determine the class from which to flip a spin one calculates the numbers

$$Q_m=\underset{im}{}P(\mathrm{\Delta }E_i),m=1,\mathrm{},M\text{and}Q_0=0,$$
(8)
which are the integrated probabilities for a spin flip within the first $`m`$ classes. Hence a class is selected by generating a random number $`0<r<Q_M`$ and taking class $`m`$ if $`Q_{m1}<r<Q_m`$. The spin to be overturned is chosen from this class with equal probabilities. Due to the flip, the spin and its interacting neighbors will change classes and correspondingly the numbers $`n(\sigma ,\mathrm{\Delta }E_i)`$ will differ from their predecessors. Finally, one has to determine the average life-time $`\tau `$ of the resulting state, i.e., one has to find out how many times the move just made would be rejected on average in the usual update scheme. The probability that the first random number would produce a flip is $`\widehat{P}=Q_M/N`$. Therefore one has for the probability that exactly $`n`$ random numbers will result in a new configuration

$$\overline{P}_n=\widehat{P}(1\widehat{P})^{n1}.$$
(9)
Thus the average life-time becomes

$$\tau =\underset{n=1}{\overset{\mathrm{}}{}}n\overline{P}_n=\underset{n=1}{\overset{\mathrm{}}{}}n\widehat{P}(1\widehat{P})^{n1}=\frac{N}{Q_M}.$$
(10)
Now, the steps which are actually carried out by the N-fold way version of the Wang-Landau method are the following:

1. choose an initial configuration and set $`H(E)=0`$, $`S(E)=0`$ for all $`E`$ and $`f_0=\mathrm{log}e0.4342944`$ and also fix $`f_{final}`$.
2. determine (update) the probabilities $`p(EE+\mathrm{\Delta }E_i)`$ and the $`Q_m`$’s of the (initial) configuration using eqs. (6), (7) and (8).
3. determine average lifetime $`\tau `$ of (initial) configuration via eq. (10).
4. increment histogram, density of states and update $`f_i`$:

$`H(E)`$ $``$ $`H(E)+\tau `$
$`S(E)`$ $``$ $`S(E)+\mathrm{\Delta }S(E)`$
$`f_i`$ $``$ $`f_{i+1}`$
with:

$`\mathrm{\Delta }S(E)`$ $`=`$ $`\{\begin{array}{cc}\hfill f_i\tau \text{if}& f_i\tau \mathrm{log}e0.4342944\hfill \\ \hfill \mathrm{log}e\text{if}& f_i\tau >\mathrm{log}e\hfill \end{array}`$ (13)
$`f_{i+1}`$ $`=`$ $`\{\begin{array}{cc}\hfill f_i\text{if}& f_i\tau \mathrm{log}e\hfill \\ \hfill \mathrm{\Delta }S(E)/\tau \text{if}& f_i\tau >\mathrm{log}e,\hfill \end{array}`$ (16)
in case of $`f_i/f_{i+1}>2`$ we set $`f_{i+1}f_i/2`$.
5. after some fixed sweeps check $`H(E)`$ and refine $`f_j`$ according to $`f_{j+1}=f_j/2`$ if $`H(E)`$ is flat. Stop if histogram is flat for a $`f_if_{final}`$.
6. determine the $`Q_m`$’s (the $`p(EE+\mathrm{\Delta }E_i)`$’s are not updated here).
7. choose and flip spin as described above.
8. go to $`2.`$

Equations (13) and (16) ensure that the increment $`\mathrm{\Delta }S(E)`$ is kept below or equal to $`\mathrm{log}e`$. In the single spin flip case this constraint is always guaranteed by the choice of the initial modification factor. Larger upper limits for $`\mathrm{\Delta }S(E)`$ result in large statistical errors in $`\widehat{g}(E)`$ as was reported in . Furthermore they can lead to $`Q_M=0`$ in the context of the N-fold way which will definitely destroy the iteration procedure described above. We have observed that adjusting $`f_i`$ according to eqs. (13) and (16) only applies to the first stage of iteration. Once the histogram is flat for the first time the increment $`\mathrm{\Delta }S(E)`$ is simply $`f_i\tau `$. Controling $`f_i`$ in other ways may also work.

## 3 Simulations

### 3.1 Two-dimensional Ising model

First we have tested our algorithm with the two-dimensional Ising model, eq. (2), on $`L\times L`$ square lattices with $`L=8,16,32,50`$ over the entire energy-range $`E/JN[2,2]`$. For these system sizes, the density of states $`\widehat{g}(E)`$ can be calculated exactly  and we will denote it by $`\widehat{g}_{ex}(E)`$. As already mentioned the method only provides a relative density of states. To obtain the absolute density of states $`\widehat{g}(E)`$, we first utilize (4) by taking

$$g(E)(g(E)+g(E))/2.$$
(17)
In order to normalize $`g(E)`$ one can now use eq. (3) which will yield the best accuracy in the density of states for low and high-lying energy levels or one uses eq. (5) which will lead to best accuracy in $`\widehat{g}(E)`$ for energies for which the corresponding density of states has large values. We have to stress here that this is the simplest method to obtain an absolute density of states, but it suffices in order to demonstrate the improvements made by combining the Wang-Landau method with the N-fold way. For further improvement of accuracy concerning the density of states, one has to resort to optimization methods. This was done succesfully in , where a least-squares method was considered with eqs. (3),(4) and (5) as constraints.

To compare our data with the tests provided in  we consider the relative error per energy level for the density of states

$$\epsilon (E)=\left|\frac{\widehat{g}(E)}{\widehat{g}_{ex}(E)}1\right|$$
(18)
and its average

$$\overline{\epsilon }=\frac{1}{N2}\underset{E}{}\epsilon (E).$$
(19)
The factor $`1/(N2)`$ is due to the fact that we have $`N1`$ different energy levels whereby the groundstate is excluded from averaging because we have imposed the exact value for the corresponding density of states, see eq. (3). For each system size we have performed $`30`$ runs, starting always with all spins up. Monte Carlo time was measured in units of sweeps, whereby one sweep was taken to be N spin flips. The results are given in table 1 and figures 2 and 3. The used parameters $`ϵ`$ and $`f_{final}`$ are stated in the corresponding captions. We also measured the tunneling-time $`\tau _t`$, i.e., the average number of sweeps which is elapsing during a transition from the groundstate to a state of highest energy or vice versa. The tunneling time is averaged over the complete run, i.e., the various stages with different modification increments. For the two-dimensional Ising model we get

$$\tau _t0.79L^{2.58}$$
(20)
from a fit which is depicted in figure 1. For comparison, the transition matrix Monte Carlo method yielded

$$\tau _t0.4L^{2.8}.$$
(21)
We have to add here that the tunneling-time is difficult to define within the Wang-Landau method because it depends on the density of states which is not constant during the simulation. Thus $`\tau _t`$ may also depend on the interplay between the parameters $`ϵ`$ and $`f_{final}`$. Thus we have chosen the typical values $`ϵ=0.85`$ and $`f_{final}1.010^6`$ for all considered system sizes.

### 3.2 Errors at right edges

As already mentioned in section 2.1, one shortcoming of the single spin flip Wang-Landau method are the relatively large errors in the density of states at right edges of an energy interval over which one wants to determine the density of states . In the N-fold way version presented here such systematic errors do not occur. In order to demonstrate this we have calculated the density of states for the first $`25`$ levels of a $`L=32`$ two-dimensional Ising model using single spin flips as well as N-fold way updates, see figure 4. This test also shows that the simulational effort is enormously reduced in regions where the acceptance rate for a spin flip is low.
In the single spin flip Wang-Landau method systematic errors occur at right edges due to a boundary effect which sets in as soon as the energy interval does not cover all the possible energies the system may have. If the random walk then hits an energy level which is outside the allowed range this move is simply rejected and the old level is counted once more analogous to the case when no boundary is involved in a transition, i.e., no difference is made between the density of states at boundaries and those away from the boundary, whereas in our N-fold way algorithm the density of states at edges is treated correctly since it is forbidden by the definition of the flip rates (eq. (7)) to choose a spin whose flipping would result in an energy outside the allowed range.

## 4 Thin Ising films

The two-dimensional Ising model is an ideal testing ground for Monte Carlo algorithms since it can be solved exactly and therefore permits one to verify simulational results directly. To show that our algorithm works also well with other more complicated Ising systems, we have considered the following Hamiltonian

$$=J\underset{i,j}{}s_is_jH\underset{i\text{bulk}}{}s_ih_1\underset{i\text{surface}\mathrm{\hspace{0.17em}1}}{}s_ih_2\underset{i\text{surface}\mathrm{\hspace{0.17em}2}}{}s_i\text{with}s_i=\pm 1,$$
(22)
which is essentially a three-dimensional $`L\times L\times D`$ Ising model confined by two walls with fields $`h_1,h_2`$ acting on them and periodic boundary conditions in the $`L\times L`$ planes. A sketch of this geometry is depicted in figure 5. With $`H=0`$ and $`h_1=h_2`$ this model exhibits the critical behavior of an interface localization/delocalization transition. For finite $`D`$ one has a single transition at $`T_c(D)`$ which presumeably belongs to the two-dimensional Ising universality class. For $`T<T_c(D)`$ the interface is bound either to the left or the right wall and for $`T_c(D)<T<T_{cb}`$ it is fluctuating delocalized around the center of the film, where $`T_{cb}`$ is the critical bulk temperature. If $`T`$ satisfies $`T>T_{cb}`$ the film is disordered with exception of the response to the surface fields $`h_1,h_2`$ near the walls. It has become clear  that such systems are extremely difficult to simulate, due to the fact that one has to deal with a much more severe slowing down than in simpler models, which results from the presence of a very large length scale $`\xi _{}`$ parallel to the surfaces in the phase with the delocalized interface which is related to the bulk correlation length $`\xi _b`$ via $`\xi _{}\mathrm{exp}(D/4\xi _b)`$ . The corresponding correlation time is then described by $`\tau \xi _{}^z`$, with a dynamic exponent of $`z2`$ . One can argue  that the asymptotic region of the two-dimensional Ising critical behavior is very narrow and that the true nature of the transition can be observed only for $`L\xi _{}`$. Consequently one is forced to study rather large linear dimensions $`L`$, depending on the thickness $`D`$ of the film, since for $`D=12`$, for example, one has to satisfy $`L\mathrm{exp}(D/4\xi _b)=26`$ ($`\xi _b0.92`$ for $`J/k_bT=0.25`$ ). In the case considered here ($`D=6`$, $`L=32`$) we have $`\mathrm{exp}(D/4\xi _b)=5.1`$ and the maximum of the specific heat is still rather broad, as can be seen from figure 7. While for standard Ising models cluster-algorithms  have proven to be a remedy for critical slowing down their utility has to be tested in detail when it comes to thin Ising films with additional surface fields. In  such a model has been studied with a ghostspin Swendsen-Wang cluster algorithm  and a significantly reduced auto-correlation time was achieved only for relatively small surface fields. Therefore we find it promising to test the applicability of our new algorithm to the model in question. For this reason we have performed $`5`$ runs of a system with size $`32\times 32\times 6`$ and competing surface fields $`h_1/J=h_2/J=0.55`$. In order to calculate canonical averages of themodynamic observables we have estimated the density of states $`\widehat{g}(E)`$ over the energy interval $`\frac{E}{JN}[2.8\overline{3},0.2]`$ which consists of already $`185974`$ different energy levels despite the moderate size of the system. The whole interval was partitioned into $`100`$ adjacent intervals with small overlaps, which are needed to join the density of states afterwards. For the normalization of the density of states we solely used the twofold degeneracy of the groundstate. The average total number of sweeps for the complete energy range stated above was $`(1.14\pm 0.18)10^7`$, the average total CPU time amounted to $`242.4h`$ on a HP v-class, whereby a single energy interval needed $`(8.7\pm 0.2)10^3s`$ CPU time on average. The simulation was carried out in a parallel fashion such that the effective average running time was $`12h`$ for a complete run. The simulation was performed with $`ϵ=0.95`$ and $`f_{final}=(7.05\pm 0.03)10^{10}`$.
In particular we have calculated the internal energy

$$U(T)=E_T=\frac{\underset{E}{}E\widehat{g}(E)e^{\beta E}}{\underset{E}{}\widehat{g}(E)e^{\beta E}},$$
(23)
the specific heat

$$C(T)=\frac{E^2_TE_T^2}{T^2},$$
(24)
as well as the Gibbs free energy

$$F(T)=kT\mathrm{ln}\left(\underset{E}{}\widehat{g}(E)e^{\beta E}\right),$$
(25)
and the entropy

$$S(T)=\frac{U(T)F(T)}{T},$$
(26)
over a range of temperatures near $`T_C(D)`$ and checked it against conventional Monte Carlo data obtained from a heat bath algorithm . We also considered canonical averages of the order parameter $`\left|m\right|=\left|\frac{1}{N}s_i\right|`$:

$$\left|m\right|_T=\frac{\underset{E}{}\left|m\right|_E\widehat{g}(E)e^{\beta E}}{\underset{E}{}\widehat{g}(E)e^{\beta E}},$$
(27)
its fourth order cumulant

$$U_4(T)=1\frac{m^4_T}{3m^2_T^2},$$
(28)
and the susceptibility

$$\chi (T)=\frac{N}{T}\left(m^2_T\left|m\right|_T^2\right),$$
(29)
by sampling microcanonical averages $`_E`$ during the last stage of the simulations. The results are shown in figure 6 \- 12. The value of the extrapolated ($`L\mathrm{}`$) critical temperature $`J/k_BT_C(D=6)=0.2655\pm 0.0002`$ estimated in ref.  is marked by an arrow. One can see from these figures that quantities like the specific heat, the magnetization, the susceptibility and the fourth order cumulant can be obtained with an accuracy that is clearly better than that of the corresponding Monte Carlo data obtained with the standard single spin flip heat bath algorithm , although the effort in computer time was comparable. Of course, for a precise analysis of the critical behaviour of the present model it is necessary to repeat the study for several choices of $`L`$ and perform a finite size scaling analysis. This will be presented elsewhere . A further merit of the present approach is that it immediately yields precise data for the entropy and free energy of the model as well (figs. 8,9), which would be accessible from the standard algorithm  only by tedious techniques of “thermodynamic integration” . The availability of the free energy is particularly advantageous in the case of first-order interface localization-delocalization transitions , where the relaxation time increases exponentially with the linear dimension $`L`$ of the system, and straightforward Monte Carlo studies would be hampered by hysteresis for large $`L`$ (and huge statistical errors for intermediate values of $`L`$ ).

## 5 Conclusions

In this paper, it has been shown that the efficiency of the flat histogram method of Wang and Landau  can be improved significantly by combining it with the so-called “N-fold way”  technique of sampling states of a lattice model via Monte Carlo methods. In particular, the problem that errors get strongly enhanced near the right edge of an energy interval that is sampled is eliminated, and performing systematic comparisons for the two-dimensional Ising model with transition matrix Monte Carlo approaches proposed by Wang and Swendsen  it is shown that the present algorithm is competitive in efficiency with techniques described in the latter study , but the merit of the present method is that it is rather straightforward to implement, and it should be immediately useful for wide classes of lattice models. Unlike cluster algorithms, there is no problem with the inclusion of magnetic fields. As a nontrivial example, we show new results of a first application to a thin Ising film with competing boundary fields.

## Acknowledgement

One of us (B. J. S.) received financial support from the Deutsche Forschungsgemeinschaft (DFG) under grant No. Bi314/16-2. We are grateful to David P. Landau for numerous stimulating discussions, and for providing preprints of ref.  prior to publication. We are grateful to Jian-Sheng Wang and Robert H. Swendsen for sending a preprint of ref. . This work also received current support from the Bundesministerium für Bildung und Forschung (BMBF) under grant No. 03N6015.


text:
# SETHi@Berkeley– A Piggyback 21-cm Sky Survey at Arecibo

## Observing Methodology

The UCB SETI searches use the 1420 MHz line feed on Carriage House 1 at the National Astronomy and Ionospheric Center’s 305 meter radio telescope in Arecibo, Puerto Rico. This unique arrangement allows observations to be conducted without interference with other uses of the telescope. This results in two main modes of observation. If the primary observers feed is stationary or stowed the beam scans across the sky at the sidereal rate. If the primary observer’s feed is tracking a position on the sky, the beam scans the sky at twice the sidereal rate. At twice the sidereal rate, the beam width corresponds to a 12 second beam transit time. (Korpela et al. 2001) Figure 1 shows the path of the telescope beam. The upper figure show the path of the beam over the course of 15 hours. The lower figure shows the sky coverage for about 2 years of observations.

The time domain data for the sky survey is recorded as follows: first, a 30 MHz band from the receiver is converted to baseband using a pair of mixers and low pass filters. The resulting complex signal is digitized and then filtered to 2.5MHz using a pair of 192 tap FIR filters in the SERENDIP IV instrument. (Werthimer et al. 1997) One bit samples are recorded on 35 GByte DLT tapes (one bit real and one bit imaginary per complex sample). These tapes are shipped to Berkeley for use in the SETI@home program. Tsys for the system is typically 75K.

The SETHi@Berkeley program analyzes these tapes to extract hydrogen spectra. The 2.5 MHz time series data are converted to raw spectra using 2048 point FFTs ($`\mathrm{\Delta }\nu `$=1220 Hz). 6144 FFTs are accumulated into a single power spectrum of 5.033 second integration time. The resulting power spectrum is corrected for 1 bit sampling effects by using the Van Vleck correction. The spectrum, its start and and end coordinates, and its time are stored in a database for future use. This database will be queried to develop spectral maps of the neutral hydrogen distribution.

## Example Spectra

Figure 2 shows a spectrum taken along a line of sight near 3C192 (l,b)=(197.7, 24.4). The hydrogen column along this line of sight is $`4.2\times 10^{20}`$ cm-2 The upper plot shows a single raw 5 second integration spectrum with 1.22 kHz resolution. The lower plot represents the sum of 5 adjacent spectra, and shows the typical SNR that will be achieved in a 0.1 degree skymap pixel.

Figure 3 shows 27 spectra accumulated over 136 seconds. During the accumulation of these spectra the telescope beam was moving at more than twice the sidereal rate ($``$1.6 degrees over this duration). Significant changes are seen in the spectra are seen on scales approximating the beam width.

## Program Status

Data accumulation for SETI@home began in December 1998. We currently have accumulated about 58 Msec of observation time, which will translate into 11.5 million spectra. The survey has covered 79.4% of the accessible sky. If mapped into pixels 1 beam width in size, the median exposure per pixel is $``$20 seconds. We anticipate that data collection will continue for at least 1 additional year. Generation of the spectral database is just beginning. We anticipate analysis of the existing data to be complete in $``$12 months.

## Acknowledgements

SETI@home is primarily funded through private donations. We especially would like to thank The Planetary Society, Cosmos Studios, The SETI Institute, Sun Microsystems, The University of California, and The Friends of SETI@home for their generous support. We would also like to thank the staff NAIC Arecibo observatory. NAIC is operated by Cornell University under a cooperative agreement with the National Science Foundation.


text:
# Cosmic Ray composition around the knee from EAS electromagnetic and muon data

## 1 Introduction and the detectors

The study of the primary composition in the Extensive Air Shower energy region requires the use of different observables in order to cross check the information and reduce the dependence on the interaction model and propagation codes used. At the National Gran Sasso Laboratories a program of exploiting the surface shower size measurements from EAS-TOP (2005 m a.s.l.), and the high energy muon measurements ($`E_\mu ^{th}=1.3`$ TeV) performed in the deep underground laboratories (MACRO) has been developed. Such muons in fact originate from the decays of mesons produced in the first interactions in the atmosphere and from a quite different rapidity region than the GeV muons usually used for such analysis ($`x_F>0.1`$ or $`0.2`$).

The two experiments operated in coincidence for a live time of $`\mathrm{\Delta }T=`$ 960.12 days between 1990 and 2000. The number of coincident events collected with the two detectors operating in the final configuration amounts to 28160, of which 3752 have shower cores inside the edges of the array (”internal events”) and shower size $`N_e>2.10^5`$, and 409 have $`N_e>10^{5.92}`$, i.e. above the observed knee at the corresponding zenith angle. We present here an analysis of the full data set in terms of the QGSJET interaction model as implemented in CORSIKA. For a comparison of different interaction models see .

The EAS-TOP array is located at Campo Imperatore (2005 m a.s.l., $`30^o`$ with respect to the vertical of the underground Gran Sasso Laboratories, corresponding to 930 gr cm-2 atmospheric depth). Its e.m. detector (to which we are mainly interested in the present analysis) is built of 35 scintillator modules 10 m2 each, including an area A $``$ 105 m2. The array is fully efficient for $`N_e>10^5`$. Its reconstruction capabilities of the EAS parameters (for internal events) are: $`\frac{\mathrm{\Delta }N_e}{N_e}10\%`$ above $`N_e10^5`$, and $`\mathrm{\Delta }\theta 0.9^o`$ for the EAS arrival direction. The array and the reconstruction procedures are fully described in .

MACRO, in the underground Gran Sasso Laboratory at 963 m a.s.l., with 3100 m w.e. of minimum rock overburden, is a large area multi-purpose apparatus designed to detect penetrating cosmic radiation. The lower part of the MACRO detector has dimensions $`76.6\times 12\times 4.8`$ m3. A detailed description of the apparatus can be found in .

In this work we consider only muon tracks which have at least 4 aligned hits in both views of the horizontal streamer tube planes over the 10 layers composing the whole detector. The standard reconstruction procedure of MACRO  has been used.

The two experiments are separated by a thickness of rock ranging from 1100 up to 1300 m, depending on the angle.

The muon energy threshold at the surface for muons reaching the MACRO depth ranges from $`E_\mu ^{th}=1.3`$ TeV to $`E_\mu ^{th}=1.8`$ TeV inside the effective area of EAS-TOP. Event coincidence is established off-line, using the absolute time given by a GPS system with an accuracy better than 1 $`\mu s`$.

Independent analyses of the two arrays are in  and .

## 2 Analysis and results

The analysis technique has to be adapted to the specific trigger requirements (both surface and underground detectors fired) with defined acceptance area for the EAS array (internal events), but undefined for the underground one. Therefore, the main experimental feature to be considered is the muon multiplicity distribution in different intervals of shower sizes. We have chosen six intervals of shower sizes around the knee:
5.20 $`<Log_{10}(Ne)`$ 5.31; 5.31 $`<Log_{10}(Ne)`$ 5.61; 5.61 $`<Log_{10}(Ne)`$ 5.92; 5.92 $`<Log_{10}(Ne)`$ 6.15; 6.15 $`<Log_{10}(Ne)`$ 6.35 and 6.35 $`<Log_{10}(Ne)`$ 6.7.

Within each size bin the muon multiplicity distribution has been fitted with a superposition of pure $`p`$ and $`Fe`$ components, or light $`(L)`$ and heavy $`(H)`$ admixtures containing equal fractions of $`p`$ and $`He`$ or $`Mg`$ and $`Fe`$ respectively. All spectra in the simulation have slope $`\gamma =2.62`$.

The fit has been performed in the quoted six windows by minimizing the following expression for each multiplicity distribution:

$$\chi ^2=\underset{i}{}\frac{(N_i^{exp}p_1N_i^pp_2N_i^{Fe})^2}{\sigma _{i,exp}^2}$$
(1)
where $`N_i^{exp}`$ in the number of observed events in the $`i`$-th bin of multiplicity, $`N^p`$ ($`N^L`$) and $`N^{Fe}`$ ($`N^H`$) are the number of simulated events in the same $`i`$-th multiplicity bin from the $`p`$ ($`L`$) and $`Fe`$ ($`H`$) components, respectively; $`p_1`$ and $`p_2`$ are the parameters (to be fitted) defining the fraction of each mass component contributing to the same multiplicity bin. The results of the fits have been normalized to reproduce the observed number of coincident events in each size bin.

To infer the mass composition evolution as a function of energy, for each size bin we take from the simulation the $`log(E)`$ distributions of contributing mass groups weighted by the parameters $`p_1`$ and $`p_2`$ and with weights $`w_k`$ representing the relative efficiency to trigger the underground apparatus. The resulting distributions from different size bins are summed together, and so we eventually obtain the simulated energy spectra of the two basic components that reproduce the experimental data.

These resulting spectra for $`p`$-$`Fe`$ and $`Light`$-$`Heavy`$ components, further corrected for the detection efficiency of the surface array, are reported in Fig. 1 together with their sum. Within such procedure, for each bin in $`log(E)`$ we calculate $`<logA>`$ from the expression:

$$<logA>=\frac{p_1M_k^pLog(A^p)+p_2w_kM_k^{Fe}Log(A^{Fe})}{M_k^p+p_2w_kM_k^{Fe}}$$
(2)
where $`M_k^p`$ and $`M_k^{Fe}`$ are the numbers of simulated events from $`p`$ and $`Fe`$ contributing to the k-th energy bin. The same expression, with obvious modifications, is used for the $`Light`$-$`Heavy`$ composition. The results are depicted in Figs. 2, the two curves delimiting the error band. The results show an increase in $`<logA>`$ through the energy corresponding to the knee position.

## 3 Conclusions

The analysis of the Ne-Nμ ($`E_\mu >1.3`$ TeV) data collected by the MACRO/EAS-TOP collaboration at the Gran Sasso Laboratories points to a primary composition becoming heavier around the knee of the primary spectrum (in the energy region $`10^{15}10^{16}`$ eV). The result is in very good agreement with the measurements of EAS-TOP alone at the surface, wich detects muons with energy $`E_\mu >1`$ GeV, using the same (QGSJET) interaction model . GeV and TeV muons are produced in different kinematic regions: in the central one and at the edges of the fragmentation region respectively. The present measurements show therefore a good consistency of the interaction model in describing the yield of secondaries over a wide rapidity region.


text:
# Probing New Scales at a 𝑒⁺𝑒⁻ Linear Collider

## I Introduction

The Lhc is expected to directly probe possible New Physics beyond the Standard Model (SM) up to a scale of a few TeV. While its data should provide answers to several of the major open questions in the present picture of elementary particle physics, it is important to start examining how this sensitivity can be further extended at a next generation of colliders. Today we have a number of indications that New Physics could be of supersymmetric nature. If this is the case, the Lhc will have a variety of signals to discover these new particles and the linear collider (LC) will be required to complement the probe of the SUSY spectrum with detailed measurements. However, beyond Supersymmetry there is a wide range of other scenarios invoking new phenomena at, and beyond, the TeV scale. They are aimed at explaininig the origin of electro-weak symmetry breaking, if there is no light elementary Higgs boson, at stabilising the SM, if SUSY is not realised in nature, or at embedding the SM in a theory of grand unification.

Many of such scenarios predict the existence of new particles that would be manifested as rather spectacular resonances in $`e^+e^{}`$ collisions, if the achievable centre-of-mass energy is sufficient. A high energy LC represents an ideal laboratory for studying this New Physics zp1 ; dominici . It also retains an indirect sensitivity, through a precision study of the virtual corrections to electro-weak observables, when their mass exceeds the available centre-of-mass energy. This paper summarises the results of a series of studies aimed at quantifying the potential of a high energy, high luminosity $`e^+e^{}`$ LC in extending to high scales the probe for New Physics. While a significant activity has already addressed the issues related to a TeV-class collider, we now review the potential of a multi-TeV LC, such as Clic.

## II Electro-weak Observables at a Multi-TeV LC

The analysis of the Lep and Slc data has provided a significant experience in the extraction of electro-weak observables, optimising their statistical sensitivity and controlling their systematic uncertainties. At larger centre-of-mass energies, the relevant $`e^+e^{}f\overline{f}`$ cross sections are significantly reduced and the experimental conditions at the interaction region need to be taken into account in validating the anticipated accuracies on the cross section $`\sigma _{f\overline{f}}`$, forward-backward asymmetries $`A_{FB}^{f\overline{f}}`$ and left-right asymmetries $`A_{LR}^{f\overline{f}}`$ determination at $`\sqrt{s}`$ = 1 TeV - 5 TeV. Since the two-fermion cross section is of the order of only 10 fb, it is imperative to achieve high luminosity by reducing the beam-spot sizes. In this regime the beam-beam effects are important and the primary $`e^+e^{}`$ collision is accompained by several $`\gamma \gamma \mathrm{hadrons}`$ interactions. Being mostly confined in the forward regions, this $`\gamma \gamma `$ background reduces the polar angle acceptance for quark flavour tagging and dilutes the jet charge separation using jet charge techniques. These experimental conditions require efficient and robust algorithms to ensure sensitivity to flavour-specific $`f\overline{f}`$ production. The statistical accuracy for the determination of $`\sigma _{f\overline{f}}`$, $`A_{FB}^{f\overline{f}}`$ and $`A_{LR}^{f\overline{f}}`$ has been studied, for $`\mu ^+\mu ^{}`$ and $`b\overline{b}`$, taking the Clic parameters at $`\sqrt{s}`$ = 3 TeV. The Simdet parametrised detector simulation has been used and the $`\gamma \gamma \mathrm{hadrons}`$ background, corresponding to 10 overlayed bunch crossings, has been added to $`e^+e^{}\mu ^+\mu ^{}`$, $`b\overline{b}`$ events. $`b\overline{b}`$ final states have been identified using an algorithm based on the sampling of the decay charged multiplicity of the highly boosted $`b`$ hadrons at Clic energies bmult . Similarly to Lep analyses, the forward-backward asymmetry has been extracted from a fit to the flow of the jet charge $`Q^{jet}`$ defined as $`Q^{jet}=\frac{_iq_i|p_iT|^k}{_i|p_iT|^k}`$, where $`q_i`$ is the particle charge, $`p_i`$ its momentum, $`T`$ the jet thrust axis and the sum is extended to all the particles in a given jet. Here the presence of additional particles, from the $`\gamma \gamma `$ background, causes a broadening of the $`Q^{jet}`$ distribution and thus a dilution of the quark charge separation. The track selection and the value of the power parameter $`k`$ needed to be optimised as a function of the number of overlayed bunch crossings. The results are summarised in terms of the relative statistical accuracies $`\delta 𝒪/𝒪`$ in Table 1. Another important issue is the accuracy on the luminosity determination, that needs to be controlled to 0.5%, or better.

## III Sensitivity to a $`Z^{}`$ Boson

There is a wide range of New Physics scenarios predicting the existence of new vector particles with masses in the TeV range. One of the simplest extensions of the SM consists in the introduction of an additional $`U(1)`$ gauge symmetry, whose breaking scale is close to the Fermi scale. This extra symmetry is for instance predicted in some grand unified theories. The extra $`Z^{}`$ associated to this symmetry naturally mixes with the SM $`Z^0`$. The mixing angle is strongly constrained by precision electroweak data to be of the order of few mrad while direct searches at Tevatron for a new $`Z^{}`$ boson set a lower mass limit around 600 GeV. Further improvements are expected by precision electro-weak data obtained at a Giga-Z facility.

The search for such an extended gauge sector offers an interesting framework for studying the sensitivity to scales beyond those directly accessible, even at a multi-TeV collider. It also raises the issue of the ability to discriminate between different models. Main classes of models with additional $`Z^{}`$ bosons include $`E_6`$ inspired models and the left-right models (LR) (for a recent review see leike ). Concerning the $`E_6`$ models, the $`Z^{}`$ fermion couplings depend on the angle $`\theta _6`$ which defines the embedding of the extra $`U(1)`$ in the $`E_6`$ group.

In this study we have considered the so-called $`\eta `$ model with $`\theta _6=\mathrm{arctan}\sqrt{5/3}`$, the $`\chi `$ model with $`\theta _6=0`$ and, as a reference model, the so called sequential SM (SSM) which has an additional $`Z^{}`$ boson with SM-like couplings.

There has been a significant interest in the Lhc and LC potential in the search for a new $`Z^{}`$ boson. At the LC, the indirect sensitivity to its mass, $`M_Z^{}`$, can be parametrised in terms of the available integrated luminosity $``$, and centre-of-mass energy, $`\sqrt{s}`$. In fact the scaling law for large $`M_Z^{}`$ can be obtained by considering the effect of the $`Z^{}\gamma `$ interference in the cross section $`\sigma `$. For $`s<<M_Z^{}^2`$ and assuming that the uncertainties $`\delta \sigma `$ are statistically dominated, we get the range of mass values giving a significant difference:

$$\frac{|\sigma ^{SM}\sigma ^{SM+Z^{}}|}{\delta \sigma }\frac{1}{M_Z^{}^2}\sqrt{sL}>\sqrt{\mathrm{\Delta }\chi ^2}$$
(1)
and the sensitivity to the $`Z^{}`$ mass scales as:

$$M_Z^{}(sL)^{1/4}$$
(2)
This relationship shows that there is a direct trade-off possible between $`\sqrt{s}`$ and $``$, which should be taken into account when optimising the parameters of a high energy LC.

The $`\sigma _{f\overline{f}}`$ and $`A_{FB}^{f\overline{f}}`$ ($`f=\mu ,b`$) values have been computed, for 1 TeV $`<\sqrt{s}<`$ 5 TeV, both in SM and including the corrections due to the presence of a $`Z^{}`$ boson with 10 TeV $`<M_Z^{}<`$ 40 TeV, with the couplings predicted by the models mentioned above. This has been obtained by implementing them in the Comphep program. The relative statistical errors on the electroweak observables are obtained by rescaling the values of Table I for different energies and luminosities. The sensitivity has been defined as the largest $`Z^{}`$ mass giving a deviation of the actual values of the observables from their SM predictions corresponding to a SM probability of less than 5%. This sensitivity has been determined, as a function of the $`\sqrt{s}`$ energy for an integrated luminosity $``$ of 1 ab-1 and 5 ab-1, and rescaled to other values of $``$ using the formula (2) and assuming the uncertainties to be dominated by statistics. Results are summarised in Figure 1. For the $`\eta `$ model the sensitivity is lower: for example to reach a sensitivity of $`M_Z^{}`$=20 TeV, 10 ab-1 of data at $`\sqrt{s}`$=5 TeV would be necessary.

## IV Sensitivity to Kaluza-Klein Excitations in Theories with Extra-Dimensions

Theories of quantum gravity have considered the existence of extra-dimensions for achieving the unification of gravity at a scale close to that of electroweak symmetry breaking. String theories have recently suggested that the SM could live on a $`3+\delta `$ brane with $`\delta `$ compactified large dimensions while gravity lives on the entire ten dimensional bulk. The corresponding models lead to new signatures for future colliders ranging from Kaluza-Klein (KK) excitations of the gravitons HLZ  to KK excitations of the SM gauge fields with masses in the TeV range antoniadis .

Among the models with extra dimensions we have considered the five dimensional extension of the SM with fermions on the boundary which predicts KK excitations of the SM gauge bosons with couplings $`\sqrt{2}`$ larger than those of the SM. KK masses are given by $`M_n=nM`$, with $`M`$ the compactification scale of the fifth dimension. Indirect limits on $`M4`$ TeV from electroweak measurements already exist ewlimits . These models predict excitations of the $`Z^0`$ and of the photon which are almost degenerate in mass.

In this analysis we have included in the cross section calculations only the exchange of the first KK excitations $`Z^{(1)}`$ and $`\gamma ^{(1)}`$, neglecting the effect of the remaining towers which give a small correction. The scaling law for the limit on $`M`$ can be obtained by considering the interference of the two new nearly degenerate gauge bosons with the photon in the cross section and taking the $`s<<M^2`$ limit. The result is the same as eq. (2). The analysis closely follows that for the $`Z^{}`$ boson discussed above. In Figure 1 we give the sensitivity contours as a function of $`\sqrt{s}`$ for different values of $`M`$. We conclude that the sensitivity achievable for the compactification scale $`M`$ for an integrated luminosity of 1 ab-1 in $`e^+e^{}`$ collisions at $`\sqrt{s}`$ = 3 TeV - 5 TeV is of the order of 40 TeV - 60 TeV. Results for a similar analysis, including all electro-weak observables are discussed in rizzosnow .

## V Sensitivity to Contact Interactions

The scenarios investigated above address specific models of New Physics beyond the SM.

Fermion compositeness or exchange of very heavy new particles can be described in all generality by four-fermion contact interactions ref:ci . These parametrise the interactions beyond the SM by means of an effective scale, $`\mathrm{\Lambda }_{ij}`$,

$$_{CI}=\underset{i,j=\mathrm{L},\mathrm{R}}{}\eta _{ij}\frac{g^2}{\mathrm{\Lambda }_{ij}^2}(\overline{\mathrm{e}}_i\gamma ^\mu \mathrm{e}_i)(\overline{f}_j\gamma ^\mu f_j).$$
(3)
The strength of the interaction is set by convention as $`g^2/4\pi =1`$ and models can be considered by choosing either $`|\eta _{ij}|=1`$ or $`|\eta _{ij}|=0`$ as detailed in Table II. The contact scale $`\mathrm{\Lambda }`$ can be interpreted as effect of new particles at a mass $`M_X`$, $`1/\mathrm{\Lambda }^2\lambda ^2/M_X^2`$.

In order to estimate the sensitivity of electro-weak observables to the contact interaction scale $`\mathrm{\Lambda }`$, the statistical accuracies discussed in Section II have been assumed for the $`\mu \mu `$ and $`b\overline{b}`$ final states. The systematics of the assumed 0.5% include the contributions from model prediction uncertainties. Results are given in terms of the lower limits on $`\mathrm{\Lambda }`$ which can be excluded at 95% C.L., in Figure 2. It has been verified that, for the channels considered in the present analysis, the bounds for the different $`\mathrm{\Lambda }_{ij}`$ are consistent. High luminosity $`e^+e^{}`$ collisions at 3 TeV can probe $`\mathrm{\Lambda }`$ at scales of 200 TeV, and beyond. For comparison, the corresponding results expected for a LC operating at 1 TeV are also shown. Beam polarisation represents an important tool in these studies. First, it improves the sensitivity to new interactions, through the introduction of the left-right asymmetries $`A_{LR}`$ and the polarised forward-backward asymmetries $`A_{FB}^{pol}`$ in the electro-weak fits. If both beams can be polarised to $`𝒫_{}`$ and $`𝒫_+`$ respectively, the relevant parameter is the effective polarisation defined as $`𝒫={\displaystyle \frac{𝒫_{}+𝒫_+}{1𝒫_{}+𝒫_+}}`$. In addition to the improved sensitivity, the uncertainty on the effective polarization, can be made smaller than the error on the individual beam polarization measurements. Secondly, in the case of a significant deviation from the SM prediction would be observed, $`e^{}`$ and $`e^+`$ polarization is greatly beneficial to determine the nature of the new interactions. This has been studied in details for a LC at 0.5–1.0 TeV ref:lc-sr  and those results also apply, qualitatively, to a multi-TeV collider.

## VI Conclusions

Extending the sensitivity to New Physics beyond the anticipated reach of the Lhc, is a prime aim of future colliders. By accurately measuring electro-weak observables, a LC able to achieve $`e^+e^{}`$ collisions at and beyond 1 TeV, with high luminosity, can indirectly probe scales extending from tens to several hundreds TeV.