text:
# 1 Introduction: exponential vs. power–law PR

## 1 Introduction: exponential vs. power–law PR

As is well known any trajectory of a bounded in phase space motion of Hamiltonian system recurs infinitely many times to some neighborhood of its initial position, for both regular (with discrete spectrum) as well as chaotic (with continuous spectrum) motion. These Poincaré recurrences (PR) do not imply a quasiperiodic motion which is still a widespread delusion (see, e.g., ). The difference between regular and chaotic motions lies in the statistics of recurrences which is usually described by the integral distribution $`P(\tau )`$ that is by the probability for a recurrence time to be larger than $`\tau `$. In a regular motion such a survival probability $`P(\tau )`$ has a strict upper bound in $`\tau `$ while for a chaotic motion $`\tau `$ can be arbitrarily long. In both cases PR characterize some fluctuations including arbitrarily large ones in chaotic motion. The PR statistics proved to be a very powerful and reliable method in the studies of chaotic dynamics due to its statistical stability.

To my knowledge, such a method was first used (implicitly) in Ref. for the study of a narrow chaotic layer along the separatrix of a nonlinear resonance. The result ($`\tau 1`$)

$$P(\tau )\frac{1}{\sqrt{\tau }}$$
$`(1)`$
was a surprise as it contradicted the bounded motion in chaotic layer. Indeed, the total sojourn time $`\tau P(\tau )`$ of a trajectory, which is prportional to the measure of the chaotic component of the motion, diverges as $`\tau \mathrm{}`$. Later , this apparent contradiction has been resolved simply by increasing $`\tau `$ which showed that the exponent of the power–law decay also increased from the initial $`\nu =1/2`$ to $`\nu 3/2`$.

It is instructive to mention that the origin of a short–time computation in Ref. was in apparently reasonable decision to avoid any rounding–off errors by enormous increase of the computation accuracy. As a result, the computation speed, and the available motion time, dropped by several orders of magnitude. Generally, for exponentially unstable (chaotic) motion such an approach is prohibited whatever the computer power. Fortunately, it is also unnecessary for calculating statistical characteristics of the motion like $`P(\tau )`$ since most of the latter are robust. True, the corresponding Anosov theorem  was (and can be so far) proved for the very simple Anosov systems only. Moreover, such a theorem is even wrong for discontinuous (discrete) perturbations like rounding–off ones (see, e.g., Refs.). Nevertheless, all the numerical experience confirms a sort of robustness of the statistical behavior of chaotic systems, at least with some minimal precautions (see, e.g., Refs. for discussion). Notice that without such an ’empirical’ robustness the numerical experiments with always approximate models would lose any physical meaning!

A power–law decay $`P(\tau )\tau ^\nu `$, whatever the exponent $`\nu `$, found in  for a bounded motion, was at variance with the exponential decay believed to be a generic case. In Ref. the former was interpreted as a characteristic of a qualitatively new structure of the motion near the chaos border in phase space. Later, it was termed the critical structure, which was described by a renormalization group  (see also review  and references therein).

Since then, the exponential decay has been considered as a property of ergodic chaotic motion without any chaos borders. However, in recent numerical experiments  with an asteroid motion a fairly long transient exponential decay was found. Moreover, it persists in the separatrix map also used, just the same map which seemed to have been well studied in many previous works  (see also  and references therein).

The main purpose of this paper is to reconsider various regimes of PR, and to formulate the conditions for their realization using two relatively simple models: separatrix and standard maps. Only bounded motion will be considered, with or without chaos borders. First, a classical problem of PR in an ergodic system will be discussed in some details in Section 2. Then, in Section 3, the analysis of various PR regimes in the separatrix map will be presented aimed to resolution of the apparent contradiction mentioned above. In Section 4, PR in the standard map in accelerator (microtron) regime will be described. The latter model presents a unique possibility for quantitative study of the global critical structure. Particularly, a new part of this structure has been found which size was surprisingly large. Finally, in Section 5, the main results of the present study are summarized. In addition, the first preliminary empirical evidence is presented for a new regime of Poincaré recurrences including the transition from exponential to exponential statistics.

## 2 PR in ergodic system: standard map

Consider, first, an elementary example of 1D homogeneous diffusion in momentum $`p`$. It can be described by a Gaussian distribution function

$$f_G(p,t)=\frac{\mathrm{exp}\left(\frac{p^2}{2tD}\right)}{\sqrt{2\pi tD}}$$
$`(2)`$
where $`D=<(\mathrm{\Delta }p)^2>/t`$ is the diffusion rate. Derivative $`f_P(p,t)=df_G/dp`$ with boundary condition $`f_P(0,t)=0`$ which obeys the same diffusion equation

$$\frac{f}{t}=\frac{D}{2}\frac{^2f}{p^2}$$
$`(3)`$
describes, then, PR to $`p=0`$. The distribution of recurrence times (1) is simply related to an auxiliary function $`f_G`$ by

$$P(\tau )=A_0^{\mathrm{}}f_P(p,\tau )𝑑p=Af_G(0,\tau )=\frac{A}{\sqrt{2\pi \tau D}}\sqrt{\frac{\tau _0}{\tau +\tau _0}}$$
$`(4)`$
Here $`A`$ is normalizing factor, and parameter $`\tau _0`$ provides a necessary truncation of the preceding diverging expression at small $`\tau `$. It characterizes the dynamical time scale of the diffusion (cf., e.g., free path in molecular diffusion). If the motion in $`p`$ is actually bounded (see below), Eq.(4) describes initial free diffusion.

### 2.1 A little of theory

Now, consider in more details another simple model – the kicked rotator – described by the so–called standard map:

$$\begin{array}{cc}\overline{p}=p+K\mathrm{sin}x\hfill & (\mathrm{mod}L)\hfill \\ \overline{x}=x+\overline{p}p_0\hfill & \end{array}$$
$`(5)`$
on a torus ($`0x<2\pi ,0p<L=2\pi n,n=1,2,\mathrm{}`$).

We seek a solution $`f(p,t)`$ of diffusion equation (3) with the boudary condition

$$f(0,t)=\mathrm{\hspace{0.17em}0}$$
$`(6)`$
which provides a loss of probability because of PR to $`p=0`$ (and to $`p=L`$). The orthogonal and normalized eigenfunctions of the diffusion equation for this problem have the form ($`k1`$ is integer)

$$g_k(p)=\sqrt{\frac{2}{L}}\mathrm{sin}\left(\frac{\pi kp}{L}\right)$$
$`(7)`$
with the corresponding eigenvalues

$$\gamma _k=\left(\frac{\pi k}{L}\right)^2\frac{D}{2}$$
$`(8)`$
which describe the decay rate of the eigenmodes (7). In Eq.(8) the diffusion rate is

$$D(K)=\frac{K^2}{2}C(K)$$
$`(9)`$
with the dynamical correlation function

$$C(K)\mathrm{\hspace{0.17em}1}\mathrm{\hspace{0.17em}2}J_2(K)+\mathrm{\hspace{0.17em}2}J_2^2(K)$$
$`(10)`$
where $`J_2(K)`$ is the Bessel function.

The set of eigenfunctions (7) and eigenvalues (8) provides a general solution of the diffusion equation with boundary condition (6) for an arbitrary initial distribution $`f_0(p)=f(p,0)`$. Peculiarity of PR statistics $`P(\tau )`$ is just in a very particular initial condition. Specifically, for a single trajectory in numerical experiments the recurrence time $`\tau `$ is determined by the two successive crossings of the exit line which is, in the model under consideration, $`p=0\mathrm{mod}L`$. Hence, the initial distribution is concentrated right here: $`|p|K`$. The condition for a trajectory with initial $`p>0`$ to cross the exit line reads: $`p+K\mathrm{sin}x<0`$. Whence, the probability of crossing is proportional to $`\mathrm{arccos}(p/K)`$, and the normalized initial distribution can be taken in the form:

$$f_0\left(\frac{p}{K}\right)=\mathrm{arccos}\left(\frac{p}{K}\right),0<pK<L$$
$`(11)`$
An example of $`f_0`$ is shown in the insert to Fig.1. It is convenient to chose $`f_0(p)`$ on the one side of the exit line which is possible due to the symmetry of eigenfunctions (7).

The difficulty with such an initial condition is in its narrow width which is always comparable with the dynamical scale (both are $`K`$, a single kick). This violates the diffusion approximation for the exact integro–differential kinetic equation. A simple remedy is well known, for example, from the theory of neutron diffusion where the dynamical scale is the transport free path $`l_n`$ (see, e.g., Refs. and , p.689). A simple correction improving the diffusion approximation amounts to a relatively small shift of the boundary condition (6) from $`p_b=0`$ to $`p_b=\alpha l`$ where $`l`$ is the dynamical scale in our problem, and $`\alpha 1`$ is unknown numerical factor to be determined below from the numerical experiments. This implies an increase of the global scale: $`L=L_0+2\alpha l`$ while the initial distribution remains unchanged as it is obtained directly from the dynamics (5). Notice the corresponding change in eigenfunctions (7).

The general solution of the diffusion problem is given by

$$f(p,t)=\underset{k=1}{\overset{\mathrm{}}{}}f_kg_k(p)\mathrm{e}^{\gamma _kt}$$
$`(12)`$
where the expansion coefficients $`f_k`$ are determined by the initial condition (11):

$$f_k=_0^K\frac{dp}{K}\sqrt{\frac{2}{L}}\mathrm{sin}\left(\frac{\pi kK}{L}\left(\frac{p}{K}+\mathrm{\Lambda }\right)\right)\mathrm{arccos}\left(\frac{p}{K}\right)=$$

$$\sqrt{\frac{\pi }{2L}}\frac{\mathrm{cos}(sk\mathrm{\Lambda })(1\mathrm{J}_0(sk))+\mathrm{sin}(sk\mathrm{\Lambda })\mathrm{H}_0(sk)}{sk}$$

$$\frac{\pi sk}{4\sqrt{2L}}\left(1+\frac{8\mathrm{\Lambda }}{\pi }\right)\left(1+𝒪(s^2k^2)\right)$$
$`(13)`$
Here $`\mathrm{H}_0`$ is the Struve function, $`\mathrm{\Lambda }=\alpha l/K1`$ (see below), and $`s=\pi K/L1`$ is a small diffusion parameter. The latter approximate expression in (13) holds true for $`sk\stackrel{<}{}\text{ }1`$.

Now, the PR statistics is described by

$$P(\tau )=_0^Lf(p,\tau )𝑑p=\underset{k=1}{\overset{\mathrm{}}{}}f_k\mathrm{e}^{\gamma _k\tau }\sqrt{\frac{2}{L}}_0^L\mathrm{sin}\left(\frac{\pi kp}{L}\right)𝑑p=$$

$$=\frac{2}{\pi }\sqrt{2L}\underset{m=1}{\overset{\mathrm{}}{}}\frac{f_k}{k}\mathrm{e}^{\gamma _k\tau }$$
$`(14)`$
with $`k=2m1`$ because only odd modes contribute to the integral.

Asymptotically, as $`\tau \mathrm{}`$, PR decay exponentially (Poisson statistics)

$$P(\tau )\frac{s}{2}\left(1+\frac{8\mathrm{\Lambda }}{\pi }\right)\mathrm{e}^{\gamma _1\tau }=F_1\mathrm{e}^{\gamma _1\tau }$$
$`(15)`$
with the characteristic time

$$\tau _1=\frac{1}{\gamma _1}=\frac{2}{\pi ^2}\frac{L^2}{D}$$
$`(16)`$
which is determined by the first (most slow) mode $`m=k=1`$, and which is of the order of the global diffusion time.

The factor $`F_1`$ in Eq.(15) characterizes the share of asymptotic exponential decay which is small in the diffusive regime due to $`s1`$. The main, initial, decay is a power–law one. Again, due to small $`s`$, the sum in Eq.(14) can be approximately replaced by the integral over $`m`$ to obtain:

$$P(\tau )\sqrt{\frac{\tau _0}{\tau }}\sqrt{\frac{\tau _0}{\tau +\tau _0}}$$
$`(17)`$
where

