text:
# Thermal measurements of stationary nonequilibrium systems: A test for generalized thermostatistics

## Abstract

We show that a gas thermometer in contact with a stationary classical system out of thermal (Boltzmann) equilibrium evolves, under very general conditions, towards a state characterized by a Lévy velocity distribution. Our approach is based on a kinetic-like equation that applies to a wide class of models for the system-thermometer interaction. The results clarify the role of non-exponential energy distributions as possible generalizations of the Boltzmann distribution for systems where the usual formulation of thermostatistics may not apply. In particular, they show that the power-law distributions derived from Tsallis’s nonextensive formalism are irrelevant to the stationary state of the thermometer, thus failing to give a consistent description of the system-thermometer equilibrium. We point out the need of a generalized thermostatistical formulation able to give a unified frame to Lévy and Maxwell distributions.

Lévy distributions  have been extensively applied to the description and modeling of a wide class of physical processes, ranging from anomalous transport in disordered media  and turbulent flows , to phase-space diffusion in dynamical systems  and polymer dynamics . They have also found application in other branches of science, such as in biology . The ubiquity of Lévy distributions in many natural phenomena is a straightforward consequence of their stable character under summation of random variables whose distributions have diverging moments, as shown by Paul Lévy in his generalization of the central limit theorem . This essential result gives Lévy distributions the same status as the Gaussian distribution in the statistical description of stochastic processes.

While most applications of Lévy distributions deal with dynamical and transport processes, we show in this Letter that they replace the equilibrium Maxwell distribution—a Gaussian in the velocity variable—in an extended version of the scenario of a thermal measurement. We analyze the asymptotic energy distribution of a thermometer in thermal contact with a classical system in a nonequilibrium stationary state. While the nature of this nonequilibrium state is not explicitly specified, we typically refer to a situation where the system is not isolated and suffers the effect of external stationary forces that maintain its energy distribution apart from the prediction of Boltzmann-Gibbs statistics. It is an essential fact of thermostatistics  that, if the system is in thermodynamical equilibrium, its Boltzmann energy distribution is “copied” by the thermometer in such a way that its average energy per degree of freedom, i.e. its temperature, becomes equal to that of the system. Now, what properties of the state of the system are detected by the thermometer in the case that the system is maintained out of equilibrium?

The motivation of this question is two-fold. First, since all physical systems are to some extent thermodynamically open and subject to external influence, the possible effect of these factors in a thermal measurement constitutes a problem of empirical relevance. In general, this problem transcends the limits of equilibrium thermostatistics and calls for a nonequilibrium, dynamical formulation. We choose a kinetic-like approach which describes a wide class of interaction models for the thermal contact between the system and a gas thermometer. Our results show that the thermometer attains a velocity distribution which is fully determined by the high-energy distribution of the system, and results to be insensitive to other details in the system state. In the case where the energy distribution of the system decays as a power law, as advanced above, the thermometer velocities approach a Lévy distribution. This distribution “copies” the singularities associated with divergent energy moments in the system, and reduces to a Maxwellian when the average energy is finite.

The second motivation for our question has to do with the possibility of generalizing thermostatistics to the description of physical systems where the usual Boltzmann-Gibbs formulation may not apply . Tsallis’s nonextensive thermostatistics, for instance, is claimed to replace the usual formulation for systems where long-range interactions lead extensivity assumptions to break down, while containing Boltzmann-Gibbs statistics as a limiting case . In Tsallis’s formalism, the Boltzmann exponential distribution is replaced by a class of power-law functions derived from a variational principle for a generalized entropy. Whether this formalism is compatible with the basic facts of equilibrium between systems in thermal contact can be decided, precisely, by studying the interaction of a thermometer and a system with such power-law energy distribution. Our results show that, since Lévy distributions are not regarded in this formalism as possible equilibrium states, Tsallis’s thermostatistics fails to give a consistent picture of thermodynamical equilibrium.

Consider a system in a stationary state which does not necessarily coincide with Boltzmann equilibrium. We assume that the state of the system is macroscopically characterized by an energy distribution $`F_0(ϵ)`$. The system interacts with a thermometer consisting of an ensemble of independent particles of mass $`m`$ moving in a $`d`$-dimensional domain. Within this domain, the distribution of particles is supposed to be spatially homogeneous and isotropic in velocities, such that the associated distribution function $`F(ϵ,t)`$ depends on the energy $`ϵ`$ and the time $`t`$ only. The interaction between system and thermometer is assumed to fulfill the following conditions. (i) Interaction events are time localized and infrequent, such that the typical time between interactions is large as compared with the relaxation time of the energy distribution of the system toward $`F_0(ϵ)`$. This insures that (a) at each event, the system is found at its stationary state and (b) any correlation between the states of system and thermometer, created by the interaction, dies out before the next event takes place. (ii) The amount of energy interchanged at each event is small as compared with the total energy of both the system and the thermometer.

Within these conditions, it is possible to describe the evolution of $`F(ϵ,t)`$ by means of a kinetic-like equation. On the other hand, as a consequence of condition (i), no evolution equation is required for the system distribution. To represent the interaction between system and thermometer we choose a very generic class of models—closely related to stochastic Maxwell models —which considerably simplifies the mathematical treatment of the kinetic problem and, at the same time, does not imply severe limitations to the physics of the interaction. In fact, these models include interactions where momentum conservation can either hold or be violated, and admit arbitrary angular dependence in the associated cross sections. They assume that the interaction rate is independent of the energy and that energy itself is conserved. The corresponding evolution equation for $`F(ϵ,t)`$ stands for the balance between events where, due to interaction with the system, the thermometer reaches or abandons its states of energy $`ϵ`$. Such events yield, respectively, positive and negative contributions to $`_tF`$. The equation is more conveniently written in the Fourier representation, for the Fourier-transformed distribution

$$\mathrm{\Phi }(u,t)\varphi (k,t)=\mathrm{exp}(i𝐤𝐯)f(v,t)𝑑v^d,$$
(1)
with $`u=k^2/2m`$. Here,

$$f(v,t)=\frac{\mathrm{\Gamma }(d/2)}{2\pi ^{d/2}}mv^{2d}F(mv^2/2,t)$$
(2)
is the (isotropic) velocity distribution associated with $`F(ϵ,t)`$ . Due to the normalization of $`F(ϵ)`$, we have $`\mathrm{\Phi }(0,t)=1`$ for all $`t`$. With an analogous definition for $`\mathrm{\Phi }_0(u)`$ in terms of $`F_0(ϵ)`$, the Fourier-transformed evolution equation reads :

$$_t\mathrm{\Phi }(u,t)=_0^1w(s)\mathrm{\Phi }_0(su)\mathrm{\Phi }[(1s)u,t]𝑑s\mathrm{\Phi }(u,t).$$
(3)
The kernel $`w(s)`$ characterizes the interaction, and satisfies the normalization $`_0^1w(s)𝑑s=1`$. When necessary, we assume that the moments

$$w_\mu =_0^1s^\mu w(s)𝑑s$$
(4)
are well defined.

The above condition (ii) imposes that $`w(s)`$ is appreciably different from zero only for $`s0`$. Taking this condition into account, we write $`\mathrm{\Phi }[(1s)u,t]\mathrm{\Phi }(u,t)su_u\mathrm{\Phi }(u,t)`$, to obtain an approximate form of Eq. (3). Its solution can be shown to asymptotically approach the stationary distribution

$$\mathrm{\Phi }(u)=\mathrm{exp}\left[\frac{1}{w_1}_0^u\frac{du^{}}{u^{}}_0^1w(s)\mathrm{ln}\mathrm{\Phi }_0(su^{})𝑑s\right].$$
(5)
Since $`w(s)0`$ except for $`s0`$ the asymptotic solution $`\mathrm{\Phi }(u)`$ is fully determined by the behavior of $`\mathrm{\Phi }_0(u)`$ close to $`u=0`$. We pay particular attention to the special case where, for $`u0`$,

$$\mathrm{\Phi }_0(u)1\alpha _0u^\nu ,$$
(6)
with $`\alpha _0>0`$ and $`0<\nu 1`$. For $`\nu 1`$, this form corresponds to an energy distribution with a power-law tail for high energies, $`F_0(ϵ)ϵ^{\nu 1}`$. While $`F_0`$ is still normalized, its first moment—namely, the average energy—diverges. In the limit $`\nu =1`$ a finite average energy is recovered and the long high-energy tail is lost. For $`\mathrm{\Phi }_0(u)`$ as in Eq. (6), the stationary Fourier distribution for the thermometer is

$$\mathrm{\Phi }(u)=\mathrm{exp}(\alpha u^\nu ),$$
(7)
with $`\alpha =\alpha _0w_\nu /\nu w_1`$. Note that, for $`0<\nu <1`$, $`\mathrm{\Phi }(u)`$ reproduces the singularity of $`\mathrm{\Phi }_0(u)`$ at $`u=0`$. The corresponding Fourier-transformed velocity distribution, $`\varphi (k)=\mathrm{\Phi }(k^2/2m)`$ \[see Eq. (1)\], reads

