text:
# Modelling a Dune Field

## 1 Introduction

A wide variety of dune shapes can be found in deserts, in coastal areas, on the sea–bottom , and even on Mars  and over one hundred different types of dunes have been classified by geomorphologists and their ancient predecessors so that some of them bear old arabic names. The shapes depend mainly on the amount of available sand and on the change in the direction of the wind over the year. If the wind blows steadily from the same direction throughout the year and there is not enough sand to cover the entire area, dunes called barchans develop having the shape of a crescent. They typically appear in huge fields of several hundred kilometers length stretching along coastlines and driven by trade-winds. A picture from a section of such a field is shown in Fig. 1.

Measurements of the relationship between height, width and length have been performed for instance in Peru or Morocco  and led to a good understanding of the barchan shape. In particular we know that the barchans move without changing their shape with a velocity which is inversely proportional to their height. However, the dynamics of dunes in a dune field has not yet been studied very much. In particular, there seems to be a contradiction between the fact that the barchans are not all of equal height and the stability of the dune field since when the dunes have different velocities they will eventually bump into one another. We do not yet understand the underlying stabilizing mechanism. Looking at the images in more detail also other questions arise. Are the relative positions of the dunes correlated? Do they form patterns? When dunes are getting closer they can interact through the non–local wind field or the sand flux between the dunes. In this work, we propose a simple model that defines the interaction between the dunes caused by the inter–dune sand flux. We neglect further effects of the wind field. The model predicts the time evolution of the entire field. It can reproduce the spatial correlations between the dunes that have been observed in the field, cf. Fig. 1. Furthermore, the model shows that the entire dune field does not broaden but stays focused. However, the stability of the dune field depends on the sand flux balance of a single dune, which is not well known so that we will have to make some simplifying assumptions.

## 2 Wind and sand flux

Sand is transported if the wind speed exceeds a certain threshold velocity. In equilibrium the amount of transported sand or the flux of sand $`q(u_{})`$ can be expressed as a function of the shear velocity $`u_{}`$ . The shear velocity $`u_{}`$ characterizes the turbulent boundary layer of the atmosphere and thus its well known logarithmic velocity profile $`v(z)=u_{}/\kappa \mathrm{ln}z/z_0`$, where $`v(z)`$ is the wind speed at height $`z`$ over the ground, $`\kappa 0.4`$ the universal von Kármán constant for turbulent flow, and $`z_0`$ the roughness length of the surface. For shear velocities well above the threshold, e.g. $`u_{}=0.5`$m s-1, the simplest sand flux relation, proposed by Bagnold , should give a reasonable prediction for the sand flux $`q`$,

$$q=C_B\frac{\rho _a}{g}u_{}^3,$$
(1)
where $`C_B2`$ is a phenomenological constant , $`\rho _a=1.225`$kg m-3 the density of air, and $`g=9.81`$m s-1 the acceleration due to gravity. Using the bulk density $`\rho _s=1650`$kg m-3 of a typical sand we can rewrite Eq. (1) into a relation for the volume flux $`\mathrm{\Phi }`$,

$$\mathrm{\Phi }=\frac{q}{\rho _s}$$
(2)
For the discussion of the dune field dynamics we chose the shear velocity over the flat ground $`u_0=0.5`$m s-1 and obtain a volume sand flux $`\mathrm{\Phi }_0`$ of 0.031 m2s-1. Furthermore, we assume that the sand velocity $`v_s`$ is proportional to the shear velocity $`u_{}`$ and choose a typical value of $`v_s=5`$m s-1.

## 3 Barchan dunes

A dune is characterized by a single scalar quantity, e.g. its width $`w`$, and its position $`(x,y)`$ in the dune field. Further properties such as the height $`h`$,

$$h=aw$$
(3)
or volume $`V`$,

$$V=bw^3$$
(4)
can be calculated using simple phenomenological formulas. Here, the parameters $`a=0.09`$ and $`b=0.05`$ have been obtained from field measurements .

The migration velocity $`v_d`$ of a dune depends on its height $`h`$ and the sand flux $`\mathrm{\Phi }_b`$ over the brink of the dune . The latter can be calculated through Eq. (1) if the shear velocity at the brink is known. However, the dune itself disturbs the wind field such that it becomes stronger above the dune than over a flat ground. The ratio between the shear velocity over a flat ground $`u_0`$ and the shear velocity at the brink of the dune is called the “speed–up” $`n`$. This speed–up depends on the shape of the dune and mainly on the average windward slope. A typical value for the speed–up over a sand dune is 1.4 . Finally, we can write for the dune velocity $`v_d`$,

$$v_d=\frac{\mathrm{\Phi }_b}{h}=\frac{\mathrm{\Phi }_0n}{aw}.$$
(5)
We want to note that Eqs. (3) and (4) implicitly assume scale invariant dune shapes, which is an approximation that is only valid for large dunes, much larger than the minimal dune size .

The growth and shrinkage of a dune depends on the balance of in– and out–going sand fluxes. The sand flux that arrives from a upwind location, for instance a beach or other dunes, feeds the dune with new sand, whereas at the tip of the horns sand leaves the dune. The change of volume $`\dot{V}_i`$ due to the in–flux $`\mathrm{\Phi }(y)`$ can be calculated by integrating over the dunes width,

$$\dot{V}_i=_0^w\mathrm{\Phi }(y)𝑑y.$$
(6)
The sand that leaves the dune at the horns gives rise to a loss of volume $`\dot{V}_o`$ that is modeled by the following expression,

$$\dot{V}_o=w\mathrm{\Phi }_0(a_0+a_1w),$$
(7)
where $`a_0`$ and $`a_1`$ are parameters that depend on the complex processes occurring on the dune’s surface and are not known a priori. Recent numerical calculations of the sand outflux from barchan dunes suggest that the outflux is nearly proportional to the dunes width and motivated the expression above .

For a stable dune, the gain and loss of volume must exactly compensate $`\dot{V}_i\dot{V}_o=0`$. The average inter–dune flux or influx $`\mathrm{\Phi }=f\mathrm{\Phi }_0`$ causes an average increase in volume,

$$\dot{V}_i=wf\mathrm{\Phi }_0,$$
(8)
where $`f`$ defines the saturation of the inter–dune flux. Eqs. (7) and (8) define finally an average dune width $`w`$,

$$w=\frac{fa_0}{a_1}.$$
(9)
Using the average dune width $`w`$ and $`\gamma =a_1w/f`$ as new parameters, instead of $`a_0`$ and $`a_1`$, we can rewrite Eq. (7),

$$\dot{V}_o=wf\mathrm{\Phi }_0\left(1+\gamma \frac{ww}{w}\right).$$
(10)
Finally, the change in volume $`\dot{V}`$ of a dune can be obtained from Eqs. (6) and (10),

$$\dot{V}=_0^w\mathrm{\Phi }(y)𝑑ywf\mathrm{\Phi }_0\left(1+\gamma \frac{ww}{w}\right)$$
(11)
Typical values for $`f`$ are of the order of 0.1, i.e. the inter–dune flux is far from being saturated.

## 4 The barchan field

A barchan field of length $`l_f`$ and width $`w_f`$ as depicted in Figure 3 is composed of many dunes that migrate in the $`x`$–direction, the direction of the wind. The center of mass ($`x_i`$, $`y_i`$) of a dune is used to localize the dune $`i`$ in the field. The dynamics of the dune field is described by the temporal change of the variables $`x_i`$, $`y_i`$, and $`w_i`$ which will be defined in the following.

From Eqs. (11) and (4) we obtain the change of the width of a dune, which depends on the inter–dune flux,

$$\dot{w}_i=b^{1/3}\left[_{y_iw_i/2}^{y_i+w_i/2}\mathrm{\Phi }(y)𝑑yw_if\mathrm{\Phi }_0\left(1+\gamma \frac{w_iw}{w}\right)\right]^{1/3}$$
(12)
The velocity of the dune in the wind direction is given in Eq. (5) and depends only on the width of the dune,

$$\dot{x}_i=v_d(w_i)=\frac{\mathrm{\Phi }_0n}{aw_i}$$
(13)
An inter–dune flux that is asymmetric in the lateral direction with respect to the center of mass $`y_i`$ of the dune gives rise to a change in lateral position $`y_i`$,

$$\dot{y}_i=\frac{V_iy_i+_{y_iw_i/2}^{y_i+w_i/2}\stackrel{~}{y}\mathrm{\Phi }(\stackrel{~}{y})𝑑\stackrel{~}{y}}{V_i+_{y_iw_i/2}^{y_i+w_i/2}\mathrm{\Phi }(\stackrel{~}{y})𝑑\stackrel{~}{y}}y_i,$$
(14)
where $`V_i`$ is the volume of the dune $`i`$.

Dunes can interact either due to the inter–dune flux or due to coalescence. The first case is defined by the equations above and the inter–dune sand flux $`\mathrm{\Phi }(x,y)`$ defined in the next section. For the latter case we define a characteristic distance of coalescence $`d_{\mathrm{col}}`$. If two dunes become closer than this distance they are merged to a single dune. The mass is conserved by this process and the new position is calculated according to their centers of mass and relative volumes.

### 4.1 The inter–dune sand flux in the field

The flux of sand in the field, i.e. between the dunes, is a complex problem and causes an interaction between the dunes. Dunes can act here as sources and sinks for the sand at the windward foot and at the horns, respectively. Furthermore, the position and width of the dunes changes in time and thus their influence on the inter–dune sand flux. We propose the following model for the inter–dune flux: First, the sand is convected with the constant velocity $`v_s`$ in the wind direction. Secondly, a lateral movement is caused by diffusion and obeys a diffusion equation,

$$\frac{\mathrm{\Phi }}{t}D\frac{^2\mathrm{\Phi }}{y^2}=0,$$
(15)
where the time scale of the diffusion process is connected to the transport in wind direction by the sand velocity $`v_s`$. Hence, we start solving the diffusion equation at the influx–border $`x=0`$ of the field, where the sand flux is spatially constant. When integrating the solution forward in time by displacing its $`x`$–position we obtain the stationary solution of the two–dimensional inter–dune flux field $`\mathrm{\Phi }(x,y)`$,

$$v_s\frac{\mathrm{\Phi }}{x}D\frac{^2\mathrm{\Phi }}{y^2}=0.$$
(16)
For each dune in the field, the flux $`\mathrm{\Phi }_i`$ is set to zero along the entire width of the dune (from $`y_iw_i/2`$ to $`y_i+w_i/2`$). Then, at the horns (i.w. the positions $`y_iw_i/2`$ and $`y_i+w_i/2`$) the outflux is imposed as two $`\delta `$-functions. Finally, the new dune mass is calculated by the in– and out–flux balance, Eq. (12).

### 4.2 Initial configuration and Boundary conditions

The field is initialized with dunes at random positions until a certain density $`\rho _d`$ is reached, which is a parameter of the model. The dune width is chosen between a minimum and a maximum width, $`w_{min}ww_{max}`$, having an average width $`w=(w_{min}+w_{max})/2`$. From the initial density $`\rho _d`$ we can define the rate $`r`$ with which new dunes enter the field at the in–flux boundary,

$$r=W\rho _dv(w),$$
(17)
where $`v(w)`$ is the dune velocity according to Eq. (5) and $`W`$ is the width of the dune field, which is chosen to be only one third of the total simulation area, cf. Fig 3 where $`w_f=3W`$. The width of the dunes that enter with the rate $`r`$ at the influx boundary is chosen randomly as in the initial configuration. If a dune crosses the outflow boundary (at the right side of the field) it is removed from the simulation.

The inter–dune flux of sand $`\mathrm{\Phi }`$ is chosen to be constant between $`w_f/3<y<2w_f/3`$ and zero outside. The constant influx is chosen to be $`\mathrm{\Phi }=f\mathrm{\Phi }_0`$ and thus corresponds to an inter–dune flux of a dune field with average dune width $`w`$, cf. Eq. (8). In the direction perpendicular to the wind ($`y`$–direction) periodic boundary conditions are used.