$$\tau _0=\frac{\pi }{16C(K)}\left(1+\frac{8\mathrm{\Lambda }}{\pi }\right)^2$$
$`(18)`$
function $`C(K)`$ is given by Eq.(10), and the approximate expression for $`f_k`$ in Eq.(13) is used. The latter is not applicable for $`\tau 0`$, so that the final expression in Eq.(17) is an approximate truncation of the preceding diverging relation (cf. Eq.(4)).

The power–law/exponential crossover time $`\tau _{cro}`$ is obtained from the comparison of Eqs. (15) and (17), and is given approximately by the relation:

$$\tau _{cro}\frac{L^2}{8\pi D},\gamma _1\tau _{cro}=\frac{\pi }{16}$$
$`(19)`$
Again, in the diffusive regime ($`L^2D`$) the intermediate power–law decay may be very long until the exponential asymptotics is reached.

### 2.2 Numerical experiments

An example of PR in ergodic case is shown in Fig.1. We use the standard map (5) on a torus of sufficiently large circumference $`LK`$ to provide a diffusive relaxation ($`s1`$, for the opposite limit of ballistic relaxation $`s\stackrel{>}{}1`$ see Section 4 below). How strange it may seem, the conditions for ergodicity even in such an apparently ’simple’ model are still unknown ! However, numerical experiments (see, e.g., Ref. indicate that, at least, for a particular value of the parameter $`K=7`$ the share of the regular domains, if any, is negligible ($`\stackrel{<}{}\text{ }10^9`$) besides the two small islets (per map’s period, see Section 4 below). Fortunately, their effect on PR is also negligible because they are related to the accelerator mode in which the momentum $`p`$ quickly moves around the torus, so that a trajectory immediately crosses the exit line $`p=0\mathrm{mod}L`$ (cf. Section 4).

In Fig.1 empirical data for a particular value of $`L=50\pi `$ are shown which corresponds to $`25`$ periods of map (5) in $`p`$. All the data were obtained from the run of a single trajectory over $`10^7`$ iterations. Transition from a power law (straight dashed line) to an exponential (dashed curve) is clearly seen.

For a quantitative comparison with the theory above (Section 2.1) we fix the dynamical parameter $`lK\sqrt{C(K)}=12.1`$ where the value $`C(7)3`$ is used which has been obtained from a special numerical experiment. It considerably differs from the value $`C(7)1.78`$ according to approximate relation (10) just because of accelerator islets mentioned above. Since our model is a map, the minimal empirical recurrence time is $`\tau _{em}=1`$ instead of $`\tau _{th}=0`$ in a continuous theory (for example, in numerical data $`P(1)1`$). The corresponding corrections are negligible except the initial dependence for $`\tau 1`$ (see below).

Numerical data in Fig.1 were fitted to Eq.(15) in the interval $`\tau =5001000`$ iterations, and the empirical values of the characteristic time $`\tau _1=125`$, and of the factor $`F_1=0.26`$ were obtained. The corresponding values of the correction parameter are $`\alpha _\tau =2.3`$, and $`\alpha _F=0.71`$. The difference in these two values of $`\alpha `$ characterizes the accuracy of the correction which is rather poor because of a very narrow initial distribution (see Eq.(11), and discussion around). Without correction ($`\alpha =0`$) the theoretical values would be: $`\tau _1=69`$, and $`F_1=0.07`$ which both are substantially underestimated.

For a more systematic study the similar numerical data were computed for a number of $`L`$ values specified by the integer $`n=L_0/2\pi `$. The results are shown in Fig.2.

Dependence $`\tau _1(n)`$ is well described by the uncorrected relation (16) for large $`n`$ as expected. In intermediate region ($`n3`$) the agreement is father improved by the correction which provides a smooth transition to the ballistic limit (see Eqs.(28) and (29) in Section 4). In other words, the correction is not very important for the asymptotic decay rate because it is determined by the first eigenfunction which is only slightly disturbed, for large $`n`$, by the shift of the boundary. This is no longer the case for the amplitude $`F_1`$ which strongly depends just on the distorted region near the boundary $`p=0`$. As a result the correction is most important for large $`n`$. The dependence $`F_1(n)`$ in the intermediate region remains unclear. For $`n\stackrel{>}{}10,s\stackrel{<}{}\text{ }1/3`$ both relations, Eqs.(15) and (16), are in a reasonable agreement with the numerical data for the same average value of the correction parameter $`<\alpha >=1.22`$.

Coming back to Fig.1, we see that the initial power–law decay is well described by a simple relation (17) with $`\tau _0=1`$ which is shown by the dashed straight line, and which would correspond to $`\alpha _0=0.66\alpha _F`$.

## 3 PR with a chaos border: separatrix map

Now we consider an opposite limit of essentially nonergodic system with a large chaos border and the critical structure. As an example we take the separatrix map which was studied in many papers (see, e.g., Refs.), and for which a new regime of PR has been recently observed . The latter was the main motivation for the present studies. We take the separatrix map in the form :

$$\begin{array}{cc}\overline{p}=p+\mathrm{sin}x\hfill & \\ \overline{x}=x\lambda \mathrm{ln}(|\overline{p}|)p_0\hfill & \end{array}$$
$`(20)`$
Here the motion is always strictly confined to the so–called chaotic layer: $`|p|p_b(x)`$. Previously, the most studied case corresponded to big parameter $`\lambda 1`$. In this limit $`p_b\lambda `$, so that the width of the layer ($`2\lambda `$) is much larger than the dynamical scale of the diffusion (a single ’kick’) which, for map (20), is unity (cf. Eq.(5)). Besides the critical structure along the two borders, the average diffusion rate within the layer is nearly constant (see Eq.(10)):

$$<D>\frac{<C(\lambda /|p|)>}{2}\frac{1}{2}$$
$`(21)`$
Hence, the initial decay of PR is a simple power law (1) which was observed, indeed, from the beginning  (Section 1). The crossover time to a diferent law is given by a simple diffusion estimate:

$$\tau _{cro}\frac{p_b^2}{D}\lambda ^2$$
$`(22)`$
Unlike the ergodic case, the asymptotics of PR in the presence of chaos border is also a power law but with a different exponent $`\nu 1.5`$. This is explained by a very specific critical structure near the border where the diffusion rate rapidly drops. As a result no trajectory can ever reach the exact border, even though it is approaching, from time to time, the border arbitrarily close (see Refs.  for details).

An example of this well known behavior is shown in Fig.3 (upper solid curve). A transition between the two different power laws (dashed straight lines) at $`\tau 100`$ is clearly seen in agreement with estimate (22). There is no sign of any exponential decay. Now, how does it appear in a similar model ?

The first observation is that in application to celestial mechanics (dynamics of asteroids) the parameter $`\lambda `$ of map (20) is typically rather small: $`\lambda 1`$ . This drastically changes the structure of the layer. First of all, the layer width is reduced down to the size of a single kick. An example is shown in Fig.4. Hence, the diffusion approximation becomes inapplicable. Instead, the so–called ballistic relaxation comes into play which is much quicker. In other words, a slow diffusive motion from the exit line to a critical structure is replaced now by rapid jumps of a trajectory over the whole layer with some probability to get into the critical structure. Since those jumps are very irregular in a chaotic layer the PR are expected to decay exponentially. This is the case indeed as an example in Fig.3 demonstrates (lower solid curve, $`\lambda =1`$ ). The exponential decay can be intermediate only as the trajectory is evetually captured into the critical structure, and the decay turns to a power law. Generally, the initial part of the power law is an approximate relation in that its exponent is not universal, and is even varying with $`\tau `$. In the latter example $`\nu 1.1`$ which is rather different from $`\nu 1.5`$ for the upper curve in Fig.3.

Another interesting and important question is how long is the intermediate exponential? For the lower curve in Fig.3 it is rather long: $`\tau _{cro}150`$ which corresponds to the PR crossover as low as $`P_{cro}10^6`$ ! However, under different conditions with the same $`\lambda =1`$ the exponential is much shorter: $`\tau _{cro}50`$, and $`P_{cro}10^2`$. The difference is in the exit line as shown in Fig.4.

In the latter case the exit line is usual: $`p=0`$. The critical structure is determined by the two big islands comparable in size with that of the whole layer. This entails a rapid capture of a trajectory into the critical structure, and a fast transition to a final power law (with the local exponent $`\nu 1.65`$). The lower curve in Fig.3 corresponds to the same $`\lambda =1`$ but to a different exit line:

$$p_{ex}=\mathrm{cos}(x)$$
$`(23)`$
It is chosen in such a way to cut through both stability islands and, thus, to suppress any sticking to their critical structure. Then, the final power law is determined by the critical structure at the layer borders which is apparently very narrow and cannot be discerned ’by eye’ in Fig.4. Nevertheless, it does exist as the asymptotic power law of PR in Fig.3 proves. Moreover, the latter even allows us to estimate the size of the critical structure: its relative area (with respect to that of the layer) is $`A_{cr}4\times 10^4`$, or the width $`(\mathrm{\Delta }p)_{cr}10^3`$ (see Section 4). This exponential transient is well fitted by the relation similar to (15) (up to $`\tau \tau _{cro}`$) with $`\tau _112`$, and $`F_10.32`$. Both values are in a surprisingly well agreement with the uncorrected theory ($`\alpha =0,L4.5`$, see Fig.4) which gives $`\tau _18.2`$, and $`F_10.35`$. Apparently, this is because the diffusion parameter $`s0.71`$ is still not large enough.

Now, we can summarize the conditions for the transient exponential in PR for a nonergodic motion: (i) fast, ballistic, relaxation, and (ii) a small measure of the regular domains. Besides, it turns out that the exponential PR allow for, at least, some estimates of that measure. A more quantitative study of this interesting relation is convenient to continue with the standard map again. This is because the latter has an infinite series of the special values of parameter $`K=K_n2\pi n`$ for which there are well studied islands of regular motion with a simple scaling and of rapidly decreasing area.

## 4 PR in microtron: the standard map again

The main advantage of this microtron model is in that it is very simple, especially for numerical experiments, and well studied already. Here we are interested primarily in the domains of regular motion which exist for an infinite series of the special values of parameter $`K=K_n2\pi n`$ where $`n>0`$ is any integer. Within these domains (islands) $`|p|`$ grows indefinitly proportional to time which is the so–called microtron acceleration. It was well studied since the celebrated paper due to Veksler in 1944 (see, e.g.,  and references therein). However, in the present paper, as well as in Ref., the main object for study is not the regular acceleration itself but rather the chaotic motion outside the microtron islands which is generally affected by the critical structure at the island borders. A picture of this scale–invariant border is shown in Fig.5a in dimensionless variables

$$x_s=(xx_0)K,p_s=(pp_0)K$$
$`(24)`$
where $`p_0`$ is a parameter of map (5), and

$$K\mathrm{sin}x_0=\mathrm{\hspace{0.17em}2}\pi n,K^2=\sigma ^2+(2\pi n)^2,\sigma =K\mathrm{cos}x_0,4<\sigma <\mathrm{\hspace{0.17em}0}$$
$`(25)`$
The latter inequalities determine the stability region around a fixed point
$`\pm x_0,p_0\mathrm{mod}2\pi `$. In Fig.5a and below $`\sigma =2`$ (the center of stability). For each integer $`n`$ there are two islets per phase space bin $`2\pi \times \mathrm{\hspace{0.17em}2}\pi `$ one of which is presented in Fig.5a. The picture shows a single trajectory of $`5000`$ iterations. During this time interval the trajectory is sticking to the critical structure very close to the exact chaos border which results, under particular conditions (see below), in an asymptotic power–law decay of PR (cf. Fig.3 above). The island relative area (with respect to that of the phase–space bin) is given also by a dimensionless relation :

$$A_nK_n^2=A_0(\sigma )\mathrm{\hspace{0.17em}0.17}$$
$`(26)`$
where the latter value corresponds to $`\sigma =2`$. This area rapidly decreases as island’s number $`n`$ grows. Yet, for any $`n\mathrm{}`$ it determines the asymptotic PR decay, as we shall see below.

In Fig.5b another, much smaller, microtron island is shown for comparison. In this case an outside, and much longer, trajectory was used which cannot ever cross the chaos border and enter the island. Its area is given by the same estimate (26) with $`A_0(3.1)0.0038`$.

The main difficulty with the microtron model for our purposes here is the rapid growth in $`|p|`$ within and around the chaos border. This destroys any long sticking of a trajectory whatever the exit line for PR (cf. Section 3, Fig.4). To overcome this difficulty we used the following method. First, we have chosen the exit line in such a way not to cut any island. It was done simply by fixing parameter $`p_0=\pi (0\mathrm{mod}2\pi )`$ in map (5) without any change in the configuration of map’s torus. Second, we compensated acceleration by adding the term $`2\pi n`$ to the first Eq.(5). This helps, of course, for one island of each pair only.

Now, we need to provide the ballistic regime of relaxation that is a sufficiently large parameter $`s=\pi K/L`$ (Section 2.1). It is convenient to take $`L=2\pi n`$, so that the parameter

$$s=\frac{\pi K_n}{2\pi n}\pi $$
$`(27)`$
is nearly independent of $`n`$ except a few small values of the latter.

Neglecting any dynamical correlations of the motion (particularly, those caused by the presence of small microtron islets including the compensation of acceleration) it is straightforward to calculate the probability $`w`$, per map’s iteration, for a trajectory to stay within the torus without crossing the exit line. As is easily verified, it is given by the relation:

$$w=\frac{dxdp}{2\pi L/4}=\frac{2}{\pi L}_0^{x_m}𝑑x(LK\mathrm{sin}x)=\frac{2}{\pi }\left(x_m\frac{s}{\pi }(1\mathrm{cos}x_m)\right)$$

$$1\frac{2}{\pi }\mathrm{\hspace{0.17em}0.363}$$
$`(28)`$
where

$$x_m=\{\begin{array}{ccc}\mathrm{arcsin}\left(\frac{\pi }{s}\right)& ,\hfill & s\pi \hfill \\ \pi /2& ,\hfill & s\pi \hfill \end{array}$$
$`(29)`$
These general relations were used in Section 2.2 (Fig.2) to draw the ballistic approximation.

The latter expression in Eq.(28) corresponds to the value $`s=\pi `$ used in numerical experiments. Without additional shift $`\mathrm{\Delta }p=2\pi n`$ discussed above the average time of the exponential decay would be

$$<\tau >=\frac{1}{\mathrm{ln}w}\mathrm{\hspace{0.17em}0.988}$$
$`(30)`$
For $`s=\pi `$ the shift increases $`w`$ and $`<\tau >`$ up to

$$\stackrel{~}{w}=\frac{w}{2}+\frac{1}{2}=\mathrm{\hspace{0.17em}1}\frac{1}{\pi }\mathrm{\hspace{0.17em}0.682}$$

$$\stackrel{~}{<\tau >}=\frac{1}{\mathrm{ln}\stackrel{~}{w}}\mathrm{\hspace{0.17em}2.61}$$
$`(31)`$
Now we can turn to numerical experiments with this microtron model.

### 4.1 PR in microtron: numerics

The main results of numerical experiments are presented in Fig.6, and in the Table below. In Fig.6 the points show numerical data computed from a single trajectory (for each $`n`$) up to $`3\times 10^{11}`$ iterations (for the largest $`n=5000`$). The straight solid line is the fitted intermediate exponential with the decay time $`<\stackrel{~}{\tau }>=2.41`$ in a good agreement with the expected theoretical value $`2.61`$ in Eq.(31). This justifies neglecting dynamical correlations assumed in the above theory in ballistic regime.

The exponential/power–law crossover time systematically increases with $`n`$ that is with the decrease of the microtron island area (see Table). The power–law tails of PR were fitted by the expression

$$P_n(\tau )=\frac{A_{cr}(n)}{\tau ^\nu }$$
$`(32)`$
Remarkably, all values of the exponent were found to be close: $`\nu 2`$. The relation of this expression to the size of the critical structure is based on the following hypothesis: dependence (32), fitted to the tail of PR, can be extrapolated back to $`\tau =1`$. If true, it allows us to interpret the parameter $`A_{cr}(n)`$ as the relative area of the whole (global) critical structure around the corresponding microtron island of area $`A_n`$.

One could expect that both areas are comparable: $`A_{cr}(n)A_n`$. Surprisingly, this is not the case (see Table, third column; the data in fourth column will be discussed below). Their ratio $`R=A_{cr}/A_n100`$ is not only very large but also slowly increasing with $`n`$ according to the following approximate empirical relation

$$R(n)\mathrm{\hspace{0.17em}50}n^{1/4}$$
$`(33)`$
The origin of this small correction to a simple scaling $`R`$ const remains unclear.

| Table. Global critical structure around microtron islets | | | |
| --- | --- | --- | --- |
| $`n`$ | $`A_n`$ | $`R`$ | $`R_{ex}`$ |
| islet number | islet area | Fig.6 | Fig.7 |
| 1 | $`4.30\times 10^3`$ | 36 | 25.0 |
| 5 | $`1.71\times 10^4`$ | 124 |  |
| 10 | $`4.31\times 10^5`$ | 65 | 13.9 |
| 100 | $`4.31\times 10^7`$ | 202 | 10.6 |
| 500 | $`1.72\times 10^8`$ | 291 |  |
| 1000 | $`4.31\times 10^9`$ | 176 | 9.8 |
| 5000 | $`1.72\times 10^{10}`$ | 461 |  |
| $`10^4`$ | $`4.31\times 10^{11}`$ |  | 8.6 |
| $`10^5`$ | $`4.31\times 10^{13}`$ |  | 8.3 |
| $`10^6`$ | $`4.31\times 10^{15}`$ |  | 8.1 |
| $`10^7`$ | $`4.31\times 10^{17}`$ |  | 7.7 |

In any event, the size of the whole critical structure seems to be much larger than expected. This main outer part of the structure looks ergodic, and forms a sort of halo around the usually narrow inner part with a typical admixture of chaotic and regular components of motion. The former reminds the ergodic critical structure around a parabolic fixed point, that is the limiting case of an island of zero size, studied in Ref.. In a sense, such a halo is some ’hidden’ critical structure, without internal chaos borders but with apparently strong correlations in the motion which keep a trajectory within this relatively small domain.

Now, the principal question to be aswered reads: is the observed halo a real physical structure or the result of a wrong interpretation of the empirical data using the above extension hypothesis?

### 4.2 Exit time statistics

To clarify this question a new series of numerical experiments was undertaken. To this end, the exit times from the halo, instead of recurrences, were measured. Such a method was recently successfully used in the studies of the critical structure in Ref.. In the problem under consideration here the measurement of exit times was organized as follows. A number (typically 100) of trajectories with the initial conditions homogeneously distributed over the circle around a microtron island (see Fig.5a) were run until they leave the interval ($`0<x<\pi `$). The dependence of the average exit time $`\tau _{ex}`$ for a series of the circles with increasing radius $`\rho _s`$ as a function of the area within a circle $`A_s=\pi \rho _s^2`$ (in scaled variables (24)) was thus computed. The minimal circle of radius $`\rho _s=3`$ touches the island, and comprises the area $`A_{min}=28.3`$ while island’s area in these units is $`A_{sn}=6.72`$, the minimal ratio being $`R_{ex}=A_s/A_{sn}=4.21`$.

The main results of this measurement are shown in Fig.7 for 8 different values of parameter $`n`$ up to $`n=10^7`$ with the island area as small as $`A_n4\times 10^{17}`$ ! This is completely out of reach for the PR method (cf. Ref.). The difference is in a rather short exit time from the halo, we are interested in, as compared to the long recurrence time on the tail where it is eventually separated from the exponential (Fig.6).

The main result revealed in Fig.7 is a transition between the two different scalings. One, for relatively large $`\tau _{ex}`$, is the standard critical scaling shown by the dashed line which is the fitting of numerical data to the relation

$$R_{ex}(\tau _{ex})\frac{R_{ex}(1)}{\tau _{ex}}$$
$`(34)`$
with $`R_{ex}(1)50`$. As expected, this part of the data does not depend on $`n`$. Moreover, scaling (34) is in a good agreement with the PR tail in Fig.6. The relation between the two is well known . Generally, the power–law PR statistics is descibed by (cf. Eq.(32)):

$$P(\tau )\frac{P(1)}{\tau ^\nu }A(\tau )\frac{<\tau >}{\tau }$$
$`(35)`$
where $`<\tau >`$ is the average PR time, and the latter expression is obtained from the ergodicity within the chaotic component of the motion. Using approximate relation $`\tau 2\tau _{ex}`$  we obtain

$$R_{ex}(1)\frac{A(2)}{A_n}\frac{R}{<\tau >\mathrm{\hspace{0.17em}2}^{\nu 1}}$$
$`(36)`$
where $`R=A_{cr}/A_nA(1)/A_n`$ (see above). For integer map’s time

$$<\tau >=\underset{\tau =1}{\overset{\mathrm{}}{}}\tau \left(\frac{1}{\tau ^\nu }\frac{1}{(\tau +1)^\nu }\right)=\zeta (\nu )$$
$`(37)`$
where $`\zeta `$ is the Riemann function. Whence, for $`\nu =2`$ the relation (36) gives (see Eq.(34))

$$R\frac{\pi ^2}{3}R_{ex}(1)\mathrm{\hspace{0.17em}160}$$
$`(38)`$
which is in a reasonable agreement with numerical data (third column in Table). However, unlike the data in Fig.6 where the actual power–law scaling is not seen under much larger exponential transient, the data in Fig.7 clearly demonstrate that the critical scaling does not reach the limit $`\tau =1`$ assumed above. Moreover, the crossover time $`\tau _{cro}`$ increases, and hence the size of the global critical structure ($`R_{ex}`$) decreases, as $`n`$ grows (Fig.7, insert). The increase of $`\tau _{cro}`$ must have an upper bound because otherwise the critical structure near the chaos border would be also destroyed in contradiction to the detailed studies of that in anomalous diffusion . Indeed, the empirical dependence $`\tau _{cro}`$ in Fig.7 (insert) can be fitted reasonably well by the expression

$$\tau _{cro}=\mathrm{\hspace{0.17em}6.69}\frac{4.66}{n^{0.172}}$$
$`(39)`$
The upper limit in $`\tau _{cro}`$ corresponds, according to Eq.(34), to the lower limit in the ratio $`R_{ex}>7.47`$. Combining Eqs. (34) and (39) we obtain approximately

$$R_{ex}(n)\frac{R_{ex}(1)}{\tau _{cro}(n)}$$
$`(40)`$
The empirical values of this dependence are given in Table (fourth column). They indicate much smaller, yet still a fairly large, size of the critical structure as compared to the limiting estimate for PR (third column). The former seems to be more reliable and realistic. A different, new and unknown, scaling in Fig.7 for $`\tau _{ex}<\tau _{cro}`$ requires further studies. What is of importance here is the termination of the critical scaling at a finite $`\tau _{ex}=\tau _{cro}`$. This determines the outer border of the critical structure.

## 5 Discussion: a new puzzle

The original motivation of these studies was the unusual exponential transient observed in PR in the presence of chaos border . However, in the course of investigating the mechanism and conditions of this phenomenon a more interesting observation has come out. It suggests the existence of a new, unknown to my knowledge, part of the critical structure surrounding, like a halo, the well–known inner part close to the chaos border. In spite of some contradictory empirical evidence the halo apparently occupies the most of the global critical structure. In any event, in the microtron model considered in this paper the area of the halo is much larger than that of the regular island inside it, even according to the minimal estimates (see Table and Fig.7).

As is well known, the scaling of the peripheral part of the critical structure is generally nonuniversal, at least quantitatively, in the sense of the corresponding power–law exponents, for example . However, it might be nevertheless typical qualitatively as it appears in our model. In this respect, it would be interesting to look at different examples of the global crtical structure. One possibility is to use the same model with a fixed parameter $`L=2\pi n`$ (Section 4) but for different values of the stability parameter $`\sigma `$ in Eq.(25). First preliminary numerical experiments have been done for 9 values of $`\sigma `$ within the whole stability interval ($`4\sigma 0;2\pi K7.45`$) including the ’quasiergodic’ case $`K=7`$ used in Section 2 for other purposes. In all cases but the latter the PR behavior was similar to that in the main series of numerical experiments (Fig.6), at least qualitatively. However, just for $`K=7`$ a sudden surprise has emerged which is presented in Fig.8.

In spite of a very long run ($`10^{11}`$ iterations) no clear sign of the expected power–law decay is seen. A small deviation from the final exponential at the end of the dependence is a typical feature due to a poor statistics (cf., e.g., Fig.6). The first exponential is close to the expected one with the fitted decay time $`<\stackrel{~}{\tau }>=2.41`$ as compared to the theoretical $`<\stackrel{~}{\tau }>=2.38`$ (see Section 4). For the second exponential the empirical decay time $`<\stackrel{~}{\tau }>=23.1`$ is about 10 times longer. This means that a trajectory is kept within (sticks to?) a certain domain but not in a way it does so in the usual critical structure. Moreover, the relative area $`A_d2\times 10^4`$ of this peculiar domain, estimated similarly to $`A_{cr}(n)`$ in Eq.(32), is small and is comparable with that of the island inside (Fig.5b): $`A_77.8\times 10^5`$. This island does have a chaos border, yet contrary to usual behavior, it does not produce any appreciable power–law decay of PR. Another preliminary remark is that a more careful inspection of Fig.5b seems to suggest a different, more regular than usual, structure of the chaos border for $`K=7`$ (cf. Fig.5a). Certainly, this ’anomaly’ deserves further investigation.

Acknowledgements. I am grateful to I.I. Shevchenko for many interesting discussions and important remarks. This work was partially supported by the Russia Foundation for Fundamental Research, grant 97–01–00865.


format: text

text:
# 1 Origin of Neutrino Mass

## 1 Origin of Neutrino Mass

Recent neutrino experiments , which are most naturally explained by neutrino oscillations, are strong indications that neutrinos are massive and mix with one another. So where do neutrino masses come from? In the minimal Standard Model, neutrinos appear only as left-handed fields in three electroweak doublets $`(\nu _i,l_i)_L`$, where $`i=1,2,3`$ is the family index. Charged leptons $`l_i`$ have right-handed components which are singlets, but not neutrinos. Hence neutrinos are massless two-component fermions, as long as there is no physics beyond the minimal Standard Model. Otherwise, there may be an effective dimension-5 operator

$$\frac{1}{\mathrm{\Lambda }}(\nu _i\varphi ^0l_i\varphi ^+)(\nu _j\varphi ^0l_j\varphi ^+),$$
(1)
where $`(\varphi ^+,\varphi ^0)`$ is the standard Higgs doublet and $`\mathrm{\Lambda }`$ is a large mass, which yields a nonvanishing Majorana neutrino mass matrix as $`\varphi ^0`$ acquires a nonzero vacuum expectation value in the spontaneous breaking of the $`SU(2)\times U(1)`$ electroweak gauge symmetry.

All models of neutrino mass with the same low-energy particle content as that of the minimal Standard Model differ only in the way the above effective operator is realized . The most well-known such model by far is the canonical seesaw model , where three right-handed neutrino singlets with large Majorana masses are added. This amounts to inserting a heavy internal fermion line between the two singlet factors of Eq. (1). The corresponding diagram can be read off to obtain the neutrino mass matrix as

$$(_\nu )_{ij}=\frac{f_if_j\varphi ^0^2}{m_N},$$
(2)
where $`f_i`$ are the Yukawa couplings linking $`\nu _i`$ with the heavy singlet $`N`$ with mass $`m_N`$. Obviously, we need three $`N`$’s to obtain three naturally small seesaw masses for $`\nu _e`$, $`\nu _\mu `$, and $`\nu _\tau `$. On the other hand, other mechanisms are available , the simplest alternative being the addition of a heavy scalar triplet . This is easily recognized if we rewrite Eq. (1) as

$$\frac{1}{\mathrm{\Lambda }}[\nu _i\nu _j\varphi ^0\varphi ^0(\nu _il_j+l_i\nu _j)\varphi ^0\varphi ^++l_il_j\varphi ^+\varphi ^+],$$
(3)
hence an insertion of the heavy scalar triplet

$$\xi =(\xi ^{++},\xi ^+,\xi ^0)$$
(4)
into the above yields a neutrino mass matrix given by

$$(_\nu )_{ij}=\frac{2f_{ij}\mu \varphi ^0^2}{m_\xi ^2},$$
(5)
where $`f_{ij}`$ are the Yukawa couplings of $`\xi `$ to the lepton doublets and $`\mu `$ is its coupling (with the dimension of mass) to the scalar doublets. Note that Eq. (5) can also be interpreted as due to $`\xi ^0`$, i.e.

$$\xi ^0=\frac{\mu \varphi ^0^2}{m_\xi ^2}<<m_\xi .$$
(6)
This shows explicitly that the vacuum expectation value of a heavy scalar field may in fact be very small.

## 2 Leptogenesis

Given that lepton number is not conserved in models of Majorana neutrino masses, the attractive possibility exists that a primordial lepton asymmetry may be created in the early Universe, which then gets converted into the present observed baryon asymmetry through the $`B+L`$ violating, but $`BL`$ conserving interactions of the electroweak sphalerons . In the canonical seesaw model, this is accomplished  by the decays of $`N_i`$. Being Majorana fermions, $`N_{1,2}`$ may decay into either $`l^{}\varphi ^+`$ with $`L=1`$ or $`l^+\varphi ^{}`$ with $`L=1`$. Consider now the amplitude for $`N_1l^+\varphi ^{}`$. It is the sum of 3 terms: the obvious tree graph, the one-loop vertex correction with

$$N_1l^{}\varphi ^+\varphi ^{}l^+$$
(7)
through the exchange of $`N_2`$, and the one-loop self-energy correction with

$$N_1l^{}\varphi ^+N_2l^+\varphi ^{}.$$
(8)
Let this amplitude be denoted by $`A+iB`$, where $`B`$ is the absorptive part, i.e. from putting the intermediate state $`l^{}\varphi ^+`$ on the mass shell. Then the lepton asymmetry from the decay of $`N_1`$ is proportional to

$$|A+iB|^2|A^{}+iB^{}|^2=4Im(A^{}B).$$
(9)
This means that $`CP`$ violation is essential and that the presence of a different $`N`$, i.e. $`N_2`$, in the loop is necessary for leptogenesis.

Since 1995, there has been a resurgence of activity  in this topic. Consider $`N_{1,2}`$ and the mass matrix linking $`(\overline{N}_{1L},\overline{N}_{2L},\overline{N}_{1R}\overline{N}_{2R})`$ with $`(N_{1L},N_{2L},N_{1R},N_{2R})`$:

$$_N=\left[\begin{array}{cc}0& \\ \stackrel{~}{}& 0\end{array}\right],$$
(10)
where

$`=\left[\begin{array}{cc}M_1+H_{11}& H_{12}\\ H_{12}& M_2+H_{22}\end{array}\right],\stackrel{~}{}=\left[\begin{array}{cc}M_1+H_{11}& \stackrel{~}{H}_{12}\\ \stackrel{~}{H}_{12}& M_2+H_{22}\end{array}\right],`$ (15)
and

$`H_{12}`$ $`=`$ $`[M_1{\displaystyle \underset{\alpha }{}}h_{\alpha 1}^{}h_{\alpha 2}+M_2{\displaystyle \underset{\alpha }{}}h_{\alpha 1}h_{\alpha 2}^{}]\left[g_{\alpha 12}^{disp}{\displaystyle \frac{i}{32\pi }}\right],`$ (16)
$`\stackrel{~}{H}_{12}`$ $`=`$ $`[M_1{\displaystyle \underset{\alpha }{}}h_{\alpha 1}h_{\alpha 2}^{}+M_2{\displaystyle \underset{\alpha }{}}h_{\alpha 1}^{}h_{\alpha 2}]\left[g_{\alpha 12}^{disp}{\displaystyle \frac{i}{32\pi }}\right],`$ (17)
$`H_{jj}`$ $`=`$ $`2M_j{\displaystyle \underset{\alpha }{}}|h_{\alpha j}|^2\left[g_{\alpha jj}^{disp}{\displaystyle \frac{i}{32\pi }}\right].`$ (18)
Note that $`\stackrel{~}{H}_{12}H_{12}^{}`$ because of the absorptive contribution, hence $`_N`$ is not hermitian. This is analogous to $`\overline{K}K`$ or $`\overline{B}B`$ mixing when the decay of the particles is also taken into account. The self-energy correction to $`_N`$ implies that mass eigenstates need not be $`CP`$ eigenstates, i.e. indirect $`CP`$ violation as exemplified by the $`ϵ`$ parameter. The vertex correction to $`_N`$ implies that $`CP`$ may be violated in the decay itself, i.e. direct $`CP`$ violation as exemplified by the $`ϵ^{}`$ parameter. $`N_1N_2`$ oscillations may also occur. However, all these things happen in an expanding Universe, i.e in a dense, hot medium which is changing with time, hence the exact details are complicated and are still being worked out .

The primordial lepton asymmetry is generated from the decay of the lightest $`N`$, i.e. $`N_1`$:

$$\delta \frac{G_F}{2\pi \sqrt{2}}\frac{1}{(m_D^{}m_D)_{11}}\underset{j=2,3}{}Im(m_D^{}m_D)_{1j}^2\frac{M_1}{M_j},$$
(19)
where $`m_D`$ is the Dirac mass matrix linking $`\nu `$ with $`N`$, and $`M_1<<M_j`$ has been assumed. This expression  can then be used to study neutrino masses and mixing from atmospheric  and solar  neutrino oscillations and to extract information  on $`N_1`$. The range $`10^910^{13}`$ GeV for $`M_1`$ is found to be consistent with $`n_B/g_{}n_\gamma 10^{10}`$. Since the self-energy contribution is proportional to $`M_1/(M_2M_1)`$, $`\delta `$ may be enhanced if $`M_2M_1<<M_1`$, but the limit $`M_2M_10`$ is not singular , as it is bounded by the decay width of $`N_1`$.

## 3 Triplet Higgs Model

If neutrino masses come from heavy triplet scalars , then the mixing of $`\xi _1`$ and $`\xi _2`$ through their absorptive parts, i.e. self-energy contributions, leads to the physical mass eigenstates $`\psi _1`$ and $`\psi _2`$ which are not $`CP`$ eigenstates. Their decay asymmetries are given by

$$\delta _i\frac{Im[\mu _1\mu _2^{}\underset{k,l}{}f_{1kl}f_{2kl}^{}]}{8\pi ^2(M_1^2M_2^2)}\frac{M_i}{\mathrm{\Gamma }_i},$$
(20)
and for $`M_{1,2}10^{13}`$ GeV, realistic neutrino masses and leptogenesis are possible.

## 4 Neutrino Masses from $`R`$ Parity Violation

I now come to my third example which is the generation of neutrino masses through $`R`$ parity violation in supersymmetry . The well-known superfield content of the Minimal Supersymmetric Standard Model (MSSM) is given by

$`Q_i=(u_i,d_i)_L(3,2,1/6),`$ (21)
$`u_i^c(3^{},1,2/3),`$ (22)
$`d_i^c(3^{},1,1/3),`$ (23)
$`L_i=(\nu _i,l_i)_L(1,2,1/2),`$ (24)
$`l_i^c(1,1,1);`$ (25)
$`H_1=(\overline{\varphi }_1^0,\varphi _1^{})(1,2,1/2),`$ (26)
$`H_2=(\varphi _2^+,\varphi _2^0)(1,2,1/2).`$ (27)
Given the above transformations under the standard $`SU(3)\times SU(2)\times U(1)`$ gauge group, the corresponding superpotential should contain in general all gauge-invariant bilinear and trilinear combinations of the superfields. However, to forbid the nonconservation of both baryon number ($`B`$) and lepton number ($`L`$), each particle is usually assigned a dicrete $`R`$ parity:

$$R(1)^{3B+L+2j},$$
(28)
which is assumed to be conserved by the allowed interactions. Hence the MSSM superpotential has only the terms $`H_1H_2`$, $`H_1L_il_j^c`$, $`H_1Q_id_j^c`$, and $`H_2Q_iu_j^c`$. Since the superfield $`\nu _i^c(1,1,0)`$ is absent, $`m_\nu =0`$ in the MSSM as in the minimal Standard Model. Neutrino oscillations  are thus unexplained.

Phenomenologically, it makes sense to require only $`B`$ conservation (to make sure that the proton is stable), but to allow $`L`$ violation (hence $`R`$ parity violation) so that the additional terms $`L_iH_2`$, $`L_iL_jl_k^c`$, and $`L_iQ_jd_k^c`$ may occur. Note that they all have $`\mathrm{\Delta }L=1`$. Neutrino masses are now possible  with Eq. (1) realized in two ways. From the bilinear terms

$$\mu H_1H_2+ϵ_iL_iH_2,$$
(29)
a $`7\times 7`$ neutralino-neutrino mass matrix is obtained:

$$\left[\begin{array}{ccccc}M_1& 0& g_1v_1& g_1v_2& g_1u_i\\ 0& M_2& g_2v_1& g_2v_2& g_2u_i\\ g_1v_1& g_2v_1& 0& \mu & 0\\ g_1v_2& g_2v_2& \mu & 0& ϵ_i\\ g_1u_i& g_2u_i& 0& ϵ_i& 0\end{array}\right],$$
(30)
where $`v_{1,2}=\varphi _{1,2}^0/2`$ and $`u_i=\stackrel{~}{\nu }_i/2`$, with $`i=e,\mu ,\tau `$. Note first the important fact that a nonzero $`ϵ_i`$ implies a nonzero $`u_i`$. Note also that even if $`u_i/ϵ_i`$ is not the same for all $`i`$, only one linear combination of the three neutrinos gets a tree-level mass . From the trilinear terms, neutrino masses are also obtained, but now as one-loop radiative corrections. Note that these occur as the result of supersymmetry breaking and are suppressed by $`m_d^2`$ or $`m_l^2`$.

## 5 $`L`$ Violation and the Universe

As noted earlier, the $`R`$ parity violating interactions have $`\mathrm{\Delta }L=1`$. Furthermore, the particles involved have masses at most equal to the supersymmetry breaking scale, i.e. a few TeV. This means that their $`L`$ violation together with the $`B+L`$ violation by sphalerons  would erase any primordial $`B`$ or $`L`$ asymmetry of the Universe . To avoid such a possibility, one may reduce the relevant Yukawa couplings to less than about $`10^7`$, but a typical minimum value of $`10^4`$ is required for realistic neutrino masses. Hence the existence of the present baryon asymmetry of the Universe is unexplained if neutrino masses originate from these $`\mathrm{\Delta }L=1`$ interactions. This is a generic problem of all models of radiative neutrino masses where the $`L`$ violation can be traced to interactions occuring at energies below $`10^{13}`$ GeV or so.

## 6 Leptogenesis from $`R`$ Parity Violation

Once the notion of $`R`$ parity violation is introduced, there are many new terms to be added in the Lagrangian. Some may be responsible for realistic neutrino masses and may even participate in the erasure of any primordial $`B`$ or $`L`$ asymmetry of the Universe, but others may be able to produce a lepton asymmetry on their own which then gets converted into the present observed baryon asymmetry of the Universe through the sphalerons. A recent proposal  shows how this may happen in a specific model.

Consider the usual $`4\times 4`$ neutralino mass matrix in the $`(\stackrel{~}{B},\stackrel{~}{W}_3,\stackrel{~}{h}_1^0,\stackrel{~}{h}_2^0)`$ basis:

$$\left[\begin{array}{cccc}M_1& 0& sm_3& sm_4\\ 0& M_2& cm_3& cm_4\\ sm_3& cm_3& 0& \mu \\ sm_4& cm_4& \mu & 0\end{array}\right],$$
(31)
where $`s=\mathrm{sin}\theta _W`$, $`c=\mathrm{cos}\theta _W`$, $`m_3=M_Z\mathrm{cos}\beta `$, $`m_4=M_Z\mathrm{sin}\beta `$, and $`\mathrm{tan}\beta =v_2/v_1`$. The above assumes that $`ϵ_i`$ and $`u_i`$ are negligible in Eq. (26), which is a good approximation because neutrino masses are so small. We now choose the special case of

$$m_3,m_4<<M_2<M_1<\mu .$$
(32)
As a result, the two higgsinos $`\stackrel{~}{h}_{1,2}^0`$ form a heavy Dirac particle of mass $`\mu `$ and the other two less heavy Majorana fermion mass eigenstates are

$`\stackrel{~}{B}^{}`$ $``$ $`\stackrel{~}{B}+{\displaystyle \frac{sc\delta r_1}{M_1M_2}}\stackrel{~}{W}_3+\mathrm{},`$ (33)
$`\stackrel{~}{W}_3^{}`$ $``$ $`\stackrel{~}{W}_3{\displaystyle \frac{sc\delta r_2}{M_1M_2}}\stackrel{~}{B}+\mathrm{},`$ (34)
where $`\delta =M_Z^2\mathrm{sin}2\beta /\mu `$, and

$$r_{1,2}=\frac{1+M_{1,2}/\mu \mathrm{sin}2\beta }{1M_{1,2}^2/\mu ^2}.$$
(35)
We now observe that whereas $`\stackrel{~}{B}`$ couples to both $`\overline{l}_L\stackrel{~}{l}_L`$ and $`\overline{l}_L^c\stackrel{~}{l}_L^c`$, $`\stackrel{~}{W}_3`$ couples only to $`\overline{l}_L\stackrel{~}{l}_L`$ because $`l_L^c`$ is trivial under $`SU(2)_L`$. On the other hand, $`R`$ parity violation implies that there is $`\stackrel{~}{l}_Lh^{}`$ mixing as well as $`\stackrel{~}{l}_L^ch^+`$ mixing. Therefore, both $`\stackrel{~}{B}^{}`$ and $`\stackrel{~}{W}_3^{}`$ decay into $`l^\pm h^{}`$ and may be the seeds of a lepton asymmetry in such a scenario.

Let the $`\stackrel{~}{l}_Lh^{}`$ mixing be very small (which is a consistent assumption for realistic neutrino masses from bilinear $`R`$ parity violation). Then $`\stackrel{~}{W}_3^{}`$ decays only through its $`\stackrel{~}{B}`$ component. Hence the decay rate of the LSP (Lightest Supersymmetric Particle), i.e. $`\stackrel{~}{W}_3^{}`$, is very much suppressed, first by $`\delta `$ and then by the $`\stackrel{~}{l}_L^ch^+`$ mixing which will be denoted by $`\xi `$. Our construction is aimed at satisfying the out-of-equilibrium condition:

$$\mathrm{\Gamma }(\stackrel{~}{W}_3^{}l^\pm h^{})<H=1.7\sqrt{g_{}}(T^2/M_{Pl})$$
(36)
at the temperature $`TM_2`$, where $`H`$ is the Hubble expansion rate of the Universe with $`g_{}`$ the effective number of massless degrees of freedom and $`M_{Pl}`$ the Planck mass. This implies

$$\left(\frac{\xi |\delta |r_2}{M_1M_2}\right)^2\frac{1}{M_2}<1.9\times 10^{14}\mathrm{GeV}^1,$$
(37)
where we have used $`g_{}=10^2`$ and $`M_{Pl}=10^{18}`$ GeV.

The lepton asymmetry generated from the decay of $`\stackrel{~}{W}_3^{}`$ has both vertex and self-energy loop contributions from the insertion of $`\stackrel{~}{B}^{}`$. However, the coupling of $`\stackrel{~}{B}^{}`$ to $`l^\pm h^{}`$ is suppressed only by $`\xi `$ and not by $`\delta `$, thus a realistic asymmetry may be established if $`\xi `$ is not too small. Let $`xM_2^2/M_1^2`$, then the decay asymmetry of $`\stackrel{~}{W}_3^{}`$ is given by

$$ϵ=\frac{\alpha \xi ^2}{2\mathrm{cos}^2\theta _W}\frac{Im\delta ^2}{|\delta |^2}\frac{\sqrt{x}g(x)}{1x},$$
(38)
where

$$g(x)=1+\frac{2(1x)}{x}\left[\left(\frac{1+x}{x}\right)\mathrm{ln}(1+x)1\right].$$
(39)
The phase of $`\delta `$ comes from the relative phase between $`M_1`$ and $`M_2`$.

To make sure that at $`TM_2`$, the $`L`$ violating processes $`l^\pm h^{}l^{}h^\pm `$ through $`\stackrel{~}{B}^{}`$ exchange do not erase $`ϵ`$, we require

$$\left(\frac{2e^2\xi ^2}{\mathrm{cos}^2\theta _W}\right)^2\frac{1}{M_1^2}\frac{T^3}{32\pi }\frac{f(x)}{(1x)^2}<H$$
(40)
at $`TM_2`$, where

$`f(x)=1+{\displaystyle \frac{2(1x)}{x^2}}[(1+3x)\mathrm{ln}(1+x)x(1+x)],`$ (41)
which implies

$$\frac{\xi ^4}{M_2}\frac{xf(x)}{(1x)^2}<2.6\times 10^{14}\mathrm{GeV}^1.$$
(42)
A sample solution is

$`M_1=3\mathrm{TeV},{\displaystyle \frac{\delta }{M_1M_2}}=8.3\times 10^4,`$
$`M_2=2\mathrm{TeV},\xi =2\times 10^3.`$ (43)
In that case,

$$ϵ=3.6\times 10^8Im\delta ^2/|\delta |^2,$$
(44)
and

$$\frac{n_B}{g_{}n_\gamma }\frac{ϵ}{3g_{}}10^{10}\frac{Im\delta ^2}{|\delta |^2}.$$
(45)
Hence realistic leptogenesis is possible if $`\xi 10^3`$ can be obtained.

The origin of $`\stackrel{~}{l}_L^ch^+`$ mixing in $`R`$ parity violation is usually the term $`H_1\stackrel{~}{L}\stackrel{~}{l}^c`$, which is very small because $`\stackrel{~}{\nu }`$ has to be very small. To obtain $`\xi 10^3`$, we need to add the nonholomorphic  term $`H_2^{}H_1\stackrel{~}{l}^c`$ which is generally unconstrained.

## 7 Conclusion

Given a mechanism for generating small Majorana neutrino masses, it is often a bonus to find that leptogenesis is possible at the same time. In the canonical seesaw and triplet Higgs models, the same new physics is responsible for both. In $`R`$ parity nonconserving supersymmetry, they may come from different sectors of the theory.

## Acknowledgments

I thank Enrico Nardi and William Ponce and the other organizers for their great hospitality at Cartagena. This work was supported in part by the U. S. Department of Energy under Grant No. DE-FG03-94ER40837.


text:
# Spin-splitting in GaAs 2D holes

## I Introduction

Zero-magnetic-field ($`B=0`$) spin-splitting of energy bands has long been a subject of considerable experimental and theoretical effort because of its role in the electronic properties of two-dimensional (2D) systems . It concerns details of energy band structure that are of fundamental interest. The eigenstates in a system that has inversion symmetry in both space and time are spin-degenerate. If the system is made spatially inversion asymmetric, this degeneracy is lifted by spin-orbit interaction for all momentum $`\stackrel{}{k}0`$, even in the absence of an external magnetic field. In 2D systems, we can create this spatial inversion asymmetry simply by making the potential which confines the carriers to 2D asymmetric.

In this work, we study the $`B=0`$ spin-splitting in GaAs 2D holes confined to a quantum well (QW) grown on the (311)A surface of a GaAs substrate. Section II provides some experimental details. In Section III, we show that by using two gates on the back and front of the sample to control the electric field applied perpendicular to the QW, the asymmetry and therefore the spin-splitting can be tuned while the carrier density is kept constant. In the process we discover that the low-$`B`$ Shubnikov-de Haas oscillations are not as simply related to the $`B=0`$ spin-splitting as has often been presumed .

Tuning the spin-splitting without changing the density allows us to isolate the effect of spin-splitting on various phenomena. In recent years there has been significant interest in the surprising metallic behavior of various high-quality 2D systems . In Section IV we show the effect of spin-splitting on the metallic behavior of the 2D holes in GaAs .

Finally, motivated by recent theoretical predictions , we have applied an in-plane $`B`$ to these samples and studied the magnetoresistance (Section V). We have found that, remarkably, the in-plane $`B`$ depopulates the upper spin-subband at a rate which depends on the relative orientations of $`B`$ and the crystal axes! This is a reflection of the anisotropic band structure and Zeeman splitting in the system. Furthermore, when the upper spin-subband is depopulated, the metallic behavior turns to insulating behavior. We also see an in-plane magnetoresistance anisotropy that depends on the relative orientations of $`B`$ and the current $`I`$; the finite layer thickness of our 2D system may be responsible for this observed anisotropy .

## II Experimental

A schematic cross section of our samples is shown in Fig. 1A. The samples are grown, via molecular beam epitaxy (MBE), on GaAs (311)A substrates. Each sample contains a 20 nm-wide GaAs QW flanked by AlGaAs barrier layers which are modulation doped with Si. On GaAs (311)A substrates, Si is normally incorporated as an acceptor, leading to a high-quality 2D hole system in the GaAs QW. We make ohmic contacts to the 2D holes using In:Zn and alloying in a reducing atmosphere. As shown in Fig. 1A, our samples have back and front gates to change the density of the 2D hole system and also tune the symmetry of the confining QW potential. The back gate is made by contacting a layer of In left over on the back of the substrate from the MBE growth, and the front gate is evaporated metal (Ti-Au).

GaAs 2D systems grown on (311)A substrates have an intrinsic mobility anisotropy due to an anisotropic interface roughness . The mobility for current along the $`[\overline{2}33]`$ direction can be up to several times the mobility for current along the $`[01\overline{1}]`$ direction. For this reason, we did our measurements on $`L`$-shaped Hall bars, which allow simultaneous measurement of the resistivity $`\rho `$ for current along these two directions (Fig. 1B).

Measurements were done in dilution and 3He refrigerators, with base temperatures of 30 mK and 0.3 K respectively. The systems are fitted in a superconducting magnet capable of magnetic fields up to 16 T. The 3He refrigerator has a tilting stage for the sample that can be rotated in-situ, so the relative angle between the magnetic field and 2D hole system plane can be varied.

Figure 2 shows an example of data for one of the samples from both mobility directions, measured simultaneously at 30 mK . The inset to A shows the low-$`B`$ data. The main part of A demonstrates the high quality of the sample. The inset to B shows raw $`T`$-dependence of $`\rho `$ data, and the main part of B shows it scaled to $`\rho _0`$, the resistivity at 30 mK.

## III Tuning Spin-splitting

We use the gates to tune the asymmetry of the QW by applying an electric field ($`E_{}`$) perpendicular to its plane. Figure 1C schematically demonstrates the procedure. We set the front gate ($`V_{FG}`$) and back gate ($`V_{BG}`$) voltages, and measure $`\rho `$ as a function of $`B`$. Then, at a small $`B`$, $`V_{FG}`$ is increased and the change in the hole density is measured from the change in the Hall coefficient. $`V_{BG}`$ is then reduced to recover the original density. This procedure changes $`E_{}`$ while maintaining the same density to within 3%, and allows calculation of the change in $`E_{}`$ from the way the gates affect the density. These steps are repeated until we have probed the range of $`V_{FG}`$ and $`V_{BG}`$ that are accessible without causing gate leakage . Increasing the magnitude of $`E_{}`$ increases the asymmetry of the sample, which increases the spin-splitting.

Figure 3A shows examples of the Shubnikov-de Haas (SdH) oscillations measured for current along the $`[01\overline{1}]`$ direction in our 20 nm QW, following the procedure described above to change $`E_{}`$ from about 5 kV/cm (pointing towards the front gate) in the top panel through to about -6 kV/cm in the bottom panel. Beating can be clearly seen in the oscillations of all traces except the center trace where $`E_{}0`$. Fast Fourier transforms (FFTs) of $`\rho `$ vs. $`B^1`$ quantify the frequencies present. Figure 3B shows the FFTs of the SdH oscillations at all of the measured sets of $`V_{FG}`$ and $`V_{BG}`$. These frequencies have long been thought to be directly proportional to the spin-subband densities ($`p_\pm `$) following :

$$p_\pm =\frac{e}{h}f_{SdH}^\pm .$$
(1)
With this interpretation it is clear in Fig. 3B that we are tuning the spin-splitting through a minimum as we change $`E_{}`$ from 5 kV/cm to -6 kV/cm: the two peaks, corresponding to two spin-subband densities, get closer together as $`E_{}`$ approaches zero, finally merging, only to separate again as $`E_{}`$ is decreased away from zero. Below, we show that Eq. 1 is only approximately true.

In order to compare the theory of spin-splitting to the experimental results, we have performed self-consistent subband calculations that have no adjustable parameters (further details are in Refs.  and ). These calculations produce, for a given 2D hole density and $`E_{}`$, both the $`B=0`$ spin-subband densities and the SdH oscillations. It is important to note that both of these results are produced by the same Hamiltonian, so they are directly comparable. These calculations yield a surprising result: FFTs of the calculated SdH oscillations show that the frequencies present, when converted to spin-subband densities using Eq. 1, do not agree with the calculated $`B=0`$ spin-subband densities .

This is highlighted in Fig. 4 which shows the calculated $`B=0`$ spin-subband densities (right axis), and the peak positions of FFTs ($`f_{SdH}`$) of both the calculated and measured SdH oscillations (left axis), as a function of $`E_{}`$ . The two $`y`$-axes of the figure are scaled by Eq. 1, so the sets of data can be directly compared. The calculated $`f_{SdH}`$ consistently underestimate the $`B=0`$ spin-splitting. The agreement between the theoretical and measured $`f_{SdH}`$, however, is very good. Especially noteworthy are the results for $`E_{}0`$. As the $`B=0`$ calculations in Fig. 4 indicate, even at $`E_{}=0`$, there should be finite spin-splitting because of the inversion asymmetry of the GaAs (zincblende) crystal structure. However, both the measured and calculated FFTs show only one peak for $`1E_{}1`$ kV/cm .

We believe that the inaccuracy of Eq. 1 is due to a breakdown of Onsager’s argument , which is based on Bohr-Sommerfeld quantization of the semiclassical motion of the Bloch carriers. The presence of spin-orbit interaction makes the system inherently quantum-mechanical, so a semiclassical picture fails. This finding is quite general in that one can expect deviations from a semiclassical picture for all systems with strong spin-orbit interaction. However, the full quantum-mechanical calculations predict that Eq. 1 can be quite accurate for some 2D systems grown on high-symmetry crystal directions, while being inaccurate for other crystal directions or other systems .

## IV Metallic behavior and spin-splitting

For many years, it was widely accepted that there can be no metallic phase in a 2D carrier system . However, recent experiments on several different high quality 2D systems have provided us with reason to re-visit this belief, as they showed that $`\rho `$ has a metallic-like temperature dependence: at very low temperatures $`\rho `$ decreases with decreasing $`T`$ . Various mechanisms  have been proposed to explain the metallic behavior, but no clear model has emerged which quantitatively describes the sizeable body of experimental data.

By tuning the spin-splitting at constant density, we isolate the effect of spin-splitting on the $`T`$-dependence of $`\rho `$ . Figure 3C shows $`\rho `$ as a function of $`T`$ for each of the measured $`E_{}`$. The traces are lined up with the corresponding FFTs in Fig. 3B. From this data it is evident that the magnitude of the change in $`\rho `$ from 30 mK to $``$0.7 K is correlated with the spin-splitting. The 30 mK, $`B=0`$ resistivity ($`\rho _0`$) for each trace is listed on the $`y`$-axis.

In order to characterize the $`T`$-dependence data in a simple way, we calculate $`\mathrm{\Delta }\rho ^T/\rho _0`$, the fractional change in resistivity from 30 mK to $``$0.7 K. This is plotted in Fig. 5 vs. the spin-subband population difference calculated at $`B=0`$ ($`\mathrm{\Delta }p_s`$). We have repeated this type of experiment on lower density samples, down to a density of $`2.5\times 10^{10}`$ cm-2. For densities from $`7.0\times 10^{10}`$ to $`3.3\times 10^{11}`$ cm-2, for which the $`\rho `$ vs. $`T`$ traces have qualitatively the same shape as in Fig. 3C, we have calculated the $`B=0`$ spin-subband densities. As shown in Fig. 5, data for these densities exhibit the same trend: $`\mathrm{\Delta }\rho ^T/\rho _0`$ is larger at larger $`\mathrm{\Delta }p_s`$, and $`\mathrm{\Delta }\rho ^T/\rho _0`$ is more sensitive to spin-splitting at the lower measured densities.

The $`T`$-dependence of $`\rho `$ for a density of $`2.5\times 10^{10}`$ cm-2 has a qualitatively different shape from the higher-density data: it exhibits a local maximum. Increasing $`\mathrm{\Delta }p_s`$ moves this local maximum to lower $`T`$. We note that previous experiments on the density-dependence of the metallic behavior show that at high densities $`\rho `$ monotonically increases with $`T`$, and that as the density is reduced and the transition to insulating behavior is approached, a local maximum appears in $`\rho `$ vs. $`T`$ . Qualitatively, in the density range we have measured, the effect of increasing spin-splitting on shape of $`\rho `$ vs. $`T`$ is similar to the effect of reducing density.

For the highest density, $`3.3\times 10^{11}`$ cm-2, we have direct experimental support for the calculated $`\mathrm{\Delta }p_s`$. At the lower densities the sample quality is typically worse and the spin-splitting is too small so that two frequencies are not resolved in the SdH oscillations. We use the calculations described in Section III to determine the expected $`\mathrm{\Delta }p_s`$ used in Fig. 5 from the $`E_{}`$. Even if a large error in $`\mathrm{\Delta }p_s`$ is allowed for, the conclusions of the previous paragraphs are still valid.

Tuning $`E_{}`$ tunes the spin-splitting, but can also affect the mobility and cause changes in $`\rho _0`$. However, a careful examination of all our data reveals that the changes in $`\rho _0`$ are not causing the changes in $`\mathrm{\Delta }\rho ^T/\rho _0`$. The variation of $`\mathrm{\Delta }\rho ^T/\rho _0`$ due to $`E_{}`$ does not correlate with the changes in $`\rho _0`$. One example of this is in Fig. 3C, where there are three traces at different $`E_{}`$ that have different $`\mathrm{\Delta }\rho ^T/\rho _0`$, but the same $`\rho _0`$ (58.3 $`\mathrm{\Omega }`$/sq.).

All of the data we have presented so far have been from the low-mobility $`[01\overline{1}]`$ arm of the $`L`$-shaped Hall bar. The SdH oscillations in the high-mobility $`[\overline{2}33]`$ direction data are very similar to those in the $`[01\overline{1}]`$ traces, and the FFTs show that the frequencies present, as expected, are the same. The $`T`$-dependence data, while qualitatively similar, are different along $`[\overline{2}33]`$. In a given measurement, $`\mathrm{\Delta }\rho ^T/\rho _0`$ is always smaller for the $`[\overline{2}33]`$ direction than for the $`[01\overline{1}]`$ direction. An example of this behavior can be seen in Fig. 2B.

In summary, we have found that in the density regime where $`\rho `$ increases monotonically with increasing $`T`$, the magnitude of the change in $`\rho `$, $`\mathrm{\Delta }\rho ^T/\rho _0`$, is correlated with the spin-splitting of the 2D hole system. As the density is reduced, $`\mathrm{\Delta }\rho ^T/\rho _0`$ becomes larger and more sensitive to increased spin-splitting. In the density regime where $`\rho `$ has a local maximum, increasing the spin-splitting moves the maximum to lower $`T`$. We also find that the direction of the current in the sample plays a surprising role: the higher-mobility direction shows a smaller $`\mathrm{\Delta }\rho ^T/\rho _0`$.

## V In-plane magnetic field

We also employed an in-plane magnetic field, using the tilting stage of the 3He refrigerator, to probe the Zeeman splitting in GaAs 2D holes . Similar magnetoresistance (MR) measurements have been recently reported for 2D systems that exhibit a low-$`T`$ metallic behavior . In our measurements we have discovered a remarkable anisotropy in the effect of $`B`$ on the spin-subbands, pointing out that the Zeeman splitting in this 2D system is very anisotropic.

The measurements were done on a sample similar to those described above, but the densities were lower still. As with previous samples, we used both front and back gates to control the density. The longitudinal and Hall MRs were first measured with a perpendicular $`B`$, and then the sample was tilted 90o and $`\rho `$ as a function of in-plane $`B`$ was measured. The measurements were made once with the sample mounted with the $`[\overline{2}33]`$ direction parallel to the tilt axis. Then the sample was warmed up, re-mounted with the $`[01\overline{1}]`$ axis parallel to the tilt axis, and the measurements were repeated. Note that this sample also had an $`L`$-shaped Hall bar geometry as shown in Fig. 1B, so for each orientation of the sample relative to the in-plane $`B`$, $`\rho `$ was measured for both current directions simultaneously.

In order to make the data easier to track, we scale the traces to the $`B=0`$, base-$`T`$ value of $`\rho `$. Figure 6 shows such data for various densities, organized by the relative orientations of $`B`$ and the crystal axes. All traces have an overall positive MR and, in addition, show a broad feature: there is an inflection point followed by a reduction in slope, followed by another inflection point beyond which the traces curve upwards again. To highlight this behavior, the arrows in Fig. 6 are placed between the two inflection points, at a $`B`$ we will refer to as $`B^{}`$. Surprisingly, for each density, $`B^{}`$ for the $`B[\overline{2}33]`$ traces is about 4 T smaller than for the $`B[01\overline{1}]`$ traces, regardless of the $`I`$ direction. Also, $`B^{}`$ becomes smaller as the density is reduced. Figure 6 reveals that the relative orientations of $`B`$ and the crystal axes play an important role in the MR features.

The existence of the MR features around $`B^{}`$ is intriguing. Similar, though sharper, features have been observed in in-plane $`B`$ measurements in systems with multiple confinement subbands when a subband is depopulated . Similarly, the features in our data may be related to a spin-subband depopulation and the resulting changes in subband mobility and inter-subband scattering as the in-plane $`B`$ is increased. In support of this hypothesis, numerical calculations similar to the ones described above show that the spin-subband depopulation happens at much lower $`B`$ for $`B[\overline{2}33]`$ than for $`B[01\overline{1}]`$ . This reflects the highly anisotropic spin-subband structure of the 2D hole systems in GaAs (311)A QWs. Our hypothesis is further supported by the observation (Fig. 6) that while $`B^{}`$ clearly depends on the orientation of $`B`$ with respect to the crystal axes, it is rather insensitive to the current direction: spin-subband depopulation should not depend on the direction of current in the sample.

At higher in-plane $`B`$, beyond the MR features around $`B^{}`$, the data in Fig. 6 are qualitatively similar. The traces for $`BI`$ have greater slope than the corresponding traces with $`BI`$, regardless of crystal axes. In this regime the magnetic confinement can become comparable to the electric confinement, and the effects due to the finite-thickness of the 2DHS may be dominant. Indeed, Ref.  predicts that MR with in-plane $`B`$ should be significantly larger for $`BI`$ than for $`BI`$, in agreement with our highest $`B`$ data.

Figure 7 shows the $`T`$-dependence of the MR at a density of $`3.9\times 10^{10}`$ cm-2, for the four measured relative orientations of $`B`$, $`I`$, and the crystal axes. For each panel, the traces exhibit a nearly $`T`$-independent magnetic field $`B_T`$ which occurs near the trace’s first inflection point. This is consistent with the data of Ref. . For $`B<B_T`$, the data show metallic behavior, and for $`B>B_T`$, insulating behavior. $`B_T`$ is different in each panel and, similar to $`B^{}`$, depends much more strongly on the orientation of the crystal axes relative to $`B`$ than on the orientation of $`I`$ relative to $`B`$. Our experiments indicate that $`B^{}`$ and $`B_T`$ depend very similarly on the parameters of our systems ($`p`$, $`E_{}`$, direction of $`B`$). Our observation, which is in agreement with the in-plane MR data of Ref. , strongly suggests that the metallic behavior is linked to the presence of two populated spin-subbands .

## VI Summary

We have demonstrated tunable $`B=0`$ spin-splitting at constant density in a GaAs 2D hole system. In the process, we have discovered that systems with a significant spin-orbit interaction show a more complicated relationship between the $`B=0`$ spin-subband populations and the frequencies present in the Shubnikov-de Haas oscillations than had previously been expected. Using the tunability of the spin-splitting to investigate its effect on the metallic behavior observed in this 2D system, we find that changing the magnitude of the spin-subband population difference changes the $`T`$-dependence of $`\rho `$.

Through the use of an in-plane magnetic field, we have measured a surprising anisotropy of the subband structure of (311)A GaAs 2D holes. When a magnetic field is applied in the plane of the 2D system parallel to the $`[\overline{2}33]`$ direction, the upper spin-subband is depopulated at a significantly lower field than if the field is applied parallel to the $`[01\overline{1}]`$ direction. Furthermore, we observe that the $`B=0`$ metallic behavior turns into insulating near $`B`$ at which the upper spin-subband depopulates.

Finally, we note that Das Sarma and Hwang have recently reported calculations aiming to explain the $`T`$-dependence of the resistivity  and the in-plane MR  of 2D systems that exhibit metallic behavior at finite $`T`$. Their calculations, which include only charged impurity scattering and the orbital motion, qualitatively reproduce some of the experimental data. We wish to point out that our results reveal the importance of the spin degree of freedom, and suggest that for an understanding of the experimental data it is important to also consider a scattering mechanism involving the spin-subbands, perhaps intersubband scattering . Also important for (311)A GaAs 2D holes is the inclusion of interface roughness scattering: both the $`T`$-dependence of $`\rho `$ at $`B=0`$ (Fig. 2b), as well as the in-plane MR data (Fig. 6), depend on the direction of the current in the crystal.

This work was supported by the ARO and the NSF.


text:
# A Comptonization Model for Phase Lag Variability in GRS 1915+105

## 1. Introduction

The galactic microquasar GRS 1915+105, discovered in 1992 by GRANAT (Castro-Tirado, Brandt & Lund castro:1992 (1992)), is a transient source with an extremely wide variety of variability modes (Greiner et al. greiner:1996 (1996); Belloni et al. belloni:1997 (1997); Muno, Morgan & Remillard muno:1999 (1999)). During the period 1996-1997 it has been extensively and repeatedly observed by the Rossi X-ray Timing Explorer (RXTE), alternating phases of dramatic variability to phases of remarkably regular behavior (Belloni et al. belloni:2000 (2000)).

GRS 1915+105 was the first galactic source to show superluminal radio expansion (Mirabel & Rodrìguez mirabel:1994 (1994)), usually interpreted in terms of emission from gas expanding at relativistic velocities. Measurements of this superluminal motion allow to estimate the source distance, which is in the range 6–12.5 kpc (Rodrìguez & Mirabel rodrig:1999 (1999); Fender et al. fender:1999 (1999)). Because of the high extinction along the line of sight, only a faint infrared counterpart has been discovered (Mirabel et al. mirabeletal:1994 (1994)). As a consequence, neither the mass function nor the orbital period are known. The source is believed to host a black hole because it shows properties (variability and superluminal motion) very similar to the other known galactic microquasar, GRO 1655-40, for which a precise dynamical measurement yields $`M7M_{}`$ (Orosz & Bailyn oro:1997 (1997)).

Quasi Periodic Oscillations (QPOs) have been observed in GRS 1915+105 in an incredibly large interval of frequencies, ranging from about 0.001 to 67 Hz (Morgan, Remillard & Greiner morgan:1997 (1997); Chen, Swank & Taam chen:1997 (1997)). In particular, in some observations a QPO with variable centroid frequency ($`\nu _{QPO}0.5`$–10 Hz), whose value and root mean square amplitude correlate with the thermal flux, is detected (Markwardt, Swank & Taam mark:1999 (1999); Trudolyubov, Churazov & Gilfanov trudo:1999 (1999)).

The X-ray spectrum of GRS 1915+105 shows all the distinctive properties of other Galactic Black Hole Candidates (BHCs). A power-law tail with variable luminosity and slope is present at high energies, while a thermal component often dominates at low energies (Belloni et al. belloni:1997 (1997); Muno, Morgan & Remillard muno:1999 (1999)). A detailed spectral classification of the different states of GRS 1915+105 has been recently presented by Belloni et al. (belloni:2000 (2000)). According to these authors, the variability of the X-ray emission is associated with the alternation of three basic states, called A, B and C (see also Markwardt, Swank & Taam mark:1999 (1999); Muno, Morgan & Remillard muno:1999 (1999)). During states A and B, corresponding to relatively high and soft flux, a strong thermal disk component dominates the emission. In state C, the 3–25 keV flux is characterized by lower luminosity and harder X-ray colors, being dominated by a hard ($`\mathrm{\Gamma }=1.5`$–2.5) power law with much softer (or even unobservable) thermal component. State C can last uninterruptedly, for a long time (weeks to months, the so called plateau intervals), in which case no soft component is observed at all. State C is always associated with the 0.5–10 Hz QPO and a flat-top component in the Power Density Spectrum (Markwart, Swank & Taam mark:1999 (1999); Trudolyubov, Churazov & Gilfanov trudo:1999 (1999); Reig et al. reig:2000 (2000)).

The analysis of a set of RXTE observations of GRS 1915+105 in this state showed that phase lags are present, both in the continuum and the QPO (Reig et al. reig:2000 (2000)). The sign of the lag is correlated with QPO frequency, count rate and X-ray colors. The lag is positive (photons in the 5–13 keV band arrive later than those in the 2–5 keV band) when the spectrum is hard and the centroid of the 0.5–10 Hz QPO is below 2 Hz; it then becomes negative (soft photons are delayed) as the spectrum softens and the QPO centroid frequency increases above 2 Hz. The smooth correlations of the lags with colors suggests that a unique mechanism is driving the observed transition from positive to negative lags.

In this letter we propose a simple comptonization model to account for the observed properties of GRS 1915+105 in state C, which is outlined in § 2. In § 3 we report the main results obtained and compare them with observations. Finally, § 4 contains a short discussion of the implications of our results.

## 2. The Model

No definite model for the accretion flow in GRS 1915+105 has been presented as yet. Here we consider a highly idealized scenario in which a Shakura-Sunyaev disk (Shakura & Sunyaev shaksuny:1973 (1973)) coexists with a hotter component. As suggested by Belloni et al. (belloni:1997 (1997)), the main cause for the large-scale time variability of this source is the onset of a thermal-viscous instability which blows off the inner, radiation-pressure supported region of the disk. In particular, according to this scenario, the hard state corresponds to the replenishment by the disk of the previously evacuated region on a viscous timescale. Further support to this picture has been provided by numerical simulations (Szuszkiewicz & Miller ewa-john:1997 (1997); ewa-john:1998 (1998)).

Fitting the observed spectra of BHCs in the hard state requires the presence, in addition to a standard disk, of a hot, less dense phase, formed e.g. by the debris of the puffed-up inner region or by the evaporation of the disk itself. The upscattering of soft disk photons by a hot cloud at $`T10`$ keV can account also for the positive lags of the continuum. In fact, Compton cloud models at constant temperature have been already proposed in connection with other sources which show (positive) time delays (Miyamoto et al. miyamo:1991 (1991); Vaughan et al. vaug:1994 (1994); see also Poutanen put:2000 (2000) for a review). Compton models suffer from the problem that soft photons cool down the cloud to the temperature of the impinging radiation, unless an efficient heating mechanism is at work. However, if cooling of ions is inefficient, like in advection-dominated flows, protons decouple from electrons and stay hot at a fraction of their virial temperature, $`T_p(r)\eta GMm_p/r`$. The electron thermal balance is fixed by the competition between ion heating and Compton cooling. A simple estimate of the equilibrium electron temperature gives

$$T(r)10(c^2r_{in}/GM)^{2/5}(1r_{in}/r)^{2/5}\mathrm{keV}$$
(1)
for a scattering depth $`1`$, $`\eta =0.01`$, and $`L_{soft}(1/r_{in}1/r)`$. As expected (see Poutanen put:2000 (2000)), the cloud temperature diminishes as the disk rim moves in and $`L_{soft}`$ increases, but the dependence is weak. In an ion-heated corona the temperature is larger at smaller radii because $`L_{soft}`$ decreases with $`r`$ so that less photons are available to cool the inner part of the cloud.

A non-homogeneous Compton cloud model may indeed explain the positive/negative lags and the decrease in the hard color observed in the state C of GRS 1915+105, provided that the corona becomes denser as the disk extends inwards. The density may increase in time if a cloud of fixed mass contracts as the result of angular momentum losses caused by viscous and/or radiation stresses. In order to illustrate the basic idea of the model, we sketch the corona as formed by two regions (referred to as “warm” and “hot”), each at uniform density and temperature, $`T_W`$ and $`T_H`$, with $`T_H>T_W`$ (see Figure A Comptonization Model for Phase Lag Variability in GRS 1915+105). The scattering depths are $`\tau _W`$ and $`\tau _H`$, $`\tau _H>\tau _W`$. When the density of the corona is low, soft photons get preferentially scattered in the inner part of the corona which has $`\tau _H1`$ (while $`\tau _W1`$), gaining energy and producing a positive lag. As density increases, the Compton $`y`$-parameter of the corona exceeds unity below a certain radius $`r_{sc}`$. Photons emitted by the inner part of the disk undergo repeated scatterings with the hot electrons and start to fill the Wien peak at $`T_H`$. The spectrum emitted by the disk plus the hot, dense part of the corona is therefore the superposition of a multitemperature blackbody and a Bose-Einstein distribution at $`T_H`$; the latter has the same number of photons originally emitted by the disk in the region $`r_{in}rr_{sc}`$. These photons have still to traverse the warm part of the corona before escaping to infinity. In doing so, they are downscattered and this produces the negative lag and the observed spectrum. A somehow similar geometry (a standard disk with variable radius and a hotter, comptonizing cloud), although in a different context, has been successfully invoked to explain the spectral properties of Nova Muscae 1991 (e.g. Misra & Melia mimel:1997 (1997); Esin, McClintock & Narayan esmcnar:1997 (1997)).

In the next section we show that quantitative estimates, based on this simple picture, are indeed in agreement with the observed colors and lags for reasonable values of the parameters. In our model the inner radius changes only by a factor $`3`$, so, for the sake of simplicity, in the following we neglect the small ($`1.5`$, see eq. ) change of $`T`$ with $`r_{in}`$.

## 3. Numerical Estimates

The calculation of the emerging spectrum proceeds in two steps and follows the standard procedure discussed by Sunyaev & Titarchuk (suntit:1980 (1980)). First, we compute the time-dependent comptonized spectrum for a given seed photon distribution by solving the Kompaneets equation in a semi-infinite homogeneous medium and then use an escape probability algorithm to obtain the time evolution of the spectrum emerging from a medium of finite depth. We introduce a dimensionless time $`u=(n_e\sigma _Tc\theta )t`$ and energy $`x=h\nu /kT`$, where $`n_e`$ is the electron density, $`\sigma _T`$ the Thomson cross-section and $`\theta =kT/m_ec^2`$. In the following all lengths are expressed in units of $`GM/c^2`$. The photon occupation number $`n(u,x)`$ is the solution of the Kompaneets equation with an appropriate initial condition $`n(0,x)`$. For $`r_{in}>r_{sc}`$ the initial photon occupation number is that of a multitemperature blackbody (calculated including GR and hardening effects, e.g. Novikov & Thorne notho:1973 (1973); Shimura & Takahara giap:1995 (1995)) truncated at $`r_{in}`$ whereas, for $`r_{in}<r_{sc}`$, it corresponds to the superposition of a Bose-Einstein distribution at $`T_H`$ and a multitemperature blackbody truncated at $`r_{sc}`$.

Since in our model the long-term evolution corresponds to a variation of the inner edge of the disk, the Kompaneets equation has been solved numerically for a sequence of values of $`r_{in}`$. Once $`n(u,x)`$ has been evaluated, the time evolution of the escaping photons is given by

$$N(u,x)=x^2n(u,x)P(u,\tau )$$
(2)
where $`N`$ is proportional to the count rate and $`P(u,\tau )`$ is the escape probability distribution. The explicit form of $`P(u,\tau )`$ depends on the geometry (Sunyaev & Titarchuk suntit:1980 (1980)). In our calculations we have used a weighted sum of the expressions for a spherical cloud with either a central or a diffused photon source.

The total count is proportional to the convolution

$$F(x)=_0^{\mathrm{}}x^2n(u,x)P(u,\tau )𝑑u.$$
(3)
We note that eq. (3) gives exactly the spectrum produced by a stationary source of photons. Therefore it will be used to compute the hardness ratios on time intervals substantially larger than typical variability timescales in the continuum. On the other hand, if one is interested in following the time evolution of the spectrum over timescales $`1/\nu _{QPO}`$ the relevant quantity to consider is $`N(u,x)`$.

To reproduce the observed time evolution of the color in a given energy band we have convolved $`N(u,x)`$ with the interstellar absorption and with the detector response function

$$C(E_1E_2)_{E_1}^{E_2}N(u,x)\mathrm{exp}[N_H\sigma (x)]A(x)𝑑x$$
(4)
where $`N_H`$ is the column density. The lags in the continuum are computed comparing the time evolution of the the soft and hard colors, $`C(25)`$ and $`C(513)`$, as obtained from eq. (4). The hardness ratios HR1 and HR2, defined as $`C(513)/C(25)`$ and $`C(1360)/C(25)`$, are calculated again starting from eq. (4) but this time replacing $`N(u,x)`$ with $`F(x)`$.

The calculations presented here have been obtained for $`T_H=15`$ keV, $`T_W=1.5`$ keV and $`r_{sc}=20`$; the black hole mass and angular momentum are $`M=6M_{}`$ and $`a=0.1M`$. The model with $`r_{in}=16`$ corresponds to the usual picture in which a hot corona ($`T=T_H`$) of moderate depth ($`\tau =\tau _H3`$) upscatters the soft photons produced by the disk; under these conditions the energy exchange in the warm region can be reasonably assumed to be negligible. The opposite situation is represented by the model with $`r_{in}=6`$. Now the disk inner region lies entirely within the hot part of the corona which, because of the increase in the density, has $`\tau _H100`$. Escaping photons are then reprocessed in the outer, warm part of the corona ($`T=T_W`$, $`\tau =\tau _W10`$). Intermediate cases are more difficult to model in the present framework, because photons scatter in both regions and this would require a treatment of comptonization in a non-uniform medium. To avoid undue complications, at this stage we assume that the effects of a temperature profile can be roughly described in terms of an “effective” temperature $`T`$, $`T_W<T<T_H`$. For these models $`T`$ and $`\tau `$ were chosen in such a way to maintain $`y1`$.

Our results are illustrated in Figure A Comptonization Model for Phase Lag Variability in GRS 1915+105. As shown in Figure A Comptonization Model for Phase Lag Variability in GRS 1915+105a, the observed hardness ratios are quite well reproduced by the model for a value of the column density, $`N_H=2\times 10^{22}`$ cm<sup>-2</sup>, consistent with observations (Belloni et al. belloni:2000 (2000)). When the inner edge of the disk is far out, the spectrum is a power-law produced by thermal comptonization of soft disk photons. As $`r_{in}`$ decreases, downscattering becomes progressively more important with respect to upscattering in the hot corona, the spectrum softens and, at the same time, the phase lag of the continuum decreases. Eventually, systematic downscattering dominates and soft (2–5 keV) photons, which underwent more scatterings, escape later than hard (5–13 keV) photons, causing an inversion in the sign of the lag. Figures A Comptonization Model for Phase Lag Variability in GRS 1915+105b and A Comptonization Model for Phase Lag Variability in GRS 1915+105c show the lag of the continuum as a function of HR1 and QPO frequency (the scattering of the computed points is partially due to the poor sampling of the parameter space). Time lags have been computed assuming that the typical size of the comptonizing corona varies with time in accordance with $`r_{in}`$ in the range $``$20–400. The relation between $`\nu _{QPO}`$ and $`r_{in}`$ was derived associating the model with $`r_{in}=6`$ (16) to the observation with the largest (smallest) QPO frequency, 6 (0.6) Hz and yields

$$r_{in}13\nu _{QPO}^{0.43}.$$
(5)
## 4. Discussion

In this Letter we have shown that a simple Compton model can account for the main observational properties of GRS 1915+105 in state C. In particular, the fundamental correlations (color-color, lag-color and lag-QPO frequency) are successfully reproduced for a reasonable choice of the model parameters. Although we have adopted here a quite definite picture for the accretion flow, it is important to stress that the basic results are largely independent of the specific choice for the geometry. The decrease of the lag in the continuum and the softening of the spectrum as the disk fills in are simply a consequence of the increase of the coronal depth and of the decrease of the temperature where radiation scatters before escaping. Thus, only three main ingredients are required to account for the observed correlations within the framework of a simple thermal comptonization model: a non-isothermal corona of varying optical depth, a truncated standard disk and a source of more energetic ($`10`$ keV) primary photons at small radii. Further assumptions (like the variable size of the corona and the existence of an “effective” comptonization temperature) are introduced only to match the observed values. It should be noted, however, that other geometrical settings for the time-varying accretion flow may work equally well.

On the other hand, albeit not directly supported by hydrodynamics computations, the flow structure adopted here exhibits a certain degree of consistency in the way in which the QPO frequency and the coronal optical depth vary with inner disk radius. In fact, without entering into the mechanism that originates the observed QPO, we have assumed that it is somehow produced at the inner edge of the disk and that its frequency varies with $`r_{in}`$. It is interesting to note that the relation between $`\nu _{QPO}`$ and $`r_{in}`$ derived in § 3 (eq. ) is similar to that found by Di Matteo & Psaltis (dima:1999 (1999)) combining the empirical relation of Psaltis, Belloni & van der Klis (psaltis:1999 (1999)) with the upper bound provided by the keplerian orbital frequency. We stress, however, that our result has been obtained in an independent way and is a direct consequence of the values of the inner disk radius needed to match the observed colors and time lags. The contraction of the corona (a factor $`20`$) required to reproduce the largest/smallest lags implies an increase in the density ($`\rho R^3`$ at constant coronal mass) consistent with the assumed ratio of the scattering depths, $`\tau _H(r_{in}=6)/\tau _H(r_{in}=16)100`$. The smallest value of the inner disk radius required by the model, $`r_{in}6M`$, is consistent with the assumed value $`a/M=0.1`$. If the black hole is rapidly rotating, as claimed by Zhang, Cui & Chen (zhang:1997 (1997)) and Cui, Zhang & Chen (cui:1998 (1998)), then in low luminosity state C the disk does not reach the innermost stable circular orbit.

A complete discussion of the behaviour of GRS 1915+105 in states A and B is outside the scope of the present investigation. Here, we simply mention that models with high coronal optical depth but hard-photon starved (because, e.g., of the final collapse of the innermost, denser part of the corona) have hardness ratios in agreement with those observed in states A and B. Although no QPO has been observed in these states as yet, the positive detection of time lags correlated with colors in states A and B may lend further support to the present Compton model.

Very recently, complex phase lag behavior in XTE J1550-564 has been reported by Cui, Zhang & Chen (czc:2000 (2000)). The power spectra are very similar to those in Reig et al. (reig:2000 (2000)). However, the QPO phase-lag dependence is different. The fundamental/first harmonic shows a negative/zero lag, increasing in magnitude as the spectrum softens. In the light of the higher complexity of the timing properties of XTE J1550-564, any attempt to extend the present model to this source is premature.