$$\varphi (k)=\mathrm{exp}(bk^\gamma ),$$
(8)
with $`b=\alpha /(2m)^\nu `$ and $`\gamma =2\nu `$ ($`0<\gamma 2`$). Note that, in the velocity representation, this corresponds to a Lévy distribution for each velocity component. In fact, the Fourier transform of, say, $`f(v_x)=𝑑v_y𝑑v_z\mathrm{}f(v)`$ is $`\varphi (k_x,0,0,\mathrm{})=\mathrm{exp}(b|k_x|^\gamma )`$. Though, in general, it is not possible to write explicitly the corresponding form of $`f(v_x)`$, it is known that for large $`|v_x|`$ and $`\gamma 2`$ it behaves as $`f(v_x)|v_x|^{\gamma 1}`$ . For $`\gamma =2`$, $`f(v_x)`$ is a Gaussian distribution.

We stress that the result (8) for the Fourier-transformed velocity distribution of the thermometer is independent of the detailed form of the distribution $`F_0(ϵ)`$ of the system at low energies. As long as condition (ii) is satisfied, $`\varphi (k)`$ is completely determined by the power-law high-energy tail of $`F_0(ϵ)`$.

Let us analyze these results in the few cases where the velocity distribution can be explicitly obtained. We first consider the situation where the average energy of the system under study is finite, $`E=ϵF_0(ϵ)𝑑ϵ<\mathrm{}`$. Under this condition we have, in Eq. (6), $`\alpha _0=2E/d`$ and $`\nu =1`$. This gives, for the Fourier-transformed velocity distribution of the thermometer, the Gaussian $`\varphi (k)=\mathrm{exp}(Ek^2/md)`$. In the velocity representation we get the Maxwellian

$$f(v)=\left(\frac{md}{4\pi E}\right)^{d/2}\mathrm{exp}\left(\frac{md}{4E}v^2\right).$$
(9)
The average value of the kinetic energy over this distribution is, as may be expected, $`mv^2/2=E`$. In other words, the stationary velocity distribution of the thermometer is a Maxwellian whose average energy per particle coincides with that of the system. This result includes of course the case where the system is itself in thermodynamical (Boltzmann) equilibrium at temperature $`T`$, for which $`ET`$.

For $`\nu =1/2`$, the energy distribution of the system decays as $`F_0(ϵ)ϵ^{3/2}`$ and, consequently, the average energy diverges. The Fourier-transformed velocity distribution of the thermometer results to be $`\varphi (k)=\mathrm{exp}(bk)`$, with $`b=\alpha _0w_{1/2}\sqrt{2/mw_1^2}`$. Taking into account its units, we write $`b=p_0/m`$, where $`p_0`$ is a characteristic linear momentum given by the energy distribution $`F_0(ϵ)`$ and the specific form of the interaction kernel $`w(s)`$. In velocity space, each component $`v_i`$ is distributed according to a Cauchy distribution,

$$f(v_i)=\frac{mp_0}{\pi }\frac{1}{p_0^2+m^2v_i^2}.$$
(10)
Note that, here, the relevant dynamical variables in the thermometer distribution are the momentum components $`mv_i`$ and, accordingly, the relevant parameter is the momentum $`p_0`$. Compare with the result for $`\nu =1`$, Eq. (9), where the distribution is a function of the energy and the relevant parameter is $`E`$.

The present results shed light on the role of non-Boltzmannian energy distributions in defining the stationary states of out-of-equilibrium interacting systems. In particular, they show that, under the above conditions (i) and (ii), the Gaussian stationary velocity distribution of an ensemble of independent particles interacting with a non-equilibrium stationary system results to be a special instance of the more general situation where the velocity components have Lévy distributions. We stress that, despite the fact this conclusion is limited by the validity of the assumptions (i) and (ii), these conditions are fully compatible with the paradigm of thermal measurement  and therefore refer to a realistic, experimentally accessible situation.

As advanced above, our conclusions are particularly significant in the evaluation of generalized formulations of thermostatistics as plausible descriptions of physical systems to which the usual Boltzmann-Gibbs theory may not apply. We focus our attention on Tsallis’s formulation, which has motivated large amounts of work over the last two decades . Within Tsallis’s statistics, canonical probability distributions are obtained from maximization of a generalized entropy

$$S_q=\frac{1_ip_i^q}{1q},$$
(11)
where $`p_i`$ is the probability of state $`i`$, and the sum runs over all the states of the system. The parameter $`q`$ quantifies the nonadditivity of $`S_q`$, and is therefore a measure of nonextensivity in the system under consideration. The maximization process is subject to the generalized constraint

$$E_q=\frac{_iϵ_ip_i^q}{_ip_i^q},$$
(12)
with $`ϵ_i`$ the energy of state $`i`$, and leads to the probability distribution

$$p_ip(ϵ_i)=Z_q^1[1(1q)\beta _q(ϵ_iE_q)]^{1/(1q)},$$
(13)
where $`\beta _q`$ is a (“renormalized” ) Lagrange multiplier, analogous to the inverse temperature, and $`Z_q`$ is a normalization constant, analogous to the partition function. Boltzmann-Gibbs formulation is recovered for $`q1`$. Note that, for systems whose density of states behave as $`\rho (ϵ)ϵ^\sigma `$—a class which includes most classical systems—the probability $`p(ϵ)`$ corresponds to a Fourier-transformed energy distribution of the form of Eq. (6), with $`\nu =\sigma +(q2)/(1q)`$.

An attractive feature of Tsallis’s generalized thermostatistics is that it preserves most of the mathematical structure of the usual theory, including the definition of thermodynamic functions, Legendre transformations, and even linear nonequilibrium properties . In fact, a large part of the literature on this topic is devoted to the formal extension of thermodynamical relations to the generalized formulation. Thus, putting aside its phenomenological applications to nonequilibrium processes such as anomalous diffusion  and turbulence , the real tour de force of Tsallis’s formalism lies in the description of thermal equilibrium for systems where Boltzmann-Gibbs theory is supposed to fail. This is in fact the original and most frequently invoked motivation of the formulation and, consequently, constitutes its genuine source of validation. Our results clarify whether Tsallis’s statistics is relevant to an aspect that—in spite of its essential role in thermodynamics—has been scarcely treated in the profuse literature on Tsallis’s theory , namely, thermal equilibrium between interacting systems.

Assume to have a system which, due to its nonextensive nature, exhibits an energy distribution of the form of Eq. (13), as predicted by Tsallis’s thermostatistics. Allow furthermore the system to interact with a thermometer, as specified above. Consider first that the nonextensivity index $`q`$ and the density of states of the system are such that $`E=_iϵ_ip_i<\mathrm{}`$. According to our results, the stationary velocity distribution of the thermometer is a Maxwellian, corresponding to a Boltzmann exponential probability for the energy. As expected for an ensemble of independent particles, the thermometer behaves as an extensive system ($`q=1`$). The value of the parameters in the thermometer distribution predicted by Tsallis’s formalism, however, results to be wrong. Abe and Rajagopal  have shown that a formal extension of equilibrium conditions for two systems in thermal contact yields, in Tsallis’s theory, $`\beta _q=\beta _q^{}`$, where $`q`$ and $`q^{}`$ are the nonextensivity indices corresponding to the two systems . Our results show that this relation is generally not satisfied. In fact, taking $`q^{}=1`$ as the nonextensivity index of the thermometer, Eq. (9) implies $`\beta _q^{}E^1`$. On the other hand, $`\beta _q`$ shows in general a more complicated dependence on $`E`$. For instance, it can be readily shown that, for a system with a density of states $`\rho (ϵ)ϵ^\sigma `$, one has $`\beta _qE^{(q1)(\sigma +1)1}`$. Let us also point out that the generalized average $`E_q`$ of Eq. (12) plays no role in the parameter that defines the equilibrium distribution of the thermometer—namely, the inverse temperature $`\beta _q^{}`$—in spite of the fact that this generalized average replaces $`E`$ in the formal extension of thermodynamics to Tsallis’s formalism.

In the case that the average $`_iϵ_ip_i`$ diverges, our results imply that the thermometer approaches an energy distribution not belonging to the class of the system distribution, Eq. (13). Thus, no well-defined index $`q`$ can be assigned to the thermometer. In this situation, the assumption that Tsallis’s thermostatistics would describe the system-thermometer equilibrium fails drastically. The only trace of the energy distribution of the system in the thermometer distribution arises from its power-law tail, which defines the Lévy exponent $`\gamma =2\nu `$ in Eq. (8). Note that, again, the generalized average $`E_q`$ is irrelevant to the determination of the equilibrium distribution of the thermometer.

