text:
# Interplay of disorder and magnetic field in the superconducting vortex state

## I Introduction

The interaction between the superconductor quasiparticles and the vortices induced by an external magnetic field has been a subject of considerable recent interest Gorkov ; Anderson ; FT . In the presence of vortices the quasiparticles feel the combined effect of the external magnetic field and of the spatially varying field of the chiral supercurrents. By performing a gauge transformation to effectively reduce the system to the one in a zero average magnetic field it was shown FT  that the natural low energy quasiparticle modes are Bloch waves rather than the Dirac Landau levels proposed in Refs. Gorkov,  and Anderson, . The results or Ref. FT  also showed that the quasiparticles besides feeling a Doppler shift caused by the moving supercurrents Volovik  (Volovik effect) also feel a quantum “Berry” like term due to a half-flux Aharonov-Bohm scattering of the quasiparticles by the vortices.

The effect of disorder on the low-energy density of states of superconductors has also been a subject of much recent activity report , particularly in the case of $`d`$-wave symmetry. Several conflicting predictions have appeared in the literature which have mainly concentrated on the effect of the presence of impurities. Some progress toward understanding the disparity of theoretical results has been achieved realising that the details of the type of disorder affect significantly the density of states report . Particularly in the case of $`d`$-wave superconductors, in contrast to conventional gapped $`s`$-wave superconductors, the presence of gapless nodes is expected to affect the transport properties. Using a field theoretic description and linearizing the spectrum around the four Dirac-like nodes it has been suggested that the system is critical. It was obtained that the density of states is of the type $`\rho (ϵ)|ϵ|^\alpha `$, where $`\alpha `$ is a non-universal exponent dependent on the disorder, and that the low energy modes are extended states (critical metal) NTW . Taking into account the effects of inter-nodal scattering (hard-scattering) it has been shown that an insulating state is obtained instead, where the density of states still vanishes at low energy but with an exponent $`\alpha =1`$ independent of disorder Senthil . The addition of time-reversal breaking creates two new classes designated spin quantum Hall effect I and II, due to their similarities to the usual quantum Hall effect, corresponding to the hard and soft scattering cases, respectively report . The proposed formation of a pairing with a symmetry of the type $`d+id`$ breaks time-inversion symmetry did  but up to now remains a theoretical possibility. On the other hand applying an external magnetic field naturally breaks time-reversal invariance and therefore it is important to study the density of states in this case.

In general, disorder is due to the presence of impurities which may either scatter the quasiparticles and/or may serve as pinning centers for the field induced vortices. The density of states of a dirty but homogeneous $`s`$-wave superconductor in a high magnetic field, where the quasiparticles scatter off scalar impurities, was considered using a Landau level basis Dukan . For small amounts of disorder it was found that $`\rho (ϵ)ϵ^2`$ but when the disorder is higher than some critical value a finite density of states is created at the Fermi surface. In the same regime of high magnetic fields, but with randomly pinned vortices and no impurities, the density of states at low energies increases significantly with respect to the lattice case suggesting a finite value at zero energy Sacramento . Refs. Ye,  and Khveshchenko,  considered the effects of random and statistically independent scalar and vector potentials on $`d`$-wave quasiparticles and it was predicted Khveshchenko  that at low energies $`\rho (ϵ)\rho _0+aϵ^2`$, where $`\rho _0B^{1/2}`$. The effect of randomly pinned discrete vortices on the spectrum of a $`d`$-wave superconductor was considered recently and a preliminary report was presented in Ref. prl,  for the isotropic case. In this work we study the density of states as a function of vortex density and consider the effects of the Dirac anisotropy. We also study the local density of states (LDOS) and the inverse participation ratio (IPR) and study the effect of disorder on the localization of the low energy states.

The nature of the low energy states in the presence of a single vortex with $`s`$-wave symmetry was solved long ago Caroli . There are localized bound states in the vortex core. The $`d`$-wave symmetry case led to some early controversy but it was eventually clearly demonstrated that the low energy states, even though strongly peaked near the vortex core, extend along the four nodal directions and are indeed delocalized as shown by the behavior of IPR dvortex . In the vortex lattice the states are naturally also extended FT ; PRB . In the clean limit it is therefore expected that the external field will increase the low-energy density of states. The expectation that these are delocalized is evidenced by the increase of the thermal conductivity with field Chiao1  in contrast to the reducing effect of conventional superconductors Bessa .

The addition of impurities in zero magnetic field has been studied using the Bogoliubov-de Gennes (BdG) equations. It was found that the $`d`$-wave superconductivity is mainly destroyed locally near a strong scatterer. The superfluid density is strongly suppressed near the impurities but only mildly affected elsewhere Atkinson . No evidence for localization of the low energy states was found and accordingly the superfluid density is indeed suppressed but less than expected Franzet ; Randeria1  and, accordingly, the decrease of the critical temperature with disorder is much slower than previously expected in accordance with experiments Ulm . Similar results of an inhomogeneous order parameter were also obtained for $`s`$-wave superconductors Randeria2 .

The question we address in this paper is the influence of the positional disorder of the vortices on the quasiparticle states of either a $`s`$\- or a $`d`$-wave superconductor in an external magnetic field and we will not consider the scattering off impurities. We obtain that the low energy states are strongly peaked near the vortex locations as evidenced by the LDOS. At higher energies the states have a very uniform distribution throughout the system. For these states the IPR scales as $`1/L^2`$ indicative of extended states. In the case of the low energy states the same quantity does not follow this scaling but does not saturate either. This is probably a consequence of finite size effects indicative that the system sizes considered are smaller than the localization length. The results indicate at best a large value for the localization length but are more consistent with extended states rather than with localized states.

In section II we describe the model and present the BdG equations to be solved and consider the effect of the vortices on the supercurrent profiles. In section III we study the effects of positional disorder on s-wave superconductors. In section IV we consider d-wave superconductors for the isotropic and anisotropic cases. Also we study the spatial structure of the low-lying quasiparticle states. We end with the conclusions.

## II Description of the model

The strong correlations of the high $`T_c`$ materials are typically modelled using a Hubbard-like Hamiltonian for the electrons. The superconductivity is considered using a BCS-like formulation which provides good agreement with experimental results. Even though the normal phase of these materials is particularly challenging it is by now accepted that the phenomenology in the superconducting phase is reasonably well described by BCS theory.

### II.1 Bogoliubov-de Gennes Hamiltonian

We will consider the lattice formulation of a disordered $`d`$-wave superconductor in a magnetic field. We start from the Bogoliubov-de Gennes equations $`\psi =ϵ\psi `$ where $`\psi ^{}(𝐫)=(u^{}(𝐫),v^{}(𝐫))`$ and where the matrix Hamiltonian is given by

$$=\left(\begin{array}{cc}\widehat{h}& \widehat{\mathrm{\Delta }}\\ \widehat{\mathrm{\Delta }}^{}& \widehat{h}^{}\end{array}\right)$$
(1)
with FT ; PRB

$$\widehat{h}=t\underset{𝜹}{}e^{\frac{ie}{\mathrm{}c}_𝐫^{𝐫+𝜹}𝐀(𝐫)𝑑𝐥}\widehat{s}_𝜹ϵ_F$$
(2)
and

$$\widehat{\mathrm{\Delta }}=\mathrm{\Delta }_0\underset{𝜹}{}e^{\frac{i}{2}\varphi (𝐫)}\widehat{\eta }_𝜹e^{\frac{i}{2}\varphi (𝐫)}.$$
(3)
The sums are over nearest neighbors ($`𝜹=\pm 𝐱,\pm 𝐲`$ on the square lattice); $`𝐀(𝐫)`$ is the vector potential associated with the uniform external magnetic field $`𝐁=\mathbf{}\times 𝐀`$, the operator $`\widehat{s}_𝜹`$ is defined through its action on space dependent functions, $`\widehat{s}_𝜹u(𝐫)=u(𝐫+𝜹)`$, and the operator $`\widehat{\eta }_𝜹`$ describes the symmetry of the order parameter.

The application of the uniform external magnetic field $`𝐁`$ generates in type II superconductors compensating vortices each carying one half of the magnetic quantum flux, $`{\displaystyle \frac{\varphi _0}{2}}`$. We take the London limit, which is valid for low magnetic field and over most of the $`HT`$ phase diagram in extreme type-II superconductors like the cuprates, assuming that the size of the vortex cores is negligible and placing each vortex core at the center of a plaquette (unit cell). The $`N_\varphi `$ vortices are distributed randomly over the $`L\times L`$ plaquettes. In this limit we can take the order parameter amplitude $`\mathrm{\Delta }_0`$ constant everywhere in space and we can factorize the phase of the order parameter as shown in Eq. 3.

At this stage, it is convenient to perform a singular gauge transformation to eliminate the phase of the off-diagonal term (3) in the matrix Hamiltonian. We consider the unitary FT gauge transformation $`𝒰^1U`$, whereFT

$$𝒰=\left(\begin{array}{cc}e^{i\varphi _A(𝐫)}& 0\\ 0& e^{i\varphi _B(𝐫)}\end{array}\right)$$
(4)
with $`\varphi _A(𝐫)+\varphi _B(𝐫)=\varphi (𝐫)`$. The phase field $`\varphi (𝐫)`$ is then decomposed at each site of the two-dimensional lattice in two components $`\varphi _A(𝐫)`$ and $`\varphi _B(𝐫)`$ which are assigned respectively to a set of vortices $`A`$, positioned at $`\{𝐫_i^A\}_{i=1,N_A}`$, and a set of vortices $`B`$, positioned at $`\{𝐫_i^B\}_{i=1,N_B}`$. The phase fields $`\varphi _{\mu =A,B}`$ are naturally defined through the equation

$$\mathbf{}\times \mathbf{}\varphi _\mu (𝐫)=2\pi 𝐳\underset{i}{}\delta (𝐫𝐫_i^\mu )$$
(5)
where the sum runs only over the $`\mu `$-type vortices. After carying out the gauge transformation (4) the Hamiltonian (1) reads

$$^{}=\left(\begin{array}{cc}t\underset{𝜹}{}e^{i𝒱_𝜹^A\left(𝐫\right)}\widehat{s}_𝜹ϵ_F& \mathrm{\Delta }_0\underset{𝜹}{}e^{i\frac{\delta \varphi }{2}}\widehat{\eta }_𝜹e^{i\frac{\delta \varphi }{2}}\\ \mathrm{\Delta }_0\underset{𝜹}{}e^{i\frac{\delta \varphi }{2}}\widehat{\eta }_𝜹^{}e^{i\frac{\delta \varphi }{2}}& t\underset{𝜹}{}e^{i𝒱_𝜹^B\left(𝐫\right)}\widehat{s}_𝜹+ϵ_F\end{array}\right).$$
(6)
The phase factors are given by PRB  $`𝒱_𝜹^\mu (𝐫)=_𝐫^{𝐫+𝜹}𝐤_s^\mu 𝑑𝐥`$ and $`\delta \varphi (𝐫)=\varphi _A(𝐫)\varphi _B(𝐫)`$, where $`\mathrm{}𝐤_s^\mu =m𝐯_s^\mu =\mathrm{}\mathbf{}\varphi _\mu \frac{e}{c}𝐀`$ is the superfluid momentum vector for the $`\mu `$-supercurrent. Physically, the vortices $`A`$ are only visible to the particles and the vortices $`B`$ are only visible to the holes. Each resulting $`\mu `$-subsystem is then in an effective magnetic field

