text:
# Variation of Fine Structure Constant from Non-Universal Gravity

(August 28, 2002)

## Abstract

We relate the reported variation in the value of the fine structure constant to a possible non-universality of the gravitational interaction with respect to different particle generations.

Recent measurement of the fine structure constant $`\alpha `$ in the remote, yet on cosmological scale recent past yielded a $`4\sigma `$ deviation $`\mathrm{\Delta }\alpha /\alpha =(0.72\pm 0.18)10^5`$ from laboratory value within the red shift domain $`z`$0.6–2 Web01 . This poses the question about the origin of such an effect Ban02 ; Dav02 , identified as either temporal (general presumption) or spatial variation. Since the experimental effect which is observed is seen within a scheme of a systematic change in absorption lines by nebula in line of sight of quasar illumination of laboratory sensors, the effect could be originating in properties of the nebula. Such spatial interpretation of the data analysis is more compatible with other limits on variation of $`\alpha `$ Web01 .

As recent flurry of activity has shown that it is quite difficult to find a self consistent interpretation of this experimental result. We propose here a quite different approach:
a) We interpret the change in $`\alpha `$ as being due to fermion masses $`m_f`$ which enter vacuum polarization loops.
b) We will argue that the effect cannot involve electrons, thus only the heavy fermions such as the muon contribute.
c) Assuming that heavy fermion masses have varied by about -0.15–0.45% this experimental observation can be understood.
d) There is no consistent interpretation assuming time variation of $`m_f`$ (or $`\alpha `$) since there is no sensitivity to the time evolution of the late Universe we live in.
e) When the effect is associated with spatially localized gravitational mass it requires non-universal gravitational coupling.
f) We show how a theory of gravity is formulated in principle which combines different ratios of gravitational to inertial mass for the three fermion generations. This is not in contradiction with the Eötvös experimental evidence for universality of the ratio of inertial to gravitational mass, and all principal results of Einstein gravity and cosmology are retained.
g) The 1000-fold enhanced strength of gravity for 2nd and/or 3rd generation, which is required to understand $`\mathrm{\Delta }\alpha /\alpha `$, also explains the dark matter in terms of the more strongly gravitating cosmic background 2nd and/or 3rd generation neutrinos.
h) We do not present here a theoretical model how a gauge universal, gravity non-universal, theory could arise.

Conventional theoretical framework: Within the current paradigm of a unified structure of all interactions, the value of $`\alpha `$, the electromagnetic coupling in the infrared large distance $`r\mathrm{}`$ limit is not fundamental, the determining value is the strength of the universal interaction at the unification scale, which in case of ‘Grand Unification’ is $`M_{\mathrm{GUT}}310^{16}`$ GeV. As energy scale diminishes, broken symmetries emerge and individual interactions separate, and their strength evolves with the scale. Electromagnetic (EM) interaction emerges when the electro-weak (EW) interactions separate at the scale of $`𝒪`$(100 GeV), and where Bur01 ; Jeg01 ,

$$\alpha (M_Z=91.2\text{GeV})=\frac{1}{128.936(46)}.$$
(1)
Unification of interactions poses a problem for an apparent time variation of $`\alpha `$, since it implies a much larger associated variation in strong interaction coupling parameters, implicating significant changes in quantities such as the nucleon mass and nuclear interactions Lan02 ; Cal02 , and indeed leading to possible and probable disagreement with latest limits on cosmological variation of e.g. the proton to electron mass Lev02 . In addition we need to resolve the question, how can sensitivity to the structure of the universe arise at a relatively late stage of its evolution, or if spatial interpretation is preferred, in localized nebula.

Introduction of relatively large compactified extra dimensions reduces the unification scale considerably, to the TeV range Die98 . Association of the time dependence of $`\alpha `$ or as we argue here $`m_f`$ with large extra dimensions has the inherent flaw that these dimensions do not evolve beyond the early, formative period of the Universe Die02 .

The case is quite strong today, in view of the many ramifications Oli02  of a time-changing strong interaction strength accompanying the change in $`\alpha `$, that it is necessary to search for other possible sources of, on cosmological scale ‘nearby’, time, or spatial, dependence of $`\alpha `$. Regarding the forthcoming discussion of the three particle generations we recall that it is not understood today how the masses of the ‘heavy’ fermion families are generated from a unified theory involving radically different scales (hierarchy problem). It is expected that in a unified theory which comprises the quantum gravity, not only the interactions, but also the three particle family generations will be naturally explained.

The fundamental generation which is at the origin of the matter around us, also provides the measuring sticks in the Universe. The unified theory than explains the mass ratios such as the muon to electron mass is $`m_\mu /m_e=206.77`$. Furthermore, if electron mass (along with other lepton masses) were to vary, in its units the nuclear world around us would evolve, and this is not safe from contradictions with precision experiments. One can more generally argue considering the impact of masses on fine-tuned properties of known interactions that all the fundamental generation masses ($`m_e,m_{\nu _e},m_d,m_u`$) must not vary. Thus only variations in the masses of the second ($`m_\mu ,m_{\nu _\mu },m_s,m_c`$) and third ($`m_\tau ,m_{\nu _\tau },m_b,m_t`$) generation could be considered to be at the origin of the variation of $`\alpha `$.

Although we make this argument with time dependence of $`\alpha `$ in mind, it applies equally for the case we explore below of a spatial variation of $`\alpha `$. When we speak of variation of mass, we also include in this the gravitational mass defect.

Relation of $`\alpha `$ to fermion masses: The scale evolution of the EM coupling strength $`\alpha \alpha (M=0)=1/137.03599976(50)`$, arises from the interaction of photons with the virtual fluctuations of charged fermi particle pairs of mass $`m_f`$ with mass thresholds within the range of energy scale considered. Quantum electrodynamics allows to evaluate this vacuum polarization effect for leptons with great precision. For our consideration it is sufficient to recall the lowest order in electromagnetic interaction, and the large $`M/m`$ behavior:

$`\mathrm{\Delta }\alpha ^1`$ $``$ $`{\displaystyle \frac{1}{\alpha (M)}}{\displaystyle \frac{1}{\alpha (0)}}`$
$`=`$ $`{\displaystyle \frac{2}{3\pi }}{\displaystyle \underset{f}{}}Q_f^2N_f\left(\mathrm{ln}{\displaystyle \frac{M}{m_f}}{\displaystyle \frac{5}{6}}+𝒪({\displaystyle \frac{m_f}{M}})\right).`$
Here $`Q_f`$ is the charge of the fermion $`f`$ and $`N_f`$ is the internal degeneracy, such as color of quarks (=3). $`M`$ is the scale at which the interaction is stable as function of time, and/or environment. In our present considerations we assume that all the variation occurs within the energy scale domain below $`M_Z`$, in which EM emerged as separate interaction. This does not restrict the validity of our principal argument.

Eq. (Variation of Fine Structure Constant from Non-Universal Gravity) can be applied to evaluate the contribution of the electron, muon and tau $`f=e,\mu ,\tau `$ leptons, and to estimate the effect of strongly interaction quarks, $`u,d,s,c,b`$ to the ‘running’ of $`\alpha (M)`$. The contribution of quarks which also interact strongly, is more precisely determined combining dispersion methods with experimental cross section results Bur01 , rather than through the use of Eq. (Variation of Fine Structure Constant from Non-Universal Gravity). The lepton contribution is $`\mathrm{\Delta }\alpha _l^1=4.3164`$ (including a very small higher order effect) while the hadron contribution is Jeg01  $`\mathrm{\Delta }\alpha _h^1=3.823\pm 0.054`$. In comparison, the reported cosmological variation is Web01  can be expressed as

$$\delta (1/\alpha )=\frac{1}{\alpha }\frac{\delta \alpha }{\alpha }0.001,$$

which is 50 times smaller than the uncertainty remaining in the understanding how the quark degrees of freedom influence the evolution of EM interaction strength, see Eq. (1).

From Eq. (Variation of Fine Structure Constant from Non-Universal Gravity) we obtain how $`\alpha (0)`$ responds to the variation of fermion masses, assuming that $`\alpha (M)`$ is constant:

$$\delta \left(\frac{1}{\alpha (0)}\right)=\frac{2}{3\pi }\underset{f}{}Q_f^2N_f\frac{\delta m_f}{m_f}$$
(3)
In consequence, a compound variation by $`\mathrm{0.72\hspace{0.17em}10}^5\alpha ^1\mathrm{\hspace{0.17em}3}\pi /2=0.46`$% in mass of (any of the) three leptons $`e,\mu ,\tau `$ and five quarks $`u,d,s,c,b`$ is sufficient to generate the observed effect. Restricting our attention to the 2nd and 3rd generation and distributing the effect uniformly, we need less than -0.1% variation in mass for each of the 5 fermions.

Within this chain of arguments we cannot expect that the variation of the masses of second and third families is found in the conventional interactions, as these effects would also alter the fundamental generation:

1. Masses, like coupling strengths, evolve with the interaction energy scale, driven by variation of the coupling constant. It is customary to write, within QED,

$$\frac{M}{m}\frac{dm}{dM}=\gamma (\alpha (M)).$$
(4)
Given a universal function $`\gamma (\alpha )`$ for leptons, it is easy to see that the variation of all lepton masses is exactly the same and thus the ratio of lepton masses is scale independent. For quark mass ratios, where the scale evolution of strong interactions is dominated by gluons, a similar argument is applicable in a very good approximation. Thus we cannot introduce a variation of mass of 2nd and 3rd generation only by manipulating the scale running of masses.
2. Considering the case that the masses of fermions are derived from the coupling to the Higgs field, the mass ratios are ratios of coupling strengths. It is less than obvious how a fraction of percent variability of mass can arise within the Higgs scheme solely for the 2nd and 3rd generation.

Gravity for 2nd and 3rd generation: To explain the observed effect, we would like to see a variation of some of the masses, while keeping the fundamental generation masses constant. Thus by necessity we need an interaction
1) which breaks the fermion generation universality, and
2) which is sensitive to the conditions prevailing in the late Universe.

We propose a novel and speculative solution that each fermion generation may have a different ratio of inertial to gravitational mass, as is expressed by gravitational coupling. We recall that the Eötvös type experiment establishes the ratio of gravitational and inertial mass solely for the “stable” first generation which surrounds us.

To be very specific we propose a natural and very slight generalization of the Einstein general relativity equations:

$$R_{\mu \nu }\frac{1}{2}g_{\mu \nu }R=\frac{\delta }{\delta g^{\mu \nu }}\left(G_1I_1+G_2I_2+G_3I_3\right)$$
(5)
Here, $`R`$ is the Ricci curvature tensor, $`G_i8\pi g_i/c^2`$ with $`g_1=G`$ the usual gravitational constant and $`G_2,G_3`$ are the gravitational constants of the 2nd and 3rd generation of particles, while $`I_i`$ is the action of the generation $`i`$ derived from the study of inertial forces. Variation with respect to the metric $`g^{\mu \nu }`$ does not yield a global inertial energy momentum tensor $`T_{\mu \nu }`$, rather we find for each generation its $`T_{\mu \nu }^i`$-component weighted with different strength of the non-universal gravitational interaction. Since each of the components, $`T_{\mu \nu }^i`$, has the same transformation properties, the right hand side of Eq. (5) transforms like the inertial energy-momentum tensor.

We have not included gauge fields in the above. A priori their universality poses a technical challenge in our proposed approach. What needs to be derived to make our approach correctly rooted in e.g. the brane theory is a brane configuration which maintains gauge universality, while yielding gravity non-universality with respect to the different generations Die02 . However, within the brane theory more common is to find a configuration leading to gauge non-universal, gravity universal theory Die02 .

The inertial matter action $`I_i`$ is as given by Weyl, Fock and Iwanenko in late 20’s:

$$I_i=d^nx\sqrt{g}[\overline{\mathrm{\Psi }}_i(\gamma ^\mu iD_\mu m_i)\mathrm{\Psi }_i].$$
(6)
Where $`\gamma ^\mu \gamma ^\nu +\gamma ^\nu \gamma ^\mu =2g^{\mu \nu }`$ defines the Dirac-$`\gamma `$ general relativistic matrices and $`D_\mu =/x^\mu +\mathrm{\Gamma }_\mu `$ with:

$$\mathrm{\Gamma }_\mu =\frac{1}{4}\gamma _\nu \left(\frac{\gamma ^\nu }{x^\mu }+C_{\lambda \mu }^\nu \gamma ^\lambda \right),$$

and the Christoffel symbol is:

$$C_{\lambda \mu }^\nu =\frac{1}{2}g^{\nu \beta }\left(\frac{g_{\lambda \beta }}{x^\mu }+\frac{g_{\beta \mu }}{x^\lambda }\frac{g_{\lambda \mu }}{x^\beta }\right).$$

Cosmological solutions of field equations arising from Eq. (5) have the usual structure. Presence of matter with different gravitational couplings (e.g. cosmic background neutrinos) does not alter the structure of the Einstein equations. In the dust model for the energy-momentum tensor the effective equations of state are different, for the case of mass less (neutrino) gas the Einstein equations including the 2nd and 3rd generation appear with effective degeneracy $`g_{2,3}^{\mathrm{eff}}=g_{2,3}G_{2,3}/G_1`$, there seems to be more gravitating matter in the Universe, just as required. Stellar objects containing 2nd and 3rd generation neutrinos would still be defined by the asymptotic Schwarzschild mass, we could not even tell which matter is within the star, short of locally probing the flavor and inertia of each component. Important for our later consideration of mass threshold is that the background Schwarzschild metric attaches to the mass of a particle as

$$m_im_i\sqrt{1\frac{2Mg_i}{r}}m_i\left(1\frac{Mg_i}{r}\right).$$
(7)
We thus see that the local gravitational potential attaches to the gravitational charge $`g_i`$ of particle ‘i’. Accordingly the mass defect is more significant for particles with greater coupling strength.

Consequences of non-universal gravity: There is lack of any direct experimental information on the gravity action on 2nd and 3rd family. A very weak upper limit of the strength of the interaction $`G_{3,2}`$ arises noting that massive particles in these families should not be beyond Planck mass limit, which requires $`G_{3,2}/G_1<10^{17}`$$`10^{19}`$. A stronger limit can be derived from loop contributions of virtual particles to physics on Earth, which begins to play a noticeable role when $`G_{3,2}/G_110^{10}`$ An experimental study of muonium $`\mu ^+e^{}`$ fall in the gravitational field of the Earth could reveal a better limit on $`G_2`$. During muons lifespan of 2$`\mu s`$ a particle falls 20 Å, thus a measurable effect requires $`G_2G_1`$. Short of approaching the Planck mass (see above) the enhanced gravitational interaction remains weak.

Turning to cosmological consequences we note that even if the inertial masses of neutrinos are at the (small) level of the mass differences seen in neutrino oscillations, $`m_i10^2`$ eV Bar02 , a 1000 fold enhancement of second and third neutrino generation gravitational coupling allows the two neutrinos in cosmological background to contribute a gravitational mass-energy equivalent of nearly 10 eVG each, and thus to describe the ‘dark’ gravitating mass of the Universe. The presence of more strongly gravitating component in the background of our Universe is likely to help the dynamical evolution of density fluctuations into the present day large inhomogeneities. The dynamics of the Universe expansion would need to be reexamined in this new context.

One should not before a thorough exploration of the family specific gravitational interaction see the above factor $`10^3`$ as an upper limit on gravitational interaction enhancement. Namely, the neutrino oscillation phenomenon combined with the large gravitational fields in early Universe can greatly alter this naive limit.

The required change of the fermion mass by -0.1–0.5% within the nebulae which we are looking for is understood to be the gravitational mass defect. Our Sun generates for normal gravity a mass defect $`\delta m_1/m_1=\mathrm{0.3\hspace{0.17em}10}^5`$. Thus coincidental with the explanation of the dark matter, the non-universal gravity with $`G_{3,2}/G_110^3`$ yields the expected -0.3% mass defect for the 2nd and 3rd generation fermions. Note that the sign of the effect is correct, which aside of checking all equations can be seen as follows: background gravity reduces the effective threshold for creation of a pair of particles, and this reduction in effective mass allows the fine structure constant to ‘run’ faster and thus from a fixed unification value toward a smaller value, as is reported experimentally to be present within the nebula.

Once we interpret the variation of $`\alpha `$ as a spatial effect, we could go back to check if it would not be possible that all fermions are subject to the mass defect. This would require that nebula appear to have a gravitational interactions 1000 stronger compared to our sun also for regular matter. This seems not to be possible.

Enhanced gravitational coupling has the capacity to resolve some riddles in stellar structure. For example, the accretion of background neutrinos and resulting concentration of more strongly gravitating neutrinos in large stellar objects Bil01  releases greatly enhanced amount of gravitational energy as e.g. is required to fuel the quasar phenomenon. Since there are two different additional gravitational couplings possible, and there is a finite supply of cosmic neutrinos it may be possible to explain both quasars and active galactic nuclei in terms of accretion of neutrinos, and to limit the range of cosmological time during which these phenomena occur. These remarks are consistent with appearance of event horizons (black holes) in such stellar objects, indeed the non-universal gravity of some of the accreted (neutrino) particles would lead to faster formation of such a gravitational singularity. Naturally, this singularity is the same for all test particles of any of the three generations.

Final remarks: We conclude that the small deviation of $`\alpha (0)`$ from the established values could be result of a somewhat larger (by factor $`645=1373\pi /2`$) -0.46% compound total mass variation of unstable second and third generation fermions. It seems that a spatial rather than temporal variation of fermi masses, due to a gravitational mass defect, is the source of this effect. The nebula which act as source of the absorption lines are required to generate gravitational fields causing this mass defects. An enhanced anomalous gravitational coupling of heavy fermi generations is our interpretational scheme. A factor 1000 enhancement in non-universal gravity allows interpretation of dark matter as being due to $`\nu _\mu ,\nu _\tau `$ cosmic background, while the gravitational mass defect of stable matter in absorbing nebula is at $`10^5`$, just a factor 5 greater than the mass defect of the sun.

Acknowledgments: I thank Keith Dienes for stimulating discussions and valuable suggestions, and the CERN theory division for hospitality while I was discovering this subject. Work supported in part by a grant from the U.S. Department of Energy, DE-FG03-95ER40937.


format: text

text:
# The origin of the rebrightening in soft X-ray transient outbursts

## 1 Introduction

The soft X-ray transients (SXTs) are close interacting binary systems consisting of a black hole or neutron star primary and a low-mass Roche-lobe-filling secondary star. The outbursts of soft X-ray transients are quite distinct from those seen in other binary systems such as dwarf novae. The duration and recurrence times are both strikingly longer, suggesting that a large fraction of the mass of the accretion disc is consumed in each outburst and that a different mechanism from the standard disc instability picture is required (comprehensive reviews of the disc instability picture are given by Cannizzo, 1993 and Lasota, 2001).

Typical ionization instability-driven dwarf nova outbursts last a few days and consume a few percent of the disc mass, leading to a recurrence time of the order of weeks. In comparison, a typical SXT outburst lasts for 100 days or more and recurs on a time-scale of several years. The profiles of the X-ray light curves of SXT outbursts vary greatly but they do share common characteristics - the initial peak is often followed by an exponential decline and, in many cases, a rebrightening several tens of days after the first maximum. Smaller rebrightenings later on in the decay from outburst have also been observed. These features can be seen quite clearly in the X-ray light curve of the well-known transient A0620-003 (Elvis et al. 1975, Whelan et al. 1977). A general review of SXT outbursts is given by Tanaka & Lewin (1996).

Several authors have proposed models to explain the form of the SXT light curves. The long recurrence time of the outbursts is a problem which has been addressed in many ways. It can be explained by the standard ionization instability model if the Shakura-Sunyaev viscosity parameter $`\alpha `$ is greatly reduced from that typical of a cataclysmic variable disc. Such an explanation is not confined to SXTs; it has been invoked to explain the outbursts of WZ Sagittae (Smak, 1993). However, there is no obvious physical reason why $`\alpha `$ should vary so much from one accretion disc to another. The more successful models involve the removal of the inner part of the disc where the viscous time-scale is short. This can been achieved by a hot, advection-dominated accretion flow (ADAF) onto the hole : Meyer & Meyer-Hofmeister (1994) showed that gas from the cool outer accretion disc in a binary system may evaporate into a hot corona and subsequently feed an ADAF onto the comapct object. Dubus et al. (2001) have performed one-dimensional simulations of an SXT including the effects of both irradiation and an ADAF. In this work we do not address the question of recurrence time, but concentrate on the observed outburst profiles. Cannizzo, Chen & Livio (1995) found that the exponential decay could be reproduced by parameterising $`\alpha `$ in terms of the aspect ratio $`\frac{\mathrm{H}}{\mathrm{R}}`$. Augusteijn, Kuulkers & Shaham (1993) sought to explain the rebrightenings by considering irradiation of the secondary as a mechanism for increasing the mass transfer rate. However, it is unlikely that the radiation can penetrate the photosphere of the companion to produce a variation in mass transfer rate of more than a few percent (King, 1989).