Summing up, our kinetic description of a thermometer in contact with a system in a stationary state different from Boltzmann equilibrium, suggests that a suitable extension of thermostatistics should yield Lévy distributions as the generalization of the Maxwellian velocity distribution, instead of the power-like functions of Eq. (13). We conclude that the energy probabilities derived from Tsallis’s thermostatistics do not play the role of equilibrium distributions for systems in thermal contact. Though this conclusion is not completely general, it applies to a wide class of interaction models in the realistic situation where a system is put in contact with a thermometer. The question remains open as to which form of the entropy and which constraints should be used to derive Lévy distributions from a variational principle.

We thank G. Abramson and I. Samengo for their critical reading of the manuscript.


format: text

text:
# Possible improvements on the mass of 𝜈_𝜏 using leptonic 𝐷^±_𝑠 decays

## 1 Introduction

The last few years have seen growing evidence for non-vanishing neutrino rest masses in the results from neutrino oscillation experiments . However, the direct bounds on neutrino masses remain rather weak. While the electron neutrino mass is known to be smaller than about 2.2 eV , the muon neutrino mass has to be smaller than 170 keV , a bound about 5 orders of magnitude worse and still as high as 30 % of the electron mass. The situation is even worse by another two orders of magnitude for $`\nu _\tau `$ which is known from ALEPH to be smaller than 18.2 MeV . In this paper we explore leptonic $`D_S^\pm `$ decays for improving on the bound of $`\nu _\tau `$.

## 2 Leptonic $`D_S^\pm `$ decays

In the standard model the leptonic branching ratios of $`D_S^\pm `$ are given as

$$BR(D_S^\pm l\overline{\nu }_l)=\frac{G_F^2}{8\pi }V_{cs}^2f_{D_S}^2\tau _{D_S}m_{D_S}m_l^2F(1\frac{m_l^2}{m_{D_S}^2})^2$$
(1)
with $`V_{cs}`$ is the corresponding CKM matrix element, $`\tau _{D_S}`$ the $`D_S^\pm `$ life-time, $`m_{D_S}`$ the mass of the $`D_S^\pm `$, $`f_{D_S}`$ the decay constant and $`m_l`$ the lepton mass, and F is a phase space factor which depends on the neutrino mass (F = 1 when neutrino mass is zero). Several quantities cancel when the ratio of two leptonic branching ratios is taken and furthermore this ratio is quite sensitive to the $`\nu _\tau `$ mass . Neglecting the muon neutrino mass to first approximation, the ratio between muonic and tauonic decays can be parametrised as a function of $`\nu _\tau `$ mass to first order in $`m_{\nu _\tau }`$

$$R=\frac{(D_S^\pm \tau \nu _\tau )}{(D_S^\pm \mu \nu _\mu )}=R_0(m_\nu =0)\times \left(1C(\frac{m_{\nu _\tau }}{m_\tau })^2\right)$$
(2)
$`R_0`$ is obtained as

$$R_0(m_\nu =0)=\frac{m_\tau ^2}{m_\mu ^2}\frac{\left(1\frac{m_\tau ^2}{m_{D_s^\pm }^2}\right)^2}{\left(1\frac{m_\mu ^2}{m_{D_s^\pm }^2}\right)^2}=9.79$$
(3)
$`C`$ is given by

$$C=\frac{3(m_\tau /m_{D_S^\pm })^41}{(1(m_\tau /m_{D_S^\pm })^2)^2}=28.77$$
(4)
The values were obtained using $`m_\tau `$=1777 MeV and $`m_{D_S^\pm }`$=1969 MeV . The ratio R as a function of $`m_{\nu _\tau }`$ is plotted in Fig.1.

As expected, the figure shows that the ratio decreases as the mass of $`\nu _\tau `$ is increased. Furthermore, it follows from Eq. 4 that to improve the bound on $`m_{\nu _\tau }`$ down to 10 MeV would correspond to a 0.1 % effect on $`R`$ only. To achieve the latter bound would be especially interesting, because it would close the window for observation of a MeV Majorana neutrino in double beta decay via atomic mass dependent effects.

We have not yet discussed the role of radiative corrections to the ratio $`R_0`$. For the case of $`\pi `$ decays these have been carefully calculated and are well-known. A complete calculation of analogous radiative corrections for $`D_S^\pm `$ decays is yet to be done. It is expected that these corrections would be in the range of 2 to 4 %. This will modify the bound on $`\nu _\tau `$ mass. To extract a meaningful limit on the mass of $`\nu _\tau `$ it will be necessary to have the radiative corrections to the value of R available.

One can also compare inclusive decay rates, rather than the exclusive modes above. The ratio

$$R_{\mu ,\tau }=\frac{\mathrm{\Gamma }(D_S^\pm \tau \nu _\tau +D_S^\pm \tau \nu _\tau \gamma )}{\mathrm{\Gamma }(D_S^\pm \mu \nu _\mu +D_S^\pm \mu \nu _\mu \gamma )}$$
(5)
is free from any dependence on energy resolution and for a point-like $`D_S^\pm `$ the expression is as given in . This gives a correction of 2.4 % in the direction of increasing the tau branching ratio. In both cases, it would be desirable to have an estimate of the remaining structure dependent corrections, although we believe that they would not be larger than the effects already included here.

Similar considerations can be applied to the decays of $`D^+`$ as well. However, the rates in that case are suppressed by the CKM suppression as well as phase space, and the branching ratios are smaller than for $`D_s`$ decays by an order of magnitude or more. This make $`D^+`$ unsuited for extracting mass limits on the tau neutrino.

## 3 Experimental status

The current status of leptonic branching ratios of interest are compiled in Tab. 1. As can be seen, all measurements still have large errors. Taking the PDG values  implies R = 12.5 $`\pm `$ 5.5, clearly not allowing any conclusions on $`m_{\nu _\tau }`$ . However recently there have been improvements in investigations of $`D_S^\pm `$ decays at LEP. Using the branching ratio values of ALEPH  a ratio of R = 8.5 $`\pm `$ 5.2 can be obtained, unfortunately still having a much too large an error for the purpose at hand. The situation might be improved by producing a clean and statistically large sample by looking at diffractive $`D_S^\pm `$ production in antineutrino-nucleon scattering at a neutrino factory  or accurate measurements at the planned CLEO-c charm factory . For the latter an accuracy on both branching ratios of 4 % is predicted. If we assume that the central values remain the same as obtained by ALEPH this would imply a value of R = 8.5 $`\pm `$ 0.5, which is significantly away from $`R_0`$. The accuracy in R, of about 6 %, would be a great improvement over the current one, but still not quite at the level needed to get improved bounds on the $`\nu _\tau `$ mass, which calls for a level of less than 1 %.

## 4 Summary

We discuss the possibility of gaining information on the mass of the tau neutrino by investigating leptonic $`D_S^\pm `$ decays. An improvement on the $`m_{\nu _\tau }`$ down to 10 MeV is possible if the ratio of muonic and tauonic $`D_S^\pm `$ decays can be measured with an accuracy of 0.1 %. While current data do not allow to draw any conclusions , this might change in future experiments especially by using a charm factory. We emphasize the importance of having a calculation of the radiative corrections available in anticipation of a future improvement in the measurement of these branching ratios.

## Acknowledgements

We are thankful to J. Link, W. Marciano and the referee for useful discussions and critical comments. K. Zuber is supported by a Heisenberg Fellowship of the Deutsche Forschungsgemeinschaft.


text:
# The evolutionary advantage of diploid sex

## I Introduction

The reasons why sexual reproduction appeared, hundreds of million years ago, are still under investigation today 1 . One of the reasons commonly used to justify sex is the higher genetic diversity it produces if compared to asexual reproduction. On the other hand, sex is closely related to aging: all sexual species present aging, which starts as soon as reproduction starts 2 . Because of its ability to address both questions in a single context, the Penna model for biological aging 3  has been extensively used in recent years to study why sex evolved (for a review see 4 ).

Sá Martins and Moss de Oliveira, for instance, have shown that the higher diversity of sexual populations can prevent extinction if a catastrophe damaging the genome of the individuals (such as exposure to radiation) occurs 5 . Diversity was also the key issue discussed by Sá Martins 6  in his simulation of the competition of sexual and asexual varieties of a single species, coevolving under the assault of parasitic infestation. The density of parasites was shown to drive a transition between asexual and sexual orders, in agreement with the so-called Red Queen Hypothesis for the evolution of sex 7 ). The above mentioned results show that the model provides sexual varieties with survival advantages, due to their larger genetic diversity, over the asexual ones, whenever they compete in a changing environment. Nevertheless, in the absence of competition and a mutant environment, the asexual populations resulting from simulations of the standard Penna model are always much larger, since they generate twice as many offspring for the same population size. This effect is not overcome by the need, in diploid individuals, for homozygose if the harmful allele is to become active in loci for which it is not dominant. Without catastrophes or parasites, if a sexual and an asexual population are allowed to evolve under the same immutable environment, competing for the same food resources, the asexual one will win and the sexual one disappears (due to the Verhulst logistic factor explained in section two). Only very recently has a new ingredient been introduced 8  into the Penna model allowing larger sexual populations to be obtained, with sizes comparable to the asexual ones.