## 5 Numerical Scheme

The numerical scheme can be divided into several steps. First, the dune field is randomly initialized with dunes as discussed above. Then, an iterative calculation, forward in time, of the following steps is performed:

* New dunes are created with rate $`r`$ at the left boundary, in the region $`w_f/3<y<2w_f/3`$.
* The new $`x`$–position of the dunes is calculated according to Eq. (13).
* If two dunes are closer than the coalescence distance $`d_{\mathrm{col}}=50`$m they are merged.
* The inter–dune sand flux is calculated from left to right on a discrete one dimensional grid in $`y`$–direction according to Eq. (16). The iteration is started at $`x=0`$ with a constant value of $`\mathrm{\Phi }(y)=f\mathrm{\Phi }_0`$ for $`w_f/3<y<2w_f/3`$ and zero outside. When the calculation forward in the wind direction reaches a dune, the flux is set to zero along the dune width and two discrete $`\delta `$–peaks are set at the horns of the dune as described above. Furthermore, the change in volume of the dune and its lateral movement according to Eqs. (12) and (14) are calculated.

The time step between the iteration cycles above has been chosen to be $`\mathrm{\Delta }T=1/12`$ year.

## 6 Results

We are mainly interested in measuring how the dunes are distributed in the field, possible correlation between the positions of the dunes and understanding the confinement of the field on a stripe. We did simulations of our model for a field having a length of $`50`$ km and a width of $`12`$ km. The dunes enter the field only in the region between km $`4`$ and km $`8`$. The initial density of dunes in the field was chosen to be $`\rho _f=10`$ dunes/km2. We set the saturation of the flux to be $`f=0.15`$ which corresponds to a flux of $`\varphi =0.0046`$ m2/sec.

In fig. 4 we see a typical simulation of the dune field. It stays confined inside a stripe. The dune density slowly decreases in direction of the wind. This can be seen quantitatively in Fig.5 where the density averaged over the $`y`$ direction is plotted over the $`x`$-axis ($`x`$ and $`y`$ directions are defined in fig. 4) in a semilogarithmic plot. We recognize for large $`x`$ a logarithmic decay of the form $`\rho lnx`$.

The situation seen in fig. 4 is a stationary state. As seen in fig. 6 the average size of the dunes fluctuates around a value that is constant in time and except for a short transient at early times also the number of dunes is essentially constant in time.

### 6.1 The confinement of the field

It is not clear up to now why real dune fields are typically confined on a stripe all the time. In our simple model, despite the fact that the diffusion of the sand tends to spread the sand to the borders, this behavior appears by itself. As seen in fig. 7, the sand flux decreases outside the central region of the field. Hence, dunes at the borders of that region receive less sand. The ones close to the center of the field receive more sand. This flux gradient acts like a “force” which pulls the dunes into the region where there is more flux, i.e. the central region where the other dunes are.

A way to see this effect is the distribution of dunes and their mean size at the end of the field. This is shown in figure (8). In the inset we plotted the number of dunes that leave the field at the right boundary. Near to the center of the field the dunes have approximately the same size, however the ones near to the borders are smaller, since they receive less sand. Since the dunes in the borders are smaller, they are faster and this explains why we have more dunes leaving the field close to the borders than in the center.

### 6.2 Spatial correlations

Looking at fig. 1 we can observe that the dunes tend to form lines, the horn of one dune pointing to the center of the back of another one. These patterns are apparently due to the coupling through the inter-dune sand flux and we will show how they are reproduced by our model. In order to quantitatively measure these structures, we will introduce two types of correlation functions which we will monitor as well in our model as also directly from a photograph of a real dune field in order to compare the patterns formed by our model to the observed ones.

Let us define the longitudinal correlation as

$$C_l(d)=\underset{\mathrm{i}}{}\underset{j>i}{}\delta (d|x_ix_j|)$$
(18)
This function measures the correlation of dune positions in the wind direction.

In fig. 9 we present the results for our model. To get better statistics we take the average over several thousand time-steps. We observe very pronounced peaks at regular distances of about 50 m showing a characteristic distance between the dunes in wind direction.

It is also interesting to calculate the angular correlations in the positions of the dunes which can be quantified by the angular correlation function

$$C_p(\theta )=\underset{\mathrm{i}}{}\underset{j>i}{}\delta \left(\theta arctan\left(\frac{y_iy_j}{x_ix_j}\right)\right)$$
(19)
where $`i`$ and $`j`$ are only chosen within the central part of the stripe (about one third of all dunes) in order to avoid boundary effects and for which the distance $`(x_ix_j)^2+(y_i+y_j)^2<(5\mathrm{k}\mathrm{m})^2`$. In fig. 10 we see our results for $`C_p(\theta )`$ as obtained after averaging over all dunes in the central part of the field and several thousand time-steps.

We see that while most directions have no particular structure in the forward direction there are two favorite directions symmetrical around the exact wind direction ($`x`$-axis) of high correlation. In the wind direction itself the correlation is particularly low. The selected direction corresponds to the arrangements where the dune in front is just shifted by half a dune width, so that the horns directly point to the center of the dune in front. Such configurations can be recognized with the naked eye on the photo of fig. 1.

## 7 Results from the experiments

We analysed aerial photos of the dune field of Laâyounne in Southern Morocco and obtained the position of 734 dunes as shown in fig. 11. In figs. 12 and 13 we show our results for both correlation functions corresponding to the correlation functions obtained for the model in figs. 9 and 10. Comparing figs. 9 and 12 one sees that the measured longitudinal correlation does not show much periodicity or characteristic distance. We believe that this lack of structure is rather due to insufficient statistics. After all the model gives the possibility to the average over time while the measurements are obtained from a single photo.

Concerning the polar correlations, fig. 13 clearly shows the same two peaks pointing in both directions with an angle of $`\theta `$ from the wind direction. But while in the model $`\theta `$ was about 80, it is about 150 for the real data.

## 8 Conclusions

In this paper we introduced a simple model for a field of Barchan dunes. The dunes are characterized by their width and position in the field. The model includes a diffusion of the sand. Assuming a simple law for the inter-dune sand flux we could show that the simulated field is stably confined and presents the same kind of spatial correlations as real fields.

An extension of this work could be to introduce variations in the wind direction. This could help to explain more complex patterns observed in deserts.

## 9 Acknowledgement

ARL acknowledges a fellowship from CNPq (Brazilian Agency) and Gerd Sauermann a fellowship from DAAD.


format: text

text:
# On the cosmological variation of the fine structure constant

(November 1, 2001)

## Abstract

A phenomenological model is proposed to explain the recent observed cosmological variation of the fine structure constant as an effect of the quantum vacuum, assuming a flat universe with cosmological constant $`\mathrm{\Lambda }`$ in the cases $`(\mathrm{\Omega }_M,\mathrm{\Omega }_\mathrm{\Lambda })`$ equal to $`(0.3,0.7)`$ and $`(1,0)`$. Because of the fourth Heisenberg relation, the lifetime of the virtual pairs of the zero-point radiation must depend on the gravitational potential $`\mathrm{\Phi }`$, so that the quantum vacuum changes its density and acquires a relative permittivity different from one. Since the matter was more concentrated in the past, the gravitational potential of all the universe was stronger and the optical density of the vacuum higher, the electron charge being therefore more renormalized and smaller than now. The model is based on a first order Newtonian approximation that is valid for the range of the observations, but not for very high redshift, the prediction being that $`\mathrm{\Delta }\alpha /\alpha `$ is proportional to $`\{\mathrm{\Omega }_M[a(t)^11]2\mathrm{\Omega }_\mathrm{\Lambda }[a(t)^21]\}`$, $`a(t)`$ being the scale factor. This agrees with the observations.

PACS numbers: 06.20.Jr, 95.30Sf, 98.80.Es

Short title: Cosmological variation of $`\alpha `$

Introduction. The observations of the absorption lines of distant quasars by Webb et al seem to indicate that the fine structure constant was smaller in the past . The idea that the fundamental constants are actually changing goes back to the “large number hypothesis” by Dirac , what stimulated a number of proposals to formulate the variation of the constants, see for instance ref.  in which $`\alpha `$ is written in terms of a scalar field, or other theories in which the variation of $`\alpha `$ is induced by a variation of the light velocity  (see also ).

This letter proposes a way to understand this phenomenon as an effect of the quantum vacuum. More precisely, it is argued here that, because of the fourth Heisenberg relation, the density of the sea of virtual particles in the quantum vacuum must change in a gravitational field, with a corresponding variation of its permittivity and permeability. As a consequence the quantum vacuum was optically denser in the past, when the universe was more contracted than now and the gravitational potential due to all the universe was stronger. The observed electron charge $`e`$ is the result of the renormalization of its bare charge $`e_{\mathrm{bare}}(>e)`$ by the quantum vacuum , so that a denser vacuum implies a smaller observed charge. It would be desirable to study this question in the frame of a rigorous quantum field theory, but the quantum vacuum is not understood well enough today in order to do that (for instance it is not known why its energy density seems to be so small). The only available alternative (other than waiting for a future theory) is to try a phenomenological approach, as is done here.

In the model here proposed, the quantum vacuum is treated as a transparent optical medium characterized by its permittivity and permeability, the change of $`\alpha `$ being a consequence of the fourth Heisenberg relation applied to the gravitational interaction of the virtual pairs in the zero point radiation with all the universe. The analysis is nonrelativistic (although with the inclusion of the rest-mass energy $`mc^2`$) and only weak Newtonian gravitational potentials $`\mathrm{\Phi }`$ are considered for which $`|\mathrm{\Phi }|/c^21`$.

The fourth Heisenberg relation and the vacuum density. Traditionally, quantum physics has stated that the sea of virtual pairs that are created and destroyed constantly in the quantum vacuum, i. e. the zero-point energy, has infinite density, as follows from the simple application of its basic principles. However, there is now evidence that this density may be finite. On the average and phenomenologically, a virtual pair created with energy $`E`$ (including rest-mass energy, kinetic energy and electromagnetic energy) will live during a time $`\tau _0=\mathrm{}/E`$, according to the fourth Heisenberg relation. A constant number density of pairs is established in this way at a certain equilibrium between the number of particles created and destroyed per unit time. This has an important consequence: virtual pairs must live longer in the gravitational field created by a mass distribution, because they have an extra negative potential energy $`E\mathrm{\Phi }/c^2`$ there ($`\mathrm{\Phi }`$ being the Newtonian potential). Indeed their lifetime must be

$$\tau _\mathrm{\Phi }=\mathrm{}/(E+E\mathrm{\Phi }/c^2)=\tau _0/(1+\mathrm{\Phi }/c^2).$$
(1)
The consequence is clear if unexpected: the number density of pairs $`𝒩`$ depends on the gravitational potential $`\mathrm{\Phi }`$ as

$$𝒩_\mathrm{\Phi }=𝒩_0/(1+\mathrm{\Phi }/c^2).$$
(2)
In other words, the density of the quantum vacuum depends on the gravitational potential because it must have there an extra number density of pairs depending on $`\mathrm{\Phi }/c^2`$. If $`\mathrm{\Phi }<0`$, the quantum vacuum becomes denser, if $`\mathrm{\Phi }>0`$, it becomes thinner. These two cases correspond to the gravity created by the mass and the cosmological constant, respectively. However, this change of the number density can produce no gravitational effects, since the gravitational mass of the virtual pairs in a volume depends on $`\mathrm{\Phi }`$ through the product of their mass $`E/c^2`$ by the number density per unit volume and unit energy $`𝒩(E)`$. This product does not change since the second factor acquires a factor $`(1+\mathrm{\Phi }/c^2)^1`$ in a gravitational field, the first a factor $`(1+\mathrm{\Phi }/c^2)`$, as we have seen. The effect that we are considering is optical and dielectric, not gravitational. There is no contradiction with the usual Lorentz invariance of the quantum vacuum, since the presence of a mass or any non uniform gravitational field breaks the Lorentz symmetry. This variation of the virtual pairs lifetime in a gravitational field is not an ad hoc hypothesis, but an unavoidable consequence of the fourth Heisenberg relation. Note also that the virtual particles considered here are not created by gravity: they are just the usual zero point particles that fill the space everywhere, according to elementary quantum physics: they live a bit longer (or shorter), that is all.