$$𝐁_{\mathrm{eff}}^\mu =\frac{mc}{e}\mathbf{}\times 𝐯_s^\mu =𝐁\varphi _0𝐳\underset{i}{}\delta ^2(𝐫𝐫_i^\mu )$$
(7)
where each vortex carries now an effective quantum magnetic flux $`\varphi _0`$. For the case of a regular vortex latticeFT ; PRB , these effective magnetic fields vanish simultaneously on average if the magnetic unit cell contains two vortices, one of each type. More generally, in the absence of spatial symmetries, as it is the case for disordered systems, these effective magnetic fields $`𝐁_{\mathrm{eff}}^{\mu =A,B}`$ vanish if the numbers of vortices of the two types $`A`$ and $`B`$ are equal, i.e. $`N_A=N_B`$, and their sum equals to the number of elementary quantum fluxes of the external magnetic field penetrating the system. As the number of vortices $`N_\varphi `$ in the two-dimensional system of size $`L\times L`$ is proportional to the quantized magnetic flux piercing through the system, we choose to parametrize the magnetic field intensity by the ratio of the number of vortices, $`N_\varphi =N_A+N_B`$, by the number of lattice sites, $`B=\frac{N_\varphi }{L\times L}`$.

We consider here the disorder induced by the random pinning of the $`N_\varphi `$ vortices over the two-dimensional $`L\times L`$ lattice. As the effective magnetic fields (7) experienced by the particles and the holes vanish on average within the gauge transformation (4) we are allowed to use periodic boundary conditions on the square lattice ($`\mathrm{\Psi }(x+nL,y)=\mathrm{\Psi }(x,y+mL)=\mathrm{\Psi }(x,y)`$ with $`n,m`$). The $`L\times L`$ original lattice becomes then a magnetic supercell where the vortices are placed at random with a mean intervortex spacing $`\mathrm{}=1/\sqrt{B}`$. The disorder in the system is then established over a length $`L`$.

In order to compute the eigenvalues and eigenvectors of the Hamiltonian (6) we seek for eigensolutions in the Bloch form $`\mathrm{\Psi }_{n𝐤}^{}(𝐫)=e^{i𝐤𝐫}(U_{n𝐤}^{},V_{n𝐤}^{})`$ where $`𝐤`$ is a point of the Brillouin zone. In the following we will label these eigensolutions with the index $`𝐧=(n,𝐤)`$. We diagonalize then the Hamiltonian $`e^{i𝐤𝐫}^{}e^{i𝐤𝐫}`$ for a large number of points $`𝐤`$ in the Brillouin zone and for many different realizations (around 100) of the random vortex pinning. The results of these computations are shown in section III for the $`s`$-wave disordered superconductor case and in section IV for the $`d`$-wave disordered superconductor case.

### II.2 Profiles of supercurrents

The $`\mu `$-superfluid wave vector $`𝐤_s^\mu (𝐫)`$ characterizes the supercurrents induced by the $`\mu `$-vortices. This vector can be calculated for an arbitrary configuration of vorticesPRB  like

$$𝐤_s^\mu (𝐫)=2\pi \frac{d^2k}{(2\pi )^2}\frac{i𝐤\times 𝐳}{k^2+\lambda ^2}\underset{i=1}{\overset{\mathrm{}}{}}e^{i𝐤(𝐫𝐫_i^\mu )}.$$
(8)
Here $`\lambda `$ is the magnetic penetration length and the sum extends over the infinite number of $`\mu `$-type vortex positions. We consider the case $`\lambda \mathrm{}`$ which is an excellent approximation in extreme type-II systems like high temperature superconductors. In this limit the repulsive interaction between vortices is not screened and therefore the vortex distribution is not strictly arbitrary – long range interactions will try to force the vortex system to take only incompressible configurations. For the purposes of this paper, we assume that the pinning centers are strong enough to overcome the vortex repulsion over the lengthscales relevant to experiments. Thus, we are considering the limit of strong pinning. We have checked that the case with a random distribution of vortices is qualitatively the same as for the case where the vortex positions are allowed to vary with a radius of a few unit cells around a regular lattice position.

The $`\mu `$-superfluid wave vector distribution is completely determined by the configuration of the $`\mu `$-vortices, as can be seen from Eq. 8 and is independent of the pairing symmetry. Since we are taking the London limit, where the vortex core size is negligible, we neglect any possible different symmetry contributions from inside the vortex core Inside  and the chiral supercurrents simply reflect the circulation of the superfluid density around the vortex.

In order to compute efficiently the $`\mu `$-superfluid wave vector we use the translational symmetry introduced by the periodic repetition of the $`L\times L`$ supercells and then we are able to use the Fourier series representation of expression (8), e.g.

$$𝐤_s^\mu (𝐫)=\frac{2i\pi }{(L𝜹)^2}\underset{𝑸\mathrm{𝟎}}{}\frac{𝐐\times 𝐳}{𝐐^2}\underset{i=1}{\overset{N_\mu }{}}e^{i𝐐\left(𝐫𝐫_i^\mu \right)}$$
(9)
where $`𝐐=\frac{2\pi }{L\delta }(n_1,n_2)`$ with $`n_1,n_2`$. We remark that the sum over vortex positions runs now only over the $`\mu `$-type vortices within a supercell and not, as in Eq. 8, over the infinite number of $`\mu `$-type vortices present in the two-dimensional system. The price to pay for this procedure is the preceding infinite sum over the reciprocal vectors $`𝐐`$. We truncate this latter sum when the convergence of each components of $`𝐤_𝐬^\mu (𝐫)`$ in each lattice point $`𝐫`$ is attained.

Figures 1 and 2 show examples of distributions of the conventional superfluid wave vector which is half the sum of the two types of superfluid wave vector,

$$𝐤_s(𝐫)=\frac{𝐤_s^A(𝐫)+𝐤_s^B(𝐫)}{2}=\frac{1}{2}\mathbf{}\varphi \frac{e}{c}𝐀.$$
(10)
In Fig. 1 we show the profile of the supercurrent velocities in the lattice case. The unit cell has two vortices, one of each type $`A`$ and $`B`$. In Figure 2 we show the supercurrents for an arbitrary configuration of the vortices. For example, the $`24\times 24`$ lattice shown in Fig. 2 can illustrate a disordered supercell corresponding to a magnetic field $`B=1/144`$ with two $`A`$-type vortices and two $`B`$-type vortices. As it should be, the labelling of the $`A`$\- and the $`B`$-type vortices is completely arbitrary. The intensity of the supervelocities for the regular and the disordered vortex lattice peak around each vortex and decrease rapidly outside the vicinity of the vortex core ($`5\delta `$).

## III Disordered $`s`$-wave superconductors

For the conventionnal $`s`$-wave case the operator characterizing the symmetry of the order parameter is constant $`\widehat{\eta }_\delta =\frac{1}{4}`$ and the off-diagonal terms of the Hamiltonian (6) are then considerably simplified

$$^{}=\left(\begin{array}{cc}t\underset{𝜹}{}e^{i𝒱_𝜹^A\left(𝐫\right)}\widehat{s}_𝜹ϵ_F& \mathrm{\Delta }_0\\ \mathrm{\Delta }_0& t\underset{𝜹}{}e^{i𝒱_𝜹^B\left(𝐫\right)}\widehat{s}_𝜹+ϵ_F\end{array}\right).$$
(11)
The $`s`$-wave Hamiltonian (11) describes coupled particles and holes evolving each in an effective magnetic field composed by the external magnetic field and the counteracting field of $`\varphi _0`$ flux carying vortices. The effective magnetic fields are characterized by the Peierls phase factors $`𝒱_𝜹^{\mu =A,B}`$. Independently of the system being disordered or not, the tails of the quasiparticles spectrum ($`ϵ\mathrm{\Delta }_0`$) are expected to be described by quantized Landau levels.

In Fig. 3 we show the density of states for a field $`B=1/18`$ in the case of no disorder. The gap is clearly seen at low energies characteristic of $`s`$-wave symmetry. The states inside the gap (for $`|ϵ|<t`$ in Fig. 3) are the typical Caroli - de Gennes - Matricon (CdGM) bound statesCaroli . At higher energies the Landau levels are clearly visible.

In Fig. 4 we show the effects of disorder for the same $`B`$ field in a lattice of $`18\times 18`$ sites ($`18`$ vortices). The gap is filled by the disorder. The reader should observe here that our disorder is in certain sense “infinite” since we place vortex positions completely at random. Thus, it appears that such full randomness in vortex positions, generated by strong pinning, suffices to close the gap in the single particle density of states. In constrast, the Landau level quantization structure clearly persists at high excitation energies. In general the increase of the magnetic field modifies the curvature of the quasiparticle density of states for the low energies. As it is shown in Fig. 5 the density of states seems to be decribed by the power-law formula

$$\rho (ϵ)ϵ^\alpha .$$
(12)
In Fig. 6 we fit the magnetic field dependence of the exponent $`\alpha `$. This exponent obeys the following law

$$\alpha c\mathrm{}d$$
(13)
where $`\mathrm{}=1/\sqrt{B}`$ is the mean intervortex spacing and $`d`$ is found close to $`1`$. In the cases of a strong magnetic field the density of vortices is high and strictly we are in a regime where the size of the vortex cores can not be neglected. In this regime the Landau level structure at high energies is clear and it extends to low energies superimposed by the effects of disorder. The lower limit for this high magnetic field regime is the value $`B0.16`$ for which $`\alpha =0`$ in Fig. 6.

## IV Disordered $`d`$-wave superconductors

For the unconventional $`d`$-wave case the operator $`\widehat{\eta }_𝜹`$ takes the form $`\widehat{\eta }_𝜹=(1)^{\delta _y}\widehat{s}_𝜹`$, we recall that $`\widehat{s}_𝜹`$ acts on spatial dependent functions as $`\widehat{s}_𝜹u(𝐫)=u(𝐫+𝜹)`$ and that $`𝜹=\pm 𝐱,\pm 𝐲`$ characterizes unit displacements (hops) on the lattice. With these definitions the $`d`$-wave Hamiltonian can be derived from the Hamiltonian (6) and reads

$$^{}=\left(\begin{array}{cc}t\underset{𝜹}{}e^{i𝒱_𝜹^A\left(𝐫\right)}\widehat{s}_𝜹ϵ_F& \mathrm{\Delta }_0\underset{𝜹}{}e^{i𝒜_𝜹(𝐫)+i\pi \delta _y}\widehat{s}_𝜹\\ \mathrm{\Delta }_0\underset{𝜹}{}e^{i𝒜_𝜹(𝐫)i\pi \delta _y}\widehat{s}_𝜹& t\underset{𝜹}{}e^{i𝒱_𝜹^B\left(𝐫\right)}\widehat{s}_𝜹+ϵ_F\end{array}\right)$$
(14)
where the phase factor $`𝒜_𝜹(𝐫)`$ has the form

$`𝒜_𝜹(𝐫)`$ $`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle _𝐫^{𝐫+𝜹}}\left(\mathbf{}\varphi _A\mathbf{}\varphi _B\right)𝑑𝐥`$ (15)
$`=`$ $`{\displaystyle \frac{1}{2}}{\displaystyle _𝐫^{𝐫+𝜹}}\left(𝐤_s^A𝐤_s^B\right)𝑑𝐥.`$
In the Hamiltonian (14) and in Eq. 15 the vector