We now use this modified Penna model to address a tantalizing issue: since a diploid sexual population has the upper hand when competing against an asexual one due do the diversity generated by the use of genetic material coming from two different parents, why does not nature enhance this effect by allowing the genome of the offspring to benefit from three different templates? Is the need for homozygose mentioned above enough to overcome the burden of using three individuals to generate one offspring? To provide answers to these questions, we compare two different kinds of sexual populations: one involving the mating of two diploid individuals as parents (normal sexual reproduction) and the other involving three triploid individuals.

The paper has the following structure: In section 2 we describe the model, first for diploids and then for triploids, in section 3 we present our results and in section 4, our conclusions.

## II The model

### II.1 Diploid population

Each individual of the population is represented by a “chronological genome”, which consists of two bit-strings of 32 bits (32 loci or positions) each, that are read in parallel. One string contains the genetic information inherited from the mother and the other, from the father. Each position of the bit-strings is associated to a period of the individual’s life, which means that each individual can live at most for 32 periods (“years”). Each step of the simulation corresponds to reading one new position of all individuals’ genomes. Diseases are represented by bits 1. If an individual has two bits 1 in the $`i`$-th position of both bit-strings (homozygote), it will start to suffer the effects of a disease at his $`i`$-th year of life. If the individual is homozygous with two bits zero, no disease appears in that age. If the individual is heterozygous in that position, he will become sick only if that locus is one for which the harmful allele is dominant. The dominant loci are randomly chosen at the beginning of the simulation and remain fixed. If the actual number of accumulated diseases reaches a threshold $`T`$, the individual dies.

In order to avoid an exponential increase of the population and to introduce a dispute for food and space, a logistic Verhulst factor is used. Every time step and for each individual, a random number between zero and one is generated. This number is compared with $`V=N(t)/N_{max}`$, where $`N(t)`$ is the actual size of the population and $`N_{max}`$ is the carrying capacity. If the random number is smaller than $`V`$, the individual dies, independent of its age or genome. For a discussion of other alternatives of implementation of this mean-field-like interaction between individuals of the population we direct the reader to Ref. cebrat .

If a female succeeds in surviving until the minimum reproduction age $`R`$, it generates $`b`$ offspring every period until death. The female chooses randomly a male to mate, with age also greater than or equal to $`R`$. The offspring genome is constructed from the parents’ ones; first, the strings of the mother are cut in two at a randomly chosen position (“crossing”) and two complementary pieces, one from each string, are recombined to generate the female gamete (one string of 32 bits). $`M`$ deleterious mutations are then randomly introduced. The same process occurs with the father’s genome and the union of the two resulting gametes forms the new genome. The sex of the baby is randomly chosen, with equal probability. Deleterious mutation means that if a bit 0 is randomly chosen in the parent’s genome, it is set to 1 in the offspring genome. However, if a bit already set to 1 is randomly chosen, it remains 1 in the offspring genome (no back mutations).

The description given above corresponds to the original sexual version of the Penna model 4 ; 5 . The new ingredient, mentioned in section 1 and introduced by Sá Martins and Stauffer 8 , consists in assuming that harmful mutation reduces the survival probability. Thus, at each iteration, or “year,” each individual survives with probability $`\mathrm{exp}(mϵ)`$ if it has a total of $`m`$ harmful mutations (taking into account dominant positions) in it’s full genome (it is killed if a random number is tossed that is smaller than the survival probability). $`ϵ`$ is a parameter of the simulation, fixed from the start. To summarize, an individual may now die for any one of three reasons: i) randomly, due to the Verhulst logistic factor; ii) if its actual number of accumulated diseases reaches the limit $`T`$; iii) due to its survival probability being too small.

### II.2 Triploid Population

In this case, we assume that mating involves three triploid individuals (two males and one female or vice-versa). The chronological genomes consist of three bit-strings that are read in parallel. Homozygous positions are those with three equal bits at homologous loci. Harmful mutations are active only if there are three bits 1 at that same position, or at a heterozygous locus at which harmful mutations are dominant. Only females generate offspring. There are random crossing and recombination to produce the offspring genome (see fig. 1), and deleterious mutations are randomly introduced in each of the three gametes. The baby is a male or a female, with equal probability.

The first question we want to investigate relates to the competition between the benefit provided by a triploid genome and the effort involved in a mating that needs three individuals, instead of two, to generate the offspring. The benefit is the fact that, for triploids, mutations that happen in loci where the harmful allele is not dominant need to appear in all three bit-strings to become active.

## III Results

The curves presented below correspond to the average of the results obtained for 20 different populations (20 different initial seeds for the random number generator), using the following parameters:

Initial population = 10,000;

Maximum population size $`N_{max}=100,000`$ individuals;

Maximum number of genetic diseases $`T=3`$;

Minimum reproduction age $`R=8`$;

Birth rate $`b=2`$;

Mutation rate $`M=1`$ per bit-string (or gamete);

Number of dominant positions $`d=6`$;

Decrease in survival probability $`ϵ=0.015`$;

Total number of Monte Carlo steps = 800,000;

In figure 2 we present the time evolution of a diploid population (upper curve) and of two triploid populations, one for which reproduction involves the mating of one male and two females (central curve) and the other corresponding to the mating of two males and one female (lower curve). In figure 3 we present the corresponding dimensionless survival rates $`S(a)`$ as a function of age $`a`$, $`S(a)=N(a+1)/N(a)`$, where $`N(a+1)`$ is the number of individuals with age $`a+1`$ and $`N(a)`$ is the number of individuals with age $`a`$.

From these figures we can see that the diploid population is not only larger than the other two, but also presents a higher survival probability. However, it may be not enough to guarantee that the diploid population is better than the triploid ones. A second important measure, as mentioned in the introduction, is the genetic diversity of the populations. It has already been shown that the survival probability of a sexual population is the same as that of a diploid asexual population that reproduces by meiotic parthenogenesis. However, the larger genetic diversity of the sexual population may prevent it from extinction if, for example, exposure to radiation occurs . The genetic diversity is calculated by measuring the Hamming distance, in this case defined by the number of different loci (bits) between the genomes, for all pairs of individuals. The probability distribution of these distances is obtained by making a histogram of the fraction of pairs, out of all possible pairs in the population, that present a given Hamming distance, normalized by its maximum possible value (64 for diploids and 96 for triploids). Figure 4 shows the resulting distributions for the diploid and triploid populations. It is clear that the diploid population presents both a larger mean distance between pairs, indicated roughly by the position of the peak of the distribution, and a larger variance, measured by the width at half the maximum height of the curves. This assures a larger diversity within its genome space for the diploid populations, and will give them the upper hand when competing against the triploids under the pressure of a rapidly mutating environment. It is worth while mentioning that the results are essentially the same if we allow a double crossing of the triploid genome during reproduction, and there is no benefit for the triploids to ensure that the offspring have their genetic material gathered from all three parents.

## IV Conclusions

We generalized the sexual version of the Penna model for biological aging to simulate also triploid populations. In these populations individuals present three sets of homologous chromosomes. A harmful mutation needs to appear in all of them to become active or to be in a locus where the harmful allele is dominant. We showed that normal diploid sexual populations have higher survival rates and are larger than the triploid ones, for the same carrying capacity of the environment. The genetic diversities of the populations, measured by the Hamming distance distribution, were also compared. The result shows that the diploid population presents a higher genetic diversity than the triploid ones. We may thus conclude that, concerning survival probabilities, population sizes and diversity, usual sex is better than that involving three individuals. Sexual diploid populations would be favored in direct competition with triploids either in a stable or a mutating environment, and we claim that this is the reason for normal sex to have been chosen by evolution as the dominant reproduction strategy. These results mean that the fact that triploids need mutations to appear in all three sets of homologous cromossomes to become effective is not enough to overcome the effort of a mating involving three individuals, and that more is not necessarily better insofar as genetic diversity is concerned.

Acknowledgments: We thank D. Stauffer and P.M.C. de Oliveira for a critical reading of the manuscript, and the Brazilian agencies CNPq, FAPERJ and CAPES for partial financial support.


text:
# Are there cosmological evolution trends on Gamma-Ray Burst features?

## 1 Introduction

The study of gamma-ray bursts (GRBs) afterglows has enable the measurement of their distances, so far GRBs are known as an explosive phenomenon occurring at cosmological distances, emitting large amount of energy mostly in the gamma-ray range (see, e.g. Piran 1999; Cheng & Lu 2001 for a review). So GRBs can provide useful information about the early epochs in the history of the universe.