We consider the effect of disc irradiation by the X-ray emission from the compact object. King & Ritter (1998) showed analytically that the assumption that the central X-ray source can efficiently irradiate the disc leads automatically to long exponential X-ray decays with the observed time-scales. Observations of SXTs suggest that the optical/X-ray flux ratio for these objects is simply too large for the optical emission to result from intrinsic viscous dissipation in the disc alone (van Paradijs & McClintock, 1994; Shahbaz & Kuulkers, 1998). Moreover, the optical brightness is strongly correlated with the irradiating X-ray flux. However, detailed one-dimensional modelling of the disc structure (e.g. Cannizzo 1994) suggests that the response of the disc to the irradiation is for the surface to become convex rather than concave, thus shielding the outer regions from the irradiating flux. The observations seem to contradict this picture, and for the purposes of this paper we assume that the disc surface remains concave.

In section 2 we describe the irradiated disc model of King and Ritter in more detail, and outline the importance of tidal forces in systems with extreme mass ratios in section 3. The numerical method used to implement the model is then described in section 4. In section 5, we show that the observed peak luminosity of the 1975 outburst of A0620-003 almost certainly implies that the irradiation radius (the outer radius of the region kept heated by the central X-rays) was larger than the accretion disc. We then show in simulations that no rebrightening occurs if the whole disc is kept in the hot viscous state by the irradiation; however, if only a very small area of the outer disc is shielded in some way (by a radiation-induced warp or self-shadowing, for example), the observed rebrightening is produced.

## 2 Disc Irradiation Model

Irradiation of an accretion disc by the central X-ray source heats the exposed disc faces. In SXTs, provided the disc faces remain concave, the surface temperature of the disc is given by the formula of de Jong, van Paradijs & Augusteijn (1996) :

$$\mathrm{T}_{\mathrm{irr}}^4(\mathrm{R})=\frac{\eta \dot{\mathrm{M}}\mathrm{c}^2(1\beta )}{4\pi \sigma \mathrm{R}^2}\left(\frac{\mathrm{H}}{\mathrm{R}}\right)^\mathrm{n}\left[\frac{\mathrm{d}\mathrm{lnH}}{\mathrm{d}\mathrm{lnR}}1\right]$$
(1)
where $`\eta `$ is the accretion efficiency for producing X-rays, Ṁ is the mass accretion rate onto the primary, H is the scale-height of the disc and $`\beta `$ is the disc albedo. The index n takes the values 1 or 2 for irradiation by a point source or by the inner disc respectively. Accretion onto a neutron star will always have n = 1, while n = 2 is possible in the absence of a hard surface on which to accrete, i.e. a black hole primary. During outburst, however, the disc will probably develop an X-ray corona, as suggested by observations of a strong power-law continuum. In this case, the appropriate index reverts to n = 1. In a large enough disc, $`\mathrm{T}_{\mathrm{irr}}`$ is quite sufficient to dominate the effective temperature and exceed the hydrogen ionization temperature ($`\mathrm{T}_\mathrm{H}6500\mathrm{K}`$), thus stabilising the disc on the hot branch of the $`\mathrm{T}_{\mathrm{eff}}\mathrm{\Sigma }`$ relation (see for example, Cannizzo, 1993). $`\mathrm{T}_{\mathrm{irr}}`$ is the key factor which determines whether a low mass X-ray binary is transient or not - in a persistent source $`\mathrm{T}_{\mathrm{irr}}`$ exceeds $`\mathrm{T}_\mathrm{H}`$ throughout the entire disc (van Paradijs, 1996).

We consider a situation in which the disc undergoes an outburst which is initiated by the standard hydrogen ionization disc instability. Subsequently, as the mass accretion rate onto the primary increases, the disc is irradiated by the enhanced X-ray emission. The disc is unable to return to a cool quiescent state until the mass accretion rate has fallen sufficiently to switch off the enhanced irradiation, allowing a cooling front to sweep through the disc. Following the analytic treatment of King & Ritter (1998), one can use equation 1 to determine the variation of $`\mathrm{R}_\mathrm{h}`$ with $`\dot{\mathrm{M}}_\mathrm{c}`$ by noting that the edge of the irradiated region is determined by the point in the disc at which $`\mathrm{T}_{\mathrm{irr}}=\mathrm{T}_\mathrm{H}`$. Then :

$$\mathrm{R}_\mathrm{h}^2=\frac{\eta \dot{\mathrm{M}}\mathrm{c}^2(1\beta )}{4\pi \sigma \mathrm{T}_{}^{4}{}_{\mathrm{H}}{}^{}}\left(\frac{\mathrm{H}}{\mathrm{R}}\right)^\mathrm{n}\left[\frac{\mathrm{d}\mathrm{lnH}}{\mathrm{d}\mathrm{lnR}}1\right]=\mathrm{B}_\mathrm{n}\dot{\mathrm{M}}$$
(2)
Using typical parameters for neutron star and black hole systems, King & Ritter showed that the mass accretion rates for a given $`\mathrm{R}_\mathrm{h}`$ are then

$$\dot{\mathrm{M}}_{1,\mathrm{crit}}=4.1\times 10^{10}\mathrm{R}_{11}^{}{}_{}{}^{2}\mathrm{M}_{}\mathrm{yr}^1,$$
(3)
and

$$\dot{\mathrm{M}}_{2,\mathrm{crit}}=2.9\times 10^9\mathrm{R}_{11}^{}{}_{}{}^{2}\mathrm{M}_{}\mathrm{yr}^1,$$
(4)
where $`\mathrm{R}_{11}`$ is $`\mathrm{R}_\mathrm{h}`$ in units of $`10^{11}\mathrm{cm}`$. King and Ritter and subsequently King (1998) extend this analysis to show that the first part of the outburst decay, with the whole disc irradiated, will be exponential, but ultimately the final part of the decay will be roughly linear due to the gradual shrinking of the irradiated region. However, they could not fully explain the rebrightening, and its origin remained unknown. A schematic of the theoretical profile is given in fig. 1.

The recurrent transient A0620-003 (V616 Mon) is perhaps the best known of the SXTs. Its 1975 outburst was observed by Ariel V (Elvis et al. 1975), during which second and third maxima were observed 55 and 200 days after the first. The decay either side of the first rebrightening is exponential and the rise-time of the first rebrightening was $`1015`$ days. A0620-003 has a short orbital period of 7.75 hours (McClintock et al., 1983), and is likely to contain a black hole primary. The primary mass has been estimated to be $`7.3\mathrm{M}_{}`$ or more (McClintock & Remillard, 1986), but a more recent work concludes that it is likely to be in the range $`3.30<\mathrm{M}_1<4.24\mathrm{M}_{}`$ (Marsh, Robinson & Wood 1994). The reader is referred to Whelan et al. (1977) or Tanaka & Lewin (1995) for the full X-ray light curve of the 1975 outburst.

The analysis of King & Ritter correlates well with the observed behaviour of A0620-003 and SXTs in general, but the precise origin of the rebrightening and the behaviour of the disc remain to be explained. In this work we incorporate King & Ritter’s analytic model into an existing numerical code and perform two-dimensional simulations of a SXT with the orbital parameters of A0620-003. We find a natural explanation for the rebrightening in terms of the physical state of the disc. In a forthcoming paper we will report results from simulations of outbursts in longer period SXTs (with correspondingly larger discs) and investigate the resulting complex outburst morphologies.

## 3 Tides

The mass ratio of A0620-003 is very well-constrained : $`\mathrm{q}=\frac{\mathrm{M}_2}{\mathrm{M}_1}=0.067\pm 0.01`$ (Marsh. Robinson & Wood, 1994). Such a low value suggests that tidal resonances will be accessible to the disc. Tidal forces will therefore be an important factor in the dynamics of the accretion disc. The radius of the 3:1 tidal resonance, which has been shown to be important in dwarf novae, is

$$\frac{\mathrm{R}_{3:1}}{\mathrm{a}}=\left(\frac{1}{3}\right)^{\frac{2}{3}}(1+\mathrm{q})^{\frac{1}{3}}$$
(5)
where a is the binary separation (see, for example, Frank, King & Raine, 2002). For $`\mathrm{q}=0.07`$, this radius is $`\mathrm{R}_{3:1}=0.47\mathrm{a}`$. This is well inside the radius at which the disc is ultimately truncated by the tides. $`\mathrm{R}_{\mathrm{tide}}`$ can be estimated from the Roche lobe radius of the primary - it is always roughly 90% of this value. Using the formula of Eggleton (1983),

$$\frac{\mathrm{R}_\mathrm{L}(1)}{\mathrm{a}}=\frac{0.49\mathrm{q}^{\frac{2}{3}}}{0.6\mathrm{q}^{\frac{2}{3}}+\mathrm{log}_\mathrm{e}(1+\mathrm{q}^{\frac{1}{3}})},$$
(6)
$`\mathrm{R}_{\mathrm{tide}}0.60\mathrm{a}`$ for q = 0.07. Therefore we expect tidal effects to become important and the disc to become eccentric when it reaches $`\mathrm{R}_{3:1}`$.

## 4 Numerical Method

We use a two dimensional smoothed particle hydrodynamics (SPH) accretion disc code. The great advantage of this method is that it incorporates the full tidal potential of the binary system without need for an approximate prescription for the removal of angular momentum at the disc edge employed in 1-D simulations. The importance of using a two- or three-dimensional scheme with the tidal effects included has been demonstrated effectively with this code in recent simulations of cataclysmic variables (Truss et al. 2001). Full details of the code are described in Truss et al. (2000) and a complete review of the SPH technique is given by Monaghan (1992). The implementation here is very closely related to the dwarf nova study of Truss et al. (2000,2001), in which we incorporate the disc instability in terms of a viscosity switch based on surface density and build up the disc by a mass transfer stream. The disadvantage of the scheme is that there is no self-consistent calculation of vertically-averaged thermal structure, as in one-dimensional codes. There is a trade-off in computational efficiency between calculating the thermodynamics in a self-consistent manner in 1-D and performing the hydrodynamics in 2-D. We choose to present a 2-D calculation in this work to examine the role of tidal forces, and in this way the scheme is complementary to modern 1-D treatments.