$$𝐚_s=\frac{1}{2}\left(𝐤_s^A𝐤_s^B\right)$$
(16)
acts as an internal gauge field independent of the external magnetic fieldPRB . The associated internal magnetic field $`𝐛=\mathbf{}\times 𝐚_s`$ consists of opposite $`AB`$ spikes fluxes carying each one half of the magnetic quantum flux $`\varphi _0`$, centered in the vortex cores and vanishing on average since the numbers of $`A`$\- and $`B`$-type vortices are the same.

### IV.1 Differences from the regular vortex lattice case

The situation where the vortices are regularly distributed in a lattice was treated before PRB . Since the average effective magnetic field vanishes it is possible to solve the BdG equations using a standard Bloch basis– the supercurrent velocities are periodic in space and there is no need to consider the magnetic Brillouin zone. Taking the continuum limit and linearizing the spectrum around each node effectively decouples the nodes. It was shown that the low-energy quasiparticles are then naturally described as Bloch waves FT  and not Dirac-Landau levels as previously proposed Gorkov ; Anderson . However, it was found that in the linearized problem different assignments of the $`A`$ and $`B`$ vortices lead to somewhat different spectra, which was unexpected PRB . It was found that defining theory on the lattice regularized this problem and indeed the system has a manifest internal gauge symmetry such that the spectrum is independent of the $`A`$-$`B`$ vortex assignments, as it should be. Moreover, the lattice formulation explicitly involves internodal contributions which are important for the properties of the density of states in the disordered case. In the vortex lattice case, however, it was found that only in special commensurate cases (for the square lattice) the inclusion of the internodal contributions is relevant since only in such cases a gap develops due to the interference terms between the various nodes, estimated to be of the order of $`\sqrt{B}`$. In a general incommensurate case the interference is not relevant leading to qualitatively similar spectra. In the $`d`$-wave case the spectrum is gapless with a linear density of states at low energy PRB . One would therefore expect that in a general disordered vortex case internodal scattering might not be relevant (particularly for high Dirac cone anisotropy $`\alpha _D=v_F/v_\mathrm{\Delta }=t/\mathrm{\Delta }_01`$).

In Fig. 7 we compare the densities of states for the lattice case and a case with disorder. At weak fields the density of states is small at low energies having a dip close to zero energy. We have checked for finite size effects on the spectrum. For system sizes larger than $`16\times 16`$ the density of states at not very low energies converges and the finite size dependence is negligible.

### IV.2 Isotropic case

In this subsection we will focus our attention on the isotropic case $`\alpha _D=1`$ in order to extract the general dependencies of the quasiparticle density of states in the magnetic field $`B`$.

In Fig. 8 we show the density of states for a system with size $`20\times 20`$ and for various magnetic fields. The density of states at small energies is finite up to quite small energies where there is a dip to a value that decreases as the magnetic field decreases. Only for quite small magnetic fields the density of states approaches zero at the origin. Neglecting the narrow region close to the origin we have fitted the density of states using the power law

$$\rho (ϵ)=\rho _0+\beta |ϵ|^\alpha .$$
(17)
In the inset of Fig. 8 we show the fits for the various values of the magnetic field. Reasonable fits are obtained taking $`\alpha 2`$ and we obtain that $`\rho _0B^{1/2}`$. The various system sizes fit in the same universal curve indicating that the finite size effects are negligible. Note that in the lattice case the density of states at low energies is linear FT ; Marinelli ; PRB ; AV  (this result differs from the behavior obtained by others for a $`d`$-wave superconductor with no disorder Volovik ; Wang , where $`\rho (ϵ0)B^{1/2}`$). The finite density of states at zero energy is therefore a consequence of finite disorder.

In Fig. 9 we focus on the narrow region close to $`ϵ=0`$ for the same set of parameters considered in Fig. 8. Except for the lowest field case $`B=1/200`$ the density of states seems to be finite at zero energy. Here the field density $`B=1/200`$ corresponds to the particular case where only two vortices pierced the $`20\times 20`$ disordered supercell. As shown in Ref. Marinelli,  in this case the spectrum is usually gapped and therefore the density of states vanishes at zero energy. Performing a fit like in Eq. 17 we obtain an exponent which is now close to $`1`$. In this regime $`\rho _{0\mathrm{d}\mathrm{i}\mathrm{p}}`$ also scales with $`\sqrt{B}`$ and the slope scales linearly with $`B`$. In this low energy regime the finite size effects are still noticeable but the dependence on the magnetic field is common to the various system sizes. At these low energies the density of states for the various system sizes appears to be of the following approximate form

$$\rho (ϵ)\frac{1}{\omega _H}\frac{1}{l^2}\left(\frac{ϵ}{\omega _H}\frac{\delta ^2}{l^2}\right),$$
(18)
where $`\omega _H\sqrt{\mathrm{\Delta }B}`$ and $``$ is a universal function. In the left panel of Fig. 9 we show $`\rho (ϵ)`$ for various fields while in the right panel we illustrate the near scaling at low energies consistent with $`(x)c_1+c_2x`$ at small $`x`$. Note that $`\beta (B0)`$ appears to be small but finite, consistent with a crossover to a “Dirac node” scaling $`\rho (ϵ)(1/\omega _H)(1/l^2)\left(ϵ/\omega _H\right)`$ at very low fields.

### IV.3 Anisotropic case

We investigate now the dependence of the quasiparticle density of state on the value of the Dirac anisotropy ratio $`\alpha _D`$. Such a study is interesting in order to compare our results with experiments; indeed for high-$`T_c`$ superconductors such as YBa<sub>2</sub>Cu<sub>3</sub>O<sub>7</sub> and Bi<sub>2</sub>Sr<sub>2</sub>CaCu<sub>2</sub>O<sub>8</sub> the Dirac anisotropy ratio isChiao1 ; Chiao2  respectively $`\alpha _D14`$ and $`\alpha _D19`$.

In the lattice case a high anisotropy increases the density of states at low energies and leads to lines of quasi-nodes PRB ; Marinelli . At high anisotropy ($`\mathrm{\Delta }_0<<t`$) the nodes are very narrow and a one-dimensional like character is evidenced.

As shown in Fig. 10, the low-energy quasiparticle density of states at constant field is filled when the anisotropy parameter $`\alpha _D`$ is increased. There is a narrow linear region close to the origin that decreases as $`\alpha _D`$ increases. Fig. 11 shows the power law fits of the data presented in Fig. 10 neglecting a very narrow region of width $`0.025t`$ around $`ϵ=0`$. It turns out that the low-energy quasiparticle density of states has the form

$$\rho (ϵ)\rho _0+\mathrm{\Gamma }\left(\frac{ϵ}{\mathrm{\Delta }_0}\right)^\alpha $$
(19)
with zero-energy density of states

$$\rho _0\sqrt{\frac{t}{\mathrm{\Delta }_0}},$$
(20)
the exponent

$$\alpha 1+\sqrt{\frac{\mathrm{\Delta }_0}{t}}$$
(21)
and the parameter

$$\mathrm{\Gamma }\frac{\mathrm{\Delta }_0}{t}.$$
(22)
The dependence of the exponent is interesting since in the isotropic case we retrieve $`\alpha =2`$ and with increasing anisotropy the exponent decreases to $`1`$ (see lower left panel of Fig. 11) characteristic of the Dirac nodes. At $`ϵ0`$ and for a high Dirac anisotropy ratio the quasiparticle density of states flattens since the coefficient $`\mathrm{\Gamma }`$ in (19) decreases with $`\mathrm{\Delta }_0`$.

### IV.4 Spatial structure of the low-lying quasiparticle states

As discussed above, the nature of the low lying states is an important point which is relevant due to the presence of disorder. A standard way to analyse the nature of the states is to study the IPR. The IPR is defined in the usual way

$$\xi (ϵ)=\frac{_{𝐧,𝐫}\left(\left|u_𝐧\left(𝐫\right)\right|^4+\left|v_𝐧\left(𝐫\right)\right|^4\right)\delta \left(ϵE_𝐧\right)}{_{𝐧,𝐫}\left(\left|u_𝐧\left(𝐫\right)\right|^2+\left|v_𝐧\left(𝐫\right)\right|^2\right)\delta \left(ϵE_𝐧\right)^2}.$$
(23)
The brackets denote the averaging over disorder configurations. The IPR $`\xi `$ is a direct measure of the spatial extend of the quasiparticle wavefunctions. It scales as $`1/L^2`$ for extended states and is constant for localized states with localization length $`\mathrm{}_c<L`$.

Fig. 12 presents the IPR $`\xi `$ as a function of the quasiparticle energy $`ϵ`$ and for different system sizes $`L`$. All the states are well extended, since there is no size dependence for the quantity $`L^2\xi (ϵ)`$, except those in the close vicinity of the Fermi surface ($`|ϵ|<0.5t`$). In the inset of Fig. 12 we present the scaling of the IPR $`\xi `$ with $`1/L^2`$. We consider an average over states at the Fermi surface with $`|ϵ|<0.05t`$ ($`\xi _0`$ in the inset) and an average over states with $`ϵt`$ ($`\xi _{1.0}`$ in the inset). The states close to $`ϵ=t`$ are clearly extended since the scaling with $`1/L^2`$ is quite accurate. The nature of the states close to $`ϵ=0`$ is however less clear. From Fig. 12 we see that the quasiparticle states close to $`ϵ=0`$ deviate from the strict $`1/L^2`$ scaling law (see the weak intercept of $`\xi _0`$ in the inset of Fig. 12), however our results show an obvious size dependence for $`\xi _0`$ and do not show a saturation of $`\xi _0`$ with the linear system size $`L`$. This latter fact indicates, in the hypothesis of an eventual localization of these low-lying states, that we are still far from the thermodynamic limit with the linear system sizes we can currently attain within our model.

To gain further insight into the nature of the low lying states we show in Fig. 13 the LDOS for the states close to $`ϵ=0`$ and those close to $`ϵ=t`$. At low energies the LDOS defined by

$$\rho (𝐫,ϵ)=\underset{𝐧}{}\left(\left|u_𝐧(𝐫)\right|^2+\left|v_𝐧(𝐫)\right|^2\right)\delta (ϵE_𝐧)$$
(24)
is strongly peaked near the vortex cores (upper left panel of Fig. 13). At higher energies ($`ϵ1`$) the LDOS is much more homogeneously spread over the system indicative of extended states; the values of the LDOS over the lattice are approximatively constant and slightly fluctuate around the expected value for extended states $`1/L^20.0017`$ (lower left panel of Fig. 13).

The right side panels of Fig. 13 show the profiles of the internal gauge field $`𝐚_s`$ (upper panel) and of the conventional superfluid wave vector $`𝐤_s`$ (lower panel) corresponding to the same configuration of vortex pinning as used in the computation of the LDOS (left side panels). For the profile of the superfluid wave vector $`𝐤_s`$ the intensity of the supercurrents is peaked around each vortices and rapidly decreases as the inverse of the square of the distance to the vortices. The profile of the internal gauge field $`𝐚_s`$ is more interesting. It depicts the two types of vortices $`A`$ and $`B`$ whirling in opposite directions (see upper right panel of Fig. 13). Due to this opposite circulation around $`A`$\- and $`B`$-type vortices, the interaction between a $`A`$\- and a $`B`$-type vortex produces significant field currents between the vortices. Fig. 13 shows a spatial correlation between these intervortex currents and the inhomogeneity of the local density of states between the vortices for the low-lying states $`ϵ0`$ (see upper left panel of Fig. 13).