In the following we will deal with the variation of the observed values of certain quantities between a spacetime point of the universe with potential $`\mathrm{\Phi }`$ and a terrestrial observatory. As they will be expressed at first order in the potential, $`(\mathrm{\Phi }\mathrm{\Phi }_{})`$ will be written instead of $`\mathrm{\Phi }`$, $`\mathrm{\Phi }_{}`$ being the present gravitational potential of all the universe here at Earth.

Vacuum permittivity and permeability in a gravitational potential. Following the previous considerations, we will admit as a phenomenological hypothesis that the quantum vacuum can be considered as a substratum, similar to an ordinary transparent optical medium and characterized by a permittivity and a permeability that depend on $`\mathrm{\Phi }`$. As $`\mathrm{\Phi }`$ becomes more negative, its density increases and the bare electron charge $`e_{\mathrm{bare}}`$ is renormalized further or, otherwise stated, the observed charge must become smaller. In such a view the permittivity and the permeability of the quantum vacuum must change to $`ϵ_\mathrm{r}ϵ_0`$ and $`\mu _\mathrm{r}\mu _0`$, the first factor expressing the effect of its thickening (or lightening). As we assume a weak field, we can express the relative permittivity and permeability at first order in the potential. Since both are equal to 1 now at Earth, their value at a spacetime point with potential $`\mathrm{\Phi }`$ can be expressed as:

$$ϵ_\mathrm{r}=1\beta (\mathrm{\Phi }\mathrm{\Phi }_{})/c^2,\mu _\mathrm{r}=1\gamma (\mathrm{\Phi }\mathrm{\Phi }_{})/c^2,$$
(3)
$`\beta `$ and $`\gamma `$ being certain coefficients, which must be positive since the quantum vacuum is dielectric but paramagnetic (its effect on the magnetic field is due to the magnetic moments of the virtual pairs). It must be stressed that eqs. (3) are a plausible hypothesis, consequence of the fourth Heisenberg relation.

Because of (3) the observed values of the electron charge and light velocity at potential $`\mathrm{\Phi }`$ must be equal to $`e/ϵ_\mathrm{r}=e(1+\beta (\mathrm{\Phi }\mathrm{\Phi }_{})/c^2)`$ and $`c/\sqrt{ϵ_\mathrm{r}\mu _\mathrm{r}}=c[1+(\beta +\gamma )(\mathrm{\Phi }\mathrm{\Phi }_{})/2c^2]`$, at first order (the more negative is $`\mathrm{\Phi }`$, the smaller the charge and the slower the light). Consequently $`\alpha =e^2/4\pi \mathrm{}cϵ_0`$ must change to $`\alpha \sqrt{\mu _\mathrm{r}/ϵ_\mathrm{r}^3}`$. This means that the change in the observed value of $`\alpha `$, as computed at Earth from the light absorbed or emitted at a spacetime point with potential $`\mathrm{\Phi }`$, must be

$$\mathrm{\Delta }\alpha /\alpha =\xi (\mathrm{\Phi }\mathrm{\Phi }_{})/c^2,$$
(4)
where $`\xi =(3\beta \gamma )/2`$ and $`\mathrm{\Phi }_{}`$ is the potential at Earth.

On the gravitational potential of all the universe. As said before, this work considers only weak gravitational potentials that verify $`|\mathrm{\Phi }|/c^21`$ and can be studied therefore in the Newtonian limit . When that condition is no longer verified, the gravitational energy of a body is comparable to its rest energy and a relativistic theory must be used instead. This point, however, needs a clarification, since there are two independent different aspects in a Newtonian approximation. One is that the force on a particle must be weak and equal to $`\mathrm{\Phi }`$, the potential $`\mathrm{\Phi }`$ being determined up to an additive constant, other is the value of the potential energy of a particle. When studying a bounded system, as the Earth-Moon pair, a galaxy or a cluster, the effect of the faraway bodies can be neglected in the study of its motion, taking only the potential $`\mathrm{\Phi }`$ caused by the nearby masses. The reason is clearly that the distant bodies contribute to $`\mathrm{\Phi }`$ with a term that is almost independent of the space coordinates, at the scale of the system studied, and has therefore a negligible effect on the forces on the system. This is correct for all practical purposes concerning the motion. However, the gravitation being a long range force, one must be careful to include all the faraway matter that could have an effect on the problem if something depends on its gravitational potential energy. For this reason, the potentials $`\mathrm{\Phi }`$ and $`\mathrm{\Phi }_{}`$ that appear in eqs. (1)-(4) will be taken in this work as being due to all the universe. Note that this potential must include the contribution of matter, either ordinary or dark, as well as the effect of the cosmological constant.

Indeed the gravitational potential energy of a body is defined as the energy needed to bring it from the infinity to its actual position, without changing its kinetic energy or any other nongravitational energy. It could be argued that a body can not be brought actually from the infinity, since this is farther away than the horizon of the visible universe. However, it is clear that less energy is needed to create a virtual pair if $`\mathrm{\Phi }`$ is negative than if it is zero, the difference playing the same role as the gravitational energy when bringing a body to a spacetime point at potential $`\mathrm{\Phi }`$. Therefore, the variation of the lifetime of the virtual pairs is due in this model to the gravitational interaction of the virtual pairs of the quantum vacuum with all the matter and radiation in the universe, and also with all the quantum vacuum itself.

It is assumed here that the universe is flat (i.e. $`k=0`$) and consists in ordinary plus dark matter (with zero pressure) and dark energy due to the cosmological constant $`\mathrm{\Lambda }`$. The densities over the critical density at present time are noted as usual $`\mathrm{\Omega }_M,\mathrm{\Omega }_\mathrm{\Lambda }`$ (with $`\mathrm{\Omega }_M+\mathrm{\Omega }_\mathrm{\Lambda }=1`$). Taking $`p_\mathrm{\Lambda }=\rho _\mathrm{\Lambda }c^2`$ as the equation of state of the quantum vacuum, the source of its gravity in the Newtonian approximation is $`\rho _\mathrm{\Lambda }+3p_\mathrm{\Lambda }/c^2=2\rho _\mathrm{\Lambda }`$. Consequently, the space average potentials at present time due to the matter and to the quantum vacuum are $`\mathrm{\Phi }_{\mathrm{av},\mathrm{M}}=\mathrm{\Omega }_M\mathrm{\Phi }_0`$ and $`\mathrm{\Phi }_{\mathrm{av},\mathrm{\Lambda }}=2\mathrm{\Omega }_\mathrm{\Lambda }\mathrm{\Phi }_0`$, respectively, where $`\mathrm{\Phi }_0/c^2=_0^{R_U}G\rho _{\mathrm{cr}}4\pi rdr/c^20.3`$ is the potential that would be created by a mass distribution with the critical density, $`R_U=3,000`$ Mpc being taken as the visible universe radius. Most of these potentials is due to the larger distances. This means that the average gravitational potential produced by the visible universe is now equal to $`\mathrm{\Phi }_{\mathrm{av}}=\mathrm{\Phi }_0(\mathrm{\Omega }_M2\mathrm{\Omega }_\mathrm{\Lambda })`$. The effect of the inhomogeneities can be neglected. For instance, the contributions to the potential of the Sun, of Earth and of the Galaxy at a terrestrial laboratory are, respectively, $`10^8c^2`$, $`7\times 10^{10}c^2`$ and $`6\times 10^7c^2`$, which are much weaker than $`\mathrm{\Phi }_{\mathrm{av}}`$.

Cosmological variation of $`\alpha `$. Webb et al measured lines absorbed by distant gas clouds at high redshift, so that $`\mathrm{\Phi }`$ in (4) must be the potential at the absorption time $`t`$, which can be taken to be approximately uniform, except for small scale inhomogeneities. In the past, when the universe was more compact and dense, the distances varied as the scale factor $`a(t)`$, the mass density as $`a(t)^1`$ while the density of the quantum vacuum was constant (except for its variation due to eq. (2) that would give a second order effect). Consequently, the space average potential at time $`t`$ was $`\mathrm{\Phi }_{\mathrm{av}}(t)=\mathrm{\Phi }_0(\mathrm{\Omega }_M/a(t)2\mathrm{\Omega }_\mathrm{\Lambda }a^2(t))`$. Introducing this expression in eq. (4) and neglecting the local inhomogeneities, it turns out that $`(\alpha _z\alpha )/\alpha `$, can be expressed as

$$\frac{\mathrm{\Delta }\alpha }{\alpha }=\xi \frac{\mathrm{\Phi }_0}{c^2}\left[\mathrm{\Omega }_M\left(\frac{1}{a(t)}1\right)2\mathrm{\Omega }_\mathrm{\Lambda }\left(a^2(t)1\right)\right].$$
(5)
Assuming a flat dust universe, the time evolution of the scale factor is

$$a(t)=\left(\frac{\mathrm{\Omega }_M}{\mathrm{\Omega }_\mathrm{\Lambda }}\right)^{1/3}\mathrm{sinh}^{2/3}\left[\frac{(3\mathrm{\Lambda })^{1/2}t}{2}\right],$$
(6)
with $`\mathrm{\Lambda }=8\pi G\mathrm{\Omega }_\mathrm{\Lambda }\rho _{\mathrm{cr}}`$ (eq. (6) reduces to $`a(t)=(t/t_0)^{2/3}`$ if $`\mathrm{\Lambda }=0`$, $`t_0`$ being the age of the universe). Equations (5)-(6) give the main result of this work: the time dependence of $`\mathrm{\Delta }\alpha /\alpha `$. Note that, since $`\xi \mathrm{\Phi }_0<0`$, $`\mathrm{\Delta }\alpha <0`$ for $`t<0`$, so that the quantum vacuum must have been optically denser in the past (the electron charge was more renormalized than now since $`\mathrm{\Phi }_{\mathrm{av}}(t)\mathrm{\Phi }_{\mathrm{av}}(t_0)<0`$).

The thick line in Figure 1 shows the relative change of the fine structure constant given by eqs. (5)-(6) in the case $`(\mathrm{\Omega }_M=0.3,\mathrm{\Omega }_\mathrm{\Lambda }=0.7)`$ as compared with the observations by Webb et al  (in units of $`10^5`$) versus the look back-time (in units of the age of the universe), with $`\xi =1.3\times 10^5`$. This is the value that gives the best fit (remind that $`\mathrm{\Phi }_0/c^20.3`$); it was obtained by minimizing $`\chi ^2`$, the minimum value being $`0.63`$ per point. The thin line shows the same result for $`(\mathrm{\Omega }_M=1,\mathrm{\Omega }_\mathrm{\Lambda }=0)`$. The best fit was obtained here for $`\xi =1.9\times 10^5`$, the value of $`\chi ^2`$ being of $`0.56`$ per point. Although the latter gives a slightly better fit, both curves fit well the observational points taking into account the large error bars, their difference being small. As explained before, the light velocity must be affected also by the quantum vacuum and be time dependent. A similar argument shows that the refractive index at time $`t`$ is $`n(t)=c/c(t)=1(\beta +\gamma )\mathrm{\Phi }_0[\mathrm{\Omega }_m(1/a1)2\mathrm{\Omega }_\mathrm{\Lambda }(a^21)]/2c^2(>1)`$, so that the light velocity is an increasing function of time.