### 4.1 Viscosity

When the local surface density at the position of any particle in the disc exceeds the critical value for ionization of hydrogen to occur as calculated by Cannizzo, Shafter & Wheeler (1988),

$$\mathrm{\Sigma }_{\mathrm{max}}=11.4\mathrm{R}_{10}^{1.05}\mathrm{M}_1^{0.35}\alpha _{\mathrm{cold}}^{0.86}\mathrm{gcm}^2,$$
(7)
the viscosity acting on that particle is increased. The change in viscosity is made on an appropriate thermal time-scale (much shorter than the local viscous time-scale) from a value $`\alpha _{\mathrm{cold}}`$ to a value $`\alpha _{\mathrm{hot}}`$, where $`\alpha `$ is defined in terms of the sound speed of the fluid, $`\mathrm{c}_\mathrm{s}`$, and the scale height of the disc, H, in the usual viscosity formulation of Shakura & Sunyaev (1973):

$$\nu =\frac{2}{3}\alpha \mathrm{c}_\mathrm{s}\mathrm{H}.$$
(8)
The local viscosity remains at $`\alpha _{\mathrm{hot}}`$ until the local surface density decreases below a critical value

$$\mathrm{\Sigma }_{\mathrm{min}}=8.25\mathrm{R}_{10}^{1.05}\mathrm{M}_1^{0.35}\alpha _{\mathrm{hot}}^{0.8}\mathrm{gcm}^2.$$
(9)
when it is returned to the low state value $`\alpha _{\mathrm{cold}}`$.

Equation 8 is employed almost directly in the SPH formulism. The SPH momentum equation for a particle i is:

$$\frac{\mathrm{d}𝐯_\mathrm{i}}{\mathrm{dt}}=\underset{\mathrm{j}}{}\mathrm{m}_\mathrm{j}\left(\frac{\mathrm{P}_\mathrm{i}}{\rho _\mathrm{i}^2}+\frac{\mathrm{P}_\mathrm{j}}{\rho _\mathrm{j}^2}+\frac{\beta \mu _{\mathrm{ij}}^2\zeta \overline{\mathrm{c}}_{\mathrm{ij}}\mu _{\mathrm{ij}}}{\overline{\rho }_{\mathrm{ij}}}\right)_\mathrm{i}\mathrm{W}_{\mathrm{ij}}$$
(10)
where $`\overline{\rho }_{\mathrm{ij}}`$ and $`\overline{\mathrm{c}}_{\mathrm{ij}}`$ are the average density and sound speed of the pair of particles i and j respectively. h is the smoothing length and $`\mathrm{W}_{\mathrm{ij}}`$ is the interpolating kernel. $`\beta `$ and $`\zeta `$ are the quadratic and linear artificial viscosity parameters and