Although a great achievements have been made about the GRB afterglows, we still know little about the origin of gamma-ray bursts, the reason is that the GRBs with known redshifts are relatively rare, now there are only about 20 GRBs with known redshifts, so it is difficult to do some statistics about GRBs features, such as their luminosity function, duration distribution, etc.. However, two important correlations have been discovered, i.e. between the degree of variability of the gamma-ray burst light curve and the GRB luminosity (Ramirez-Ruiz & Fenimore 1999; Feminore & Ramirez-Ruiz 2001), and between the differential time lags for the arrival of burst pulses at different energies and the GRB luminosity (Norris, Marani & Bonnell 2000), although these correlations are still tentative, they offer the possibility to derive independent estimates of the redshifts of GRBs.

Recently Amati et al. (2002) have reported the spectral and energetic properties of several GRBs with known redshifts, these bursts were all detected by BeppoSAX and have good-quality time-integrated spectra. In addition, Lloyd-Ronning & Ramirez-Ruiz (2002) have found that bursts with highly variable light curves have greater $`\nu F_\nu `$ spectral peak energies in their cosmological rest frames. These results reinforce the validity of the redshift estimates derived from the luminosity - variability relation. Here we will discuss the properties of two group GRBs, one group includes the bursts with known redshifts and well-defined spectra detected by BeppoSAX, and another group consists of bursts whose redshifts are derived from the luminosity \- variability relation. We will show that the properties of these two group GRBs are all correlated with redshift, which suggests that the luminosity - variability relation may be reliable, and furthermore the GRBs’ features are redshift dependent.

## 2 The Properties of Two Group GRBs

Amati et al. (2002) have analysed the spectral properties of the X-ray and gamma-ray emission from GRBs with known redshifts. Those bursts were all detected by BeppoSAX satellite and have good-quality spectra. The extension of the spectral analysis to the X-ray energy band allows a better determination of the continuum spectrum, reducing the bias in the measurement of the spectral slope below the peak energy $`E_p`$ of $`\nu F_\nu `$ spectra. In their sample total 12 gamma-ray bursts were included, however among them, there are three bursts (GRB980326, GRB980329 and GRB000214) whose redshifts are not determined accurately, only redshift intervals are given, so in our one group we ignore these three bursts, and contains the other 9 bursts with firm redshifts.

GRB temporal profiles are so complicated that, at first sight, their behavior obeys no simple rule. However, several authors have suggested that there may be correlation between the properties of burst time structure and burst luminosity (e.g. Feminore & Ramirez-Ruiz 2001; Reichart et al. 2001; Norris, Marani & Bonnell 2000). In particular, Fenimore & Ramirez-Ruiz (2001) explored the possibility of using the ”spikiness” of the time structure, combined with the observed flux, to obtain the GRB redshifts. They have analysed several hundred long and bright bursts and derived the redshifts and luminosities for 220 BATSE bursts. Lloyd-Ronning & Ramirez-Ruiz (2002) used total 159 bursts from the above 220 GRBs to investigate the dependence of the burst spectra on variability. The observed spectra of these 159 bursts can be well characterized by the Band function (Band et al. 1993), defined by a low energy spectral index $`\alpha `$, a high energy spectral index $`\beta `$, and a peak energy of the $`\nu F_\nu `$ spectra $`E_p`$. So our second group consists of these 159 bursts with known peak energy, and the redshifts of them are derived from the luminosity - variability relation (Feminore & Ramirez-Ruiz 2001).

Fig.1 gives the GRB duration in their cosmological rest frame, $`T^{}=T/(1+z)`$, versus the redshift, where $`T^{}`$ is the intrinsic duration and $`T`$ is the observed duration of a burst. The filled circles are bursts with secure redshifts estimates, while the empty triangles are bursts in which the redshifts are derived using the luminosity - variability distance indicator. Fig.2 and Fig.3 show the luminosity and radiated energy of GRBs in their cosmological rest frame versus the redshift, Fig.4 shows the relation between the peak energy $`E_p`$ of the $`\nu F_\nu `$ spectra in their cosmological rest frame and the redshift.

From Fig.1 we see that although the distribution of the burst duration is somewhat scatter, there is still a clear trend that the intrinsic duration decreases with the redshift, if we fit a power law to the data, we find $`T^{}(1+z)^{0.85\pm 0.08}`$. Fig.2 and Fig.3 show that the isotropic luminosity and radiated energy have a positive correlation with the redshift, excluding the 4 bursts with smallest redshifts, we have $`L(1+z)^{2.5\pm 0.1}`$, $`E(1+z)^{1.62\pm 0.14}`$. Fig.4 shows that the peak energy of $`\nu F_\nu `$ spectra also has a positive correlation with the redshift, the power law fit is $`E_p(1+z)^{0.76\pm 0.07}`$.

However, we know that the BATSE GRB sample is flux truncated, i.e. only those bursts whose flux exceeding the threshold flux can be detected, this effect can give rise to several apparent correlations. For example, the burst duration (e.g. $`T_{90}`$) are defined by the time interval over which $`90\%`$ of the fluence is detected. High-z bursts are fainter, therefore some faint components of an intrinsic long burst could have been cut out by the flux threshold effect, so that their apparent durations are shortened. For another example, the correlation trend in both Fog.2 and Fig.3 are also contaminated by the flux truncation effect, as has been demonstrated by Lloyd-Ronning et al.(2002). So we need some way to estimate the true correlations between parameters that suffering from flux truncation.

Here we use the Monte Carlo simulation method. We simulate the burst sample, assuming their luminosity satisfy $`L=L_0(1+z)^\alpha `$, and their redshifts are distributed uniformly between 0 and 10. In our simulation, we take the value $`L_0=3\times 10^{51}`$ erg s-1, the flux threshold is $`10^7`$ erg cm-2 s-1, and the Hubble constant $`H_0=65`$ km s-1 Mpc-1, $`\mathrm{\Omega }_m=0.3`$, $`\mathrm{\Omega }_\mathrm{\Lambda }=0.7`$. Under this circumstance, we produce 300 bursts in original sample, and find that 100 bursts of them with flux lower than the threshold, so there are 200 bursts in our ”observed” sample. Fig.5 shows our simulation results, where the circles are the 300 bursts produced by simulation, and the solid line is the flux threshold for BATSE. We find that, when $`\alpha =1.7\pm 0.5`$, there is the relation $`L(1+z)^{2.50\pm 0.08}`$, which is consistent with the observed value. So we see that the flux truncation has great effect on the correlation coefficient.

In order to discuss this flux truncation effect on the burst duration, one needs to know the intrinsic pulse profile. For simplicity, here we assume the intrinsic pulse profile $`L=L_\mathrm{p}e^{t^2/t_0^2}`$, since only the component with $`LL_{\mathrm{lim}}`$ can be observed (where $`L_{\mathrm{lim}}`$ is the limit luminosity corresponding to the flux threshold), so we can get the observed pulse duration $`T_{}=t_0(\mathrm{ln}\frac{L_\mathrm{p}}{L_{\mathrm{lim}}})^{1/2}`$. The intrinsic pulse duration should be proportional to $`t_0`$, which is assumed to have the form $`t_0(1+z)^\beta `$, $`L_\mathrm{p}`$ is assumed to take the form described above, $`L_\mathrm{p}(1+z)^{1.7\pm 0.5}`$. As above, we produce a GRB sample with 300 bursts, 100 bursts of them are not ”observed” due to their low flux, so the remaining 200 bursts consists of our ”observed” sample. We find that when we take the value of $`\beta =0.45\pm 0.20`$, the ”observed” sample has the relation $`T_{}(1+z)^{0.84\pm 0.14}`$, which is compatible with the observed value.

Fig.6 shows the relation between the intrinsic peak energy of the $`\nu F_\nu `$ spectra and the luminosity for group one, the solid line represents the relation $`E_pL^{1/2}`$. Fig.7 gives the intrinsic peak energy of $`\nu F_\nu `$ spectra versus the luminosity for group two, the two solid lines represent the relation $`E_pL^{1/2}`$. We see that the relation $`E_pL^{1/2}`$ can account for the observed data quite well.

## 3 Discussion and conclusion

In previous section we have shown that the properties of our two group GRBs, one group includes 9 GRBs with secure redshifts derived from the afterglow observation, another group consists of 159 GRBs with known peak energy and their redshifts are derived from the luminosity - variability relation, are all correlated with the redshifts, which reinforces the validity of the luminosity - variability distance indicator. Our results further support the conclusion obtained by Lloyd-Ronning & Ramirez-Ruiz (2002), they found that there is positive correlation between the peak energy in the cosmological rest frame and the variability for gamma-ray bursts whose redshifts are derived either from optical spectral features or from the luminosity - variability distance indicator.