Validity of the model. Being based on a first order approximation in $`\mathrm{\Phi }`$, this model is valid only if $`\mathrm{\Omega }_M|\mathrm{\Phi }_0|/a(t)c^21`$, what means that the gravitational potential energy of a particle is much smaller than its rest energy. This condition is verified at present time, approximately at least, but looses progressively its validity when $`a(t)`$ decreases so that the model can not be extended arbitrarily towards the past (in its present nonrelativistic formulation).

It must be stressed that the singularity of (5) at $`t=0`$ has no physical meaning since it is outside of the applicability of the model. Neither it can be said that $`\alpha `$ is predicted to change too rapidly or to vanish near time zero. One can not apply these ideas, therefore, to the recombination era, which happened at a very high redshift. But the interval of validity does cover the range of the observations by Webb et al, at least as a good approximation. To extend this model farther away in the past, it must be reformulated in relativistic terms.

This is not to say that the Newtonian approximation can not be used farther away in the past, for higher redshift. It can certainly be used whenever it is not necessary to consider the gravitational interaction of a system with all the universe.

On the Oklo and other data. Damour and Dyson  analyzed the data from the natural reactor which operated 1.8 billion years ago at Oklo (Gabon) and concluded that the relative change of $`\alpha `$ from then to now is in the interval $`(0.9\times 10^7,1.2\times 10^7)`$ (assuming that other constants like the Fermi constant do not vary). Equation (5) gives for that time the relative change $`\mathrm{\Delta }\alpha /\alpha 1.3\times 10^6`$ (resp. $`5.4\times 10^7`$) if $`\mathrm{\Omega }_M=0.3,\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ (resp. $`\mathrm{\Omega }_M=1,\mathrm{\Omega }_\mathrm{\Lambda }=0`$), values which are outside but not far from that interval.

The Oklo study is thought to give the most powerful method to determine the variation of $`\alpha `$ with geochemical data but, according to Uzan  “one has to understand and to model carefully the correlations of the variations of $`\alpha _\mathrm{w}`$ and $`g_\mathrm{s}`$ as well as the effect of $`\mu `$ ($`=m_e/m_p`$). This difficult but necessary task remains to be done”. There are several others studies which set bounds on the variation of $`|\mathrm{\Delta }\alpha |`$, using a number of different data . The results of the present model are compatible with these other bounds, except for the one set by Olive et al using data of $`\beta `$-decay in meteorites , $`|\mathrm{\Delta }\alpha /\alpha |<3\times 10^7`$ during the past 4.6 Gyr (redshift about 0.45), while this model gives $`3\times 10^6`$ (resp. $`1.7\times 10^6`$) if $`\mathrm{\Omega }_M=0.3,\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ (resp. $`\mathrm{\Omega }_M=1,\mathrm{\Omega }_\mathrm{\Lambda }=0`$). However, this bound by Olive et al “can also be altered if the neutrinos are massive”, according to Uzan .

In any case, this model gives a fair account of the observations by Webb et al. Taking everything into account, it seems worth of consideration to explore its consequences.

Comparison with the gravitational redshift. The effect described here produces a frequency shift which may seem similar to the gravitational redshift $`\mathrm{\Delta }\omega /\omega =\mathrm{\Delta }\mathrm{\Phi }/c^2`$. The two effects are different, however. The one observed by Webb et al is due to the variation of the fine structure constant. The principal part of the frequencies emitted by an atom is proportional to $`(e^2/4\pi \mathrm{}ϵ_0)^2`$, not exactly to $`\alpha ^2`$ while the relative separation between lines in the relativistic fine structure splitting does depend on $`\alpha ^2=(e^2/4\pi \mathrm{}cϵ_0)^2`$. On the other hand, the gravitational shift affects equally to all the lines in a multiplet. In experiments that disregard the width of the multiplets, the contribution of the change of $`\alpha `$ to the frequency shift is obtained from (4), but with $`2\beta `$ instead of $`\xi `$, i.e. $`\mathrm{\Delta }\omega /\omega =4\beta \mathrm{\Delta }\mathrm{\Phi }/c^2`$. In that case, the observed shift would be $`\mathrm{\Delta }\omega /\omega =(1+4\beta )\mathrm{\Delta }\mathrm{\Phi }/c^2`$, i.e. the addition of the two effects. The best confirmations of the gravitational redshift, those by Pound, Rebka and Snider, agree with the prediction of General Relativity up to about 1 $`\%`$ , but they refer to nuclear levels in which the electromagnetism plays only a part. This means that the results of this work are certainly compatible with the experiments on the gravitational redshift if $`\beta 2.5\times 10^3`$, an inequality that this work assumes to be satisfied. However, there is a problem: as $`\xi <3\beta /2`$, a bound on $`\beta `$ is also a bound on $`\xi `$ but the converse is not true. In any case, a necessary condition for the compatibility of the results of this work and the gravitational redshift experiments is that $`\xi <4\times 10^3`$, while this model predicts that $`\xi `$ is equal to about $`1.3\times 10^5`$ and $`1.9\times 10^5`$ in the two cases considered.

Two final comments. First, note that the dressed electron charge and the fine structure constant are not treated here as universal constants but as the result of the interaction of point charges with the quantum vacuum, their variations being an indirect consequence of the universal expansion, through the modification of their renormalization.

Second, because of the conservation of the charge conjugation, the pairs are created with $`L=S=0`$, so that their energy in a magnetic field is $`E2\mu _BB`$ with probability 1/2, $`E`$ being the non-magnetic energy. For magnetic field close to $`B_0=mc^2/\mu _B8.8\times 10^{13}`$ gauss, some pairs would have zero total energy, their lifetime being infinite according to the fourth Heisenberg relation). If $`B`$ approaches $`B_0`$ from below, the creation of the pairs would begin to be energetically free. There would be a threshold to some peculiar phenomenon, although it is not clear what this can be (note that the inequality $`BB_0`$ is a necessary condition for the validity of the Euler-Heisenberg Lagrangian to study in QED the interactions of the quantum vacuum with a magnetic field ). Intriguingly, it turns out that $`B_0`$ is close to the highest magnetic field measured in magnetars, see Kouveliotou et al . Commenting ref. , Kulkarni and Thompson  say, “our failure to detect radio pulsars with magnetic fields greater than $`2\times 10^{13}`$ gauss is because, at such high fields, quantum electrodynamic effects help to damp the radio emission” this being “direct evidence of magnetic fields strong enough to perturb the very structure of the vacuum”. This gives support to the present work: the upper limit for the magnetic field of magnetars could be just a manifestation of the fourth Heisenberg relation similar to the one considered here.