If we compare now the results for the LDOS with those for the IPR $`\xi `$ for the low-lying states, we can argue the following: as the participation ratio $`\xi _0^1`$ physically counts the number of lattice sites occupied by the quasiparticle wave function $`\mathrm{\Psi }(𝐫)`$, and as the low energy states are mainly located around the vortices, the participation ratio $`\xi _0^1`$ should scale for a fixed magnetic field as the number of vortices in a supercell $`\xi _0^1N_\varphi BL^2`$. This argument explain the fact that $`\xi _0`$ in the inset of Fig. 12 seems to follow the $`1/L^2`$ without saturation, the weak intercept being then negligible. The low-lying quasiparticle states appear then to be delocalized although strongly peaked around the vortex cores.

We present now a spatial scanning of the LDOS in the vicinity of a vortex core for the case of a regular vortex lattice (Fig. 14) and for the case of a disordered vortex lattice (Fig. 15). For both cases the LDOS away from the vortex core are qualitatively the same and are comparable to that of a $`d`$-wave superconductor in a zero-magnetic field. Also for both cases the low-lying states are predominant in the close vicinity of the vortex core but their respective spectra are different. For the regular vortex lattice case (Fig. 14) we remark a double peak structure around the vortex core. This result is in qualitative agreement with the double peak in conductance observed in YBa<sub>2</sub>Cu<sub>3</sub>O<sub>7-δ</sub> by scanning tunneling microscopy (STM) Maggio . For the disordered vortex case (Fig. 15) zero energy peaks (ZEP) appear in the close vicinity of the vortex core (closest neighbor sites) and rapidily disappear (typically over $`3\delta `$) when moving away from the vortex core. We note also that close to the vortex cores the coherence peaks dissapear in both the regular and the disordered cases.

## V Conclusion

In summary, we have calculated the density of states of a disordered superconductor in a pinned fully random vortex array. Both the disorder and the magnetic field fill the density of states at low energies. In the $`s`$-wave case the density of states behaves as a power law with an exponent that scales with $`1/\sqrt{B}`$. In general we find a finite density of states at zero energy for the $`d`$-wave case except in the limit of very small magnetic fields. The density of states deviates from the zero energy value by a power law. The zero energy density of states scales with the inverse of the magnetic length ($`\sqrt{B}`$). In the $`d`$-wave case the Dirac anisotropy further increases the weight of the density of states at low energies. Also it affects the exponent of the power law. In the isotropic case the exponent is $`1`$ at very low energies and around $`2`$ neglecting this narrow region. In the linear regime the density of states is of the form of a scaling function of the energy and the magnetic field. As the anisotropy increases this narrow regime shrinks considerably and, neglecting this region, the exponent of the density of states interpolates to $`1`$ which is the Dirac limit. This limit is also obtained in the zero field limit in the isotropic case. Except for the zero energy finite value the energy dependence of the density of states in the case with disorder is similar to the lattice case. This suggests that the vortex disorder does not dramatically affect the density of states at low energies. An analysis of the IPR and the LDOS shows that the lowest lying states are delocalized, even though strongly peaked at the vortex cores. In the gapped $`s`$-wave case however the disorder introduces states in the gap thereby changing qualitatively the low energy density of states, as in the high field limit Sacramento . We found a power law behavior with an exponent that scales with the magnetic length.

###### Acknowledgements.

This work was supported in part by NSF grant DMR00-94981 (ZT) and by FCT Fellowship SFRH/BPD/5602/2001 (JL).


format: text

text:
# 1 𝑑𝐼/𝑑𝑉 obtained from 𝑅_𝑐(𝑇) of [1] (shown in the insert) assuming Joule heating origin of I(V) nonlinearities.

Intrinsic tunneling or Joule heating?

In a recent letter  Yurgens et al. reported superconducting $`\mathrm{\Delta }_s`$ and normal $`\mathrm{\Delta }_p`$ (pseudo)gaps in the c-axis differential conductance of Bi2201-La mesas and ruled out a precursor Cooper pair scenario as an explanation of $`\mathrm{\Delta }_p`$, while mechanisms based on the Van Hove singularity in the density of states, or on the resonant tunnelling between $`CuO_2`$ planes, were outlined as possible explanations of the normal state pseudo-gap. I will show that a Joule heating of mesas may cause $`I(V)`$ nonlinearities similar to those observed in Ref..

This approach provides a natural explanation of another finding of Yurgens et al. , the surprising similarity between the ‘intrinsic Josephson effect’ (IJE) characteristics of Bi2201 and Bi2212 despite of a threefold difference in their critical temperatures. Moreover, similar (to those applied to Bi2212) levels of dissipation (estimated as $`VI/A`$) are required to achieve the characteristic features of IJE-spectra of Bi2201 at the same temperature, 4.2-4.5K (Fig.1,2a ) namely, 0.03-0.1, 0.6-3.9, and 2.5-15.6$`kW/cm^2`$ for the end of ‘brush’, $`\mathrm{\Delta }_s`$ and $`\mathrm{\Delta }_p`$ respectively (the scatter reflects an uncertainty in the mesa area, $`A`$, in ). The relevance of heating in Bi2212 is experimentally confirmed by different methods; in particular, significant overheating of the mesa’s surface was reported by  even in the ‘brush’-like part of I(V). I believe that the character of Joule heating of a sample with typical c-axis R(T) dependence, may be responsible for the similarities reported in .

Here I will show that the ‘tunnelling characteristics’ from Fig.4  could be reproduced qualitatively and quantitatively using the experimental out-of-plane normal state resistance $`R_c(T)`$ (inset in Fig.4 of Ref.) and assuming that the heating of the mesa caused by the Joule dissipation is the only reason for effects observed at high bias. Using Newton’s Law of Cooling (1701), the temperature of a thin mesa is given by

$$T=T_0+IV/(Ah),$$
(1)
where $`T_0`$ is the temperature of a coolant medium (liquid or gas) and $`h`$ is the heat transfer coefficient. According to Eq.(1), monitoring I(V) at certain base temperature $`T_0`$ results in a sample temperature rise that entails non-linearity in $`I=V/R_c(T)`$. Thus constructed $`dI/dV`$ curves resemble those presented in fig.4 . This similarity allows for an estimate of the inverse heat transfer coefficient $`(hA)^162.5K/mW`$ and overheating at high bias, $`80K`$ for $`T_0=200K`$; the set of curves accounting for this coefficient is shown in Fig.1. As it is clearly seen from Fig.1, the ‘heating’ spectra taken at $`T_0<T^{}`$ reveal a ‘pseudo-gap’ which disappears entirely as soon as $`T_0>T^{}`$. This result shows that Eq.(1) provides a natural explanation of some of the puzzling findings of .

To conclude, the nonlinear I(V) characteristics observed by  in Bi2201 mesas are related to the temperature dependence of the normal state c-axis resistance, rather than to the ‘Josephson tunnelling’ between the planes. As far as the unusual $`R_c(T)`$, itself is concerned, a few models were proposed. In particular, the $`R_c(T)`$ explanation in the framework of the bipolaron model of cuprates  was supported experimentally . In this model the temperature-independent-pseudogap is half of the binding energy of bipolarons. Although only the normal state data were considered here, there is no doubt that heating plays an even more important role at low temperatures, so that heating issues have to be accounted for in the analysis of the low temperature data also.

The financial support of the Leverhulme Trust (F/00261/H) is gratefully acknowledged.

V.N.Zavaritsky

Department of Physics, Loughborough University, Loughborough LE11 3TU, United Kingdom.

 A.Yurgens et al., Phys. Rev. Lett. 90, 147005,(2003).

 C.E.Gough et al., cond-mat/0001365(2000).

 V.N.Zavaritsky, J Superconductivity (US) 15, 567 (2002).

 A.S. Alexandrov et al., Phys. Rev. Lett. 77, 4796 (1996).

 V.N.Zavaritsky et al., Europhys. Lett, 51, 334, (2000);

V.N. Zverev and D.V. Shovkun, JETP Lett. 72, 73 (2000).


text:
# OVAL: the CMS Testing Robot

## I Introduction

New agile methologies XPROG  put the test programs at the center of the software development: they should be written before the source code, and they should be executed whenever the code is modified. This optimistic goal can only be approached if this task is really easy and as automated as possible. Oval is a tool which tries to fulfill this need.

Oval OVAL  was initially developed in the context of CMS CMS  Software Quality Improvement to support the high level validation of physics algorithms. We then felt the tool could also be used for low level unitary tests, and it appears to be actually its first use.

As a CMS tool, Oval was primarily developed for the Unix world, and given a look-and-feel close to Scram SCRAM , the CMS tool for configuration and build. That means it has a simple command line interface, similar to the CVS one CVS , the Oval configuration files are XML-like XML , and the whole is written in PERL.

Oval is mainly a robot which runs the tests for you. It does not provide any helper code, therefore it is quite independant from your programming or scripting language, but you must write your own tests. Another key point is that Oval does not analyze some dedicated log files or output channels, but the standard Unix output. This makes it only usable with programs which use this way of communication. On the other hand, this allows to scan any output from any third party software, even when this software has not been especially instrumented for Oval.

While Oval was appearing useful and usable outside of CMS, it has been progressively isolated from Scram and extended with many configuration possibilities. Today, it can be used either by a single programmer with few C++ files and a simple makefile, or by a large collaboration with complex software processes and tools.

In this paper, we will first browse the main features of Oval: how it must be configured by the users, what are the main commands and what one will find in the log files. At last, we will discuss the many ways which can be used by the administrators so to customize Oval.

## II Setup of a new test directory

Given a directory with test programs, one must create a special file called "OvalFile", where one specifies all the information needed by Oval about the local tests. The first thing to do is to list the programs which must be taken into consideration. For example:

```
<program name="Clusters.cpp">
<program name="Electrons.cpp"
<program name="EnergyFlow.cpp">

```

If your directory is dedicated to the test, you will surely list here all the programs of the directory. Then, you can ask for the compilation of the programs by typing "oval build". In this case, Oval does not do much, but simply delegates the compilation to the tool which has been specified by your Oval administrator. Yet, the output of the compilation is stored in a log file for each program. In our example, three files will be created ("Clusters.log", "Electrons.log" and "EnergyFlow.log"). The content of "Electrons.log" could be this:

```
[oval build] ===========================
[oval build] make Electrons
[oval build] ===========================

c++ -O2 -o Electrons Electrons.o
  -L/home/chamont/ORCA_5/lib/Linux__2.2
  -L/opt/cms/Releases/ORCA_5/lib/Linux__2.2
  -L/cern/2001/lib
  -lElectronFacilities -lEgammaH4Support
  -lz -lnsl -ldl -lg2c -lm

```

Some header lines (starting with "\[oval build\]") recall the actual compilation command which has been used, and the lines after are the compilation output (reformatted by Oval so to improve legibility).

## III Runtime conditions

Once all the programs compile successfully, one must care about their runtime conditions. Again, "OvalFile" is the place where to specify this information. Here is the new content we could give to the configuration file:

```
<var name="FEDERATION" value="cmsuf01">
<environment name="pt15">
  <var name="DATASET" value="eg_ele_pt15">
  <program name="Clusters.cpp">
  <program name="Electrons.cpp"
     args="-geo detailed">
</environment>
<environment name="flow">
  <var name="DATASET" value="jm_minbias">
  <program name="EnergyFlow.cpp">
</environment>

```

The second program, "Electrons", will be executed with the command-line arguments "-geo detailed". Also, before running any of the programs, Oval will set the Unix variable "FEDERATION" to the value "cmsuf01". On the contrary, the value of "DATASET" depends on the programs. The first two programs will see a value of "eg\_ele\_pt15", although the third one will see a value of "jm\_minbias". The XML tags are self-explanatory and do not require much more comments.

On top of Unix variables, on can also specify some files to be created on the fly, before running the programs. We call them ”auxiliary files”. It is rarely relevant to specify by this way the input data (often huge and/or binary), but this applies well for the configuration files of the software to be tested, especially if you want to specify different file contents for the different programs. Here is an example of such a specification, to be added to our "OvalFile":

```
<file name=".orcarc">
  GoPersistent = 1
  MaxEvents = 500
  Random:Seeds = 0 3
</file>

```

If only a subset of the file depends on the program, one can specify the common part at the top level of the "OvalFile", and each specific part in the relevant environments. When Oval prepares an auxiliary, it first searches for the different parts and concatenates them.

You can now run all your programs with the command "oval run". You can also run one program, for example "oval run Electrons". When executing each program, Oval adds a new section to the log file of the program:

```
[oval run] =============================
[oval run] LD_LIBRARY_PATH =
[oval run]   /home/chamont/ORCA_5/lib
[oval run]  :/opt/cms/Releases/ORCA_5/lib
[oval run]  :/opt/Objectivity/5.2.1/lib
[oval run] FEDERATION = cmsuf01
[oval run] DATASET = eg_ele_pt15
[oval run] .orcarc:
[oval run]   GoPersistent = 1
[oval run]   MaxEvents = 500
[oval run]   Random:Seeds = 0 3
[oval run] =============================

                Welcome to COBRA
...

```

The header lines (starting with "\[oval run\]") recall the Unix variables which have been defined and the content of the auxiliary files, followed by the real output of the program.

## IV Analysis of the programs output

If the software system to be tested is a little verbose, you will never get exactly the same output from one software release to the other, and Oval must be told what to compare, otherwise it will always warn you about irrelevant differences. As you can guess, this is done by adding tags to "OvalFile", and this can differ from one program to the others if you enclose the comparison tags and the programs in different environments. Let’s comment a trivial example of comparison tags:

```
<diffline expr="^OVAL:">
<diffnumber expr="^energy: (.*)$"
            tolerance="5%">

```

The first tag "\<diffline\>" is the default behavior of Oval: it specifies to compare all the lines which start with "OVAL:". It is not enough if you want to scan output from a third party piece of software, where you cannot insert "OVAL:" at the beginning of relevant lines. Then you can put additional tags \<diffline\>, where the value of the attribute expr is a PERL regular expression describing some relevant lines.

In the example above, the second tag is more interesting if you want to track the value of some important physics quantities. This tag defines some numbers to be checked: the regular expression specifies which lines to take into consideration, with the parenthesis specifing where the number within the line is, and the second attribute specifying how much the value can differ from the reference.

Once you have specified your comparison strategy, you must eventually modify your programs so to display the expected outputs, and run them until their outputs can be considered as a good reference. Then, the command "oval validate" will register the current outputs as the references. Basically, Oval simply copies each file "\<program\>.log" as "\<program\>.ref". For example, if we validate the first two programs with the command "oval validate Clusters Electrons", the content of the directory will be:

```
cms038:ElectronPhoton> ls
  Makefile
  Clusters.cpp
  Clusters.log
  Clusters.ref
  Electrons.cpp
  Electrons.log
  Electrons.ref
  EnergyFlow.cpp
  EnergyFlow.log
  OvalFile

```

## V Validation sessions

At this stage, your test directory is quite ready to be used. Each time the software to be tested is modified, you can move to the directory and type the top level command "oval prod". For each program which has a reference, Oval will call "oval build", then "oval run", and finally "oval diff" which will perform the comparison between the log file and the ref file. The output of Oval could look like this:

```
cms038:ElectronPhoton> oval prod
  Clusters: build, run, diff.
  Electrons: build, run, diff (DIFFS).

```

Above, Oval has noticed a difference between the Electrons output and the reference. If you want more details, you can run separately the command "oval diff Electrons" or look at the new section which has been added at the end of "Electrons.log":

```
[oval diff] ===========================
[oval diff] diff line: /^OVAL:/
[oval diff] diff number: /^energy: (.*)$/ ~5%
[oval diff] ===========================

ref#1452 != log#2053
ref: OVAL: 12 electrons
---
log: OVAL: 11 electrons

ref#1972 !~ log#2592 (>5%)
ref: energy: 29.7275
---
log: energy: 27.4728

```

In this log file, the header lines recall the comparison tags which have been used, followed by a description of all the differences, with a format inspired from the Unix diff command.

## VI Other Features

Since it has been developed and used for several years, Oval has accumulated many other features. Let’s list the important ones:

1. Oval is not restrained to programs. For complex testing tasks (for example when databases are involved), one can also run executable scripts or ready binaries which do not need to be built.
2. While a test is getting more complex and is given runtime options, one could want to run this test several times with different runtime conditions (commandline arguments, Unix variables, auxiliary files). An "OvalFile" allows several occurrences of the same program. Each occurence will have its own log file.
3. If there are many tests, they should be organized into a hierarchy of Unix directories. One can insert an "OvalFile" at any level: the specifications will be recursively propagated to the subdirectories.
4. The results of a validation session can be mailed to the watchers specified in the configuration. This is the seed of a future interface with more elaborated bug tracking tools.

## VII Site Customization

Oval is potentially collaborating with many other actors: the operating system, the configuration tool, the build tool, the compiler, the run tool, the bug tracking system, etc. No two teams are using the same tools and apply the same software process. This is why it is especially important to make Oval as flexible as possible. So to demonstrate it, this section focuses on what can be done by a site administrator. We will begin with the description of the versions and flavors control, then have a look at the interfaces with the external tools, and finally discuss how to add new local commands.

### VII.1 Versions and Flavors

All the Oval versions of a site are expected to be installed under a common directory, which we will call "OVAL\_DIR". When an oval executable is invoked, it first compares its internal version with the one required by the user, and if it is not the same, it is able to find the correct executable under "OVAL\_DIR" and execute this one instead. The user can require a specific version by setting the variable "OVAL\_DIR", or thanks to the tag "\<oval version=...\>" within an "OvalFile". The latter allows to attach each set of test programs to a given oval version, so that it will still work when the site administrator will install a new oval version and make it the default (provided this administrator does not erase the old oval versions…). This is not only useful, but absolutely necessary for a site with many users.

On top of the version, the users can also select a given ”flavor”, thanks to the variable "OVAL\_FLAVOR". Depending on the flavor, oval will use different default configurations, so to fill the needs of different user groups. For example, if the value of "OVAL\_FLAVOR" is "salty", oval will take its default values from "OVAL\_DIR/site/OvalFile.salty". The definition of the flavors and the corresponding configurations is up to the site administrators.

Actually, whenever Oval needs a file which can be customized (for example "OvalFile"), it first scans for this file in the "OVAL\_DIR/site" directory, where it will select in priority the most specialized file it can found, for a given version and flavor (for example "OVAL\_DIR/site/version/flavor/OvalFile" or "OVAL\_DIR/site/OvalFile.version.flavor"). This lets the site administrators provide default configurations for any combination of version and flavor, if they need to.

These site files are typically used for the definition of instructions which are platform-dependant. They can also be used to enforce options on the commands. Here is an example:

```
<options command="expr" value="-v">
<config name="mail instruction"
        value="/bin/mail -s %">
<config name="custom url"
        value="http://www.site.fr/oval">

```

### VII.2 Access to the External Tools

The use of external tools, such as scram or make, is done through perl wrappers. The directory "OVAL\_DIR/OVAL\_VERSION/share/Interfaces" contains the wrappers delivered with Oval. Since the features of the various tools are somehow fuzzy and overlapping, we made the list of all the functions that the external tools should implement, and we splitted them into few groups, that we will call ”interfaces”. There is a global map which associates a concrete tool to each interface. Whenever Oval needs a given external function, which is part of a given interface, Oval looks which tool is attached to the interface, and invokes the function from this tool. Currently, we have defined two interfaces: "build" which is attached by default to the tool "make", and "run" which is attached by default to the tool "oval" (the pseudo tool "oval" provides a default implementation for all the defined interfaces). If the users of a site are using Scram, the administrator can declare it in "OVAL\_DIR/site/OvalFile":

```
<config name="build tool" value="scram">
<config name="run tool" value="scram">

```

As we described above, this configuration can depend on the version and/or the flavor. If an administrator want to use a tool unknown to Oval (for which no wrapper is provided), he can copy one of the provided wrapper from "OVAL\_DIR/OVAL\_VERSION/share/Interfaces" to "OVAL\_DIR/site/Interfaces" and modify it appropriately. One can also reimplement an existing wrapper: the site wrapper will hide the default Oval one. This can also depend on the version and/or flavor. For example, a tool wrapper called "OVAL\_DIR/site/salty/Interfaces/cmt.pm" would only be visible to the users who set "OVAL\_FLAVOR" to "salty", and the tool wrapper "OVAL\_DIR/site/3\_0\_0/Interfaces/make.pm" will hide the default one only for the users who set "OVAL\_VERSION" to "3\_0\_0".

### VII.3 Commands

When a need is so specific that the Oval team can hardly include it in its standard packaging, the site administrators can provide a specific site command. If their file respects the Oval commands requirements, and is placed in a directory called "Commands" ( for example "OVAL\_DIR/site/Commands/mycom.pm" ), the new command will be detected by Oval and made available to the users. As one can guess, this requires a rather deep understanding of Oval internals.

Even more ”nasty”, the administrators can provide their own implementation of the existing commands.

## VIII Conclusion

As hopefully demonstrated, Oval can be used for simple projects as well as configured for large collaborations which use other external tools and a strict software process. The only aspect which could prove a little painful is the design of the regular expressions.

Any piece of software, written in any language, can be automatically tested all along its life (provided it is talking to the output channel). Unhapilly, the programmers still have to write the tests and maintain the reference outputs, but we are not convinced any tool will ever avoid this.

The CMS collaboration currently uses Oval for its regression tests, and the CERN LCG SPI group now recommends the use of Oval together with CppUnit CPPUNIT  which provide some C++ helper code. Oval is not only useful for the punctual validation of software releases: in the spirit of agile methods, it can be used daily so to detect the regressions as soon as they appear.

For the near future, the main task will be to add a tutorial and to react to the users feedback. We also want to better document the internals, use the CERN LCG savannah server and open the development of Oval. At last, we plan to redesign the output and log files of Oval, so to ease the life of librarians who handle huge hierarchies of tests, and enable to trace the status of all these tests.


text:
# 1 Introduction

## 1 Introduction