Since the BATSE GRB sample is flux truncated, i.e. only those bursts with flux exceeding the threshold flux can be detected, so it is necessary to discuss this selection effects on the statistical results. Here we use the Monte Carlo simulation method. We find that, after correcting the selection effects, the cosmological evolution trends are still exist, although they are much shallower than the trends found in the face values of the data. For example, if ignored the flux truncation effect, we obtained the relation $`L(1+z)^{2.5\pm 0.1}`$, while when the selection effects are taken into account, we got the relation $`L(1+z)^{1.7\pm 0.5}`$, which is consistent with the value $`L(1+z)^{1.4\pm 0.5}`$ obtained by Lloyd-Ronning et al. (2002).

Fig.5 illustrates our simulated results, where the circles are the 300 bursts produced by the simulation, the solid line is the flux threshold for BATSE, and the dashed line is the flux threshold for Swift. From this figure it is obvious that the property of the ”observed” sample should depend on the adopted flux threshold, for different flux threshold their statistical properties are different. For example, when taking the BATSE flux threshold, there are 100 bursts with flux lower than the threshold, while when taking the Swift flux threshold, there are only 7 bursts with flux lower than the threshold, so in this case the observed $`Lz`$ relation should close to the true $`Lz`$ relation. Therefore we expect that the luminosity - redshift relation for bursts observed by Swift should be shallower than that for bursts observed by BATSE.

Fig.6 and Fig.7 show that there is a good correlation between the peak energy and luminosity for both group GRBs, and the relation $`E_pL^{1/2}`$ can account for the observed data quite well. Up to now the location of the GRB emission site is still unsettled, although the internal shock model is thought to be more reasonable than the external shock model. Zhang & Meszaros (2002b) analysed various fireball models within a unified picture and investigated the $`E_p`$ predictions of different models. It is known that for internal shock model, if the GRB bulk Lorentz factors are not correlated with the luminosities, then there is the relation $`E_pL^{1/2}`$, while for external shock model $`E_p\mathrm{\Gamma }^4`$, where $`\mathrm{\Gamma }`$ is the shock Lorentz factor. So our results suggest that the gamma-ray burst emission are more likely from the internal shock.

Frail et al. have discussed the afterglow properties of several GRBs with known redshifts, they assumed that the breaks in the afterglow light curves are caused by the sideways expansion of the jet, and then they concluded that the GRB emission energy is nearly a constant, $`E5\times 10^{50}`$ ergs (Frail et al. 2001). However, Fig.3 shows that the isotropic radiated energy increases with redshift, so if the conclusion of Frail et al. is true, then we notice that the jet opening angle must decrease with the redshift, as shown by Fig.8, where the data points are taken from the paper of Frail et al. (2001). This point is very interesting, since it can put constraints on the central engines of GRBs. Of course, this phenomena can also be explained within the framework of a structured universal jet model (Zhang & Meszaros 2002a; Rossi et al. 2002). In this model, an observer closer to the jet axis would detect a higher luminosity, thus at higher redshift, smaller viewing angle detections are preferred due to luminosity selection effect.

It should be noted that in our analysis we have not considered the errors coming from the luminosity - variability relation. We know that the relation between luminosity and variability is somewhat scatter, the correlation coefficient $`n`$ ($`LV^n`$) can be range from 2.2 to 5.8 (the best value is 3.3, see Fenimore & Ramirez-Ruiz 2001), so it is natural that the redshifts and luminosities inferred from this luminosity - variability relation should have large errors, and these errors should be somehow transferred to the final errors in the correlation indices. However this effect is very complicated, we hope that this effect can be taken into account in the future work.

In summary, in this paper we discuss the properties of two group GRBs, one group with known redshift from afterglow observation, and another group with redshift derived from the luminosity - variability relation. We find that the properties of these two group GRBs are all correlated with the redshifts, which reinforces the validity of the redshift estimates derived from the luminosity \- variability relation. If this is true, then we see that the burst features, such as their intrinsic duration, luminosity, radiated energy and peak energy $`E_p`$, are all redshift dependent, which means that there are cosmological evolution effects on gamma-ray bursts features, and this can provide an interesting clue to the nature of GRBs.

## 4 acknowledgements

We are very grateful to the referee for several important comments that greatly improved this paper. This work is supported by the National Natural Science Foundation (grants 10073022, 10233010 and 10225314) and the National 973 Project on Fundamental Researches of China (NKBRSF G19990754).


text:
# Mass Limits to Primordial Star Formation from Protostellar Feedback

## 1 The Collapse of Primordial Gas Clouds

Recent numerical studies have followed the gravitational collapse of perturbations from cosmological to almost stellar dimensionsAbel, Bryan, and Norman (2002); Bromm, Coppi, and Larson (1999). Baryon-dominated clouds, cooled to about 200-300 K by trace amounts of molecular hydrogen, form at the centers of dark matter halos. For $`n_\mathrm{H}>10^4\mathrm{cm}^3`$ the cooling rate becomes independent of density, and so the dissipation of gravitational energy in the densest regions gradually raises the temperature. In the simulation of Abel et al.Abel, Bryan, and Norman (2002), the gas cloud is centrally-concentrated ($`\rho r^{k_\rho };k_\rho 2.2`$) and is contracting quasi-hydrostatically with infall speeds about one third of the sound speed. In addition to thermal support, the cloud is filled with a turbulent cascade of weak shocks (T. Abel, private comm.). The structure can be described by an approximately hydrostatic polytrope with $`\gamma _p=1+1/n=k_P/k_\rho =2(11/k_\rho )=1.1`$, where $`Pr^{k_P}`$. The contraction is akin to the maximally sub-sonic HunterHunter (1977) settling solution, for which the accretion rate is a factor $`\varphi _{}2.6`$ greater than the classic ShuShu (1977) solution. This accretion rate can be expressed in terms of the entropy parameter of the polytrope, $`K=P/\rho ^{\gamma _p}`$ and the collapsed mass, $`Mm_{}`$,

$$\dot{m}_{}=\frac{8\varphi _{}}{\sqrt{3}}\left[\frac{(3k_\rho )k_P^3K^3}{2(2\pi )^{53\gamma _p}G^{3\gamma _p1}}\right]^{\frac{1}{2(43\gamma _p)}}M^j0.026K^{15/7}\left(\frac{M}{M_{}}\right)^{3/7}M_{}\mathrm{yr}^1,$$
(1)
where $`j3(1\gamma _p)/(43\gamma _p)`$. The numerical evaluation assumes $`\varphi _{}=2.6`$, since $`\gamma _p`$ is not too different from one, and $`K=1.88\times 10^{12}(T/300\mathrm{K})n_{\mathrm{H},4}^{0.1}\mathrm{cgs}1.88\times 10^{12}K^{}\mathrm{cgs}`$. We set the temperature normalization a factor $`4/3`$ higher than is seen in simulationsAbel, Bryan, and Norman (2002), to allow for partial pressure support from sonic and isotropic turbulent motions; i.e. $`T`$ is an effective temperature. The small ratio of turbulent to thermal support is in marked contrast to contemporary massive star formationMcKee & Tan (2002). In primordial clouds it is the microphysics of $`\mathrm{H}_2`$ cooling that determines both the evolution (via $`\gamma _p`$) and normalization (via $`K^{}`$ and $`\varphi _{}`$) of the accretion rate. Collapse assuming constant $`K`$ \- “isentropic accretion” - agrees with 1-D numerical studiesOmukai and Nishi (1998); Ripamonti et al. (2002) (Fig. 1a).

The mass-averaged rotation speed is a fraction, $`f_{\mathrm{Kep}}0.5`$, of Keplerian, approximately independent of radiusAbel, Bryan, and Norman (2002). Assuming angular momentum is conserved inside the sonic point, $`r_{\mathrm{sp}}`$, leads to a disk size $`r_d=f_{\mathrm{Kep}}^2r_{\mathrm{sp}}=3.4(f_{\mathrm{Kep}}/0.5)^2(M/M_{})^{9/7}K^{10/7}\mathrm{AU}`$. Matter falls onto this disk at all radii $`rr_d`$, as well as directly to the star. We follow UlrichUlrich (1976) in describing the density distribution of the rotating, freely-falling envelope.

## 2 Protostellar Evolution