$$\mu _{\mathrm{ij}}=\{\begin{array}{cc}\frac{𝐯_{\mathrm{𝐢𝐣}}𝐫_{\mathrm{𝐢𝐣}}}{𝐫_{\mathrm{𝐢𝐣}}^{}{}_{}{}^{2}+\eta ^2}\hfill & 𝐯_{\mathrm{ij}}𝐫_{\mathrm{ij}}0\hfill \\ 0\hfill & 𝐯_{\mathrm{ij}}𝐫_{\mathrm{ij}}>0\hfill \end{array}$$
(11)
where $`\eta ^2=0.01\mathrm{h}^2`$ avoids singularities.

Following Murray (1996), we set $`\beta =0`$ and in the continuum limit (i.e. $`\mathrm{h}0`$ and number of particles $`\mathrm{N}\mathrm{}`$) the acceleration due to the artificial viscosity term can be written

$$𝐚_\mathrm{v}=\frac{\zeta \mathrm{H}\kappa }{2\mathrm{\Sigma }}((\mathrm{c}_\mathrm{s}\mathrm{\Sigma }\mathrm{S})+(\mathrm{c}_\mathrm{s}\mathrm{\Sigma }𝐯)),$$
(12)
where $`\kappa `$ is a constant determined by the kernel W. In this case we use a 2D cubic spline kernel, for which $`\kappa =\frac{1}{4}`$. $`\mathrm{S}_{\mathrm{ij}}`$ is the deformation tensor, where

$$\mathrm{S}_{\mathrm{ij}}=\frac{\mathrm{v}_\mathrm{i}}{\mathrm{x}_\mathrm{j}}+\frac{\mathrm{v}_\mathrm{j}}{\mathrm{x}_\mathrm{i}}.$$
(13)
If the density and sound speed vary on lengthscales much longer than the velocity, and if we assume that the fluid is divergence-free, then

$$𝐚_\mathrm{v}=\frac{\zeta \mathrm{Hc}_\mathrm{s}}{8}^2𝐯.$$
(14)
So, in the continuum limit, the artificial viscosity introduces a shear viscosity

$$\nu =\frac{1}{8}\zeta \mathrm{c}_\mathrm{s}\mathrm{H}$$
(15)
In fact, Murray (1996) used the SPH smoothing length, h, as the lengthscale in his prescription for viscosity. We do not impose such a constraint here : the viscous energy is dissipated over the disc scale height, H.

In the simulations, we take the sound speed to be constant everywhere and vary $`\zeta `$ to change the viscosity of a particle in the way described above. With such an isothermal scheme, the net viscosity in the simulations has a slightly different radial dependence to an analytic alpha-disc. In the alpha model, the variation of sound speed and scale height with radius ($`\mathrm{c}_\mathrm{s}\mathrm{R}^{\frac{3}{8}}`$ and $`\mathrm{H}\mathrm{R}^{\frac{9}{8}}`$) lead to

$$\nu _\alpha =1.8\times 10^{14}\alpha ^{\frac{4}{5}}\dot{\mathrm{M}}_{16}^{\frac{3}{10}}\mathrm{M}_1^{\frac{1}{4}}\mathrm{R}_{10}^{\frac{3}{4}}\mathrm{cm}^2\mathrm{s}^1.$$
(16)
where $`\dot{\mathrm{M}}_{16}`$ is the mass accretion rate in units of $`10^{16}\mathrm{gs}^1`$. In the SPH simualtions, the sound speed is a constant and the scale height is a linear function of radius:

$$\nu _{\mathrm{sph}}=6.4\times 10^{14}\zeta \left(\frac{\mathrm{c}_\mathrm{s}}{5\times 10^6\mathrm{cms}^1}\right)^2\mathrm{R}_{10}\mathrm{cm}^2\mathrm{s}^1.$$
(17)
where we usually use $`\zeta _{\mathrm{cold}}=0.1`$ and $`\zeta _{\mathrm{hot}}=1.0`$. A comparison of these two viscosity variations is shown in fig. 2. A full discussion of how the results of this scheme can be interpreted in relation to observations, and its limitations, is given in Section 6.

### 4.2 Radiation

The irradiation radius is calculated by direct implementation of equation 3 at each timestep . All particles within this radius are instantaneously given the high-state SPH viscosity $`\zeta _{\mathrm{hot}}`$. Since we know the number of particles being accreted onto the primary in the simulation, we can readily calculate the mass accretion rate at each timestep, and hence the irradiation radius. We are not concerned with calculating accurate optical light curves in these simulations; instead we concentrate on the total X-ray emission, which we take to be proportional to the mass accretion rate onto the primary. This scheme is ideal for investigating how variations in the X-rays arise from waves of mass travelling through the disc.

Our model is distinct from those employed by Huang & Wheeler (1989) and Mineshige & Wheeler (1989) in that we follow the evolution of the whole disc in two dimensions and there is no need to approximate the binary potential. In particular, for numerical reasons, Mineshige & Wheeler used a very small disc in their simulation of A0620-003, using an outer radius of 3.16 $`\times 10^{10}\mathrm{cm}`$, well inside even the circularisation radius. We simulate the entire disc and show that the outermost regions can play a pivotal rôle in determining the form of the observed X-ray light curve.

## 5 Results

### 5.1 Irradiation efficiency

Given a distance of 870 pc (Oke, 1977), the peak X-ray luminosity of the 1975 outburst of A0620-003 was $`\mathrm{L}_\mathrm{x}1.3\mathrm{x10}^{38}\mathrm{erg}.\mathrm{s}^1`$. The total accretion luminosity can be written

$$\mathrm{L}_\mathrm{x}=\eta \dot{\mathrm{M}}\mathrm{c}^2,$$
(18)
where $`\eta `$ is the rest mass to energy accretion efficiency. We can use this equation with King & Ritter’s analytic expressions for $`\mathrm{R}_\mathrm{h}`$ to estimate the maximum value of $`\mathrm{R}_\mathrm{h}`$ reached at the peak of the 1975 outburst. Equations (3) and (4) were derived using specific values of $`\eta `$; (3) assumes that the efficiency for accretion onto a neutron star is $`\eta =0.15`$, while for the black hole case (4) assumes $`\eta =0.10`$. If we relax this constraint and rewrite the equations in terms of a general $`\eta `$ consistent with (18), we have :

$$\dot{\mathrm{M}}_{\mathrm{crit}}=\frac{6.2\mathrm{x10}^{11}}{\eta }\mathrm{R}_{11}^2\mathrm{M}_{}\mathrm{yr}^1(\mathrm{n}=1),$$
(19)
and

$$\dot{\mathrm{M}}_{\mathrm{crit}}=\frac{2.9\mathrm{x10}^{10}}{\eta }\mathrm{R}_{11}^2\mathrm{M}_{}\mathrm{yr}^1(\mathrm{n}=2).$$
(20)
Now, if we infer the mass accretion rate from the observed luminosity in (18), we see from (19) and (20) that the estimates of $`\mathrm{R}_\mathrm{h}`$ are independent of $`\eta `$. For n=1 we infer that $`\mathrm{R}_\mathrm{h}=6.1\times 10^{11}\mathrm{cm}`$ and for n=2, $`\mathrm{R}_\mathrm{h}=2.8\times 10^{11}\mathrm{cm}`$. For a total system mass of $`\mathrm{M}_1+\mathrm{M}_2=4\mathrm{M}_{}`$ and orbital period 7.75 hours, the binary separation of the system is $`2.2\times 10^{11}\mathrm{cm}`$. Therefore, the maximum value reached by $`\mathrm{R}_\mathrm{h}`$ during an outburst of A0620-003 is outside the disc, and quite possibly much larger than the disc. The parameter space for varying masses is shown in figure 3. Given a mass ratio, q, the disc size can be estimated using the formula for the volume radius of the Roche lobe of the primary, equation 6. This is shown in the lower curve of fig. 3, for q = 0.07.

### 5.2 Outbursts

We perform separate simulations of A0620-003, allowing the irradiation radius to take different values for the same mass accretion rate. In this way we are effectively varying the parameter $`\mathrm{B}_\mathrm{n}`$ in equation 2. We aim to discover the effects of allowing $`\mathrm{R}_\mathrm{h}`$ to become larger than the accretion disc or not. In these simulations we take $`\mathrm{P}_{\mathrm{orb}}=7.75\mathrm{h}`$, $`\mathrm{M}_1=3.7\mathrm{M}_{}`$ and $`\mathrm{M}_2=0.25\mathrm{M}_{}`$, consistent with the measurements of Marsh, Robinson & Wood (1994) and $`\dot{\mathrm{M}}_2=3\mathrm{x}10^{11}\mathrm{M}_{}\mathrm{yr}^1`$ after McClintock et al. (1983) and Cannizzo (1998). The sound speed is $`\mathrm{c}_\mathrm{s}=0.05\mathrm{a}\mathrm{\Omega }_\mathrm{b}=24.5\mathrm{kms}^1`$ and the high and low-state SPH viscosities are $`\zeta _{\mathrm{hot}}=1`$ and $`\zeta _{\mathrm{cold}}`$ = 0.1.

Figure 4 shows the variations in mass accretion rate onto the primary and heating radius $`\mathrm{R}_\mathrm{h}`$ for the simulation of A0620-003. The entire disc is irradiated and there is no large rebrightening. There is an exponential decay as the irradiated gas is accreted followed by a faster decline as the accretion rate is controlled by $`\mathrm{R}_\mathrm{h}`$ gradually retreating inside the disc. There is a small variation at 70 days after outburst maximum. This arises from the hot / cold boundary passing through structure near the disc edge associated with the spiral arms. The effect is too small to be responsible for the larger rebrightening in A0620-003, but may be responsible for the later, smaller glitches observed in the light curve.

In the second simulation, shown in fig. 5, the efficiency, $`\eta `$, is reduced such that $`\mathrm{R}_\mathrm{h}`$ is limited to a maximum value $`\mathrm{R}_{\mathrm{h},\mathrm{max}}=0.35\mathrm{R}_{\mathrm{disc}}`$. This allows us to mimic the effect of a self-shadowed outer disc. In this case, when the mass accretion rate starts to rise, only the inner part of the disc is irradiated. The accretion of gas in the irradiated region gives rise to the first peak in fig. 5, which decays exponentially on a time-scale appropriate to the inner part of the disc. Meanwhile, the disc encounters the 3:1 tidal resonance and a two-armed spiral wave is launched in the disc. Gas in the outer, unirradiated regions is accumulated in the spiral arms and the surface density rises. Eventually, $`\mathrm{\Sigma }`$ exceeds the critical limit for the viscosity to be increased. When this gas subsequently arrives at the primary, there is a rebrightening in the X-ray light curve. The second stage of the decay takes place on a longer time-scale, and the slope is correspondingly more gradual. The final, rapid drop to quiescence only takes place when $`\mathrm{R}_\mathrm{h}<0.35\mathrm{R}_{\mathrm{disc}}`$.

The two distinct regions, one irradiated and one unirradiated but subject to the disc instability, give rise to the rebrightening. The light curve looks more like that of A0620-003, but differs in one important detail - the slope of the decay either side of the rebrightening is not the same. The observed similarity in the gradients of the decay either side of the rebrightening in A0620-003 suggests that the waves of mass which are producing these parts of the outburst originate in the same region of the disc.

It is clear from the discussion in section 4.1 that $`\mathrm{R}_\mathrm{h}`$ becomes much larger than the outer edge of the disc during an outburst of A0620-003. However, if a small part of the outer disc is shadowed from the radiation in some way, the rebrightening and the decay slopes could be reproduced even when $`\mathrm{R}_\mathrm{h}>>\mathrm{R}_{\mathrm{disc}}`$. In order to demonstrate this, a simulation was performed in which the efficiency limits $`\mathrm{R}_\mathrm{h}`$ to a value near, but not quite as large as $`\mathrm{R}_{\mathrm{disc}}`$ (Run A). It is not our intention to suggest that the accretion efficiency is extremely low, but to perform a simulation that mimics the shadowing of a small part of the outer disc during outburst. The resultant light curve is shown in fig. 6, and this time the rebrightening and slopes are closer to the observations. Again, the final drop to quiescence only occurs when $`\mathrm{R}_\mathrm{h}<0.85\mathrm{R}_{\mathrm{disc}}`$.

### 5.3 Varying the simulation parameters

In order to assess the robustness of these results we performed run A using different SPH viscosities and different particle numbers. A summary of the parameters in each run is given in table 1.

Firstly, the simulation was repeated with the same viscosity as before, $`\zeta _{\mathrm{cold}}=0.1,\zeta _{\mathrm{hot}}=1.0`$, but with fewer SPH particles (runs B and C). Figure 7 shows the results with roughly one half and one quarter the number of particles of the orginal calculation. The greatest difference as particle number is decreased is the resolution of the primary maximum. With fewer particles, the spatial resolution in the inner part of the disc where the density is lower suffers to a greater extent than the outer part. Consequently, the accretion of particles from the inner disc leads to an under-resolved maximum in the accretion rate. Therefore, we expect the primary maximum to be even better resolved for a higher particle number, and the convergence between runs A and B is encouraging. The particle masses are scaled so that the total mass of the disc is the same in each simulation, so in run C the large drop in mass transfer rate before the rebrightening is caused by the resolution being very poor indeed - the accretion of just a few particles corresponds to a large response in accretion rate. The accuracy of this run should be judged accordingly.

In runs D and E we keep the particle number the same as in run A, but vary the magnitude of $`\zeta _{hot}`$. The resulting light curves are shown in fig. 8. In both runs, the rebrightening starts at the same time as it does in run A - between 15 and 20 days. However, the rise-time of the rebrightening is seen to increse as $`\zeta _{\mathrm{hot}}`$ is decreased (fig. 9). This suggests very strongly that the the rise-time of the rebrightening is governed by the viscosity, but the time delay from the start of the outburst to the rebrightening is independent of the viscosity.

## 6 Discussion

In this section we discuss the applicability and limitations of these results to interpreting observations of A0620-003 and suggest a possible mechanism for the rebrightening.

### 6.1 Thermal structure

We have already seen in section 4.1 that the SPH viscosity used in this isothermal scheme has a slightly different functional form to that of a standard $`\alpha `$-disc. Fig. 2 shows that in general the viscosity is slightly overestimated at small radii and slightly underestimated at larger radii. It should also be kept in mind that in the simulations, the irradiated region is indistinguishable from an unirradiated region in the high-viscosity state. The question thus arises as to how good an approximation this is to a real irradiated disc, i.e. does the radiation alter the thermal structure? The magnitude of the optical depth, $`\tau `$, and the origin of opacity in the disc is the key issue. $`\mathrm{T}_{\mathrm{irr}}`$ can be quite sufficient to dominate the effective temperature of the disc

$$\mathrm{T}_{\mathrm{eff}}=\left(\frac{3\mathrm{G}\mathrm{M}_1\dot{\mathrm{M}}(\mathrm{R})}{8\pi \sigma \mathrm{R}^3}\right)^{\frac{1}{4}},$$
(21)
but the vertical thermal structure of an irradiated disc will only be altered if $`\mathrm{T}_{\mathrm{irr}}`$ also exceeds the mid-plane temperature of the disc, given by $`\mathrm{T}_{\mathrm{m}}^{}{}_{}{}^{4}\tau \mathrm{T}_{\mathrm{eff}}^{}{}_{}{}^{4}`$. King (1998) has argued that $`\tau `$ is of order 10 and consequently the vertical structure of the outer disc is nearly isothermal, with the sound speed dependent on $`\mathrm{T}_{\mathrm{irr}}`$ and not $`\mathrm{T}_\mathrm{m}`$. The viscosity is therefore independent of the local value of $`\mathrm{\Sigma }`$ and the diffusion equation

$$\frac{\mathrm{\Sigma }}{\mathrm{t}}=\frac{3}{\mathrm{R}}\frac{}{\mathrm{R}}\left\{\mathrm{R}^{1/2}\frac{}{\mathrm{R}}\left[\nu \mathrm{\Sigma }\mathrm{R}^{1/2}\right]\right\},$$
(22)
leads to

$$\frac{\mathrm{\Sigma }}{\mathrm{t}}\mathrm{\Sigma }$$
(23)
and the observed exponential decays.

However, simulations by Dubus et al. (1999), in which the vertical structure is calculated self-consistently including the effects of irradiation by a point source in the orbital plane, suggest that $`\tau `$ is of the order $`10^210^3`$ and the structure is close to that of an $`\alpha `$-disc. Indeed, for an $`\alpha `$-disc with a Kramer’s opacity law,

$$\tau =1200\left(\frac{\alpha }{0.1}\right)^{\frac{4}{5}}\dot{\mathrm{M}}_{16}^{\frac{1}{5}}$$
(24)
where $`\dot{\mathrm{M}}_{16}`$ is the mass transfer rate in units of $`10^{16}\mathrm{gs}^1`$. In the more familiar case of an unirradiated accretion disc, the equation for vertical radiative energy transport is:

$$\frac{4\sigma \mathrm{T}_\mathrm{m}^4}{3\kappa \mathrm{\Sigma }}=\frac{9}{8}\nu \mathrm{\Sigma }\frac{\mathrm{GM}}{\mathrm{R}^3}$$
(25)
and we can see straight away that although the opacity $`\kappa `$ varies only slowly with $`\mathrm{\Sigma }`$ and $`\mathrm{T}_\mathrm{m}`$ in much of the disc, the viscosity $`\nu `$ is not independent of $`\mathrm{\Sigma }`$. With the alpha prescription,

$$\frac{\mathrm{\Sigma }}{\mathrm{t}}\mathrm{\Sigma }^{\frac{5}{3}}$$
(26)
gives a power-law decay.

Clearly, the magnitude of $`\tau `$ is crucial to determining any change in the vertical thermal structure of the disc. It remains unclear how the energy is deposited in the vertical direction in an accretion disc, and whether the $`\alpha `$-model is a valid description of the real processes involved. If the $`\alpha `$-model breaks down we have no way of calculating $`\tau `$. However, it remains entirely possible that $`\mathrm{T}_\mathrm{m}>\mathrm{T}_{\mathrm{irr}}>\mathrm{T}_\mathrm{H}`$ and the irradiation ionizes the disc, bringing about a change in viscosity, but does not alter the vertical thermal structure. The assumption that we make in the code is that the viscosity increases when $`\mathrm{T}_{\mathrm{irr}}>\mathrm{T}_\mathrm{H}`$ and there is no change in $`\mathrm{T}_\mathrm{m}`$ (i.e. the sound speed remains constant).

We do not claim that our approach is a more accurate representation than a 1D model that is thermodynamically self-consistent within the framework of the $`\alpha `$-model. Rather, we use the advantages of our scheme, particularly the two-dimensional nature of the code, to show that tidal forces (which are calculated self-consistently) can have a significant effect on the accretion flow in SXT discs and can account for the rebrightening observed in the X-ray light curves.

### 6.2 Rebrightening time-scales

The rebrightening in all the simulations begins 15-20 days after the onset of the outburst, while the peak of the second maximum occurs at different times according to the viscosity. This suggests that the rise-time of each rebrightening is controlled by viscous effects, but the time at which the rebrightening occurs is not. In each case, the accreted gas that is responsible for the rebrightening in the simulation originates from the outer (initially unirradiated) part of the disc.

The viscous time-scale for the outer part of a Shakura-Sunyaev $`\alpha `$-disc is longer than 50 days. The viscous time-scale is

$$\mathrm{t}_\mathrm{v}\frac{\mathrm{R}^2}{\alpha \mathrm{c}_\mathrm{s}\mathrm{H}},$$
(27)
while the $`\alpha `$-disc solution with a Kramer’s opacity law gives

$$\mathrm{t}_\mathrm{v}390\left(\frac{\alpha }{0.1}\right)^{\frac{4}{5}}\dot{\mathrm{M}}_{16}^{\frac{3}{10}}\mathrm{M}_1^{\frac{1}{4}}\mathrm{R}_{11}^{\frac{5}{4}}\mathrm{d}.$$
(28)
where $`\dot{\mathrm{M}}_{16}`$ is the mass accretion rate in units of $`10^{16}\mathrm{gs}^1`$, $`\mathrm{M}_1`$ is the primary mass in solar masses and $`\mathrm{R}_{11}`$ is the radius in units of $`10^{11}\mathrm{cm}`$. This is plotted in fig. 10 for the parameters of A0620-003, along with the SPH viscous time-scale.

Clearly, for the parameters of A0620-003 in outburst, the viscous time in the outer part of an $`\alpha `$-disc is longer than 50 days. We have seen in section 3 that for such an extreme mass ratio ($`\mathrm{q}=0.07`$), the 3:1 tidal resonance is easily accessible to the disc, and it is here that tidal effects start to affect the outburst. The tidal instability gives rise to superhumps observed in the V-band light curve of CVs, and we have previously shown for a CV with $`\mathrm{q}=0.15`$ that tidal effects prolong the outburst and drive eccentricity growth in the disc, with the launch of a two-armed spiral wave (Truss et al., 2001). Exactly the same situation applies to the transient system studied here. Fig. 13 shows the evolution of the disc over run A - the spiral wave and the eccentricity are clearly evident during the outburst.

The disc comes into contact with the 3:1 resonance soon after the onset of the outburst (fig. 13, top panel), but it takes some time for eccentricity to develop in the disc and the spiral waves to penetrate right down into the inner disc regions. Only then does the rebrightening occur (fig. 13, centre panel). By this time the disc has developed a significant eccentricity and gas near the disc edge is spread over a range of different radii. So, what is the viscous time-scale near the disc edge? Consider the disc shown in the centre panel of fig. 13. Gas at the edge of the disc nearest the secondary is at a radius $`\mathrm{R}0.3\mathrm{a}`$, while gas at the far edge is at $`\mathrm{R}0.6\mathrm{a}`$. This factor of two in radius translates to a factor of two in viscous time-scale.

In a real disc there is an additional factor which could drive accretion on a shorter-than-expected time-scale and produce the rebrightening. Consider an idealised accretion disc with a Shakura-Sunyaev $`\alpha `$ model structure. The equation of energy balance in the unirradiated part of the disc, 25, is now:

$$\frac{4\sigma \mathrm{T}_{\mathrm{m}}^{}{}_{}{}^{4}}{3\kappa \mathrm{\Sigma }}=\frac{9}{8}\nu \mathrm{\Sigma }\frac{\mathrm{GM}}{\mathrm{R}^3}+\mathrm{Q}_{\mathrm{tide}}.$$
(29)
where $`\mathrm{Q}_{\mathrm{tide}}`$ is the energy dissipation due to tides. This can be written in the form given by Papaloizou & Pringle (1977):

$$\mathrm{Q}_{\mathrm{tide}}=(\mathrm{\Omega }_{\mathrm{Kepler}}\mathrm{\Omega }_{\mathrm{orb}})\mathrm{T}_{\mathrm{tide}}$$
(30)
where $`\mathrm{T}_{\mathrm{tide}}`$ is the tidal torque per unit area. In 1-D treatments of disc instabiliies, this is usually included as

$$\mathrm{T}_{\mathrm{tide}}=\mathrm{CR}\nu \mathrm{\Sigma }\mathrm{\Omega }_{\mathrm{orb}}\left(\frac{\mathrm{R}}{\mathrm{a}}\right)^\mathrm{n}$$
(31)
(Smak 1984; Hameury et al. 1998). There are two uncertainties in this equation: the value of C, which is a numerically-determined constant, and the value of n, which determines how much of the disc is tidally heated. Nevertheless, we can use this to demonstrate the effect of tides on the mid-plane temperature of the disc. This is shown schematically for n=5 and an arbitrary value of C in fig. 11. Clearly, in order to significantly impact on the mid-plane temperature, the tidal torque on the outer disc has to be both large and sustained.

Figure 12 shows the energy dissipation rate in the disc averaged over one orbital period near the peak of the simulated outburst in run A, along with a similar plot for a simulated outburst of the dwarf nova SS Cygni (with the same sound speed). SS Cygni is not tidally unstable, since the disc never encounters the 3:1 radius - the dissipation rate falls off with radius to the edge of the disc and rises in a well-defined narrow region associated with the stream-impact hot-spot. However, in A0620-003, tidal forces produce the dominant contribution to dissipation over a surprisingly large area of the disc. The tidal heating is spread over a range of radii because the disc is eccentric - over an orbital period the disc edge is not confined to a single radius. The dissipation rate in the outer regions of the disc is consequently about an order of magnitude higher than expected. We stress that this enhanced dissipation does not have any effect on the viscous time-scale in the simulations : the viscosity is fixed at the value set in equation 17.

Figures 13 and 14 show the evolution of surface density and particle position in the disc. The disc edge very quickly encounters the 3:1 resonant radius and a two-armed spiral wave is launched (fig. 14, top panel). At the peak of the outburst, material accumulates in the spiral arms in the unirradiated outer regions of disc (solid curve in fig. 13), while gas is the inner region is accreted. The surface density in the outer ’spike’ of gas eventually surpasses the upper critical limit to trigger the disc instability. As the outburst progresses, the disc eccentricity grows. The rebrightening is associated with the arrival of the initially unirradiated gas at the black hole. This is not viscously spread on a time-scale that would apparently be appropriate to the outer disc - the spiral wave launched by the tidal instability penetrates the inner part of the disc and carries with it some of the gas in the outer disc (centre plot in fig. 14). Hence the rebrightening is caused by the accretion of a genuinely short pulse of gas. This penetration of the spiral arms into the inner regions of a hot disc has also been found in grid-based calculations performed by Stehle (1999). Successive waves of material are accreted in this way as the outburst progresses. We find that $`60\%`$ of the disc mass is accreted in the primary maximum. This figure would be higher in a calculation in which a greater fraction of the disc is exposed to the radiation.

The enhanced dissipation in the outer regions of the disc is consistent with the findings of Blondin (2000), who used an Eulerian scheme to show that tidal forces alone can produce a viscosity of the order of $`\alpha 0.1`$ in the outer regions of an accretion disc. The effective viscosity decreased rapidly with radius, and only made a significant contribution in the outer half to two-thirds of the disc: exactly what we find with SPH in fig. 12. It is still not absolutely clear that this will be sufficient to account for the 10 to 20 day time-scale for the rebrightening in A0620-003. Fig. 10 shows that the viscous time-scale at the inner edge of the region of enhanced dissipation ($`\mathrm{R}0.5\mathrm{R}_\mathrm{L}(1)6.5\times 10^{10}\mathrm{cm}`$) is $`85`$ days for a pure $`\alpha `$-disc. However, Blondin (2000) and Spruit (1987) note that the propagation of spiral waves in a real disc is more likely than in an isothermal disc, because in a real disc the sound speed increases with decreasing radius ($`\mathrm{c}_\mathrm{s}\mathrm{R}^{\frac{3}{8}}`$). Hence, with a more realistic equation of state, the effective $`\alpha `$ may be even higher in the inner disc. Furthermore, there is an additional effect at work in the SPH simulations which was not found by Blondin - the development of a significant eccentricity in the disc. These two effects may combine to make the short time-scale of the rebrightening achievable.

In the SPH simulations, both sound speed and viscosity are fixed, so the rebrightening is caused by the arrival at the primary of a density wave associated with the spiral waves penetrating deep into the inner part of a highly eccentic disc. We perform one further simulation to show that the rebrightening is indeed a tidal effect. The same input parameters as run A are used, with the exception that the mass ratio, $`\mathrm{q}=0.6`$. At this mass ratio, the 3:1 resonance should not be accessible to the disc and therefore it cannot become tidally unstable. The result is shown in fig. 15. As expected, no significant eccentricity develops and there is no rebrightening: just a simple decay with a long tail.

In summary, then, there are two time-scales associated with the rebrightening: the growth time-scale of the tidal instability in the disc, giving the time delay from the onset of the outburst to the start of the rebrightening, and the viscous time-scale of this tidally-stressed gas. In the simulations, the eccentricity of the disc causes gas initially near the disc edge to lie at much smaller radii. This explains why the rise-time of the rebrightening is so short : the gas that produces it is viscously spread on a time-scale much shorter than would be expected for gas originating near the edge of the disc. Undoubtedly, with a larger number of particles, the spatial resolution of the spiral shocks would improve and the simulated rebrightening would be even narrower. We have also proposed that in a real disc (but not the simulations) there is an additional factor that can increase the viscosity and make the rise-time of the rebrightening much shorter than expected. This is the tidal heating of the disc, raising the mid-plane temperature.

## 7 Conclusions

We have simulated complete outbursts of soft X-ray transient systems, for the first time in two dimensions. The rebrightenings seen in the X-ray light curves of several systems can be explained by competition between the ionization instability, the effects of disc irradiation and tidal forces from the secondary star. Their appearance in the light curve is ultimately determined by the irradiated area of the disc, the location within the disc where the disc instability is first triggered and the growth time-scale of the tidal instability. The main findings of this work can be summarized as:

* A rebrightening can occur if a small portion of the disc is shielded from the central X-rays in some way. This complements the conclusions of Dubus et al. (1999), who showed that a planar disc cannot be kept heated by irradiation from the central source, in contradiction to observations.
* Tidal torques can be a dominant factor in the energy balance over a large area of the disc. Such an enahnced dissipation in a disc with the classical $`\alpha `$-model structure will drive accretion on a shorter than expected time scale.
* The time delay between the onset of the outburst and the start of the rebrightening is determined by the growth time-scale of the eccentricity and spiral waves in the disc associated with the tidal instability.
* The rise-time of the rebrightening is governed by the disc viscosity - but tidal eccentricity effects cause this time-scale to be much shorter than expected for the arrival of gas from the disc edge.
* No rebrightening is reproduced if the entire disc is irradiated.

The main limitation of this approach is that the viscous evolution of the disc is governed by the local surface density and not a local mid-plane temperature. The ultimate goal is to merge the current 1D thermodynamic state-of-the-art with our 2D approach. The eventual extension to three dimensions is also particularly desirable in examining the beahviour of the tidal forces and the dynamics of the spiral arms.

There already exists a wealth of X-ray observations of transient sources, and we stress that although some aspects of the X-ray light curves are found to be common to many SXTs, they are all different to some degree. In this work we have presented a mechanism for the outbursts which explains the origin of the rebrightenings observed in the X-ray light curves of some of these systems.

Although one would not expect to observe superhumps in the X-ray light curve, superhumps have been seen in optical light curves of SXT outbursts (O’Donoghue & Charles 1996). The superhump mechanism must be slightly different to that operating in cataclysmic variables, as the optical emission from the disc is completely dominated by the reprocessed X-rays. Rather, the superhump modulation can be attributed to the change in visible area of an eccentric disc over an orbital period (Haswell et al. 2001). This is exactly the scenario which we find in this work, and is a clear demonstration of the importance of using a two- (or preferably a three-) dimensional code able to take the tidal effects into account in a natural way.

## Acknowledgments

Research in theoretical astrophysics at the University of Leicester is supported by a PPARC rolling grant. Some of the simulations were performed on the UK Astrophysical Fluids Facility (UKAFF). Others were performed on the Theoretical Astrophysics Group’s Linux Cluster which is supported by Advanced Micro Devices (AMD). The authors are grateful to the referee for extremely constructive comments and suggestions. MRT acknowledges a PPARC studentship and the support of the William Edwards Educational Charity.


text:
# Quantum Necking in Stressed Metallic Nanowires

## Abstract

When a macroscopic metallic wire is subject to tensile stress, it necks down smoothly as it elongates. We show that nanowires with radii comparable to the Fermi wavelength display remarkably different behavior. Using concepts from fluid dynamics, a PDE for nanowire shape evolution is derived from a semiclassical energy functional that includes electron-shell effects. A rich dynamics involving movement and interaction of kinks connecting locally stable radii is found, and a new class of universal equilibrium shapes is predicted.

non-linear dynamics, metallic nanowires

It has recently become possible to image metal nanowires with subatomic resolution, and record their structural evolution in real time, using transmission-electron microscopy Kondo and Takayanagi (1997); Kizuka (1998); Kondo and Takayanagi (2000); Rodrigues and Ugarte (2001, 2002), thus opening the door to a world of nonlinear dynamical phenomena at the nanoscale. Heretofore, most information about the structure and dynamics of nanowires was obtained indirectly, from transport and cohesive properties Rubio et al. (1996); Stalder and Dürig (1996); Untiedt et al. (1997); Yanson et al. (1999, 2000). Classical molecular dynamics simulations Landman et al. (1990); Todorov and Sutton (1996); Sørensen et al. (1998); Gülseren et al. (1998) provide important insights into atomic structure and dynamics of nanowires, but cannot explain the observed stability Kondo and Takayanagi (1997, 2000); Rodrigues and Ugarte (2002) of long nanowires, which would be unstable classically under surface tension Chandrasekhar (1981). Quantum molecular dynamics can in principle study this stability, but has so far been limited to systems too small to address such issues (e.g. single or double chains Hakkinen et al. (2000); Bahn and Jacobsen (2001)).

By contrast, a nanoscale free-electron model Stafford et al. (1997, 1999); Kassubek et al. (2001); Zhang et al.  has proven successful in explaining the quantum suppression of the Rayleigh instability Kassubek et al. (2001). Here we extend this approach by developing a dynamics for a free-electron nanowire under the assumption that surface diffusion is a dominant process. This dynamics, of instrinsic interest in nonlinear science, explicitly includes quantum-size effects and allows the study of a broad range of dynamical processes (see Figs. 13): Approach to equilibrium of a nanowire, propagation of an instability, and evolution of a wire under elongation/compression. Figure 1 shows how an initially random nanowire necks down to a universal shape, consisting of a central cylinder connected to contacts having the form of unduloids of revolution, providing a dynamical mechanism for the nanofabrication technique invented by Kondo and Takayanagi Kondo and Takayanagi (1997). Figure 3 compares well to early experiments on gold nanocontacts Rubio et al. (1996); Stalder and Dürig (1996), showing perfect correlation between force and conductance with steep conductance steps, and a hysteresis between elongation and compression.

The model of a nanowire consists of free electrons confined to an axisymmetric wire by a hard-wall potential Stafford et al. (1997). The wire radius is $`R(z,t)`$ in the interval $`[0,L]`$ along the $`z`$-axis, and periodic boundary conditions com (a) are used to extend $`R(z)`$. In contrast to the classical case Chandrasekhar (1981) for which all nonaxisymmetric distortions raise the energy, Jahn-Teller distortions that break axisymmetry of nanowires do lower the energy in a few cases. However, the deepest energetic minima occur for axisymmetric wires Urban et al. , and axisymmetry, once present, is preserved by the diffusion equation, so that Jahn-Teller distortions are most likely suppressed dynamically.

In the spirit of the Born-Oppenheimer approximation, the total energy of the nanowire is taken to be the electronic energy. Since we are dealing with an open system of electrons, the grand-canonical potential is used, and can be separated into Weyl and mesoscopic contributions Strutinsky (1968); Stafford et al. (1999); bra . Dropping the volume contribution, which is assumed to be constant, the energy functional is

$$\mathrm{\Omega }[T,R(z)]=\sigma (T)𝒮\left[R(z)\right]+_0^LV(T,R(z))\text{d}z,$$
(1)
where $`\sigma (T)`$ is the surface tension com (b), $`𝒮`$ is the surface area of the wire, and $`V`$ is a mesoscopic electron-shell potential whose integral gives the so-called “shell correction” to the energy (as used, for instance, in the study of clusters Yannouleas and Landman (1995)). Higher-order terms Stafford et al. (1999); bra  proportional to the mean curvature, etc., are neglected. $`V`$ can be expressed in terms of a Gutzwiller-type trace formula

$$V(T,R)=\frac{2\epsilon _F}{\pi }\underset{w=1}{\overset{\mathrm{}}{}}\underset{v=2w}{\overset{\mathrm{}}{}}a_{vw}(T)\frac{f_{vw}\mathrm{cos}\theta _{vw}}{v^2L_{vw}},$$
(2)
where the sum includes all classical periodic orbits $`(v,w)`$ in a disk billiard bra , characterized by their number of vertices $`v`$ and winding number $`w`$, $`L_{vw}=2vR\mathrm{sin}(\pi w/v)`$ is the length of orbit $`(v,w)`$, and $`\theta _{vw}=k_FL_{vw}3v\pi /2`$. The factor $`f_{vw}=1\text{ for }v=2w,2`$ otherwise, accounts for the invariance under time-reversal symmetry of some orbits, and $`a_{vw}(T)=\tau _{vw}/\mathrm{sinh}\tau _{vw}`$ ($`\tau _{vw}=\pi k_FL_{vw}T/2T_F`$) is a temperature-dependent damping factor. A similar trace formula was previously derived Kassubek et al. (2001) for small axisymmetric perturbations of a cylinder using semiclassical perturbation theory Ullmo et al. (1996); Creagh (1996); bra . Here we point out that semiclassical perturbation theory remains valid for large deformations, as long as new classes of non-planar orbits can be neglected (adiabatic approximation).

Eqs. (1) and (2) define an energy functional which is general and simple enough to solve for nontrivial nanowire geometries with a wide range of radii. The multiple deep minima of $`V(R)`$, shown in Fig. 1(a), favor certain magic radii Yanson et al. (1999) and suggest a multistable field theory analogous to the sine-Gordon model. We therefore anticipate solutions consisting of cylindrical segments connected by kinks.

Rather than directly minimizing the energy functional (1), we derive an equation describing the nanowire surface dynamics which yields the approach to equilibrium as well as the onset of instability. The model has two main assumptions. First, in the thin wires considered, the majority of atoms are at the surface and thus surface self-diffusion is the mechanism of ionic motion. Second, under the Born-Oppenheimer approximation, the electronic energy (1) acts as a potential for the ions, and the dynamics derive from ionic mass conservation:

$$\frac{R(z,t)}{t}=\frac{v_a}{R(z,t)}\frac{}{z}\left[R(z,t)J_z(z,t)\right],$$
(3)
where $`v_a=3\pi ^2/k_F^3`$ is the volume of an atom, and the $`z`$-component of the surface current is given by Fick’s law:

$$J_z=\frac{\rho _SD_S}{k_BT}\frac{1}{\sqrt{1+(R^{})^2}}\frac{\mu }{z},$$
(4)
where $`R^{}=R/z`$, and $`\rho _S`$ and $`D_S`$ are the surface density of ions and the surface self-diffusion coefficient, respectively. The precise value of $`D_S`$ for alkali metals is not known, but it can be removed from the evolution equation by rescaling time to the dimensionless variable $`\tau =(\rho _SD_ST_F/T)t`$. The chemical potential $`\mu `$ is obtained by calculating the change in the energy (1) with the addition of an atom at point $`z_0`$, $`\mu (z_0)\mathrm{\Omega }[T,R(z)+A\delta (zz_0)]\mathrm{\Omega }[T,R(z)]`$, where $`A=v_a/2\pi R`$ is chosen so that the volume of an atom is added:

$$\begin{array}{c}\mu (z)=\frac{2\epsilon _F}{5}+\frac{3\pi \sigma }{k_F^3R(z)\sqrt{1+R_{}^{}{}_{}{}^{2}}}\left(1\frac{RR^{\prime \prime }}{1+R_{}^{}{}_{}{}^{2}}\right)\hfill \\ \hfill \frac{3\epsilon _F}{k_F^2R^2}\underset{w,v}{}\frac{a_{vw}f_{vw}}{v^2}\left[\mathrm{sin}\theta _{vw}+b_{vw}\frac{\mathrm{cos}\theta _{vw}}{k_FL_{vw}}\right],\end{array}$$
(5)
where $`b_{vw}(T)=a_{vw}\mathrm{cosh}\tau _{vw}`$. Eq. (3) is then solved numerically using an implicit scheme. Our non-linear dynamical model, Eqs. (35), differs from previous studies of axisymmetric surface self-diffusion Coleman et al. (1995); Eggers (1998); Bernoff et al. (1998) by the inclusion of electron-shell effects \[2nd line of Eq. (5)\], which fundamentally alter the dynamics.

To explore the equilibrium shapes of nanowires, we considered several random initial configurations of various lengths and widths, one of which is shown in Fig. 1(b). Under surface self-diffusion, the short-wavelength fluctuations are rapidly smoothed out \[Fig. 1(c)\], and several cylindrical segments connected by “kinks” are visible. On a much longer time-scale, the wires neck down to a universal equilibrium shape consisting of a nearly perfect cylindrical region with a radius near one of the minima of $`V(R)`$, connected to thicker leads \[Fig. 1(d)\]. In all fifteen cases studied, the leads closely approximate Delaunay unduloids of revolution. The unduloid is a stationary configuration, unstable under surface tension alone Bernoff et al. (1998), but stabilized in the present model by electron-shell effects in the connecting cylinder. The unduloids arise as superpositions of large numbers of closely spaced kinks. Equilibrium wires with conductance $`G/G_0=1,3,6,12,17,23,\mathrm{}`$ were obtained, where $`G_0=2e^2/h`$ is the conductance quantum. Intermediate conductance values may occur more rarely.

Other cases of interest are the propagation of an instability Zhang et al. (1998) and the homogeneous-inhomogeneous transition in a nanowire Zhang et al. . In the former case, starting with a cylinder of an unstable radius, a maximum in $`V(R)`$ \[Fig. 1(a)\], and adding a localized Gaussian deformation of small amplitude ($`0.001R`$), we find that the perturbation decreases during an incubation period and then the instability sets in with essentially a single Fourier component corresponding to the most unstable mode Zhang et al. , saturates at a finite amplitude, and propagates with constant velocity (see Fig. 2). A similar instability is observed in the pearling of membranes Bar-Ziv and Moses (1994). The classical Rayleigh instability itself Zhang et al. (1998) propagates in much the same way, but does not saturate, instead leading to breakup of the cylinder. Both the wavelength Zhang et al.  and front velocity of the instability are strongly dependent on the initial radius. Once the instability is fully developed, kink-antikink pairs start annihilating, eventually leading to a shape consisting of a series of cylinders of different radii \[corresponding to neighboring minima in $`V(R)`$, Fig. 1(a)\], connected by kinks. Starting from a more realistic initial condition with random perturbations, the instability occurs at several places more or less simultaneously.

The dynamical model also allows for the study of the evolution of a nanowire under strain. The length $`L`$ of the wire is changed at an average rate $`dL/dt`$ in a sequence of infinitesimal steps. In each step, the wire deforms elastically, i.e. the change of length is distributed along the wire as if each of its slices were Hookean, with spring constant proportional to its cross-sectional area. The radius of each slice is simultaneously changed so as to conserve its volume. Between steps, the wire evolves according to Eq. (3). We start with an equilibrated unduloid of outer radius $`k_FR=99`$, representing a nanocontact between two quasi-macroscopic leads \[see Fig. 3 (inset)\]. We then strain the wire by elongating it at a constant speed $`k_FdL/d\tau =0.1`$ until it necks down to a long, atomically thin wire \[Fig. 3(b-g)\], at which point we let it partially equilibrate \[Fig. 3(h)\]. Subsequently, we compress the wire at the same rate until it reaches its initial length \[Fig. 3(i-n)\]. During elongation, the necking of the wire occurs by nucleation of kink-antikink pairs at the center with subsequent motion of these kinks toward the leads, leaving a long, cylindrical wire in the middle. During compression, the remaining kink-antikink pairs annihilate, but new kinks do not emerge from the leads. Similar necking was observed for various elongation rates $`k_FdL/d\tau 1`$. For much faster rates, structural relaxation due to surface diffusion is suppressed.

These different processes lead to a substantial hysteresis between elongation and compression, as can be seen by computing the electrical conductance of the contact and the force applied to it. The force is obtained from the change in the energy (1) with elongation, $`F=\mathrm{\Omega }/L`$. The conductance is computed from the Landauer formula, using the adiabatic and WKB approximations Stafford et al. (1997) to compute the transmission probabilities. Both quantities are shown in the central inset of Fig. 3. Compared to our previous results obtained without structural relaxation Stafford et al. (1997), the conductance steps and relaxations of the force are much steeper, in better agreement with experiments Rubio et al. (1996), and a substantial hysteresis is observed in the conductance. It has been argued Rubio et al. (1996); Untiedt et al. (1997) that the abruptness of the conductance steps, and their perfect correlation with the jumps in the force, rule out an electronic mechanism. However, we have shown that this behavior arises naturally in a model which takes proper account of electronic quantum-size effects on the structure, even when atomistic effects are neglected.

A final interesting feature is the spikes in the force that occur at the opening and closing of a channel, corresponding to the rapid energy relaxation during the creation or annihilation of a kink/antikink pair. These spikes are suppressed, and the hysteresis decreases, when the speed of deformation is increased, suggesting that the spikes might be observable experimentally by decreasing the elongation rate. A failure of this model is the force during compression: contrary to experiment Rubio et al. (1996); Untiedt et al. (1997), the force is attractive, suggesting that a resistance to compression due to the ions needs to be added to the free-electron model. Indeed, an interesting question raised by this work is how discrete positive ions would arrange themselves to accomodate the predicted shapes favored by the conduction electrons.

###### Acknowledgements.

This research was supported in part by NSF grants DMR0072703 and DMR0312028 (CAS & JB), and by NSF grant DMR9812526 (REG). JB also acknowledges support from the Swiss National Science Foundation and grant SFB276 of the Deutsche Forschungsgemeinschaft.


text:
# Evolution of Brown Dwarf Atmospheres: Investigating Physical Parameters from Near-IR Spectra.

## 1. Introduction

In order to study the initial mass function in young open clusters at the low mass end, one needs to have a fast and reliable way of separating pre-MS members from foreground dwarfs and background giants. These objects occupy same place on the color-magnitude diagram, but have different surface gravities. Using this fact, we decided to explore the gravity-sensitive features in the near-infrared spectra of the objects of spectral types M–L, corresponding to the range of masses near sub-stellar boundary.

## 2. Observations

Using FSpec IR spectrograph on the 6.5m MMT telescope (Mt. Hopkins, Arizona), we obtained J & K spectra of $``$ 30 brown dwarf candidates, SpT M4–L5, drawn from the field and open clusters and associations ($`\rho `$ Oph, Tau, IC 348, Upper Sco, Pleiades) to sample a range of gravities $`logg`$ $``$ 3 to 5 dex (Fig.1), as well as 3 late M giants ($`logg`$$``$0). The resolution of our smoothed spectra is $``$700 and an average S/N in flux is 40–80 (depending on the feature).

## 3. Results

We calculated gravities using masses from the theoretical $`T_{eff}`$$`L`$ diagram, assuming dwarf colors and having distances to most of our objects. The correlation of some features measured from our spectra with gravity is shown in Fig.1.

We confirm strong temperature dependence of H2O and CO bands (eg. Reid et al. 2001) and infer that gravity dependence of these features is insignificant compared to that of FeH, KI and NaI lines. While the weakness of the latter features in cluster compared to field objects was known before (eg. Martin et al. 1996, Luhman et al. 1998, Bejar et al. 1999, Lucas et al. 2001), no systematic investigation of this effect has been carried out so far, especially in the J band. IR is the only option to work with cool, low mass and embedded objects, and in the J band the atomic features are especially prominent.

We also compared our spectra to recent atmospheric calculations by Allard et al. 2001 with the conclusion that while they predict correct trends, the systematic offsets between observations and theory remain (Fig.1).

Using our results we managed to confirm low surface gravity of two brown dwarfs from the TW Hya and Taurus associations, whose membership was originally proposed based on their optical spectra (2MASS1139-3159 and CFHT-BD-Tau4).

## References

Allard, F., Hauschildt, P. H., Alexander, D. R., Tamanai, A., & Schweitzer, A., 2001, ApJ, 556, 357

Bejar, V. J. S., Zapatero Osorio, M. R., & Rebolo, R., 1999, ApJ, 521, 671

Lucas, P. W., Roche, P. F., Allard, F., & Hauschildt P. H., 2001, MNRAS, 326, 695

Luhman, K. L., Rieke, G. H., Lada, C. J., & Lada, E. A. 1998, ApJ, 508, 347

Martin, E. L., Rebolo, R., & Zapatero Osorio, M. R. 1996, ApJ, 469, 706

Reid, I. N., Burgasser, A.J., Cruz, K. L., Kirkpatrick, J. D., & Gizis J. E. 2001, AJ, 121, 1710