It is hard to deny that the Search for ExtraTerrestrial Intelligence (SETI) is one of the major scientific adventures in the history of humankind. At the beginning of XXI century it remains the oldest and perhaps the most fascinating scientific problem. However, the field is still largely qualitative and thus often not taken seriously enough. One of the attempts to overcome this circumstance is encapsulated in the famous Drake equation, developed by Frank Drake for the first SETI symposium in 1961.

The first problem any student of SETI faces is that there is no canonical form of the Drake equation. Various authors quote various forms of the equation, and it is in a sense dependent on what is the desired result of the analysis. We shall investigate the following form (e.g. Shklovskii and Sagan 1966; Walters, Hoover, and Kotra 1980):

$$N=R_{}f_gf_pn_ef_lf_if_cL,$$
(1)
while keeping in mind that other equivalent forms exist as well. In this expression, the symbols have the following meanings:

$`N`$ = the number of Galactic civilizations with whom communication is possible.

$`R_{}`$ = mean rate of star formation in the Galaxy,

$`f_g`$ = fraction of stars suitable for supporting life,

$`f_p`$ = fraction of stars with planetary systems,

$`n_e`$ = number of planets per planetary system with conditions ecologically suitable for the origin and evolution of life,

$`f_l`$ = fraction of suitable planets where life originates and evolves into more complex forms,

$`f_i`$ = fraction of planets bearing life with intelligence,

$`f_c`$ = fraction of planets with intelligence that develops a technological phase during which there is the capability for an interest in interstellar communication,

$`L`$ = mean lifetime of a technological civilization.

Almost all authors agree on the general meanings of various $`f`$-parameters and $`n_e`$ (though wildly differing in the values they ascribe to each of them!); on the other hand, the product $`R_{}L`$ is sometimes written in the form

$$R_{}L=n_{}\frac{L}{t_0},$$
(2)
where $`n_{}`$ is the current number of stars in the Galaxy, and $`t_0`$ is the age of our stellar system (currently thought to be $`t_012`$ Gyr; e.g. Krauss and Chaboyer 2003). This is useful since (i) $`R_{}`$ is not a directly measurable quantity while $`n_{}`$ and $`t_0`$ are, at least in principle, and (ii) it enables direct comparison of two characteristic timescales, cosmological ($`t_0`$) and ”astrosociological” ($`L`$). There is a catch in (2), however, since the star formation rate is not uniform throughout the history of the Galaxy, and thus in general $`R_{}n_{}/t_0`$. While this particular problem is not acute from the SETI point of view, due to the metallicity effects (early epochs of intense star formation are characterized by low metallicity), it points in the direction of similar difficulties following from unwarranted assumptions of uniformity. We argue below that the main shortcoming of the Drake equation is its lack of temporal structure, i.e., it fails to take into account various evolutionary processes which form a pre-requisite for anything quantified by $`f`$-factors and $`n_e`$.

It is important to understand that we are criticizing the Drake equation not as an expression per se, but as a guiding line for a rather specific set of programs, procedures and (in the final analysis) investments, known overall as SETI. In SETI research proposals, Eq. (1) figures very prominently. Both supporters and opponents of SETI invoke the same simple numerical relationship in order to promote their respective views. However, arguments for both side are to be seriously doubted if the underlying relationship has serious deficiencies for any practical application, e.g. for estimating the timescale for sustained SETI effort by which we might expect to detect extraterrestrial intelligent signals (or artifacts).

Probably the most telling argument as to the inadequacy of the Eq. (1) is that almost four decades of SETI projects have not given results, in spite of the prevailing ”contact optimism” of 1960s and 1970s, motivated largely by uncritical acceptance of the Drake equation. Conventional estimates of that period spoke about $`10^610^9`$ (!) advanced societies in the Milky Way forming the ”Galactic Club” (Bracewell 1975). Nowadays, even the greatest SETI optimists have abandoned such fanciful numbers, and settled to a view that advanced extraterrestrial societies are much rarer than previously thought. One of the important factors in this downsizing of SETI expectations has certainly been demonstrations by ”contact pessimists”, especially Michael Hart and Frank Tipler, that the colonization—or at least visit—of all stellar systems in the Milky Way by means of self-reproducing von Neumann probes is feasible within a minuscule fraction of the Galactic age (Hart 1975; Tipler 1980). In this light, Fermi’s legendary question: Where are they? becomes disturbingly pertinent. In addition, Carter (1983) suggested an independent and powerful anthropic argument for the uniqueness of intelligent life on Earth in the Galactic context. Today, it is generally recognized that ”contact pessimists” have a strong position. How then, one is tempted to ask, does the discrepancy with our best analyses of Eq. (1) arise?

Now we shall show that there are two serious problems which make the equation (1) much less practical from the SETI point of view than has been conventionally thought. The two have the opposite effect on $`N`$, and may well partially cancel one another out; still, by a careful consideration those effects could be disentangled. Some other criticisms of the Drake equation, from different points of view, can be found in Walters et al. (1980), Wilson (1984), Ward and Brownlee (2000), and Walker and Ćirković (2003, preprint).

## 2 Upper limit on civilization’s age

In principle, the parameter $`L`$ in Eq. (1) could be arbitrarily large, thus offsetting any exceptionally small value among different $`f`$-parameters. Historically, that was the conventional assumption of ”contact optimists” like Sagan, Shklovskii, or Drake in the earlier decades (1960s and 1970s) of SETI efforts. It is reasonable to assume that after a technological civilization overcomes its ”childhood troubles” (like the threat of destruction in a nuclear war or through the misuse of nanotechnology) and starts colonizing space, it has very bright prospects for survival on timescales of millions (if not billions) years. Since it was intuitively clear (although quantified only recently; see below) that most of the inhabitable planets in the Milky Way are older than Earth, it was hypothesized that civilizations to be found through SETI projects will be significantly older than our civilization. However, it is a leap of faith from a reasonable estimate of the temporal distribution of civilizations to the assumption that we would be able to communicate with them, or that they would express any interest in communicating with us using our primitive communication means. Even worse, people of the ”contact optimism” camp have been expressing hope that we would be able to intercept communications between such very advanced societies, which seems still less plausible.111For a profound and poignantly ironic literary account of these issues see Lem (1984, 1987).

Obviously, from (1) we have $`lim_L\mathrm{}N=\mathrm{}`$, which is senseless, for the finite spatial and temporal region of spacetime we are considering in practical SETI. And still, remarkably, it is not senseless to contemplate upon the possibility that very advanced civilizations can exist indefinitely in an open universe (e.g. Dyson 1979), i.e. that the limit $`L\mathrm{}`$ makes sense. Whether an advanced technological society can exist indefinitely—in accordance with the so-called Final Anthropic Principle of Barrow and Tipler (1986) or the Final Anthropic Hypothesis of Ćirković and Bostrom (2000)—is still an open question in the nascent astrophysical discipline of physical eschatology (Adams and Laughlin 1997; Ćirković 2003). Any results from it, albeit very exciting and interesting in their own right, are unimportant to SETI due to the large disparity of the timescales involved.

According to a recent study by Lineweaver (2001), Earth-like planets around other stars in the Galactic habitable zone are, on average, 1.8 ($`\pm 0.9`$) Gyr older than our planet. His calculations are based on chemical enrichment as the basic precondition for the existence of terrestrial planets. Applying the Copernican assumption naively, we would expect that correspondingly complex lifeforms on those others to be, on average, 1.8 Gyr older. Intelligent societies, therefore, should also be older than ours by the same amount. In fact, the situation is even worse, since this is just the average value, and it is reasonable to assume that there will be, somewhere in the Galaxy, an inhabitable planet (say) 3 Gyr older than Earth. Since the set of intelligent societies is likely to be dominated by the small number of oldest and most advanced members (for an ingenious discussion in somewhat different context, see Olum 2003), we are likely to encounter a civilization actually more ancient than 1.8 Gyr (and probably significantly more).

It seems preposterous to even remotely contemplate any possibility of communication between us and Gyr-older supercivilizations. It is enough to remember that 1 Gyr ago, the appearance of even the simplest multicellular creatures on Earth lay in the distant future. Thus, the set of the civilizations interesting from the point of view of SETI is not open in the temporal sense, but instead forms a ”communication window”, which begins at the moment the required technology is developed (factor $`f_c`$ in the Drake equation), and is terminated either through extinction of the civilization or through it passing into the realm of ”supercivilizations” unreachable by our primitive SETI means. Formally, this could be quantified by adding a term to the Drake equation corresponding to the ratio of the duration of the ”communication window” and $`L`$. Let us call this ratio $`\xi `$; we are, thus, justified in substituting $`L`$ in Eq. (1) with $`\xi L`$. Since $`\xi `$ is by definition smaller than unity (and perhaps much smaller, if the present human advances in communication are taken as a yardstick), the net effect would be to drastically reduce the value of $`N`$. Fortunately (from the SETI point of view) this is not the only evolutionary effect hidden in the Drake equation.

## 3 Simplicity of uniformitarianism

A still more important shortcoming of Eq. (1) as a guideline to SETI consists of its uniform treatment of the physical and chemical history of our Galaxy. It is tacitly assumed that the history of the Galaxy is uniform with respect to the emergence and capacities of technological societies. This is particularly clear from the form (2), as mentioned above. If, on the contrary, we assume more or less sharply bounded temporal phases of the Galactic history as far as individual terms in Eq. (1) are concerned, and take into account our own existence at this particular epoch of this history, we are likely to significantly underestimate the value of $`N`$. We shall consider such a toy model below.

Uniformitarianism has not shone as a brilliant guiding principle in astrophysics and cosmology. It is well-known, for instance, how the strictly uniformitarian (and from many points of view methodologically superior) steady-state theory of the universe of Bondi and Gold (1948) and Hoyle (1948) has, after the ”great controversy” of 1950s and early 1960s, succumbed to the rival evolutionary models, now known as the standard (”Big Bang”) cosmology (Kragh 1996). Balashov (1994) has especially stressed this aspect of the controversy by showing how deeply justified was the introduction of events and epochs never seen or experienced by the Big Bang cosmologists. Similar arguments are applicable in the nascent discipline of astrobiology, which might be considered to be in an analogous state today as cosmology was half a century ago.

The arguments of Lineweaver (2001) are crucial in this regard, too. Obviously, the history of the Galaxy divides into at least two periods (or phases): before and after sufficient metallicity for the formation of Earth-like planets has been built up by global chemical evolution. But this reflects only the most fundamental division. It is entirely plausible that the history of the Galaxy is divided still finer into several distinct periods with radically different conditions for life. In that case, only weighted relative durations are relevant, not the overall age.

Exactly such a picture is presented by a class of phase transition models (Clarke 1981; Annis 1999; see also Norris 2000), which assume a global regulation mechanism for preventing the formation of complex life forms and technological societies early in the history of the Galaxy. Such a global mechanism could have the physical form of $`\gamma `$-ray bursts, if it can be shown that they exhibit sufficient lethality to cause mass biological extinctions over a large part of the volume of the Galactic habitable zone (Scalo and Wheeler 2002). If, as maintained in these models, continuous habitability is just a myth, the validity of the Drake equation (and the spirit in which it was constructed and used) is seriously undermined.

For illustration, let us assume that the parameter $`f_l`$ has the following evolutionary behavior:

$$f_l=\{\begin{array}{cc}\hfill 10^6:& 0<tt_p\hfill \\ \hfill 0.9:& t_p<tt_0\hfill \end{array}.$$
(3)
(we put the zero of time at the epoch of the Milky Way formation). Here, $`t_p`$ is the epoch of global ”phase transition” (Annis 1999), i.e. the epoch in which the lethal Galactic processes became rare enough for sufficiently complex lifeforms to emerge. Let us take $`t_0=12`$ Gyr and $`t_p=11`$ Gyr. Naive uniformitarian application of the Drake equation would require us to find the average $`f_l`$ in particular example $`f_l=0.072`$; if we assume $`n_e=1`$, other $`f`$-factors all equal to 0.1 (rather conservative assumption), and $`R_{}=5`$ yr-1, we obtain that $`N=3.6\times 10^5\xi L`$, where $`L`$ is measured in years, and $`\xi `$ is the relative duration of the communication window discussed above. In fact, the true result is rather $`N=4.5\times 10^4\xi L`$, more than an order of magnitude higher. Such a big difference is of obvious relevance to SETI; if $`\xi L10^5`$ yrs or less, it might as well be the difference between sense and nonsense in the entire endeavor. The discrepancy increases if the epoch of the phase transition moves closer to the present time. The latter is desirable if one wishes to efficiently resolve Fermi’s ”paradox” through phase-transition models.

More realistically, we would expect several of the $`f`$-factors, as well as $`n_e`$, to exhibit secular increase during the course of Galactic history in a more complicated manner to be elaborated by future detailed astrobiological models. Yet, steps similar to the one in (3) seem inescapable at some point if we wish to retain the essence of the phase transition idea. Barring this, the only fully consistent and meaningful idea for both explanation of the ”Great Silence” and retaining the Copernican assumption on Earth’s non-special position is the ”Interdict Hypothesis” of Fogg (1987), as the generalized ”Zoo Hypothesis” (Ball 1973), which still seems inferior, since it explicitly invokes non-physical (e.g., sociological) elements.

Intuitively, it is clear that in such phase transition models it is a very sensible policy for humanity to engage in serious SETI efforts: we expect practically all ETI societies to be roughly of the same age as ours, and to be our competitors for Hart-Tiplerian colonization of the Milky Way. The price to be paid for bringing the arguments of ”optimists” and ”pessimists” into accord is, obviously, the assumption that we are living in a rather special epoch in Galactic history—i.e. the epoch of ”phase transition”. That such an assumption is entirely justifiable in the astrobiological context will be argued in a subsequent study. (Parenthetically, this is entirely in accord with the tenets of the currently much-discussed ”rare Earth” hypothesis; see Ward and Brownlee 2000.)

Note that in this case, the overall average age of a civilization ($`L`$) would give an entirely false picture at the outcome of the Drake equation. In the toy model above, any hypothetical civilization age of (say) 10 Gyr is obviously irrelevant (although possibly sociologically allowed). This conclusion is valid even if the width of the communication window is very large, and in fact spans most of the lifetime of a civilization, as SETI pioneers claimed ($`\xi 1`$). Thus, we obtain a physically more desirable picture for the explanation of Fermi’s ”paradox” in which sociological influences are much less relevant.

## 4 Conclusions

We conclude that the Drake equation, as conventionally presented, is not the best guide for both operational SETI and future policy-making in this field. The reason for this is its lack of temporal structure and appreciation of the importance of evolutionary effects, so pertinent in the modern astrobiological discourse. If we wish to go beyond the ”zeroth-order” approximation encapsulated by (1), we will need to account for evolutionary effects, such as metallicity build-up and ”catastrophic” regulation of habitability. Notably, the non-uniform history of the Galaxy as conceived in the phase transition models can accommodate both the arguments of ”contact pessimists” and the justification for SETI projects, which have been deemed incompatible in the literature so far. Future detailed modeling will show in which way we can best accommodate our knowledge on the history of the Galaxy in the overall astrobiological picture.

Acknowledgements. The author wholeheartedly thanks Alan Robertson, Saša Nedeljković, Vesna Milošević-Zdjelar, Srdjan Samurović, and Mark A. Walker for invaluable technical help. Useful discussions with Nick Bostrom, Robert J. Bradbury, Zoran Živković, Ivana Dragićević, and David Brin are also acknowledged.

## References

Adams, F. C. and Laughlin, G. 1997, Rev. Mod. Phys. 69, 337.

Annis, J. 1999, J. Brit. Interplan. Soc. 52, 19 (preprint astro-ph/9901322).

Balashov, Yu. V. 1994, Studies in History and Philosophy of Science 25B, 933.

Ball, J. A. 1973, Icarus 19, 347.

Barrow, J. D. and Tipler, F. J. 1986, The Anthropic Cosmological Principle (Oxford University Press, New York).

Bondi, H. and Gold, T. 1948, Mon. Not. R. astr. Soc. 108, 252.

Bracewell, R. N. 1975 The Galactic Club: Intelligent Life in Outer Space (W. H. Freeman, San Francisco).

Carter, B. 1983, Philos. Trans. R. Soc. London A 310, 347.

Ćirković, M. M. 2003, Am. J. Phys. 71, 122.

Ćirković, M. M. and Bostrom, N. 2000 Astrophys. Space Sci. 274, 675.

Clarke, J. N. 1981, Icarus 46, 94.

Dyson, F. 1979, Rev. Mod. Phys. 51, 447.

Fogg, M. J. 1987, Icarus 69, 370.

Hart, M. H. 1975, Q. Jl. R. astr. Soc. 16, 128.

Hoyle, F. 1948, Mon. Not. R. astr. Soc. 108, 372.

Kragh, H. 1996, Cosmology and Controversy (Princeton University Press, Princeton).

Krauss, L. M. and Chaboyer, B. 2003, Science 299, 65.

Lem, S. 1984, His Master’s Voice (Harvest Books, Fort Washington).

Lem, S. 1987, Fiasco (Harcourt, New York).

Lineweaver, C. H. 2001, Icarus 151, 307.

Norris, R. P. 2000, in When SETI Succeeds: The impact of high-information Contact, ed. A. Tough (Foundation for the Future, Washington DC).

Olum, K. D. 2003, ”Conflict between anthropic reasoning and observation” (preprint gr-qc/0303070).

Scalo, J. and Wheeler, J. C. 2002, Astrophys. J. 566, 723.

Shklovskii, I. S. and Sagan, C. 1966, Intelligent Life in the Universe (Holden-Day, San Francisco).

Tipler, F. J. 1980, Q. Jl. R. astr. Soc. 21, 267.

Walker, M. A. and Ćirković, M. M. 2003, ”The Fermi Equation,” preprint.

Walters, C., Hoover, R. A., and Kotra, R. K. 1980, Icarus 41, 193-197.

Ward, P. D. and Brownlee, D. 2000, Rare Earth: Why Complex Life Is Uncommon in the Universe (Springer, New York).

Wilson, T. L. 1984, Q. Jl. R. astr. Soc. 25, 435-448.


text:
# Jet and Di-Jet production in Photon-Photon collisions

## 1 Introduction

During the past two years different analyses  of jet and di-jet production from real photons scattering at LEP have been presented. A unique computation it is available in the literature  at NLO for this process. This computation is based on the slicing method to cancel analytically virtual and infrared111We use the word “infrared” both for the properly soft effects and for the collinear ones. divergences. This technique relies on an approximate treatment of the phase space and of the matrix elements, which induces cancellations between large numbers. Hence in particular regions of the phase space it is possible that these cancellations are difficult to achieve numerically. To avoid this problem (it could be substantial for specific exclusive observables which have been measured ), it seems appropriate to use an exact technique to handle the infrared divergences. We use the subtraction method in the formulation of ref. . This method treats exactly both the phase space and the matrix elements, thus ensuring a better behaved numerical result. We present a complete computation of jet and di-jet cross-sections and we compare our theoretical predictions with previous results and with experimental data.

## 2 Some remarks about subtraction method

To compute a generic physical observable in hadronic physics it is necessary to implement the factorization theorem in a computer code. Thus we have to deal with unphysical large oscillations (due to the cancellation of the infrared divergences following the factorization prescription) and the computer precision. The solution of this problem it is to remove divergences analitycally and then to use a numerical method (typically a MC algorithm) to integrate the finite reminders. In this section we briefly describe a simplified 1-D model  to outline the differences between the subtraction and the slicing method.

### 2.1 Application of the slicing and the subtraction method to a simple case

Let us consider the integral,

$`F={\displaystyle \frac{1}{2ϵ}}{\displaystyle _0^1}𝑑x\delta (1x)F(x)+`$
$`{\displaystyle _0^1}𝑑x(1x)^{12ϵ}F(x),`$ (1)
in the limit of $`ϵ0`$, where $`F`$ is any regular function. This is a simplified situation of what happens to cancel out the soft emission (the second term in (1)) and the virtual correction (the first term) divergent terms when we compute a generic infrared safe observable $`F`$; $`ϵ`$ represents the (dimensional) infrared regulator. The first integral is trivial, the second one is not; in general, we can not compute it analytically. To remove the $`ϵ`$ pole we need, at least, to express explicitly the divergent part of the second integral. If we use the slicing strategy we rewrite the second term of (1) as follows:

$`{\displaystyle _0^1}𝑑x{\displaystyle \frac{F(x)}{(1x)^{1+2ϵ}}}={\displaystyle _0^{1\delta }}𝑑x{\displaystyle \frac{F(x)}{(1x)^{1+2ϵ}}}+`$
$`{\displaystyle _{1\delta }^1}𝑑x{\displaystyle \frac{F(x)}{(1x)^{1+2ϵ}}},\text{ }0<\delta 1.`$ (2)
It is clear that the only divergent term is in the second piece on the RHS, and if we perform a Taylor expansion around $`\delta =0`$, we obtain the $`ϵ`$ pole we need. Thus singular terms in $`ϵ`$ in Eq. (1) cancel out, and we obtain

$`F={\displaystyle _0^{1\delta }}𝑑x{\displaystyle \frac{F(x)}{1x}}+F(1)\mathrm{log}\delta +O(\delta ,ϵ).`$ (3)
This formula is obviously finite if we set $`ϵ=0`$ whatever $`\delta `$ we choose, thus we can perform the numerical integration. Formally the $`\delta `$ dependence of (3) is due only to the Taylor approximation we have done and it is of first order because the log$`\delta `$ term is canceled by an analogous factor coming from the upper integration limit in the (3). Then it is clear that the smaller $`\delta `$, the less accurate the numerical cancellations.
On the other hand if we use the subtraction technique we can rewrite (2) as

$`{\displaystyle _0^1}𝑑x{\displaystyle \frac{F(x)F(1)\theta (x1+x_c)}{(1x)^{(1+2ϵ)}}}+`$
$`F(1){\displaystyle _0^1}𝑑x{\displaystyle \frac{\theta (x1+x_c)}{(1x)^{(1+2ϵ)}}},\text{ }0<x_c1.`$ (4)
This time we have added and removed the same term, essentially the singular behavior of the integrand, in such a way that the divergent piece is isolated in the second factor as before. But this time we can perform exactly the integration of the divergent part. After removing the $`ϵ`$ pole and we obtain

$`F={\displaystyle _0^1}𝑑x{\displaystyle \frac{F(x)F(1)\theta (x1+x_c)}{1x}}+`$
$`F(1)\mathrm{log}x_c.`$ (5)
The logarithm term the equation (5) is formally identical to the one in (3). However, here no large number cancellations are present because we are not forced to set $`x_c1`$. Thus the results from (5) are exact and free of large fluctuations, but a good numerical integration is still necessary to obtain accurate results.

## 3 Structure of the computation

As is well known, in real photon-photon scattering photons can interact directly or hadronically. Thus the interaction is conventionally divided in three classes of Feynman diagrams: double resolved class, where the two photons interact by their hadronic components, the single resolved class, where only one photon interact by its hadronic component, and the direct class, where the photons themselves are the initial states of the hard process. As we will explicitly show each contribution has no physical meaning but only the sum correspond to an observable. For the first and second class we have modified the code  used to compute photon-hadron jet observables. As far as direct contribution is concerned we have recalculated and checked with  all the relevant matrix elements, which appear in the computation at NLO ($`\alpha ^2\alpha _S`$ order). Finally we have implemented the third class of diagrams in a computer code, analogue to the previous two codes. By summing the results from the three codes, we can obtain, using the subtraction method, any infrared finite jet and di-jet observable at NLO in real photon collisions. The details of the computation can be found in .

### 3.1 Consistency checks

We start by presenting the comparison of our theoretical predictions with the LEP OPAL data . In the lower panel of Figure 1 the differential inclusive one jet cross-section is shown. As one can see our computation (we use in Figure 1 $`\overline{\mathrm{MS}}`$ scheme, $`\mu _0=\mu _R=\mu _F=E_T^{Tot}/2`$,222Note that the half of total transverse energy of the event ($`E_T^{Tot}`$) is equal to $`E_T`$-jet at LO. GRV-HO photon density functions, $`\mathrm{\Lambda }_{QCD}=130`$ MeV; see  for the jet definition, for the cuts and for the kinematic) is in decent agreement with experimental data. In the central panel the ratios of the double resolved (dots), single resolved (dashed), direct (solid) contributions and of the sum are plotted. Our computation agrees with the previous one . Let us study the consistency of our computation, for example by changing the factorization scale ($`\mu _F=2\mu _0,0.5\mu _0`$). In fact if we change the factorization scale, upper panel of Figure 1, we obtain the same theoretical prediction (circles) but the double, single resolved and direct contributions change. Because each class of Feynman diagrams is separately divergent while the sum is finite, in order to perform the computation we had subtracted divergences analytically from each class, adding (or removing) an infinite term. Such a procedure depends on the regularization we use, i.e. on the factorization scale and scheme. Hence changing regularization the partial results change, the total does not333Up to uncontrolled higher-order terms, which appear to be negligible in this case.. In conclusion Figure 1 represents a good consistency check of our computation and shows that the separation of Feynman graphs at NLO is physically meaningless.
To distinguish between the two types of interactions of the photons, it is customary to introduce the couple of variables

$$x_\gamma ^\pm =\frac{_{Jet=1,2}(E^J\pm p_z^J)}{_{Hadrons}(E^H\pm p_z^H)},$$
(6)
where the sum in the numerator runs over the two hardest jets of the event. The meaning of the definition (6) it is clear: in fact if we consider the Born process in the direct case, we have only two back-to-back jets and no hadronic remnants, which implies $`x_\gamma ^\pm =1`$. Thus generally the closer $`x_\gamma ^\pm `$ to one, the more the cross section is dominated by the direct contribution. Elsewhere the process is dominated by resolved interactions. We are able to distinguish the behavior of each photon (coming from the electron or the positron), and then we can study the single resolved case too. It is clear that we need a $`x_\gamma ^{sep}`$ value to separate (in an arbitrary manner) between the hadronic and the direct contribution: we will use $`x_\gamma ^{sep}=0.75`$, the same value used by OPAL, suggested by MonteCarlo simulations . Therefore we consider three regions on the $`x_\gamma ^\pm `$-plane, corresponding to the three types of interactions mentioned before. We expect that the cross-section in the “hadronic” region to be the most sensitive to the PDF choice. Indeed this is the case. Moreover we can predict the differential cross-sections with respect to $`x_\gamma ^\pm `$. To this aim it is interesting to check what happens if we change again the factorization scale. In Figure 2 we plot $`d\sigma /dx_\gamma ^\pm `$ for three different scales in the region $`x_\gamma ^\pm >0.75`$. For $`x_\gamma ^\pm <0.80`$ the curves are indistinguishable, but it is evident that close to one we obtain unstable results, strongly depending on the factorization scale choice.

## 4 IR-Sensitive regions

We can deduce physical observations from these consistency tests. In fact if we look at the hadronization corrections444Here, the hadronization correction is defined as the ratio between the MC cross-sections at the parton and hadron levels. taken from  to $`d\sigma /dx_\gamma ^\pm `$ in the same region we observe that they are between 40% and 260%. Hence the conclusion is that the perturbative computation at NLO fails in this region. To understand this fact let us consider the relationship between $`x_\gamma ^\pm `$ and the fraction of energy $`(1x_\gamma ^{1,2})`$ of the emitter parton carried out by a soft gluon emitted in the final state of the direct process. It is simple to show that in a three-jet configuration

$$x_\gamma ^{1,2}1x_\gamma ^\pm 1.$$
(7)
Thus the “dangerous region” showed in Figure 2 corresponds to a soft gluon emission kinematic, where the exponentiation and the resummation to all order of the large logarithms are necessary to obtain sensible predictions. In particular the dominant term of the perturbative differential cross-section in the limit $`x_\gamma ^\pm 1`$ is roughly

$$\frac{d\sigma }{dx_\gamma ^\pm }\frac{1}{1x_\gamma ^\pm }\frac{\theta (x_\gamma ^\pm 1)}{1x_\gamma ^\pm }.$$
(8)
The expression (8) explicitly indicates the finiteness (i.e. the infrared safeness) of the cross-section we expect from jet definitions, but because infrared divergence is canceled badly (exactly for $`x_\gamma ^\pm =1`$) large logarithms appear. Thus the comparison of the data in this region with NLO predictions cannot be satisfactory. In other words one has to renounce to compare the exclusive observable $`d\sigma /dx_\gamma ^\pm `$ in this region with fixed order computations. To avoid this problem one can use more inclusive quantities , obtained by integrating $`d\sigma /dx_\gamma ^\pm `$ over $`x_\gamma ^\pm 0.8`$. This improves the agreement between data and predictions.
We also point out that there is another region in the phase space which can introduce large logarithmic corrections to NLO computations. It is known that di-jet observables are sensible to the choice of the energy cuts on the phase space, in particular it is shown  that symmetric cuts over the two jets transverse energy imply large logarithmic corrections to certain observables. Thus to avoid this effect OPAL collaboration has chosen  asymmetric cuts on the $`E_T^{Jet}`$-plane defined as

$`\overline{E}_T=(E_T^1+E_T^2)/2\text{ 5 GeV,}`$
$`0\varphi ={\displaystyle \frac{\left|E_T^1E_T^2\right|}{E_T^1+E_T^2}}<{\displaystyle \frac{1}{4}}.`$ (9)
Even if these cuts are free from the troubles we cited before, we point out that other problems could arise for NLO predictions. Indeed, NLO computation imposes $`\varphi 1/3`$, and setting the $`\varphi `$-limit near this value could be problematic. Furthermore if the $`\varphi `$-limit goes to zero infrared effects arise. In fact one should expect that the cross-section goes to zero if the phase space becomes small, but as Figure 3 shows this is not the case at NLO.

The explanation of this is simple: if $`\varphi 0`$ then there is not enough phase space to emit a soft gluon to cancel out properly virtual divergences. Therefore large logarithms of $`\varphi `$  become dominant when $`\varphi 0`$. Thus even if in principle OPAL choice  could be dangerous in practice it is far from the “dangerous regions”.

## 5 Phenomenological results

In this section we compare our results to recent data. In Figure 4 our theoretical predictions are compared to OPAL measurements  (see that paper for the kinematics and the cuts) for di-jet cross-sections are plotted . The three series of curves correspond respectively to the whole $`x_\gamma ^\pm `$ region, either $`x_\gamma ^+>0.75`$ or $`x_\gamma ^{}>0.75`$, and $`x_\gamma ^\pm <0.75`$. We use for the central theoretical prediction GRV-HO in DIS scheme with $`\mu _0=\mu _R=\mu _F=E_T^{Tot}/2`$ and $`\mathrm{\Lambda }_{QCD}=130`$ MeV. The NLO predictions, for all the observables we present, are divided by the factor $`(1+\delta _H)`$, the hadronization corrections taken from . These terms are less than 10% for $`\overline{E}_T`$ grater than 10 GeV, but rise to about 20% towards small $`\overline{E}_T`$. Plots show a very good agreement between QCD predictions and data for the most inclusive observable and in the single resolved region, whereas the agreement does not seem satisfactory in the double resolved region. Even if we consider the theoretical errors due to the scale and the PDF choice or $`\mathrm{\Lambda }_{QCD}`$, as shown in Figure 4, we are not able to adsorb discrepancy in the theoretical uncertainties. Furthermore the situation is worse if we consider $`d\sigma /dx_\gamma ^\pm `$ in this region, NLO prediction cannot satisfactory describe the experimental data. The same discrepancies affect the comparison between data and the previous theoretical computation . Thence we conclude that more experimental and theoretical studies are needed.

If we compare (Figure 5) $`d\sigma /dx_\gamma ^\pm `$ with the data for $`0x_\gamma ^\pm 1`$ and $`5<\overline{E}_T<7`$ GeV (where the “resolved” component is dominant) we have only a partial agreement. We note that the correspondence we find for $`x_\gamma ^\pm >0.8`$ is not particularly meaningful. In fact, as we have discussed in Section 3 the result is not stable under variation of the factorization scale. On the other hand in the energy range $`7<\overline{E}_T<11`$ GeV data agree well with the theoretical curves outside the “dangerous region”, i.e. for $`x_\gamma ^\pm <0.8`$. Other quantities measured by OPAL have been computed and good agreement has been found with the theoretical predictions.

Recently, an L3 analysis has been presented , for the inclusive single jet $`d\sigma /dE_T`$ on a very extended $`E_T`$ range ($`E_T^{max}`$ 60 GeV). In Figure 6 we present the comparison between our prediction and L3 preliminary data. As one can see it is possible to divide the momentum range in two regions. For momenta less than 15 GeV the theoretical curve has the right shape to fit the data but it has the wrong normalization (approximately it is a factor two higher). For large momenta data and theoretical prediction are clearly in disagreement. The discrepancy is at places of one order of magnitude, which cannot be accounted for by the theoretical uncertainties. In the upper panel of Figure 6 the scale dependence of the theoretical prediction (circles) is shown, and similarly to the case of Figure 1 the variation is very small. Dotted curve re-shows the scale dependence of the “hadronic” contribution.

## 6 Conclusion

We have implemented a computer code which allows to compute jet and di-jet observables at NLO using the subtraction method. Our results are in agreement with the previous theoretical predictions obtained with the slicing method . The comparison of our results with OPAL data  is satisfactory for inclusive observables, whereas the comparison with L3 data  display sizable disagreement in the large $`E_T`$ region.
Finally, we have analyzed some theoretical difficulties to describe particular exclusive observables measured in . The theoretical uncertainties for inclusive quantity are under control.