At densities $`n_\mathrm{H}10^{16}\mathrm{cm}^3`$ an optically thick protostellar core formsOmukai and Nishi (1998), bounded by an accretion shock. The size of the protostar then changes as it accretes matter and radiates energy. For spherical geometry, the high accretion rates typical of primordial star formation lead to optically thick conditions above the accretion shockStahler, Palla, and Salpeter (1986); Omukai and Palla (2001). Accretion energy is advected, which swells the star. However, for collapse with angular momentum, a disk forms and we expect photospheric conditions over much of the stellar surface. We employ a polytropic stellar structure and extend the energy equation model of Nakano et al.Nakano (2000), that includes gravitational, ionization-dissociation, and D-burning energies, to allow for optical depth in front of the accretion shock. We also allow for some fraction of disk accretion, depending on $`f_{\mathrm{Kep}}`$. We model the “luminosity wave” expansion featureStahler, Palla, and Salpeter (1986); Omukai and Palla (2001) as a relaxation to a more compact ($`n=3`$) state after the star is older than its Kelvin time: to conserve energy the outer radius of the star grows by a factor, which we estimate from Omukai & PallaOmukai and Palla (2001). We set $`n=1.5`$ or 3, depending on convective stability. During optically thick accretion, before the luminosity wave, we derive an effective value for $`n`$ by comparing to the results of Stahler et al.Stahler, Palla, and Salpeter (1986). Figure 1b shows this case and isentropic models with $`f_{\mathrm{Kep}}=0`$ and $`0.5`$. The two spherical models undergo similar evolution, while the rotating model forms a smaller protostar once the photosphere retreats to the stellar surface. The protostar is then supported by D core and shell burning, as in present-day star formation. The large accretion rate causes the star to join the main sequence only at relatively high masses. The photospheric temperature is much hotter than in the spherical models, leading to significant fluxes of ionizing and FUV radiation. These form the inputs for the feedback model.

## 3 Feedback versus Accretion

The lifetimes of primordial massive stars converge to about 2 MyrSchaerer (2002) so eq.(1) implies an upper limit to stellar masses of $`2000M_{}`$. However, other feedback processes are likely to intervene before this. Once the flux of ionizing photons from the protostar is greater than that of neutral H to its surface, an HII region forms. Approximating sectors as independent, we calculate the extent of the ionized region. Accretion may be suppressed if the HII region expands to distances greater than $`r_g`$, where the escape speed equals the ionized gas sound speed, $`10\mathrm{km}\mathrm{s}^1`$. For the fiducial model this occurs at the poles when $`m_{}=90M_{}`$ (Fig. 1c) and at the equator when $`m_{}=140M_{}`$. In this calculation we assumed a free-fall density distribution. In reality the ionizing radiation forceOmukai and Inutsuka (2002) decelerates and deflects the flow. For collapse with angular momentum, most streamlines do not come too close to the star, so the effect discussed by Omukai & InutsukaOmukai and Inutsuka (2002) of HII region quenching due to enhanced densities is greatly weakened. In fact, we anticipate that deflection is more important in reducing the concentration of inflowing gas near the star, so that the HII region becomes larger.

A second feedback effect is radiation pressure from Ly-$`\alpha `$ photons created in the HII region and FUV photons emitted by the star. These photons are trapped by the Lyman series damping wings of the neutral gas infalling towards the HII region. The energy density builds up until the escape rate, set by diffusion in frequency as well as in spaceAdams (1972), equals the input rate. The resulting pressure acts against the infall ram pressure. We use the results of NeufeldNeufeld (1991) to aid our numerical calculations. Radiation pressure becomes greater than twice the ram pressure at $`m_{}20M_{}`$ for typical $`f_{\mathrm{Kep}}`$. The enhancement of radiation pressure above the optically thin limit is by a factor $`10^3`$. Infall is first reversed at the poles, which would allow photons to leak out and reduce the pressure acting in other directions. We shall examine this scenario in a future study.

We thank T. Abel, V. Bromm, B. Draine, J. Goodman, D. McLaughlin, and J. Ostriker for helpful discussions. JCT is supported by a Spitzer-Cotsen fellowship from Princeton University and NASA grant NAG5-10811. The research of CFM is supported by NSF grant AST-0098365 and by a NASA grant funding the Center for Star Formation Studies.


text:
# Symplectic Resolutions for Coverings of Nilpotent Orbits

## 1 Introduction

Let $`G`$ be a simply connected semi-simple complex Lie group with Lie algebra $`𝔤`$. Let $`𝒪`$ be a nilpotent adjoint orbit in $`𝔤`$. The fundamental group $`\pi _1(𝒪)`$ of $`𝒪`$ is finite and in general non-zero (see corollary 6.1.6 and section 8.4 of \[CM\]). Let $`p:M𝒪`$ be a $`G`$-homogeneous covering of degree $`d`$. We denote by $`X`$ the normalization of $`\overline{𝒪}`$ in the function field of $`M`$, which will be called a covering of the nilpotent orbit $`\overline{𝒪}`$. Then proposition 1.2 of \[BK\] says that $`X`$ contains $`M`$ as a Zariski open subset and the map $`p`$ extends to a finite surjective $`G`$-equivariant morphism $`\overline{p}:X\overline{𝒪}`$. Furthermore $`G`$ has finitely many orbits on $`X`$ and $`XM`$ is of codimension at least 2 in $`X`$. Recall that a regular 2-form on a smooth algebraic variety is symplectic if it is closed and non-degenerate at every point. The Kostant-Kirillov form on the orbit $`𝒪`$ gives a symplectic form $`\omega `$ on $`M`$.

###### Definition 1.

A symplectic resolution for $`X`$ is a projective resolution of singularities $`\pi :ZX`$ such that the 2-form $`\pi ^{}(\omega )`$ defined a priori on $`\pi ^1(M)`$ extends to a symplectic form $`\mathrm{\Omega }`$ on the whole of $`Z`$.

As shown in \[Fu\], a resolution is symplectic if and only if it is crepant. The purpose of this note is to consider the existence of symplectic resolutions, i.e. for which coverings of a nilpotent adjoint orbit there exists a symplectic resolution. It turns out that this is a difficult problem. The following proposition gives some examples of symplectic resolutions.

###### Proposition 1.

Let $`P`$ be a parabolic subgroup of $`G`$ and $`M`$ the unique open $`G`$-orbit in $`T^{}(G/P)`$. Let $`𝒪`$ be the orbit of a Richardson element in the nilradical $`𝔲`$ of the Lie algebra of $`P`$. Then the $`G`$-equivariant desingularization $`T^{}(G/P)X`$ is a symplectic resolution for $`X`$.

For the proof, see proposition 7.4 of \[BK\]. It should be pointed out that there exist some symplectic resolutions for coverings of nilpotent orbits which are not of the above form. An example is the double covering of the minimal nilpotent orbit in $`𝔰𝔭(2n,)`$ which is $`X=^{2n}\overline{𝒪}_{min}`$. However this is the only example where $`X`$ is non-singular (see theorem 4.6 \[BK\]). The purpose of this note is to prove the converse of this proposition under some additional hypotheses.

Recall that $`𝒪`$ and $`\overline{𝒪}`$ are both stable under the scaling action of $`^{}`$ on $`𝔤`$. This induces the Euler action of $`^{}`$ on $`𝒪`$ and $`\overline{𝒪}`$. Unfortunately this action does not lift to $`X`$ in general. An example is the double covering of the minimal nilpotent orbit $`𝒪_{min}`$ in $`𝔰𝔭(2n,)`$. However we have (see lemma 1.3 and proposition 1.4 \[BK\])

###### Proposition 2.

For any covering $`X`$ of $`\overline{𝒪}`$, there exists a $`^{}`$-action on $`X`$ which lifts the square of the Euler action of $`^{}`$ on $`\overline{𝒪}`$, i.e. for any $`\lambda ^{}`$ and $`xX`$, $`\overline{p}(\lambda x)=\lambda ^2\overline{p}(x)`$. Furthermore there exists a unique point $`oX`$ such that $`\overline{p}(o)=0`$. This $`o`$ is the unique $`G`$-fixed point in $`X`$ and also the unique $`^{}`$-fixed point in $`X`$.

For this $`^{}`$-action, we have $`\lambda ^{}\omega =\lambda ^2\omega `$ for any $`\lambda ^{}`$, which is different to our situation in \[Fu\], where the $`^{}`$-action satisfies $`\lambda ^{}\omega =\lambda \omega `$. This makes the situation more complicated here. However, under the hypothesis that the degree of the covering is odd, we can prove the following

###### Theorem 1.

Let $`𝔤`$ be a semisimple complex Lie algebra and $`X`$ a covering of a nilpotent adjoint orbit $`𝒪`$ in $`𝔤`$. Suppose that the degree $`d`$ of the covering $`\overline{p}:X\overline{𝒪}`$ is odd. Then for any symplectic resolution $`\pi :ZX`$, $`Z`$ is isomorphic to $`T^{}(G/P)`$ for some parabolic subgroup $`P`$ of $`G`$. Furthermore, under this isomorphism, the map $`\overline{p}\pi `$ becomes

$$T^{}(G/P)G\times ^P𝔲𝔤,(g,u)Ad(g)u,$$

where $`𝔲`$ is the nilradical of $`𝔭=Lie(P)`$.

This theorem generalizes our main theorem in \[Fu\], where we considered the case $`d=1`$.

## 2 Outline of the proof

###### Lemma 1.

The Euler action of $`^{}`$ on $`\overline{𝒪}`$ lifts to $`X`$ and for this action $`\lambda ^{}\omega =\lambda \omega `$.