Summary and conclusions. In the phenomenological model presented here, the cosmological variation of the fine structure constant is due to the combined effect of the fourth Heisenberg relation and the gravitational interaction of the virtual pairs in the zero-point radiation with all the universe. The problem is studied in the Newtonian approximation. The quantum vacuum is treated as an optical medium characterized by its relative permittivity and permeability that depend on the average gravitational potential of the universe. The model predicts that $`\mathrm{\Delta }\alpha /\alpha `$ is proportional to $`\{\mathrm{\Omega }_M[a(t)^11]2\mathrm{\Omega }_\mathrm{\Lambda }[a(t)^21]\}`$ ($`a(t)`$ being the scale factor). The results for the cases $`\mathrm{\Omega }_M=0.3,\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$ and $`\mathrm{\Omega }_M=1,\mathrm{\Omega }_\mathrm{\Lambda }=0`$ agree with the observations by Webb et al  on the cosmological variation of the fine structure constant, as is seen in figure 1. The argument goes as follows.

(i) Because of the fourth Heisenberg relation, the lifetime of the virtual pairs of the zero-point radiation depends on the gravitational potential $`\mathrm{\Phi }`$. This causes the permittivity and the permeability of the quantum vacuum, the observable electron charge and the light velocity to depend also on $`\mathrm{\Phi }`$. The consequent change of the observed fine structure constant is expressed at first order in terms of the average gravitational potential of the universe, and of a parameter $`\xi `$ related to the renormalization effects of the quantum vacuum.

(ii) As the universe was more compact in the past, its average gravitational potential was more negative (or less positive). Consequently, the lines of the spectra of distant quasars were absorbed with a more renormalized value of the electron charge than now. As a result, $`\mathrm{\Delta }\alpha /\alpha `$ is given by eq. (5), which is plotted in fig. 1. The agreement seems good. Note that, in this model, the optical density of the quantum vacuum increases towards the past and decreases along the universe history.

(iii) The model can not be extended arbitrarily towards the past, since it is based on a first order Newtonian approximation that is valid for the recent past, including the range of the observations by Webb et al, but is no longer applicable to higher redshift. In particular, it can never be applied to the recombination era. A relativistic approach should be followed in order to extend the model to the remote past. Note finally that, according to this model, the light is also affected by the gravitational potential so that it was slower in the past.

The conclusion of this paper is that the combined effect of the gravitation of all the universe and of the fourth Heisenberg relation on the density of the zero-point radiation and the corresponding cosmological variations of $`\alpha `$, $`e`$ and $`c`$ proposed here should be further investigated.

I am grateful to Profs. Claudio Aroca, Ana Inés G. de Castro, Eloísa López, Juan M. Usón and José L. Trueba for discussions.


text:
# Parton Production Via Vacuum Polarization

## References


text:
# How the X-ray spectrum of a Narrow-Line Seyfert 1 galaxy may be reflection dominated

## 1 Introduction

The basic model for the active nucleus of a typical Seyfert 1 galaxy consists of a black hole surrounded by an ultraviolet and soft X-ray emitting thin accretion disc above which a patchy active corona emits hard X-rays. The model is supported by strong X-ray reflection signatures due to Compton backscattering and fluorescence of the harder coronal X-rays by the dense disk (Pounds et al., 1990; Nandra & Pounds, 1994; Petrucci et al., 2000). In some objects an extreme red wing is found in the iron line of the reflection spectrum (Tanaka et al., 1995; Nandra et al., 1999; Fabian et al., 2000), showing that the disk can extend very close to the black hole where large relativistic effects occur. In general the strength of the reflection spectrum indicates that the disk surface is approximately flat, although disk instabilities could well make it ribbed or clumpy, particularly if the disk is radiation-pressure dominated (e.g., Lightman & Eardley, 1974; Guilbert & Rees, 1988; Krolik, 1998; Turner, Stone & Sano, 2001). Such a situation is most likely to occur when the accretion rate is high and close to the Eddington limit, which is often thought to explain the unusual properties of Narrow-Line Seyfert 1 (NLS1) galaxies (e.g., Pounds, Done & Osborne, 1995).

Recently, Boller et al. (2002) reported on a “sharp spectral feature” at about 7.1 keV in the XMM-Newton spectrum of the NLS1 1H 0707–495 ($`z=0.0411`$). This deep drop (over a factor of two in flux) in the spectrum occurs at almost the exact energy of the neutral iron edge, so absorption models (in particular, a partial covering model) were favoured to explain the feature. However, these models suffered from a very soft intrinsic power-law ($`\mathrm{\Gamma }3.5`$) and an unreasonably large value of the iron abundance ($``$ 35$`\times `$ solar). Moreover, it was difficult for the partial covering model to account for the absence of a narrow Fe K$`\alpha `$ line and to explain the rapid variability of the source over the entire waveband. The sharp drop in the spectrum could also be due to the blue wing of a relativistic iron emission line, but, as the authors noted, to explain the depth of the drop at $``$7 keV requires “invoking a very extreme Fe abundance and/or reflection fraction.” Here, we consider the possibility of a reflection dominated X-ray spectrum in detail.

The basic idea we wish to explore is that within a system which is accreting at close to the Eddington rate, disc instabilities funnel denser material into many deep rings, between which the hard X-rays are emitted by some coronal process (Figure 1). To an off-axis observer, the X-ray emission will consist mostly of a reflection component. The rapid variability seen in the X-ray emission from 1H 0707–495 could then in part be due to changes in the geometric obscuration of one sheet on another, explaining why little spectral variability is seen (Boller et al., 2002). This reflection model requires a less extreme input spectrum and hidden flux than the partial-covering model.

While this is clearly an idealized picture, recent numerical simulations of the magneto-rotational instability in radiation-pressure dominated discs have found that the turbulence can cause density contrasts up to a factor of 200 (Turner, Stone & Sano, 2001). Therefore, it does seem plausible that a radiation-pressure dominated disc will exhibit strong density inhomogeneities.

In the next section, we calculate what the hard X-ray power-law may look like in our idealized picture and show that a steep spectrum naturally arises from this situation. In Section 3 we fit reflection dominated spectra to the broadband XMM-Newton data of 1H 0707–495. Finally, we discuss the implications of this model in Section 4.

## 2 Calculation of the X-ray power-law

Here we give a simple analytic estimate of the slope of the X-ray continuum by assuming that in the inner radiation-pressure dominated part of the accretion disc the dense and cold sheet-like part of the flow is embedded in a tenuous, hot plasma (cf., Guilbert & Rees, 1988; Celotti, Fabian & Rees, 1992; Kuncic, Blackman & Rees, 1996; Kuncic, Celotti & Rees, 1997). The accretion energy is dissipated, probably via magnetic reconnection, in the hot phase only, with the cold sheets reprocessing (reflecting and thermalizing) the incident hard flux.

For ease of computation, the picture in Fig. 1 was adapted so that the cold sheets have a flat disc-like geometry with radius $`R`$. We consider two such discs facing each other and separated by a distance $`d`$. The hot phase is sandwiched by the two cold discs and therefore has a cylindrical slab geometry. From a point in the mid-plane of the cylinder, at a distance $`r`$ from its centre, the covering fraction of the two cold sheets is

$$C(r)=\left[1\frac{1}{2}\left(\frac{1}{\sqrt{1+\frac{4(Rr)^2}{d^2}}}+\frac{1}{\sqrt{1+\frac{4(R+r)^2}{d^2}}}\right)\right].$$
(1)
The average total covering fraction is then given by

$$C=\frac{1}{R}_0^RC(r)𝑑r=1\frac{d}{4R}\left[\mathrm{log}\left(\frac{4R}{d}+\sqrt{1+\frac{16R^2}{d^2}}\right)\right].$$
(2)
With the hard luminosity of the hot phase being $`L_\mathrm{H}`$, the reprocessed luminosity is given by

$$L_{\mathrm{rep}}=CL_\mathrm{H}(1a),$$
(3)
where $`a`$ is the albedo of the cold sheets. We make use of the results of Malzac, Beloborodov & Poutanen (2001) who used Monte Carlo simulations of Comptonizing coronae above dense cold material and showed that the albedo depends on the spectral index of the illuminating radiation, being smaller for steep spectra. The thermal reprocessed soft photon flux from the cold sheets is assumed to be a black-body spectrum with temperature $`T_{\mathrm{rep}}`$, and is calculated by the method outlined by Merloni & Fabian (2001). The hard luminosity is assumed to be due only to inverse Compton scattered photons and can be approximated as a sum of a cut-off power-law and Wien component (Wardziński & Zdziarski, 2000).

The emerging spectral index of the power-law is given by

$$\mathrm{\Gamma }1=\alpha =\frac{\mathrm{ln}P_{\mathrm{sc}}}{\mathrm{ln}(1+4\mathrm{\Theta }+16\mathrm{\Theta }^2)},$$
(4)
where $`\mathrm{\Theta }=kT_\mathrm{e}/m_\mathrm{e}c^2`$ is the electron temperature of the hot phase, and $`P_{\mathrm{sc}}`$ is the scattering probability averaged over the source volume and depends only on the coronal optical depth $`\tau `$. In a slab geometry $`P_{\mathrm{sc}}`$ can be approximated as (Zdziarski, 1994)

$$P_{\mathrm{sc}}=1+\frac{\mathrm{exp}(\tau )}{2}\left(\frac{1}{\tau }1\right)\frac{1}{2\tau }+\frac{\tau }{2}E_1(\tau ),$$
(5)
where $`E_1`$ is the exponential integral.

Once the geometry is fixed (by fixing the aspect ratio $`R/d`$), the temperature of the hot phase, and consequently the slope of the Comptonized continuum, can be calculated self-consistently by solving $`L_\mathrm{H}=L_\mathrm{C}=_{2.78\mathrm{\Theta }_{\mathrm{rep}}}^{\mathrm{}}L_\mathrm{C}(x)𝑑x`$. We have calculated the spectral index $`\mathrm{\Gamma }`$ for different values of the ratio $`R/d`$ and in Figure 2 we plot the spectral index for three different values of the optical depth of the hot phase. Clearly, as the covering fraction increases, the spectrum softens significantly. Also, the spectrum is steeper for a denser hot medium (larger $`\tau `$). Finally, if there is accretion energy dissipated in the cold phase, then the spectrum will be softer for a given $`R/d`$ (c.f., Malzac, 2001).

Figure 3 shows an equivalent plot where the abscissa is now the half opening angle between the two reprocessors $`\theta /2`$.

These plots show that spectral indices typical of NLS1s (i.e., $`\mathrm{\Gamma }>`$ 2.1–2.2) are a natural consequence of this scenario if the half opening angle between the two reprocessors is less than about 20 degrees. In such a case, it is likely that a observer may only detect the reflected emission arising from the walls of the reprocessors.

## 3 Application to the NLS1 1H 0707–495

In this section, we attempt to fit the broadband X-ray spectrum of the NLS1 1H 0707–495 with models of reflection dominated spectra, which might be expected if the accretion disc has broken into a number of dense segments (Fig. 1). We consider the EPIC-pn spectrum of 1H 0707–495 between 0.5 and 11 keV and fit with the constant density ionized disc models of Ross & Fabian (1993) (see also Ross, Fabian & Young 1999). The parameters of the model are the ionization parameter $`\xi =4\pi F_\mathrm{X}/n_\mathrm{H}`$, where $`F_\mathrm{X}`$ is the illuminating flux between 10 eV and 100 keV, the photon index of the power-law spectrum that strikes the slab, and a normalization constant. Relativistic blurring was applied to the model during fitting using the kernel of Laor (1991). The emissivity index was fixed at $`3`$, but the inner and outer radii, along with the inclination angle, were left free. Absorption due to the Galactic column of $`5.79\times 10^{20}`$ cm<sup>-2</sup> was included in all fits.

If one tries to fit the data with a single reflection spectrum, then the best fit ($`\chi ^2`$/d.o.f.=506/309) is obtained with a model with a three times solar Fe abundance. Recall that no power-law component has been added to these models, giving them, in effect, an infinite reflection fraction. We find that $`\mathrm{\Gamma }=2.44`$ and $`\mathrm{log}\xi =2.040`$, which results in a strong neutral Fe K$`\alpha `$ line. Extreme Kerr blurring ($`r_{\mathrm{in}}=1.24`$ $`r_g`$) was needed to fit the large red wing of the line between 2.5 and 6 keV. However, as indicated by the large $`\chi ^2`$, there were significant residuals to the fit, particularly below 1 keV where the model predicts strong emission lines from Fe and O (Ross & Fabian, 1993). To dilute the effect of these lines, a power-law was added to the model with the value of $`\mathrm{\Gamma }`$ fixed to be the same as the reflected component. Reflection still dominated the spectrum with a reflection fraction $`>`$ 10. The power-law improved the best fit to $`\chi ^2`$/d.o.f.=479/308, but required an Fe abundance greater than 5 times solar. In this fit, $`\mathrm{\Gamma }`$ softened slightly to 2.52, but there was little change in the ionization parameter, the inner and outer radii, or the inclination angle ($`23`$ degrees). There are still residuals below 1 keV, but, interestingly, there is now a clear line-like residual just above 6 keV. Adding a Gaussian to the fit drops the $`\chi ^2`$ by 38 with the addition of three degrees of freedom – significant at $`>99.99`$ per cent, according to the F-test. Here, the best-fit model has a seven times overabundance of Fe, and $`\mathrm{\Gamma }=2.56`$. The line had an energy of 6.62 keV and a width $`\sigma =0.235`$ keV, indicating that it probably also arises from somewhere in or on the disk.

To investigate this further, the power-law and Gaussian components of the model were replace by a second reflector. The second reflection component was drawn from the same group of models as the first, and was also subject to relativistic blurring. We fixed the photon-index and inclination angle to be the same for both components. The addition of this second reflector greatly improved the fit ($`\chi ^2`$/d.o.f.=342/305), as the combined emission removed most of the residuals from the soft X-ray band, and now only a five times overabundance of Fe was required. The ionization parameter of the second component is $`\mathrm{log}\xi =3.397`$ resulting in a strong He-like Fe K$`\alpha `$ line at 6.7 keV and an inner radius of 5.4 $`r_g`$. The other reflector remained relatively neutral with $`\mathrm{log}\xi =2.09`$ and strongly blurred ($`r_{\mathrm{in}}=1.24`$ $`r_g`$). The photon-index flattened slightly to $`\mathrm{\Gamma }=2.41`$. While this fit is almost acceptable there is still a line-like residual above 6 keV. Adding a Gaussian component to account for this feature significantly improves the fit ($`\chi ^2`$/d.o.f.=292/302). This new line has an energy of 6.45 keV and a width of 0.358 keV, demonstrating that the line arises from weakly ionized Fe in the accretion disc.

The two-reflector plus Gaussian emission line model gives a very good fit to the EPIC-pn data between 0.5 and 11 keV, however it is interesting to consider replacing the line component with a third reflection spectrum. This model also results in an excellent fit to the data ($`\chi ^2`$/d.o.f.=287/303; Figure 4), but requires a 7 times solar abundance

of Fe. The three components all have different ionization parameters: $`\mathrm{log}\xi =1.986`$ (with $`r_{\mathrm{in}}=2.08`$ $`r_g`$ and $`r_{\mathrm{out}}=7.25`$ $`r_g`$), 3.847 (with $`r_{\mathrm{in}}=4.18`$ $`r_g`$), and 1.722 (with $`r_{\mathrm{in}}=19.8`$ $`r_g`$). The last two reflectors had their outer radii frozen at 50 $`r_g`$. The photon index and the inclination angle for all three components were fixed to have the same values. The best fit values were $`\mathrm{\Gamma }=2.35`$ and 18.2 degrees, respectively. The model is shown in Figure 5.

## 4 Discussion

We have shown that the X-ray spectrum of 1H 0707–495 is well fitted by a reflection-dominated spectrum. The incident power-law is mildly steep with a photon index of 2.35. Furthermore, the model can completely account for the strong soft excess in this object by a combination of blurred emission lines and thermal emission from the irradiated blobs. The geometry envisaged for the source consists of a set of deep dense sheets or rings, orbiting the central black hole, between which hard X-rays are emitted. The fact that at least two reflectors are needed to provide a good fit to the data is not inconsistent with the model. There should be numerous rings or sheets in the centre of the disk with (as in the standard coronal picture) multiple hard X-ray sources active at any one time. The system is observed at sufficient inclination that the hard X-ray sources are not directly visible and only X-rays reflected from the sides of the dense matter are seen.

In fact, the model we have presented is a simplified picture of an intrinsically complex situation. In this respect, the interpretation of the fit parameters should be consider with some caution. In a more realistic situation, the reflection features will be broadened due to both relativistic effects and disc turbulence, so that parameters like $`r_{\mathrm{in}}`$ and $`r_{\mathrm{out}}`$ should not be taken at face value, but may be interpreted as the signature of the presence of different reprocessors at different locations. Also, a range of ionization parameters is to be expected from the complex configuration we envisage, as a reflector would likely have an ionization gradient on its surface.

A consistent model for both the photon index of the input power-law and the reflection component is obtained if the opening angle of the gaps between the sheets or rings is less than 40 degrees and the whole disk inclination is about 20 degrees (Fig. 3 and Sect. 3). The metallicity is still required to be high (about 5–7 times the Solar value) but less than for other models. It is possible that the metallicity can be lowered if multiple reflection is included (this will have an approximately similar effect to squaring the reflection spectrum).

The spectral-independent variability of 1H 0707–495 (Boller et al., 2002) would be mostly due to changes in the brightness of the hard X-ray source causing the observed reflection components to vary. Additional variability may also arise from changes in the geometry of the system, with sheets or rings obscuring each other as they orbit the black hole. The geometry may also account for larger scale variations. When observed by XMM-Newton, 1H 0707–495 was about ten times fainter than in previous observations (Boller et al., 2002). This could be due to slight changes in the geometry and not necessarily a large change in intrinsic luminosity. If the hard X-ray source was previously more directly visible (say the opening angle of the sheets or rings was larger, or the hard X-ray emitting corona extended nearer to the disk surface) then the observed flux would have been larger. Indeed, the spectrum observed by ASCA when the flux was higher ($`\mathrm{\Gamma }2.27`$, Leighly 1999), is consistent with being the unobscured continuum emission peeking through the cold sheets.

It is possible that the model we have developed here for 1H 0707–495 has a wider relevance to Seyfert 1 galaxies. The reflection fraction can exceed unity if the disk is in the form of deep rings or sheets, and leads to an intrinsic hard power-law steeper than $`\mathrm{\Gamma }=2`$. MCG–6-30-15 has sometimes been observed with a spectral index of $`2.1`$ or steeper (Vaughan & Edelson, 2001; Shih, Fabian & Iwasawa, 2001) and its strong iron line suffers far less variability than its continuum (Vaughan & Edelson, 2001). This could be due to corrugations in its inner disk which may cause the reflection components to be stronger and more long-lasting than if the accretion disc was perfectly flat.

## Acknowledgments

ACF thanks the Royal Society for support. DRB acknowledges financial support from the Commonwealth Scholarship and Fellowship Plan and the Natural Sciences and Engineering Research Council of Canada. Based on observations obtained with XMM-Newton, an ESA science mission with instruments and contributions directly funded by ESA Member States and the USA (NASA).


text:
# References

In building models with supersymmetry (SUSY), we must take into account that the sfermion mass spectrum of the first and the second generations is severely constrained from the flavor changing neutral current (FCNC) processes such as $`K^0\overline{K}^0`$ mixing etc. Mainly, three approaches to address this problem (SUSY flavor problem) have been discussed in the literature; (1) the degeneracy , (2) the alignment  and (3) the decoupling .

In this letter, we consider the decoupling solution to SUSY flavor problem in the context of higher dimensional theories. The basic idea is very simple. In extra dimensions, it is well known that if the matter wave functions are localized at different points in extra dimensions, Yukawa hierarchy can be obtained by the suppression factor of the overlap of wave functions . Since the fermion mass hierarchy is $`m_1<m_2<m_3`$ where $`m_i`$ is a fermion mass of the i-th generation, the matter of the third generation is localized close to the Higgs fields and the first generation is localized most distant from the Higgs fields. Introducing SUSY breaking brane in which the chiral superfield with nonvanishing F-term is localized correlates the sfermion masses with the fermion masses. If the SUSY breaking brane is put close to the first generation matter fields, the sfermion mass hierarchy is inverted, $`\stackrel{~}{m}_1>\stackrel{~}{m}_2>\stackrel{~}{m}_3`$ where $`\stackrel{~}{m}_i`$ is the sfermion mass of the i-th generation. Therefore, we expect that the decoupling solution can be a natural solution. Namely, the sfermion masses of the first and the second generation is the order of 10 TeV and the sfermion mass of the third generation is the order of 100 GeV for naturalness.

Let us discuss the model in detail. We consider an $`𝒩=1`$ supersymmetric theory in five dimensions. We introduce two 3-branes at $`y=0`$ and $`y=L`$, where $`y`$ denotes the fifth coordinate in five dimensional space-time. The gauge supermultiplets of the Standard Model (SM) gauge groups lives in the bulk and its zero mode wave functions are flat in the fifth dimension. The matter fields also lives in the bulk and its zero mode wave functions are assumed to be Gaussian.111For readers interested in the localization mechanism of the chiral superfields, see Appendix of Ref. . Higgs doublets are assumed to be localized on the brane at $`y=0`$, we refer to this brane as “H-brane”. Further, extra chiral superfields $`X,\mathrm{\Phi }^{}`$ and $`\overline{\mathrm{\Phi }}^{}`$ localized on the brane at $`y=L`$ are introduced. $`X`$ is a chiral superfield with nonvanishing F-term (i.e. $`X=\theta ^2F`$). $`\mathrm{\Phi }^{}`$ and $`\overline{\mathrm{\Phi }}^{}`$ are vector-like superfields with a mass $`M`$. A pair of vector-like superfields are introduced for each matter chiral superfields, namely $`Q^{},\overline{Q}^{}`$ for $`Q`$, and $`U^{},\overline{U}^{}`$ for $`\overline{U}`$ and $`L^{},\overline{L}^{}`$ for $`L`$ etc. We refer to the brane at $`y=L`$ as “SUSY breaking brane”.

Naive expectation is that the gaugino mass is generated from the term

$`\delta (yL){\displaystyle d^2\theta \frac{X(x)}{M_{}^2}W^\alpha (x,y)W_\alpha (x,y)}M_\lambda ={\displaystyle \frac{F}{M_{}^2L_c}},`$ (1)
where $`x`$ is a coordinate of the four dimensional space-time, $`W_\alpha `$ is the field strength tensor superfield living in the bulk, $`L_c`$ is the width of the thick wall which should be considered as the compactification length in our framework and $`M_{}`$ is the Planck scale in five dimensions. Note that the gaugino masses receive only the volume suppression factor. The soft breaking mass of the $`i`$-th generation chiral superfield $`\mathrm{\Phi }_i`$ is naively generated from the term

$`\delta (yL){\displaystyle d^4\theta \frac{X^{}(x)X(x)}{M_{}^3}\mathrm{\Phi }_i^{}(x,y)\mathrm{\Phi }_j(x,y)}`$
$`\stackrel{~}{m}_{ij}^2={\displaystyle \frac{|F|^2}{M_{}^3L_c}}\mathrm{exp}[M_{}^2(Ly_i)^2M_{}^2(Ly_j)^2],`$ (2)
where the form of the zero mode wave function of the matter fields is assumed as $`\mathrm{\Phi }_i^{(0)}\mathrm{exp}[M_{}^2(yy_i)^2]`$. The sfermion masses receives not only the volume suppression but also the exponential suppression, therefore are negligibly small compared to the gaugino masses. This spectrum is similar to the gaugino mediation scenario . Although the spectrum of the gaugino mediation is phenomenologically interesting, this is not the subject in this paper.

As is clear from the above argument, the modification is needed to obtain the spectrum of the decoupling solution. The wayout is to replace $`M_{}`$ with some scale $`M<M_{}`$ such that the enhancement by $`M`$ compensates for the exponential suppression. We consider here the following superpotential on the SUSY breaking brane

$`W`$ $`=`$ $`\delta (yL){\displaystyle 𝑑y[\lambda X(x)\mathrm{\Phi }_i(x,y)\overline{\mathrm{\Phi }}^{}(x)+M\mathrm{\Phi }^{}(x)\overline{\mathrm{\Phi }}^{}(x)]},`$ (3)
$`=`$ $`\lambda ϵ_iX(x)\mathrm{\Phi }_i(x)\overline{\mathrm{\Phi }}^{}(x)+M\mathrm{\Phi }^{}(x)\overline{\mathrm{\Phi }}^{}(x),`$ (4)
where $`\lambda `$ is a dimensionless constant of order unity. The second expression is obtained by integrating out the fifth dimensional degree of freedom. $`ϵ_i`$ is a suppression factor coming from the zero mode wave function of i-th generation. At the scale below $`M`$, we can integrate out the massive fields $`\mathrm{\Phi }^{},\overline{\mathrm{\Phi }}^{}`$. Then, the effective superpotential of Eq.(4) vanishes, and the effective Kähler potential receives the correction at tree level such as

$$\delta K=\frac{ϵ_i^2}{M^2}X^{}X\mathrm{\Phi }_i^{}\mathrm{\Phi }_i.$$
(5)
As mentioned above, it is necessary to introduce a pair of vector-like fields for each matter chiral superfield to obtain the above Kähler potential. Otherwise, the Kähler potential with the flavor mixing will arise in general. We does not consider this possibility. Sfermion masses coming from Eq. (5) becomes

$`\stackrel{~}{m}_i^2=ϵ_i^2{\displaystyle \frac{|F|^2}{M^2}}.`$ (6)
This result has the desirable features that the suppression factor is not only replaced with $`M<M_{}`$, but also proportional to $`ϵ(\mathrm{Yukawa})^1`$. Note that this argument holds in the case $`F<M^2`$. Also, we assumed here that the overall sign of the Kähler potential is positive.

In our scenario, since the information of the location where the matter fields are localized is necessary to derive the sfermion mass spectum, we briefly discuss Yukawa hierarchy. We consider the up-type Yukawa coupling for example

$`W=\delta (y){\displaystyle 𝑑yQ_i(x,y)\overline{U}_j(x,y)H(x)},`$ (7)
where the order one coefficient is implicit. Integrating out the fifth dimensional degree of freedom, we obtain the effective Yukawa coupling in four dimensions at the compactification scale $`L_c^1`$ as

$`y_{\mathrm{eff}}\mathrm{exp}[M_{}^2(y_{Q_i}^2y_{\overline{U}_j}^2)].`$ (8)
In order to realize the Yukawa hierarchy

$`y_t𝒪(1),y_c𝒪(10^2),y_u𝒪(10^5),`$ (10)
$`y_b𝒪(10^2),y_s𝒪(10^4),y_d𝒪(10^5),`$
the location of the matter fields are determined,222To be correct, the effective Yukawa couplings (8) have to be evolved down to the weak scale by the renormalization group equation (RGE) and matched to Eqs. (10), (10) after diagonalizing Yukawa matrix. We simply neglect this RGE effects and the mixing angles. for instance,

$`y_Hy_{\overline{H}}y_{Q_3}y_{U_3}0,|y_{Q_2}||y_{U_2}|\sqrt{\mathrm{ln10}}M_{}^1,|y_{D_3}|\sqrt{2\mathrm{l}\mathrm{n}10}M_{}^1,`$ (11)
$`|y_{Q_1}||y_{U_1}||y_{D_1}|\sqrt{{\displaystyle \frac{5}{2}}\mathrm{ln10}}M_{}^1,y_{D_2}\sqrt{3\mathrm{l}\mathrm{n}10}M_{}^1.`$ (12)
Now, we turn to the sfermion mass spectrum. The decoupling solution requires that the masses of the first and the second generations should be heavier than $`𝒪`$(10) TeV and that of the third generation be around $`𝒪`$(100) GeV for naturalness. We consider two cases of the decoupling constraints according to Eqs.(11) and (12). The first case is $`\stackrel{~}{m}_{Q_2}=\stackrel{~}{m}_{U_2}10`$ TeV (solution 1), and the second one is $`\stackrel{~}{m}_{D_2}10`$ TeV (solution 2).

Let us consider the solution 1 at first. This case imposes

$`10^2ϵ_2/ϵ_3\mathrm{exp}[(M_{}L)^2(M_{}L\sqrt{\mathrm{ln10}})^2],`$ (13)
which leads to

$`M_{}L{\displaystyle \frac{3}{2}}\sqrt{\mathrm{ln10}}.`$ (14)
This implies that SUSY breaking brane is located between the first and the second gengerations. This is the interesting feature of our model. In the conventional SUSY breaking models in extra dimensions, the visible brane and the SUSY breaking brane are separated in the extra dimensional spaces and the dangerous flavor violating sfermion masses are suppressed by the locality. In our model, SUSY breaking brane are not separate from the visible brane (more correctly, the visible wall) but is predicted to be located within the visible wall to yield the decoupling sfermion mass spectrum.

Requiring that the gaugino mass should be around 100 GeV,

$`M_\lambda {\displaystyle \frac{F}{M_{}}}{\displaystyle \frac{1}{M_{}L_c}}100\mathrm{GeV},`$ (15)
we obtain

$`F(M_{}L_c)\times 10^2M_{}.`$ (16)
Since $`\stackrel{~}{m}_{Q_2}=\stackrel{~}{m}_{U_2}`$ 10 TeV,

$`\stackrel{~}{m}_{Q_2}=\stackrel{~}{m}_{U_2}`$ $``$ $`{\displaystyle \frac{F}{M}}\mathrm{exp}[(M_{}L\sqrt{\mathrm{ln10}})^2],`$ (17)
$``$ $`{\displaystyle \frac{10^2(M_{}L_c)(M_{}L)}{ML}}\mathrm{exp}[({\displaystyle \frac{3}{2}}\sqrt{\mathrm{ln10}}\sqrt{\mathrm{ln10}})^2]10\mathrm{TeV},`$ (18)
$``$ $`ML(M_{}L_c)\times 10^2\times {\displaystyle \frac{3}{2}}\sqrt{\mathrm{ln10}}\mathrm{exp}[\mathrm{ln10}/4].`$ (19)
Using Eqs. (14, 16, 19), we obtain the sfermion masses for other generations

$`\stackrel{~}{m}_{Q_1}=\stackrel{~}{m}_{U_1}=\stackrel{~}{m}_{D_1}{\displaystyle \frac{F}{M}}\mathrm{exp}[(M_{}L\sqrt{5\mathrm{l}\mathrm{n}10/2})^2]18\mathrm{TeV},`$ (20)
$`\stackrel{~}{m}_{Q_3}=\stackrel{~}{m}_{U_3}{\displaystyle \frac{F}{M}}\mathrm{exp}[9\mathrm{l}\mathrm{n}10/4]100\mathrm{GeV},`$ (21)
$`\stackrel{~}{m}_{D_2}{\displaystyle \frac{F}{M}}\mathrm{exp}[(3\sqrt{\mathrm{ln10}}/2\sqrt{3\mathrm{l}\mathrm{n}10})^2]15.7\mathrm{TeV},`$ (22)
$`\stackrel{~}{m}_{D_3}{\displaystyle \frac{F}{M}}\mathrm{exp}[(3\sqrt{\mathrm{ln10}}/2\sqrt{2\mathrm{l}\mathrm{n}10})^2]17.5\mathrm{TeV}.`$ (23)
We note that these sfermion masses are generated at the compactification scale $`L_c^1`$.

Let us comment on scales in our model. Assuming $`M_{}10^{18}\mathrm{GeV}`$ and $`M_{}L_c10`$,

$`F10^{21}\mathrm{GeV}^2,M5.6\times 10^{16}\mathrm{GeV},L_c^110^{17}\mathrm{GeV}`$ (24)
are obtained. As mentioned above, $`F<M^2`$ is satisfied. It is also interesting that the mass of the additional vector-like superfields is close to the compactification scale.333Since the mass of the vector-like fields $`M`$ is a parameter in four dimensional theory, it has nothing to do with the fundamental Planck scale $`M_{}`$.

Next we will show the case of the solution 2. Instead of the condition (13), we can impose

$`10^2ϵ_2/ϵ_3\mathrm{exp}[(M_{}L)^2(M_{}L\sqrt{3\mathrm{l}\mathrm{n}10})^2],`$ (25)
which leads to

$`M_{}L{\displaystyle \frac{5}{6}}\sqrt{3\mathrm{l}\mathrm{n}10}.`$ (26)
From the gaugino mass constraint, we obtain

$`F(M_{}L_c)\times 10^2M_{}.`$ (27)
On the other hand, $`\stackrel{~}{m}_{D_2}`$10 TeV leads to

$`ML(M_{}L_c)(ML)\times 10^2\mathrm{exp}[\mathrm{ln10}/12].`$ (28)
Other sfermion masses at the compactification scale can be estimated as the same anaylses of the previous case as,

$`\stackrel{~}{m}_{Q_3}=\stackrel{~}{m}_{U_3}100\mathrm{GeV},\stackrel{~}{m}_{D_3}12.1\mathrm{TeV},`$ (29)
$`\stackrel{~}{m}_{Q_2}=\stackrel{~}{m}_{U_2}7.7\mathrm{TeV},\stackrel{~}{m}_{Q_1}=\stackrel{~}{m}_{U_1}=\stackrel{~}{m}_{D_1}11.6\mathrm{TeV}.`$ (30)
As for scales, it is almost the same as the solution 1,

$`F10^{21}\mathrm{GeV}^2,M8.25\times 10^{16}\mathrm{GeV},L_c^110^{17}\mathrm{GeV}.`$ (31)
We give some comments here. First, in our scenario, $`\stackrel{~}{m}_{D_3}`$ is the order of $`𝒪`$(10) TeV, which is somewhat large from the viewpoint of the decoupling solution. This means that the large tan$`\beta `$ case is preferable. Second, it is known that the decoupling scenario generically suffers from a problem that the third generation sfermion mass squareds are driven to be negative through the two-loop RGE effects of the heavy first-two generation sfermion masses . Since the sfermion masses we obtained are generated at the compactification scale, the third generation sfermion mass squareds might be negative at low energy. To avoid this, we have to add the extra fields with negative SUSY breaking mass squareds as discussed in Ref. . We do not discuss this point in detail in this paper.

In summary, we have proposed the decoupling scenario in the fat brane approach. In this approach, Yukawa hierarchy is determined by the overlap of wave functions of the matter fields localized at the different points in extra dimensions. The lighter matter fileds are localized closer to the point where Higgs fields are localized. We introduced the SUSY breaking brane in which the chiral superfield with nonzero F-term is localized and discussed whether the spectrum of the decoupling solution is possible or not, namely, the firsr-two generation sfermion masses are the order of 10 TeV, and the third generation sfermion masses are of order 100 GeV. Naively, sfermion masses generated from the Kähler potential suppressed by the fundamental scale are negligibly small compared to the gaugino masses since sfermion masses receive the additional exponential suppression by the overlap between the wave functions of the matter firlds and of the chiral superfield on the SUSY breaking brane. In order to obtain the spectrum of the decoupling scenario, we have introduced the extra vector-like superfields and Yukawa interaction between the vector-like superfields, MSSM superfields and the chiral superfield with nonzero F-term. Integrating out the vector-like fields leads to the Kähler potential proportional to the suppression factor with the inverse of Yukawa hiearchy and suppressed by the mass of the vector-like fields, which is smaller than the fundamental scale. This becomes the dominant source of the sfermion masses. Two explicit examples of the spectrum have been shown in this paper. It has turned out that the large tan$`\beta `$ case is preferable because $`\stackrel{~}{m}_{D_3}`$ is somewhat large. It is interesting that SUSY breaking brane is predicted to be localized within the thick visible wall unlike the conventional SUSY breaking scenario in extra dimensions.

Finally, we have touched on the negative sfermion mass squareds problem. Since the sfermion masses we have obtained is generated at the compactification scale, a detailed RGE analysis is necessary to see whether the third generation sfermion mass squareds are indeed positive at low energy. We leave this subject for future work.

Acknowledgments

N.H. is supported by the Grant-in-Aid for Scientific Research, Ministry of Education, Science and Culture, Japan (No.12740146) and N.M. is supported by the Japan Society for the Promotion of Science for Young Scientists (No.08557).


text:
# Spin-dependent Transparency of Ferromagnet/Superconductor Interfaces

## Abstract

Because the physical interpretation of the spin-polarization of a ferromagnet determined by point-contact Andreev reflection (PCAR) is non-trivial, we have carried out parameter-free calculations of PCAR spectra based upon a scattering-theory formulation of Andreev reflection generalized to spin-polarized systems and a tight-binding linear muffin tin orbital method for calculating the corresponding scattering matrices. PCAR is found to measure the spin-dependent interface transparency rather than the bulk polarization of the ferromagnet which is strongly overestimated by free electron model fitting.

\]

The spin-polarization of a ferromagnet at the Fermi energy (not to be confused with the total magnetic moment) is a parameter of crucial importance entering all models of magneto-electronic phenomena such as tunneling magnetoresistance (TMR). Although Moodera et al.’s TMR experiments  could be interpreted in terms of bulk spin-polarizations using a very simple model due to Julliere , subsequent experiments on a wider variety of materials cast doubt on the general validity of the model . They also led to the need to find simpler ways of measuring spin-polarizations than the benchmark determinations by Tedrow and Mersevey  using tunneling spectroscopy at ferromagnet$`|`$insulator$`|`$superconductor (S$`|`$I$`|`$F) junctions. These are difficult to make and characterize, especially for complex materials. Andreev reflection spectroscopy at ferromagnet-superconductor (F$`|`$S) interfaces in ballistic point contacts (PCAR) has been proposed as a simpler and complementary method to measure spin polarizations . When the diameter of the point contact is smaller than the mean free path but larger than the Fermi wavelength, transport is governed by the Sharvin resistance , appropriately modified to take into account the interface “transparency”. In the Blonder-Tinkham-Klapwijk (BTK) theory , this is done by representing the interface as a planar $`\delta `$-function barrier between free-electron materials. The Bogoliubov-de Gennes equation  for a normal-superconducting (N$`|`$S) interface is readily solved for this model and the essential physics expressed in terms of a dimensionless barrier strength $`Z`$ and the superconducting energy gap $`\mathrm{\Delta }`$. At the interface to a ferromagnetic material, the Andreev reflection is suppressed by the difference between majority and minority electronic band structures . Using the parabolic-band Stoner-Wohlfarth model for ferromagnetic materials, BTK theory can easily be generalized to spin-polarized systems. With only a single additional parameter $`P`$, the spin polarization of the ferromagnet, a good fit to experimental data is apparently achieved .

In spite of this, doubts have been raised about the definition and physical significance of the parameter $`P`$ derived from different experiments for all but 100% polarized ferromagnets . Apart from concerns about the relevant Fermi surface averaged quantity to be used to describe different transport regimes, transport through an interface should not depend solely on parameters characteristic of bulk materials. The representation of an interface between materials with quite different electronic structures as a $`\delta `$-function barrier can only be justified in terms of analytical simplicity. The free-electron Stoner model does not describe itinerant ferromagnetism in transition metals in an internally consistent manner and the probability of Andreev reflection is overestimated by such models; complex transition metal Fermi surfaces lead to a reduction of the phase space available to an electron to pair with a spin-flipped retro-reflected hole.

In this letter we present a generalized scattering formalism of Andreev reflection  to spin-polarized systems and calculate the PCAR spectra with and without interface disorder. By basing our study on ab-initio calculations we can treat the full complexity of the transition metal band structure without introducing any arbitrary parameters and make contact with experiments for a number of systems of current experimental interest. We will be forced to conclude that the generalized BTK model omits an important aspect of the problem, basically because the bimodal distribution of the transmission probabilities of real interfaces are inaccurately modelled using a $`\delta `$-function potential.

We start with the Bogoliubov-de Gennes  equation for the electron and hole wave functions of a superconductor

$$\left(\begin{array}{cc}H_{0\sigma }& \mathrm{\Delta }\\ \mathrm{\Delta }^{}& H_{0\sigma }^{}\end{array}\right)\left(\begin{array}{c}C_{e\sigma }\\ C_{h\sigma }\end{array}\right)=\epsilon \left(\begin{array}{c}C_{e\sigma }\\ C_{h\sigma }\end{array}\right),$$
(1)
where $`\sigma =\pm 1`$. $`H_{0\sigma }`$ is the single-electron Hamiltonian matrix for majority $`(\sigma =1)`$ and minority $`(\sigma =1)`$ spins. $`C_{e\sigma (h\sigma )}`$ are the coefficient-vectors of the wave functions of electrons (holes) in some convenient basis. The excitation energy $`\epsilon `$ is measured relative to the Fermi energy. Following Beenakker , wave-function matching at the F$`|`$S interface is achieved by inserting between F and S a fictitious region in which S assumes its normal state (disregarding the proximity effect). At the F$`|`$N interface, the scattering of the Bloch states for electrons and holes can be written

$$(\begin{array}{c}C_{e\sigma }^{}(F)\\ C_{e\sigma }^+(N)\\ C_{h\sigma }^+(F)\\ C_{h\sigma }^{}(N)\end{array})=(\begin{array}{cc}S^\sigma (\epsilon )& 0\\ 0& S^\sigma (\epsilon )\end{array})\left(\begin{array}{c}C_{e\sigma }^+(F)\\ C_{e\sigma }^{}(N)\\ C_{h\sigma }^{}(F)\\ C_{h\sigma }^+(N)\end{array}\right)$$
(2)
where $`+()`$ denotes right(left) going waves and

$$S^\sigma (\epsilon )=\left(\begin{array}{cc}r_{11}^\sigma (\epsilon )& t_{12}^\sigma (\epsilon )\\ t_{21}^\sigma (\epsilon )& r_{22}^\sigma (\epsilon )\end{array}\right)$$
(3)
is the normal-state scattering matrix. The subscript 1 refers to F, 2 to N. At the (fictitious) N$`|`$S interface there is only Andreev scattering:

$`C_{e\sigma }^{}(N)`$ $`=`$ $`\alpha C_{h\sigma }^{}(N)e^{i\varphi },`$ (4)
$`C_{h\sigma }^+(N)`$ $`=`$ $`\alpha ^{}C_{e\sigma }^+(N)e^{i\varphi },`$ (5)
where for $`|\epsilon |<\mathrm{\Delta }_0`$, $`\alpha =\mathrm{exp}\left[i\mathrm{arccos}(\epsilon /\mathrm{\Delta }_0)\right]`$ describes the phase shift due to the penetration of the wave function into the superconductor and for $`|\epsilon |>\mathrm{\Delta }_0`$, $`\alpha =(\epsilon \mathrm{sgn}\left(\epsilon \right)\sqrt{\epsilon ^2\mathrm{\Delta }_0^2})/\mathrm{\Delta }_0`$, $`|\alpha |<1`$. The incoming and reflected waves are related on the F side by

$$\left(\begin{array}{c}C_{e\sigma }^{}(F)\\ C_{h\sigma }^+(F)\end{array}\right)=\left(\begin{array}{cc}R_{ee}^\sigma & R_{eh}^\sigma \\ R_{he}^\sigma & R_{hh}^\sigma \end{array}\right)\left(\begin{array}{c}C_{e\sigma }^+(F)\\ C_{h\sigma }^{}(F)\end{array}\right).$$
(6)
$`R_{ee}^\sigma `$ $`=`$ $`r_{11}^\sigma (\epsilon )`$ (7)
$`+`$ $`\alpha ^2t_{12}^\sigma (\epsilon )r_{22}^\sigma (\epsilon ){\displaystyle \frac{1}{1\alpha ^2r_{22}^\sigma (\epsilon )r_{22}^\sigma (\epsilon )}}t_{21}^\sigma (\epsilon )`$ (8)
$$R_{he}^\sigma =\alpha ^{}e^{i\varphi }t_{12}^\sigma (\epsilon )\frac{1}{1\alpha ^2r_{22}^\sigma (\epsilon )r_{22}^\sigma (\epsilon )}t_{21}^\sigma (\epsilon )$$
(9)
are the reflection coefficients for the total system. In terms of the spectral conductance

$$G_{FS}(\epsilon )=\frac{e^2}{h}\underset{\sigma =\pm 1}{}Tr(1R_{ee}^\sigma R_{ee}^\sigma +R_{he}^\sigma R_{he}^\sigma ),$$
(10)
and the Fermi-Dirac distribution function $`f`$, the current-voltage relation is $`eI(V)=𝑑\epsilon [f(\epsilon )f(\epsilon +eV)]G_{FS}(\epsilon )`$. The charge transport in the presence of Andreev reflection is thus expressed in terms of the scattering matrix for the normal state. This expression is the spin polarized generalization of Beenakker’s formula  and reduces to the BTK equivalent if we use the scattering matrix of a $`\delta `$-function potential and the parabolic band model. Note that the expressions also hold for energies above the superconducting gap.

The parameter-free calculation of the reflection matrices  and the conductance (10) is based on the surface Green’s function method  implemented with a tight-binding linear muffin tin orbital basis. Because a minimal basis set is used, we are able to carry out calculations for large lateral supercells and model disorder very flexibly within such supercells without using any adjustable parameters. The electronic structure is determined self-consistently within the local spin density approximation. For disordered layers the potentials are calculated using the layer CPA approximation .

We focus our attention here on the Cu/Pb, Co/Pb and Ni/Pb systems investigated experimentally in  which were grown with an fcc(111) orientation in point contacts. We encounter the practical difficulty that the lattice constants of Pb are quite different from those of Cu, Co and Ni. Lattice matching is modelled using lateral supercells containing $`4\times 4`$ Cu (Co or Ni) atoms on one side of the interface and $`3\times 3`$ Pb atoms/monolayer on the other side or $`7\times 7`$ Cu, Co or Ni and $`5\times 5`$ Pb atoms/monolayer which gives better matching but is computationally more expensive. Disorder is taken into account by introducing a two monolayer-thick 50% interface alloy . The results are not very sensitive to moderate variations in the alloy concentration in this range. Experiments  have been carried out at slightly elevated temperatures to suppress a possible proximity effect. Our computations are therefore carried out using the corresponding values of kT in the Fermi function: 4.2 K for Cu/Pb, 1.4 K for Co/Pb and 2.5 K for Ni/Pb. The only parameter in the calculations is the superconducting energy gap $`\mathrm{\Delta }_{\mathrm{Pb}}`$ which is taken from experiment and is very close to the bulk value.

For an ideal N$`|`$S interface Andreev reflection leads to a doubling of the conductance. In Table I the zero bias values $`G_{FS}(0)/G_{FN}(0)`$, calculated for clean and alloy interfaces, are compared with the experimental values of Ref.. Suppression of the Andreev reflection is apparent in all three cases but is smallest for non-ferromagnetic Cu. Even there, however, bandstructure mismatch between Cu and Pb gives rise to a significant reduction of the interface transparency and suppression of AR. The agreement with experiment is remarkably good and only in the case of Co/Pb is there any indication of a discrepancy. This picture changes when we compare the calculated and experimental spectra in Fig. 1 shown in terms of $`g(V)\left[G_{FS}(V)G_{FN}(0)\right]/G_{FN}(0)`$ where $`G_{FN}(0)`$ is the normal state conductance and $`G_{FS}(V)`$ is the differential conductance for an F$`|`$S interface at finite temperature and finite bias.

Let us first focus on the paramagnetic Cu/Pb point contact, Fig. 1(a). Whereas results for a specular interface deviate significantly from the measured data, the decreased conductance of a rough interface leads to essentially perfect agreement with experiment. Aspects of the problem not taken into account, like spin-orbit scattering and strong-coupling superconductivity in Pb are apparently not very important.

The agreement with the PCAR spectra for the ferromagnetic systems (Figs. 1(b) and (c)) is not satisfactory close to the band edge of quasi-particle excitations, where the theoretical results appear to strongly exaggerate the Andreev reflection . Since we have no parameters with which we could improve the agreement the deviations should be sought in deficiencies of the model. For example, the close proximity of the ferromagnet might smear out the peaked density of states of the superconductor at the band edges . Whereas this point deserves attention, we focus in the following on the zero-bias conductance.

The spin-polarization of the interface conductance, $`P=(G_{maj}G_{min})/(G_{maj}+G_{min})`$, will in general differ from that of the bulk ferromagnetic material. For example, for Ni/Pb, $`P=0.5\%`$ for a specular interface and $`1.4\%`$ for a rough interface. This should be compared to the BTK fit obtained in Ref.  of $`P_{BTK}=32\%`$. Our prediction of a small spin-dependence of the interface conductance could be tested directly by measuring the magnetoresistance of a Ni/Pb/Ni trilayer (or multilayer) in the current-perpendicular-to-the-plane configuration with the layer thickness chosen sufficiently thin that the total resistance is dominated by the interface contribution. Using the expression given by Schep et al. for diffuse systems, we predict interface resistances of $`R_{maj}=0.97f\mathrm{\Omega }m^2`$ and $`R_{min}=1.27f\mathrm{\Omega }m^2`$ for a rough Ni/Pb (111) interface.

We emphasize the difference between the single-channel BTK model and the multichannel scattering theory generalization by focusing on the zero-voltage, zero-temperature conductance (10) for a non-magnetic system:

$$G_{NS}=\frac{4e^2}{h}\underset{n}{}\frac{T_{n}^{}{}_{}{}^{2}}{(2T_n)^2}=\frac{4e^2}{h}𝑑T\rho (T)\frac{T^2}{(2T)^2}$$
(11)
is expressed in terms of the eigenvalues of the transmission matrix, $`T_n`$, and the sum is over the propagating channels . For the Cu/Pb interface, the distribution function of transmission matrix eigenvalues, $`\rho (T)`$, is shown in Fig. 1(d). It turns out to be a bi-modal function with a clear dependence on interface roughness and on the constituent materials. Obviously, $`\rho (T)`$ strongly affects the Andreev reflection probability with the main contribution coming from the highly transmitting part, whereas the average conductance is much more affected by the low transmission channels. Using the Andreev conductance to derive an average transmission is therefore misleading. For a $`\delta `$-function potential and free electrons the transmission probability as a function of the in-plane wave vector has the form $`T_k_{}=1/\left(1+\frac{Z^2}{1(k_{}/k_F)^2}\right)`$. The distribution function for $`Z=0.3`$, the value which must be taken to reproduce the average conductance, is plotted as a dashed curve in Fig. 1(d) and is seen to be qualitatively incorrect compared to the realistic distributions.

For the F$`|`$S interface the situation is more complicated. Since Cooper pairs are phase sensitive, the conductance is governed by the full scattering matrix, ı.e. by the scattering amplitudes and not just the transmission probabilities. We therefore cannot use a single distribution function to characterize the interface because the Andreev reflection probability contains terms like $`1\alpha ^2r_{22}^\sigma (\epsilon )r_{22}^\sigma (\epsilon )`$ which cannot be expressed by $`T_{}`$ and $`T_{}`$. Still, it should be clear that the deficiencies of the $`\delta `$ function barrier potential will not vanish.

In summary we calculate the PCAR spectra of Co/Pb, Ni/Pb, and Cu/Pb point contacts and the spin-polarization of Co/Pb, Ni/Pb interfaces from first principles which allows us to remove what is a dubious approximation made when applying the BTK model to transition metals. We find good agreement with experimental results by the Cornell group for the full spectrum of Cu/Pb and for the zero bias conductance in the ferromagnetic systems. Qualitative differences are revealed by comparing our results with that of a simple BTK model. The physical significance of parameters obtained by fitting experiments with this simple model is doubtful.

This work is supported by the Stichting FOM, the NEDO International Joint Research Grant Program “Nano-magnetoelectronics”, the Grant Agency of the Czech Republic (202/00/0122), the European Commission’s TMR Research Network on “Interface Magnetism” (contract No. FMRX-CT96-0089) and RT Network “Computational Magnetoelectronics” (contract No. HPRN-CT-2000-00143). We acknowledge fruitful discussions with W. Belzig, A. Brataas, Yu. Nazarov, and J. Kudrnovsky and thank R.A. Buhrman for providing us with the experimental data shown in Figure 1.