###### Proof.

Let $`R`$ be the regular functions ring of $`X`$. For $`k`$, let $`R[k]=\{\varphi R|\lambda \varphi =\lambda ^k\varphi ,\lambda ^{}\}`$. Proposition 1.4 of \[BK\] implies that if the degree of the covering $`X\overline{𝒪}`$ is odd, then $`R[k]=0`$ for $`k`$ odd, i.e. $`R=_0^{\mathrm{}}R[2k]`$. This gives that the Euler action of $`^{}`$ (not only its square!) on $`\overline{𝒪}`$ lifts to a $`^{}`$-action on $`X`$. For this action, we have $`\lambda ^{}\omega =\lambda \omega `$. ∎

###### Lemma 2.

The action of $`G\times ^{}`$ on $`X`$ lifts to $`Z`$ such that $`\pi :ZX`$ and $`\overline{p}\pi :Z\overline{𝒪}`$ is $`G\times ^{}`$-equivariant. For the $`^{}`$-action on $`Z`$, we have $`\lambda ^{}\mathrm{\Omega }=\lambda \mathrm{\Omega }`$ for any $`\lambda ^{}`$.

Note that the action of $`G`$ (resp. $`^{}`$) on $`M`$ lifts to $`\pi ^1(M)X`$. Now the proof goes along the same line as proposition 3.1 of \[Fu\]. Using the two lemmas, we can apply our analysis in \[Fu\] to complete the proof of theorem 3. For reader’s convenience, we give an outline of our method.

Let $`Z^{^{}}`$ be the fixed points subvariety in $`Z`$. Since $`\overline{p}\pi `$ is proper and the $`^{}`$-action on $`\overline{𝒪}`$ extends to a $``$-action, the valuative criterion of properness shows that there exists an attraction $`q:ZZ^{^{}}`$. As show by lemma 3.5 \[Fu\], the $`^{}`$-action on the smooth connected component $`Z_0Z^{^{}}`$ containing $`q(\pi ^1(M))`$ is definite, i.e. for $`zZ_0,T_z^kZ=\{vT_zZ|\lambda _{}v=\lambda ^kv\}`$ is zero if $`k<0`$. Now the equation $`\lambda ^{}\mathrm{\Omega }=\lambda \mathrm{\Omega }`$ implies a duality between $`T_z^kZ`$ and $`T_z^{1k}Z`$. Thus for $`zZ_0`$, we have $`T_zZ_0T_z^0ZT_z^1Z`$ and $`T_z^kZ=0`$ if $`k0,1`$. This shows that $`Z_0`$ is Lagrangian in $`Z`$. By a classical result of Bialynicki-Birula, the attraction $`q:q^1(Z_0)Z_0`$ is $`^{}`$-equivariantly isomorphic to $`T^{}Z_0Z_0`$.

Now we study the $`G`$-action on $`Z_0`$, which has an open-dense orbit, namely $`q(\pi ^1(M))`$. In fact this orbit is the whole of $`Z_0`$. Note that $`Z_0Z^{^{}}\pi ^1(o)`$ is projective, so $`Z_0=G/P`$ for some parabolic subgroup $`P`$ of $`G`$. To complete the proof, we need to show that $`Z_0=Z^{^{}}`$. Note that $`\pi ^1(o)`$ is connected, so we need only to show that $`Z_0`$ is a connected component of $`\pi ^1(o)`$ or equivalently a connected component of $`(\overline{p}\pi )^1(0)`$, this is proved by using the explicit formula for the map $`T^{}(G/P)\overline{𝒪},(g,u)Ad(g)u`$.

## 3 Some corollaries and examples

###### Corollary 1.

Let $`X`$ be a covering of odd degree of a nilpotent orbit $`\overline{𝒪}`$. If $`𝒪`$ is not a Richardson orbit, then $`X`$ admits no symplectic resolution.

This follows directly from the above theorem and proposition 2. Consider the nilpotent orbits $`2A_2,2A_2+A_1,A_5,E_6(a_1),E_6`$ in the exceptional Lie algebra of type $`E_6`$, here we use notations from section 8.4 of \[CM\]. All these nilpotent orbits have $`/3`$ as fundamental group, but none of them is a Richardson orbit. The above corollary implies that the three-fold covering of any of the above orbits admits no symplectic resolution.

###### Corollary 2.

A covering $`X`$ of odd degree $`d>1`$ of a nilpotent orbit $`𝒪`$ in $`𝔰𝔩(n,)`$ admits no symplectic resolution.

This follows from the above theorem and theorem 3.3 \[Hes\], which says that for any polarization $`P`$ of $`𝒪𝔰𝔩(n,)`$, the morphism $`T^{}(G/P)\overline{𝒪}`$ is birational. Recall that every nilpotent orbit in $`𝔰𝔩(n,)`$ admits a symplectic resolution (\[Fu\]), however the above corollary shows that their odd higher degree coverings admit no symplectic resolution.

As an example, let $`𝒪`$ be the principal nilpotent orbit in $`𝔰𝔩(3,)`$, then $`\pi _1(𝒪)=/3`$. Let $`X`$ be the three-fold covering of $`\overline{𝒪}`$. Then by the above corollary, $`X`$ does not admit any symplectic resolution. In fact, $`X`$ is a so-called shared orbit in \[BK\], more precisely $`X`$ is nothing but the closure of the minimal nilpotent $`\overline{𝒪}_{min}`$ in the exceptional Lie algebra $`𝔤_2`$, which is known to admit no symplectic resolution \[Fu\].

In the list of \[BK\], shared orbits do not admit any symplectic resolution except in the following case: let $`𝒪`$ be the nilpotent orbit in $`𝔰𝔭(2n,)`$ corresponding to the partition $`[2,2,1,\mathrm{},1]`$. Then $`\pi _1(𝒪)=/2`$. By results of \[Fu\], the nilpotent orbit $`\overline{𝒪}`$ does not admit any symplectic resolution. Let $`X`$ be the double covering of $`\overline{𝒪}`$. Then $`X`$ is isomorphic to $`\overline{𝒪}_{min}`$ in $`𝔰𝔩(2n,)`$, which admits a symplectic resolution.

As shown by the above examples, there exist some nilpotent orbits which admit some symplectic resolutions, but not their coverings, and there exist some nilpotent orbits which do not admit any symplectic resolution, while some of their coverings do admit some symplectic resolutions. This indicates that the problem to determine which covering admits a symplectic resolution might be difficult and interesting.

## 4 A conjecture

From now on, we suppose that $`𝔤`$ is a simple complex Lie algebra. Let $`𝒪`$ be a nilpotent adjoint orbit in $`𝔤`$ and $`X`$ a covering of $`\overline{𝒪}`$. Recall that $`R[k]=\{\varphi R|\lambda \varphi =\lambda ^k\varphi ,\lambda ^{}\}`$. Let $`𝔤^{}=R[2]`$, which is also a simple Lie algebra (see theorem 4.2 of \[BK\]). Let $`G^{}`$ be the simply connected Lie group with Lie algebra $`𝔤^{}`$, which is the maximal connected Lie group of holomorphic symplectic automorphisms of $`X`$ (see corollary 7.3 \[BK\]). By theorem 3.1 \[BK\], $`X`$ contains a Zariski open subset $`M^{}`$ such that $`MM^{}`$, and there exists a finite $`G^{}`$-covering $`M^{}𝒪^{}`$ to a nilpotent orbit $`𝒪^{}`$ in $`𝔤^{}`$. We propose the following conjecture on symplectic resolutions for $`X`$:

###### Conjecture 1.

Let $`\overline{p}:X\overline{𝒪}`$ be a covering of a nilpotent orbit in a simple Lie algebra $`𝔤`$ such that $`R[1]=0`$. Suppose that $`\pi :ZX`$ is a symplectic resolution. Then there exists some parabolic sub-group $`P^{}`$ of $`G^{}`$ such that $`Z`$ is isomorphic to $`T^{}(G^{}/P^{})`$ and the map $`\pi `$ becomes

$$T^{}(G^{}/P^{})X.$$

###### Remark 1.

It is proved in \[BK\] (theorem 4.6) that $`R[1]0`$ if and only if $`X`$ is the double covering of $`\overline{𝒪}_{min}𝔰𝔭(2n,)`$. In this case, $`X=^{2n}`$.

###### Remark 2.

By theorem 7.5 \[BK\], if $`T^{}(G/P)X`$ is a symplectic resolution for $`X`$, then there exists some parabolic subgroup $`P^{}`$ of $`G^{}`$ such that $`G/P=G^{}/P^{}`$.

###### Remark 3.

If the covering $`\overline{p}:X\overline{𝒪}`$ is of odd degree, then by our theorem 1, the conjecture is true.

Acknowledgements: I want to thank A. Beauville and M. Brion for helpful discussions and suggestions. I am grateful to the referee for the careful reading and pertinent remarks.