text:
# The formation of close binary systems by dynamical interactions and orbital decay

## 1 Introduction

The process of star formation preferentially produces binary stellar systems (e.g. Duquennoy & Mayor 1991). The favoured mechanism for explaining this high frequency of binaries is the collapse and fragmentation of molecular cloud cores (e.g. Boss & Bodenheimer 1979; Boss 1986; Bonnell et al. 1991; Nelson & Papaloizou 1993; Burkert & Bodenheimer 1993; Bate, Bonnell & Price 1995). However, while fragmentation can readily create wide binary systems (separations $`\text{ }>10`$ AU), there are severe difficulties with fragmentation producing close binaries directly. This is a significant deficiency since approximately 20% of solar-type stars have main-sequence companions that orbit closer than 10 AU (Duquennoy & Mayor 1991), and the frequency of massive spectroscopic binaries appears to be even higher (Garmany, Conti & Massey 1980; Abt et al. 1990; Morrell & Levato 1991; Mason et al. 1998).

As a molecular cloud core begins to collapse, the formation of wide binaries through fragmentation is possible because the gas easily radiates away the gavitational potential energy that is released. The gas remains approximately isothermal and, thus, the Jeans mass decreases with density as $`\rho ^{1/2}`$. However, at densities of $`\text{ }>10^{13}`$ g cm-3 or $`n(\mathrm{H}_2)\text{ }>10^{10}`$ cm-3 (Larson 1969; Masunaga & Inutsuka 2000) the rate of heating from dynamical collapse exceeds the rate at which the gas can cool. The gas heats up, and the Jeans mass begins to increase so that a Jeans-unstable region of gas becomes Jeans-stable. This results in the formation of a pressure-supported fragment with a mass of several Jupiter-masses and a radius of $`5`$ AU (Larson 1969). Fragmentation on smaller scales is inhibited by thermal pressure. Therefore, initial binary separations must be $`\text{ }>10`$ AU. The formation of such pressure-supported fragments is frequently refered to as the opacity limit for fragmentation (Low & Lynden-Bell 1976; Rees 1976) and may set a lower limit to the mass of brown dwarfs (Boss 1988; Bate, Bonnell & Bromm 2002).

A possibility for fragmentation at higher densities (hence on smaller length scales) exists when the pressure-supported fragment has accreted enough material for its central temperature to exceed 2000 K. At this temperature, molecular hydrogen begins to dissociate, which provides a way for the release of gravitational energy to be absorbed without significantly increasing the temperature of the gas. Thus, a nearly isothermal second collapse occurs within the fragment that ultimately results in the formation of a stellar core with radius $`1R_{}`$ \[Larson 1969\]. Several studies have investigated the possibility that fragmentation during this second collapse forms close binary systems directly (Boss 1989; Bonnell & Bate 1994; Bate 1998, 2002). Boss \[Boss 1989\] found that fragmentation was possible during this second collapse, but that the fragments spiralled together due to gravitational torques and did not survive. Bonnell & Bate \[Bonnell & Bate 1994\] found that fragmentation to form close binaries and multiple systems could occur in a disc that forms around the stellar core. However, both these studies began with somewhat arbitrary initial conditions for the pressure-supported fragment. Bate \[Bate 1998\] performed the first three-dimensional calculations to follow the collapse of a molecular cloud core through the formation of the pressure-supported fragment, the second collapse, and the formation of the stellar core and its surrounding disc. In these and subsequent calculations (Bate 2002), Bate found that the second collapse did not result in sub-fragmentation due to the high thermal pressure and angular momentum transport via gravitational torques.

In this paper, we present results from the first hydrodynamical star formation calculation to produce dozens of stars and brown dwarfs while simultaneously resolving beyond the opacity limit for fragmentation. Despite the fact that no close binaries (separations $`\text{ }<10`$ AU) are formed by direct fragmentation, we find that the calculation eventually produces several close binary systems through a combination of dynamical interactions in multiple systems, and the orbital decay of wide binaries via gas accretion and their interactions with circumbinary discs.

The paper is structured as follows. In section 2, we briefly describe the numerical method and the initial conditions for our calculation. In section 3, we present results from our calculation and compare them with observations. Finally, in section 4, we give our conclusions.

## 2 Computational method and initial conditions

The calculation presented here was performed using a three-dimensional, smoothed particle hydrodynamics (SPH) code. The SPH code is based on a version originally developed by Benz (Benz 1990; Benz et al. 1990). The smoothing lengths of particles are variable in time and space, subject to the constraint that the number of neighbours for each particle must remain approximately constant at $`N_{\mathrm{neigh}}=50`$. We use the standard form of artificial viscosity (Monaghan & Gingold 1983) with strength parameters $`\alpha __\mathrm{v}=1`$ and $`\beta _\mathrm{v}=2`$. Further details can be found in Bate et al. \[Bate, Bonnell, & Price 1995\]. The code has been parallelised by M. Bate using OpenMP.

### 2.1 Opacity limit for fragmentation and the equation of state

To model the opacity limit for fragmentation without performing full radiative transfer, we use a barotropic equation of state for the thermal pressure of the gas $`p=K\rho ^\eta `$, where $`K`$ is a measure of the entropy of the gas. The value of the effective polytropic exponent $`\eta `$, varies with density as

$$\eta =\{\begin{array}{cc}\begin{array}{cc}\hfill 1,& \rho 10^{13}\mathrm{g}\mathrm{cm}^3,\hfill \\ \hfill 7/5,& \rho >10^{13}\mathrm{g}\mathrm{cm}^3.\hfill \end{array}\hfill & \end{array}$$
(1)
We take the mean molecular weight of the gas to be $`\mu =2.46`$. The value of $`K`$ is defined such that when the gas is isothermal $`K=c_\mathrm{s}^2`$, with the sound speed $`c_\mathrm{s}=1.84\times 10^4`$ cm s-1 at 10 K, and the pressure is continuous when the value of $`\eta `$ changes.

This equation of state reproduces the temperature-density relation of molecular gas during spherically-symmetric collapse (as calculated with frequency-dependent radiative transfer) to an accuracy of better than 20% in the non-isothermal regime up to densities of $`10^8\mathrm{g}\mathrm{cm}^3`$ \[Masunaga & Inutsuka 2000\]. Thus, our equation of state should model collapsing regions well, but may not model the equation of state in protostellar discs particularly accurately due to their departure from spherical symmetry.

### 2.2 Sink particles

The opacity limit results in the formation of distinct pressure-supported fragments in the calculation. As these fragments accrete, their central density increases, and it becomes computationally impractical to follow their internal evolution until they undergo the second collapse to form stellar cores because of the short dynamical time-scales involved. Therefore, when the central density of a pressure-supported fragment exceeds $`\rho _\mathrm{s}=10^{11}\mathrm{g}\mathrm{cm}^3`$, we insert a sink particle into the calculation (Bate et al. 1995). The gas within radius $`r_{\mathrm{acc}}=5`$ AU of the centre of the fragment (i.e. the location of the SPH particle with the highest density) is replaced by a point mass with the same mass and momentum. Any gas that later falls within this radius is accreted by the point mass if it is bound and its specific angular momentum is less than that required to form a circular orbit at radius $`r_{\mathrm{acc}}`$ from the sink particle. Thus, gaseous discs around sink particles can only be resolved if they have radii $`\text{ }>10`$ AU. Sink particles interact with the gas only via gravity and accretion.

Since all sink particles are created from pressure-supported fragments, their initial masses are $`10`$ Jupiter-masses (MJ), as given by the opacity limit for fragmentation \[Boss 1988\]. Subsequently, they may accrete large amounts of material to become higher-mass brown dwarfs ($`\text{ }<75`$ MJ) or stars ($`\text{ }>75`$ MJ), but all the stars and brown dwarfs begin as these low-mass pressure-supported fragments.

The gravitational acceleration between two sink particles is Newtonian for $`r4`$ AU, but is softened within this radius using spline softening \[Benz 1990\]. The maximum acceleration occurs at a distance of $`1`$ AU; therefore, this is the minimum separation that a binary can have even if, in reality, the binary’s orbit would have been hardened.

Replacing the pressure-supported fragments with sink particles is necessary in order to perform the calculation. However, it is not without a degree of risk. If it were possible to follow the fragments all the way to stellar densities (as done by Bate 1998) and still follow the evolution of the large-scale cloud over its dynamical time-scale, we might find that a few of the objects that we replace with sink particles merge together or are disrupted by dynamical interactions. We have tried to minimise the degree to which this may occur by insisting that the central density of the pressure-supported fragments exceeds $`\rho _\mathrm{s}`$ before a sink particle is created. This is two orders of magnitude higher than the density at which the gas is heated and ensures that the fragment is self-gravitating, centrally-condensed and, in practice, roughly spherical before it is replaced by a sink particle. In theory, it would be possible for a long collapsing filament to exceed this density over a large distance, thus making the creation of one or more sink particles ambiguous. However, the structure of the collapsing gas that results from the turbulence prohibits this from occuring; no long roughly uniform-density filaments with densities $`\rho _\mathrm{s}`$ form during the calculation. Furthermore, each pressure-supported fragment must undergo a period of accretion before its central density exceeds $`\rho _\mathrm{s}`$ and it is replaced by a sink particle. For example, it is common in the calculation to be able to follow a pressure-supported fragment that forms via gravitational instability in a disc for roughly half an orbital period before it is replaced. Thus, the fragments do have some time in which they may be disrupted or merge. Only occasionally during the calculation are low-mass pressure-supported fragments disrupted; most are eventually replaced by sink particles.

### 2.3 Initial conditions

The initial conditions consist of a large-scale, turbulent molecular cloud. The cloud is spherical and uniform in density with a mass of 50 M and a diameter of $`0.375`$ pc (77400 AU). At the temperature of 10 K, the mean thermal Jeans mass is 1 M (i.e. the cloud contains 50 thermal Jeans masses). The free-fall time of the cloud is $`t_{\mathrm{ff}}=6.0\times 10^{12}`$ s or $`1.90\times 10^5`$ years.

Although the cloud is uniform in density, we impose an initial supersonic turbulent velocity field on it in the same manner as Ostriker, Stone & Gammie \[Ostriker, Stone, & Gammie 2001\]. We generate a divergence-free random Gaussian velocity field with a power spectrum $`P(k)k^4`$, where $`k`$ is the wavenumber. In three dimensions, this results in a velocity dispersion that varies with distance, $`\lambda `$, as $`\sigma (\lambda )\lambda ^{1/2}`$ in agreement with the observed Larson scaling relations for molecular clouds \[Larson 1981\]. The velocity field is normalised so that the kinetic energy of the turbulence equals the magnitude of the gravitational potential energy of the cloud. The initial root-mean-square Mach number of the turbulence is $`=6.4`$.

### 2.4 Resolution

The local Jeans mass must be resolved throughout the calculation to model fragmentation correctly (Bate & Burkert 1997; Truelove et al. 1997; Whitworth 1998; Boss et al. 2000). Bate & Burkert \[Bate & Burkert 1997\] found that this requires $`\text{ }>2N_{\mathrm{neigh}}`$ SPH particles per Jeans mass; $`N_{\mathrm{neigh}}`$ is insufficient. We have repeated their calculation using different numbers of particles and find that $`1.5N_{\mathrm{neigh}}=75`$ particles is also sufficient to resolve fragmentation (Bate, Bonnell & Bromm, in preparation). The minimum Jeans mass in the calculation presented here occurs at the maximum density during the isothermal phase of the collapse, $`\rho =10^{13}`$ g cm-3, and is $`0.0011`$ M (1.1 MJ). Thus, we use $`3.5\times 10^6`$ particles to model the 50-M cloud. In fact, a gas clump with the above mass and density could not collapse because as soon as it was compressed it would heat up (equation 1) and would no longer contain a Jeans mass. Therefore, in practice, any collapsing gas clump in the simulation contains many more SPH particles than the ‘minimum’ Jeans mass above. This SPH calculation is one of the largest ever performed. It required approximately 95000 CPU hours on the SGI Origin 3800 of the United Kingdom Astrophysical Fluids Facility (UKAFF).

## 3 Results

### 3.1 Evolution of the cloud

The hydrodynamical evolution of the cloud produces shocks which decrease the turbulent kinetic energy initially supporting the cloud. In parts of the cloud, gravity begins to dominate and dense self-gravitating cores form and collapse. These dense cores are the sites where the formation of stars and brown dwarfs occurs. Although the cloud initially contains more than enough turbulent energy to support itself against gravity, this turbulence decays on the dynamical time-scale of the cloud and star formation begins after just one global free-fall time at $`t=1.04t_{\mathrm{ff}}`$ (i.e. $`t=1.97\times 10^5`$ yrs). This rapid decay of the turbulence is consistent with other numerical studies of turbulence in molecular clouds (e.g. MacLow et al. 1998; Stone, Ostriker & Gammie 1998; Ostriker et al. 2001). We followed the calculation for $`69000`$ years after the star formation began until $`t=1.40t_{\mathrm{ff}}`$ ($`t=2.66\times 10^5`$ yrs). At this point we had exhausted our allocation of computer time and the calculation was stopped.

Figure 1 illustrates the state of the cloud at the end of the calculation. The star formation occurs in three main dense cores within the cloud, resulting in three stellar groups. These groups are composed of single, binary and higher-order systems that were produced via a combination of the fragmentation of collapsing gas filaments (whose collapse is halted when the gas heats up at high densities; Inutsuka & Miyama 1992), the fragmentation of massive circumstellar discs (e.g. Bonnell 1994; Whitworth et al. 1995; Burkert, Bate, & Bodenheimer 1997), and star-disc capture (Larson 1990; Clarke & Pringle 1991a,b). At any particular time, the largest of the stellar groups contains no more than $`20`$ objects. Thus, the groups dissolve quickly. In fact, the time-scales for star formation within the dense cores and the dissolution of the groups are both $`2\times 10^4`$ years. Thus, the stellar groups undergo chaotic evolution with stars and brown dwarfs being ejected from the cloud even as new objects are forming. If the calculation were followed until well after star formation had ceased, all of the groups would dissolve on time-scales of a few tens of thousands of years.

When the calculation is stopped, the cloud has produced 23 stars and 18 brown dwarfs. An additional 9 objects have substellar masses, but are still accreting. Three of these would probably end up with substellar masses if the calculation were continued, but the other six are likely to become stars. The formation mechanism of the brown dwarfs in this calculation has been discussed by Bate et al. \[Bate, Bonnell, & Bromm 2002\]. The evolution of the cloud and the properties of the stars and brown dwarfs will be discussed in detail in subsequent papers. In this paper, we concentrate on the mechanism by which the close binary systems form.

### 3.2 The formation of close binary systems

When the calculation is stopped, amongst the 50 stars and brown dwarfs, there are 7 close binary systems with separations less than $`10`$ AU (Table 1 and Figures 1 and 2). Six of these close binaries are members of unstable multiple systems, but one has been ejected from the cloud on its own and will evolve no further (Figures 1 and 2).

As described in the introduction, it seems that binaries with separations $`\text{ }<10`$ AU cannot form by direct fragmentation because fragmentation is halted by the formation of pressure-supported fragments with radii of $`5`$ AU. Our calculation models this opacity limit for fragmentation and, as expected, none of the fragments in our calculation form closer than $`10`$ AU from each other. The three smallest separations between any existing object and a forming fragment are 9, 21, and 22 AU and only the last of these ends up in a close binary system. What mechanism(s) produce the 7 close binary systems in our calculation?

In Figure 2, we trace the histories of the stars and brown dwarfs that play a role in producing each of the 7 close binary systems. We find that the close binary systems are produced as a natural consequence of the evolution of wide binaries and multiple systems in a gas-rich environment. Three processes are involved: the accretion of gas, the interaction of a binary or triple with its circumbinary or circumtriple disc, and dynamical interactions between objects.

#### 3.2.1 Accretion

Accretion aids in the production of close binaries in two main ways. First, accretion onto a binary can decrease its separation (e.g. Figure 3; Figure 2d, steps 5–7). This effect of gas accretion on the orbital separation of a binary has been studied extensively (Artymowicz 1983; Bate 1997; Bate & Bonnell 1997; Bate 2000). Bate & Bonnell \[Bate & Bonnell 1997\] show that accretion onto a binary from a gaseous envelope decreases the orbital separation of a binary unless the specific angular momentum of the accreted gas is significantly greater than that of the binary. If a binary accretes a large amount of material compared to its initial mass (e.g. a factor $`\text{ }>10`$) this can easily reduce the separation of the binary by 1–2 orders of magnitude \[Bate 2000\].

Second, accretion can similarly destabilise stable hierarchical multiple systems by reducing the separation of wide components (e.g. Figure 4; Figure 2b, step 11; Figure 2c, steps 3 and 5; Figure 2d, step 2). This and other effects of accretion on the stability and mass ratios of triple systems were studied by Smith, Bonnell & Bate (1997). The destabilisation of multiple systems contributes to the production of close binaries because it forces the system to undergo dynamical interactions (section 3.2.3).

#### 3.2.2 Circumbinary and circumtriple discs

If a binary is surrounded by a circumbinary disc, gravitational torques from the binary transfer angular momentum from the binary’s orbit into the disc, causing the binary’s components to spiral together (Artymowicz et al. 1991; Bate & Bonnell 1997). Thus, a wide binary can be hardened into a close binary (e.g. Figure 5; Figure 2b, steps 1 and 6; Figure 2d, step 7). Similarly, a disc around a hierarchical triple can reduce the separation of the wide companion, destabilising the system (e.g. Figures 4 and 6; Figure 2b, step 11; Figure 2c, steps 3 and 5; Figure 2d, step 2).

Such disc interactions are very efficient; even relatively low-mass discs can have a significant effect over time \[Pringle 1991\]. The circumbinary and circumtriple discs in this calculation frequently contain a significant fraction of the system’s mass as can be observed by the spiral density waves caused by gravitational instabilities within them (Figures 1, 3, 5, and 6). This makes it difficult to separate the effects of accretion and gravitational torques, but, in any case, these discs play an important role in reducing the separations of binaries and triples.

#### 3.2.3 Dynamical interactions

Dynamical interactions can lead to the orbital evolution of a binary in several ways. If the orbital velocity of a binary is greater than the velocity of an incoming object while it is still at a great distance (i.e. the binary is ‘hard’), the binary will survive the encounter \[Hut & Bahcall 1983\]. However, several outcomes are possible. The binary may simply be hardened by the encounter, with the single object removing energy and angular momentum. Alternately, if the encounter is sufficiently close, an unstable multiple system will be formed. Its chaotic evolution will usually lead to the ejection of the object with the lowest mass. If the ejected object was a component of the original binary, the net effect is an exchange interaction.

The suggestion that dynamical encounters play an important role in the formation of close binaries has been made by Tokovinin (1997, 2000), who observed that many close binaries have wider companions (see below). In the calculation discussed here, dynamical encounters play the dominant role in producing close binaries (e.g. Figure 6 and the exchange interactions, hardening by fly-bys, and binary-binary interactions described in Figure 2). However, it is important to note that N-body dynamical interactions alone cannot produce the high frequency of close binaries that we obtain. Indeed, Kroupa & Burkert \[Kroupa & Burkert 2001\] find that N-body calculations which begin with star clusters (100 to 1000 stars) consisting entirely of binaries with periods $`4.5<\mathrm{log}(P/\mathrm{days})<5.5`$ produce almost no binaries with periods $`\mathrm{log}(P/\mathrm{days})<4`$. Similarly, the dissolution of small-N clusters typically results in binaries with separations only an order of magnitude smaller than the size of the initial cluster (Sterzik & Durisen 1998).

The key difference in this calculation is the effect of the gas-rich environment. Along with the effects from accretion and circumbinary/circumtriple discs discussed above, the presence of gas allows dynamical encounters to be dissipative and to transport angular momentum. Such dissipative encounters include star-disc encounters (Larson 1990; Clarke & Pringle 1991a,b; McDonald & Clarke 1995; Hall, Clarke & Pringle 1996), where the dissipation of kinetic energy allows the formation of bound systems from objects that would otherwise be unbound, and other tidal interactions (Larson 2002).

### 3.3 The resulting properties of close binary systems

As mentioned above, dynamical exchange interactions typically lead to the ejection of the least massive component. Thus, if a binary encounters a star whose mass is in between the masses of the binary’s components, the binary’s secondary will usually be replaced and the binary’s mass ratio will be equalised. In addition, gas with high specific angular momentum that accretes from a circumbinary disc or a gaseous envelope onto a binary, is preferentially captured by the secondary and, thus, also drives the mass ratio towards unity (Artymowicz 1983; Bonnell & Bastien 1992; Whitworth et al. 1995; Bate & Bonnell 1997; Bate 2000). Even if accretion comes from a relatively slowly rotating cloud, the long-term effect of this accretion is usually to equalise the binary’s components (Bate 2000). Thus, each of the mechanisms involved in producing close binaries favours the production of equal-mass systems. This is reflected in the mass ratios of the seven close binary systems that form in our calculation (Table 1), all of which have values $`q\text{ }>0.3`$ and most of which have $`q>1/2`$. Thus, the observational result that close binaries (periods $`\text{ }<10`$ years) tend to have higher mass ratios than wider binaries (Mazeh et al. 1992; Halbwachs, Mayor & Udry 1998; Tokovinin 2000) is explained naturally if these processes are an integral part of the formation mechanism for close binaries.

Overall, the calculation produces 7 close binaries among 50 stars and brown dwarfs, giving a close-binary fraction of $`7/43=16`$%. Given the small numbers, this is in good agreement with the observation that $`20`$% of solar-type stars have a less massive main-sequence companion within 10 AU \[Duquennoy & Mayor 1991\]. Thus, dynamical encounters and orbital decay can produce the observed frequency of close binaries.

However, because dynamical exchange interactions preferentially select the two most massive objects, another consequence of the way in which close binaries form is that the frequency of close binaries increases with the mass of the primary. This can be understood by considering a binary and an incoming object whose mass is greater than either of the two components of the binary. In this case, the typical result is that the incoming object and the primary from the original binary form a new binary. After several such exchange interactions, the mass of the primary in the new binary can be much greater than the original primary’s mass. In the calculation presented here, the dependence of the frequency of close binaries on stellar mass is quite strong. Of $`20`$ brown dwarfs, there is only one binary brown dwarf system (and if the calculation were followed further, at least one of the components would probably accrete into the stellar regime), while of the 11 stars with masses $`>0.2`$ M, 5 are members of close binary systems. While it is difficult to extrapolate from our low-mass star-forming cloud to larger star clusters and more massive stars, this trend of an increasing frequency of close binaries with stellar mass is supported by observational surveys (Garmany et al. 1980; Abt et al. 1990; Morrell & Levato 1991; Mason et al. 1998).

Finally, we note that at the end of the calculation, only one of our close binaries has been ejected from the cloud and is on its own. The others remain members of larger-scale bound groups and three are members of hierarchical triple systems (Table 1; Figure 1). This large number of wider companions is yet another indication of the importance of multiple systems in producing close binaries. Even allowing for the eventual dissolution of the bound groups, it seems likely that some of the hierarchical triple systems will survive. Although the true frequency of wide companions to close binaries is not yet well known, many close binaries do have wider components (e.g. Mayor & Mazeh 1987; Tokovinin 1997, 2000). Indeed, it was this observation that led Tokovinin (1997) to propose that dynamical interactions in multiple systems may play an important role in the formation of close binary systems. Further surveys to determine the true frequency of wide companions to close binary systems would be invaluable.

## 4 Conclusions

We have presented results from a hydrodynamic calculation of the collapse and fragmentation of a turbulent molecular cloud to form 50 stars and brown dwarfs. The calculation mimics the opacity limit for fragmentation by using a barotropic equation of state to model the heating of collapsing gas at high densities. This results in a minimum mass of $`10`$ Jupiter masses for the lowest-mass brown dwarfs (see Bate et al. 2002) and prevents fragmentation on scales smaller than $`10`$ AU. Despite this lower limit on the initial minimum separation between fragments, we find that as the stellar groups and unstable multiple systems evolve in the gas-rich environment, a high frequency of close binary systems (separations $`\text{ }<10`$ AU) is produced.

Examining the history of these close binary systems, we find that they are formed through a combination of dynamical interactions in unstable multiple systems, and orbital decay due to accretion and/or the interaction of binary and triple systems with circumbinary and circumtriple discs. These formation mechanisms allow realistic numbers of close binary systems to be produced without the need for fragmentation on length scales $`<10`$ AU. This avoids the difficulties associated with the fragmentation of optically-thick gas during the collapse initiated by the dissociation of molecular hydrogen (Boss 1989; Bonnell & Bate 1994; Bate 1998, 2002).

As a consequence of the dependence of close binary formation on dynamical exchange interactions and the accretion of material with high specific angular momentum, we find that close binaries tend not to have extreme mass ratios. All of our systems have mass ratios $`q\text{ }>0.3`$. Furthermore, the frequency of close binaries is dependent on mass in that massive stars are more likely to have close companions than lower mass stars. These properties are in good agreement with the results of observational surveys. At the end of our calculation, many of the close binaries are members of hierarchical triple systems. Although these systems may not yet have finished evolving, the implication is that many close binaries ought to have wider companions. Recent observations support this hypothesis, but larger surveys to determine the frequency of wide companions to close binary systems are necessary to demonstrate it conclusively.

## Acknowledgments

The computations reported here were performed using the UK Astrophysical Fluids Facility (UKAFF). VB acknowledges support by the European Community’s Research Training Network under contract HPRN-CT-2000-0155, Young Stellar Clusters.


format: text

text:
# Coreless vortex formation in a spinor Bose-Einstein condensate

## Abstract

Coreless vortices were phase-imprinted in a spinor Bose-Einstein condensate. The three-component order parameter of $`F=1`$ sodium condensates held in a Ioffe-Pritchard magnetic trap was manipulated by adiabatically reducing the magnetic bias field along the trap axis to zero. This distributed the condensate population across its three spin states and created a spin texture. Each spin state acquired a different phase winding which caused the spin components to separate radially.

Spin textures play a central role in describing the physics of elementary particles Skyrme (1961), liquid 3He-A Mermin and Ho (1976); Anderson and Toulouse (1977); Blaauwgeers et al. (2000), the quantum Hall effect Lee and Kane (1990), and gaseous Bose-Einstein condensates Ho (1998); Ohmi and Machida (1998); Yip (1999); Khawaja and Stoof (2001a, b). Topological defects vary between superfluid systems described by scalar and vector order parameters. In spin-less or spin-polarized condensates, line defects such as vortices have cores where the density of condensed particles is necessarily zero to keep the order parameter single-valued Yarmchuk et al. (1979); Madison et al. (2000); Abo-Shaeer et al. (2001). However, in condensates with an internal, spin degree-of-freedom, coreless vortices exist as spin textures Blaauwgeers et al. (2000); Matthews et al. (1999). Such structures are referred to as skyrmions (Anderson-Toulouse vortices Anderson and Toulouse (1977)) or merons (half-skyrmions, Mermin-Ho vortices Mermin and Ho (1976)) depending on the boundary conditions of the system.

In this Letter, we study spin textures in a Bose-Einstein condensate (BEC). Coreless vortices were created in $`F=1`$ spinor condensates held in a Ioffe-Pritchard magnetic trap by adiabatically reducing the magnetic bias field along the trap axis to zero. This continuously transformed the initially spin-polarized condensate into a coherent superposition of three spin states, each with a different phase winding. The resulting angular momentum per particle varied between spin states and the condensate evolved such that states with more angular momentum per particle circulated around states with less angular momentum per particle. Thus, the condensate had a net axial magnetization that varied with radial position. Previous work on vortices in a two-component system used laser, microwave, and radio frequency fields to spatially and temporally control the interconversion between components Matthews et al. (1999). However, without these applied fields the two-components evolved independently as distinguishable fluids. In our work, the spin states can freely interconvert at all points in space and time such that the spin texture would continually heal itself even in the presence of state-dependent losses.

In cylindrical coordinates, the spin$`F`$ condensate wavefunction can be written as $`|\mathrm{\Psi }(r,\varphi ,z)=\sqrt{n(r,\varphi ,z)}|\zeta (r,\varphi ,z)`$, where $`n`$ is the atomic number density and the $`2F+1`$ component spinor $`|\zeta =_{m_z=F}^F\zeta _{m_z}|F,m_z,|\zeta |\zeta |^2=1`$ describes a spin texture. A Ioffe-Pritchard magnetic trap consists of an axial bias field (with curvature) and a two-dimensional quadrupole field in the orthogonal plane Gott et al. (1962); Pritchard (1983):

$$\stackrel{}{B}(r,\varphi ,z)=B_z\widehat{z}+B^{}r(\mathrm{cos}(2\varphi )\widehat{r}\mathrm{sin}(2\varphi )\widehat{\varphi }),$$
(1)
where $`B^{}`$ is the radial magnetic field gradient and quadratic terms have been neglected. For a condensate of radial extent $`R`$ confined in a Ioffe-Pritchard magnetic trap with $`B_zB^{}R>0`$, $`|\zeta =|F,m_z=m_F`$, where $`m_z`$ and $`m_F`$ are the projection of the atomic spin along the z-axis and local magnetic field direction, respectively. Adiabatically ramping $`B_z`$ from $`B_zB^{}R>0`$ to zero rotates the atomic spin about the position-dependent axis $`\widehat{n}(\varphi )=\mathrm{sin}\varphi \widehat{x}+\mathrm{cos}\varphi \widehat{y}`$, and drives the transition $`|F,m_z=m_F_{m_z=|m_F|}^{|m_F|}\zeta _{m_z}\mathrm{exp}(i(m_zm_F)\varphi )|F,m_z`$ Isoshima et al. (2000); Leanhardt et al. (2002a). Thus, the condensate population is distributed across $`2|m_F|+1`$ spin states with each acquiring a different topological phase factor and angular momentum per particle due to the variation of Berry’s phase with magnetic quantum number Berry (1984).

The condensate remains in the state $`|F,m_F`$ with respect to the local magnetic field provided the local Zeeman energy, $`g_F\mu _B(B_z^2+(B^{}r)^2)^{1/2}`$, dominates the local kinetic energy associated with the spin texture, $`\mathrm{}^2/mr^2`$, where $`g_F`$ is the Landé g-factor, $`\mu _B`$ is the Bohr magneton, and $`m`$ is the atomic mass. For $`B_z=0`$, atomic spins aligned with the quadrupole magnetic field produce the planar spin texture in Fig. 1(a). However, the infinite kinetic energy associated with the wavefunction singularity at $`r=0`$ creates a non-planar spin texture over a disc of radius $`(\mathrm{}^2/mg_F\mu _BB^{})^{1/3}`$, with higher angular momentum spin states residing outside those with lower angular momentum.

Bose-Einstein condensates containing over $`10^7`$ 23Na atoms were created in the $`|1,1`$ state in a magnetic trap, captured in the focus of an optical tweezers laser beam, and transferred into an auxiliary “science” chamber as described in Ref. Gustavson et al. (2002). In the science chamber, the condensate was loaded into a microfabricated Ioffe-Pritchard magnetic trap formed by a Z-shaped 50 $`\mu `$$`\times `$ 10 $`\mu `$m electroplated copper wire carrying current $`I`$ and an external magnetic bias field, $`B_{}`$, as described in Ref. Leanhardt et al. (2002b). Typical wiretrap parameters were $`I=720`$ mA, $`B_{}=5.3`$ G, and $`B_z=1.3`$ G, resulting in a radial magnetic field gradient $`B^{}=180`$ G/cm. This produced axial and radial trap frequencies $`\omega _z=2\pi \times 4`$ Hz and $`\omega _{}=2\pi \times 250`$ Hz, respectively. Condensates in the wiretrap had $`2\times 10^6`$ atoms, a Thomas-Fermi radius of $`5\mu `$m, and a lifetime $`25`$ s Leanhardt et al. (2002c).

Coreless vortices imprinted onto the condensate wavefunction by adiabatically ramping $`B_z0`$ are shown in Fig. 2. To observe the nature of the spin texture, an axial bias field was switched on non-adiabatically along either the negative \[Fig. 2(a)-(d)\] or positive \[Fig. 2(e)-(h)\] z-axis. Switching the axial bias field on suddenly “froze” the atomic spins and effectively “projected” the condensate wavefunction onto a basis quantized with respect to the local (axial) magnetic field. This allowed the spin states to be separated by a magnetic field gradient applied during ballistic expansion. Switching the direction of the axial projection field exchanged the roles of the $`|1,1`$ and $`|1,+1`$ states. Figures 2(a) and 2(e) show the coreless nature of the vortices, while Figs. 2(b) and 2(f) show the concentric cylinder structure resulting from the competition between the atomic Zeeman energy and the kinetic energy of the rotating spin states. We assume that the two-dimensional ($`\omega _z\omega _{}`$) ballistic expansion process simply magnifies the condensate wavefunction, as it does in the expansion of a single-component condensate with vortices Lundh et al. (1998); Dalfovo and Modugno (2000).

To imprint the coreless vortices, the axial bias field was ramped linearly from $`B_z=1.3`$ G to $`B_z0`$ in 10 ms. Along the wiretrap axis, the magnetic field was

$$\stackrel{}{B}(r=0,\varphi ,z)=(B_z+\frac{1}{2}B^{\prime \prime }z^2)\widehat{z},$$
(2)
where quadratic terms neglected in Eq. (1) are included. The presence of axial magnetic field curvature, $`B^{\prime \prime }5`$ G/cm2, implies that the spin texture had a slight axial dependence. Ramping $`B_z0`$ quickly (10 ms) compared to the axial trap period (250 ms) compressed the condensate radially but left the axial extent of the condensate unchanged. Thus, the axial magnetic field variation, $`\mathrm{\Delta }B_z`$, was less than the radial magnetic field variation, $`\mathrm{\Delta }B_r`$, and the axial dependence of the spin texture was minimal. With $`B_z=1.3`$ G ($`B_z0`$ G), the condensate chemical potential was $`\mu (\mu _B/2)\times 3`$ mG ($`\mu (\mu _B/2)\times 27`$ mG) yielding $`\mathrm{\Delta }B_z3`$ mG ($`\mathrm{\Delta }B_r27`$ mG). The images shown in Fig. 2 integrated the atomic number density along the z-axis and therefore averaged over the minor axial variation to the spin texture.

To project the condensate wavefunction onto a basis quantized with respect to the local magnetic field, an axial bias field was switched on at a rate of $`\dot{B}_z=2\times 10^5`$ G/s along either the negative or positive z-axis. For $`0B_r\mathrm{\Delta }B_r`$, the Landau-Zener non-adiabatic transition probability, $`\mathrm{exp}(\pi \mu _BB_r^2/\mathrm{}\dot{B}_z)0.9`$, was sufficiently close to unity that the atomic spins remained “frozen” during the sudden application of the axial bias field and the spin texture could be accurately diagnosed. The magnetic trap was switched off 100 $`\mu `$s after applying the axial bias field allowing the atoms to expand ballistically. While the total condensate density monotonically decreased as a function of radial position \[Figs. 2(a) and 2(e)\], the density of each spin component peaked at a different radius signifying a variation in the angular momentum per particle between spin states \[Figs. 2(b) and 2(f)\].

With the application of a projection field along the positive z-axis, additional rings of atoms appeared in the $`|1,1`$ and $`|1,0`$ states \[Fig. 2(f)\]. Possible technical origins include the non-instantaneous turn on of the projection field and the ballistic expansion process. Also, ramping $`B_z0`$ in 10 ms was non-adiabatic for atoms near $`r=0`$ causing an atom loss of $`50\%`$ as observed in previous work Leanhardt et al. (2002a). If these atoms had not left the condensate before the projection field was applied, they may have contributed to the images displayed in Fig. 2(f). The additional rings of atoms may also correspond to a low energy, radial spin-wave excitation Ho (1998); Ohmi and Machida (1998). However, we could not identify any asymmetry between applying the projection field along the positive versus negative z-axis that would account for the presence of the extra rings in Fig. 2(f) but not in Fig. 2(b).

Engineering topological states in a Bose-Einstein condensate has received much theoretical attention Isoshima et al. (2000); Williams and Holland (1999); Ruostekoski (2000); Marzlin et al. (2000); Leonhardt and Volovik (2000); Ruostekoski and Anglin (2001). The evolution of a condensate confined in a Ioffe-Pritchard magnetic trap while ramping $`B_z0`$ is described by a position-dependent spin rotation about the $`\widehat{n}(\varphi )`$ axis through an angle $`\beta (r)`$ Leanhardt et al. (2002a)

$$|\zeta (r,\varphi ,z)=e^{i(\stackrel{}{}/\mathrm{})\widehat{n}(\varphi )\beta (r)}|\zeta _0,$$
(3)
where $`\stackrel{}{}`$ is the spin operator and $`|\zeta _0=|F,m_z=m_F`$ is a polarized spinor. For $`|\zeta _0=|1,1`$, Eq. (3) gives the condensate spinor in the laboratory frame as

$$|\zeta (r,\varphi ,z)=\begin{array}{c}\mathrm{cos}^2(\beta (r)/2)|1,1+\hfill \\ \frac{1}{\sqrt{2}}\mathrm{sin}(\beta (r))e^{i\varphi }|1,0+\hfill \\ \mathrm{sin}^2(\beta (r)/2)e^{i2\varphi }|1,+1.\hfill \end{array}$$
(4)
$`\beta (0)=0`$ and skyrmions (merons) are described by the boundary condition $`\beta (\mathrm{})=\pi `$ ($`\beta (\mathrm{})=\pi /2`$).

The radial dependence of $`\beta (r)`$ is determined by requiring the Gross-Pitaevskii energy functional,

$$E=d^3\stackrel{}{r}n\left(\frac{\mathrm{}^2}{2m}\zeta |\zeta +V+\left(\frac{c_0}{2}+\frac{c_2}{2}\frac{|\stackrel{}{F}|^2}{\mathrm{}^2}\right)n\right),$$
(5)
be stationary with respect to variations in $`\beta `$. In writing Eq. (5) we have made the Thomas-Fermi approximation, $`\zeta |\zeta =_{m_z=F}^F\zeta _{m_z}^{}\zeta _{m_z}`$, $`V=g_F\mu _B\stackrel{}{F}\stackrel{}{B}/\mathrm{}`$, $`\stackrel{}{F}=\zeta |\stackrel{}{}|\zeta `$, $`c_0=4\pi \mathrm{}^2\overline{a}/m`$, and $`c_2=4\pi \mathrm{}^2\mathrm{\Delta }a/m`$. Here $`\overline{a}=(2a_0+a_2)/3`$ and $`\mathrm{\Delta }a=(a_2a_0)/3`$ characterize two-body interactions, where $`a_0`$ and $`a_2`$ are scattering lengths for collisions with total angular momentum $`F=0`$ and $`F=2`$, respectively Ho (1998). For $`a_2>0`$ the atomic interactions are anti-ferromagnetic (polar), as in 23Na Stenger et al. (1998), while for $`a_2<0`$ they are ferromagnetic, as in 87Rb Klausen et al. (2001).

Using Eq. (4) as a trial spinor, we find that $`\beta (r)`$ satisfies the meron boundary conditions and varies from $`\beta (0)=0`$, due to the kinetic energy of the spin texture, to $`\beta (\mathrm{})=\pi /2`$, due to the atomic Zeeman energy, over a characteristic length scale given by the larger of $`(\mathrm{}^2/mg_F\mu _BB^{})^{1/3}`$ and $`|B_z|/B^{}`$, in reasonable agreement with the results presented in Fig. 2. We observed that the boundary condition $`\beta (0)=0`$ was maintained regardless of the sign of $`B_z`$, i.e. atomic spins along the trap axis always remained in the initial state. This emphasizes the non-trivial topology of the spin texture. If the atomic spins had simply followed the local magnetic field, $`\beta (r)`$ would satisfy $`\mathrm{tan}\beta (r)=B^{}r/B_z`$. Thus, scanning $`B_z`$ from slightly positive to slightly negative would instantaneous change the boundary condition at the origin from $`\beta (0)=0`$ to $`\beta (0)=\pi `$ and flip the atomic spins along the trap axis.

The trial spinor in Eq. (4) does not have the most general form since it was derived by rotating a polarized spinor \[Eq. (3)\] and is inherently ferromagnetic, $`|\stackrel{}{F}|=\mathrm{}`$. Accordingly, the spin-dependent interaction term in Eq. (5) does not vary with $`\beta `$ and therefore does not contribute to the determination of $`\beta (r)`$. This restriction is not severe since the Zeeman energy dominates the spin-dependent interaction energy in our experiment.

It should be possible to overlap the condensate with an optical dipole trap so that, after ramping $`B_z0`$, relaxing the radial magnetic field gradient would allow for the spin-dependent interactions to determine the evolution of the condensate. It is predicted that skyrmions/merons in condensates with anti-ferromagnetic (ferromagnetic) interactions are unstable (stable) Khawaja and Stoof (2001a, b); Battye et al. (2002); Mizushima et al. (2002a); Isoshima and Machida (2002); Zhang et al. (2002); Martikainen et al. (2002); Mizushima et al. (2002b). At present, excitations created during the field ramping process have prevented us from characterizing the stability of coreless vortices in an anti-ferromagnetic 23Na condensate. However, we were able to imprint the spin texture in the presence of the optical dipole trap, as well as reproduce the results of Ref. Leanhardt et al. (2002a) by fully inverting the axial bias field to produce vortices with a $`4\pi `$ phase winding in an optical dipole trap. At zero magnetic field, it may then be possible to observe multiply charged vortices in a spinor condensate “unwind” themselves as predicted in Ref. Ho (1998). The unwinding process is precisely the position-dependent spin rotation demonstrated in previous work Leanhardt et al. (2002a), with the spin texture investigated here as an intermediate state Martikainen et al. (2002).

In conclusion, we have demonstrated a robust technique for creating coreless vortices in a Bose-Einstein condensate. Our technique can be extended to generate spin textures with arbitrary winding number and variable angular momentum per particle by using higher-order, axisymmetric multipole magnetic fields and condensates with different spin. This work opens up the opportunity to study the stability of topological defects in spinor Bose-Einstein condensates.

We are grateful to K. Machida for bringing the spin texture studied here to our attention. We thank T. Pasquini for experimental assistance and J.R. Anglin and M. Crescimanno for valuable discussions. This work was funded by ONR, NSF, ARO, NASA, and the David and Lucile Packard Foundation.


text:
# Resonant spin current in nanotube double junctions

## Abstract

Zero bias conductance per spin of nanotube double junction (NTDJ) is investigated theoretically using the tight binding model, unrestricted Hartree-Fock approximation and non-equilibrium Green’s functions. NTDJ consists of two metallic nanotubes joined by a piece of semiconducting nanotube, with the transition between the nanotubes made up of sets of 5 and 7 member carbon rings. A quantum well forms in the central semiconducting NT region, bounded by Schottky barriers. Spin current occurs when Coulomb interactions raise the spin degeneracy of resonant levels in the quantum well. As long as an appropriate semiconducting NT length is chosen, spin direction can be controlled by gate voltage, i.e., NTDJ functions as a nano spin filter.

The combination of nanotechnology based on carbon nanotubes (NTs)  and spintronics  is thought to be extremely promising for future technological innovations. Tsukagoshi et al. have demonstrated efficient spin injection from cobalt electrodes in NTs. However, diffusion of magnetic atoms from the electrode into the NT results in an increase in uncontrollable spin flips due to the large spin-orbit interactions of ferromagnetic atoms. Since long spin relaxation times are needed, spin filters composed of atoms with small spin-orbit interactions are desirable. In this Letter, a spin filter composed entirely of carbon atoms, arranged as a semiconducting NT between two metallic NTs is proposed. Schottky barriers at the interfaces form a quantum well within the semiconducting NT, while the metallic NTs function as leads. When Coulomb interactions lift spin degeneracy in the quantum well, resonant current becomes spin polarized. To the best of the Authors’ knowledge, this is the first study into such resonant spin current. A possible method of constructing the NT junction is to allow open ends of $`(3m,0)`$ zigzag NT and $`(3m\pm 1,0)`$ NT to react; the NTs then become connected by adjacent regions of 5 and 7 member carbon rings, while still maintaining an sp2 network. This structure is referred to as the NT double junction (NTDJ). In order to examine spin current in such devices, zero bias conductance per spin is calculated for (9,0)-(8,0)-(9,0) NTDJ. Results are closely related to the fourfold degeneracy of the central semiconducting NT.

Using the unrestricted Hartree-Fock approximation, both the Hamiltonian $`H^\sigma `$ and the density matrix $`\rho ^\sigma `$ are dependent on spin $`\sigma =,`$, and can be expressed using the tight-binding approximation limited to the $`\pi `$ orbitals as

$$H_{i,j}^\sigma =t_{i,j}\stackrel{~}{U}_{i,j}\rho _{i,j}^\sigma +\delta _{i,j}\underset{k,\sigma }{}\stackrel{~}{U}_{i,k}(\rho _{k,k}^\sigma (1/2)),$$
(1)
where $`t_{i,j}=t<0`$ when $`i`$ and $`j`$ are nearest neighbors, $`t_{i,j}=0`$ otherwise. NTDJ is divided into the two asymptotic ’lead’ regions $`s`$, $`d`$, and a scattering region $`C`$, as shown in Figure 1, with Hamiltonian matrix elements assigned to one of these regions. The lead regions are approximated as having a periodic structure based on the charges and Hamiltonian elements in the NT unit cells at the edges of region $`C`$, as indicated by rectangles in Fig.1. The charges contribute to diagonal terms of $`H_c`$, while Hamiltonian elements determine $`H_{s,d}`$ and $`\tau _{s,d}`$. For simplicity, only nearest neighbor elements are treated for off diagonal elements of $`H_{s,d}`$ and $`\tau _{s,d}`$. $`s`$ region, $`d`$ region and gate electrode voltages are referred to as $`V_s`$, $`V_d`$ and $`V_g`$, respectively. Since the gate electrode is a planar perfect conductor, and $`V_g=0`$ in our model, Coulomb interaction can be expressed as $`\stackrel{~}{U}_{i,j}=U(|\stackrel{}{r}_i\stackrel{}{r}_j|)U(|\stackrel{}{r}_i\stackrel{}{R}_j|),`$ where $`\stackrel{}{R}_j`$ is the image of $`\stackrel{}{r}_j`$. The ’bare’ Coulomb kernel without a gate is assumed to be

$$U(r)=2.5t/\sqrt{1+16(r/a)^2}.$$
(2)
Voltages can be expressed in another way, $`\stackrel{~}{V}_s=(V_sV_d)/2`$, $`\stackrel{~}{V}_d=\stackrel{~}{V}_s`$ and $`\stackrel{~}{V}_g=(V_s+V_d)/2`$, since the difference between, not the absolute value of, voltages is of importance. In this letter, relative voltage representation is used and negative $`\stackrel{~}{V}_g`$ causing electrostatic hole doping is considered.

The density matrix $`\rho `$ is calculated from the retarded Green’s function $`G`$, self energy $`\mathrm{\Lambda }_p^\sigma \tau _p^\sigma (E+i\delta H_p^\sigma )^1\tau _p^\sigma `$ and its imaginary component $`\mathrm{\Gamma }_p^\sigma =i(\mathrm{\Lambda }_p^\sigma \mathrm{\Lambda }_p^\sigma )`$, giving

$`G^\sigma =(EH_c^\sigma \mathrm{\Lambda }_s^\sigma \mathrm{\Lambda }_d^\sigma )^1,`$ (3)
$`\rho ^\sigma ={\displaystyle \frac{1}{\pi }}{\displaystyle _{\mathrm{}}^{\mathrm{}}}𝑑EG^\sigma (\mathrm{\Gamma }_s^\sigma +\mathrm{\Gamma }_d^\sigma )G^\sigma f_{\mathrm{FD}},`$ (4)
where $`f_{\mathrm{FD}}=1/(\mathrm{exp}(\beta (EE_F))+1)`$ and $`k_BT=1/\beta =0.01t`$.  To calculate zero bias conductance, $`E_F=eV_s=eV_d=|e|\stackrel{~}{V}_g`$. Using Eqs.(1) and (4), $`H`$ and $`\rho `$ can be calculated self-consistently for a given $`\stackrel{~}{V}_g`$. When self-consistency has been achieved, zero bias conductance for each spin $`dI^\sigma /dV_{sd}`$ is calculated from

$$\frac{dI^\sigma }{dV_{sd}}=\frac{\beta e^2}{h}_{\mathrm{}}^{\mathrm{}}𝑑E\mathrm{Tr}[\mathrm{\Gamma }_s^\sigma G^\sigma \mathrm{\Gamma }_d^\sigma G^\sigma ]f_{\mathrm{FD}}(1f_{\mathrm{FD}}).$$
(5)
For small values of $`|\stackrel{~}{V}_g|`$, motion of holes into (8,0) NT is suppressed by the band-gap, resulting in more positively charged (9,0) NT than (8,0) NT and inducing Schottky barriers. Fig.2 shows potential energy as a function of position along the NTDJ for different spins and different values of $`\stackrel{~}{V}_g`$. It can be seen that as $`\stackrel{~}{V}_g`$ becomes more negative, the (8,0) NT region becomes electrostatically doped, lowering the Schottky barrier. Neighboring peaks and dips seen in the transition region between (9,0) NT and (8,0) NT are caused by an excess of electrons in the 5 member carbon rings, and a deficiency in the 7 member carbon rings . Note that while the gate electrode covers both $`C`$ and lead regions, and is much larger than a nanometer, it still exerts control over the atomistic spatial profile of the barrier.

Conductance per spin as a function of gate voltage $`\stackrel{~}{V}_g`$ is shown in Figs.3, 4 and 5 for $`L=`$13, 7 and 6 respectively, where $`L`$ is the length of (8,0) NT in units of $`3a`$. The first resonant conductance peak appears at the gate threshold, $`\stackrel{~}{V}_{gth}`$, with only a small tunneling current present when $`\stackrel{~}{V}_g>\stackrel{~}{V}_{gth}`$. For smaller values of $`L`$, enhanced tunneling effects increase hole density in the (8,0) NT region, making $`\stackrel{~}{V}_{gth}`$ more negative.

Fig.3 also shows the four eigenvalues of $`H_c`$ near the valence band edge, denoted by $`ϵ_{i,\sigma }`$, where $`i=\pm `$ is $`K`$ corner point degeneracy. These can be represented as

$`ϵ_{i,\sigma }`$ $`=`$ $`ϵ_i^{(0)}+J_{(+,)}n_{i,\sigma }`$ (7)
$`+(J_{(+,)}K)n_{i,\sigma }+J_{(i,i)}n_{i,\sigma },`$
where $`J`$, $`K`$, $`n_{i,\sigma }`$ and $`ϵ_{i,\sigma }^{(0)}`$ are the Coulomb integral, exchange integral, occupation number and eigenvalue with no Coulomb interaction, respectively. Unlike finite systems or single electron tunneling, $`n_{i,\sigma }`$ can take continuous values over the range 0 to 1. From Fig.3, it can be seen that large splitting between the levels of the different spin states only occurs when $`ϵ_{i,\sigma }`$ crosses the Fermi level $`E_F=|e|\stackrel{~}{V}_g`$, i.e., when occupation numbers change. This indicates that splitting is a function of $`Jn_{i,\sigma }`$ and $`Kn_{i,\sigma }`$, not of $`ϵ^{(0)}`$. Since the transmission rate $`T(E)=\mathrm{Tr}[\mathrm{\Gamma }_s^\sigma G^\sigma \mathrm{\Gamma }_d^\sigma G^\sigma ]`$ exhibits resonant peaks when $`Eϵ_{i,\sigma }`$, as $`|\stackrel{~}{V}_g|`$ increases two spin down peaks appear first, followed by two spin up peaks. The pairing of peaks with common spin is due to Hund’s rule, which holds when $`K>ϵ_+^{(0)}ϵ_{}^{(0)}`$.

For a given $`L`$, barrier thickness $`D`$ decreases as hole number $`N=4_{i,\sigma }n_{i,\sigma }`$ increases, as seen in Fig.2, broadening the energy of the resonant state, $`W`$. At the same time, $`J_{(i,j)}`$ decreases due to an increase in spatial width of the quantum well, $`L2D`$. The net effect is an increase in $`WJ_{(i,j)}`$ resulting in a decrease in spin polarization. To examine spin polarization effects, total energy $`E_{tot}`$ is approximated by

$`E_{tot}`$ $`=`$ $`{\displaystyle \underset{i,\sigma }{}}\left(ϵ_i^{(0)}{\displaystyle \frac{W}{2}}+{\displaystyle \frac{W}{2}}n_{i,\sigma }\right)n_{i,\sigma }+{\displaystyle \underset{i}{}}J_{(i,i)}n_{i,}n_{i,}`$ (9)
$`{\displaystyle \underset{\sigma }{}}Kn_{+,\sigma }n_{,\sigma }+{\displaystyle \underset{\sigma ,\sigma ^{}}{}}J_{(+,)}n_{+,\sigma }n_{,\sigma ^{}}`$
Equation (9) is derived from

$$\mathrm{\Delta }E_{tot}=\mathrm{\Delta }n_{i,\sigma }(W(n_{i,\sigma }1/2)+ϵ_{i,\sigma }),$$
(10)
where the energy region $`W/2<Eϵ_{i,\sigma }<+W(n_{i,\sigma }1/2)`$ is occupied so that $`n_{i,\sigma }`$ changes at the energy $`E=ϵ_{i,\sigma }+W(n_{i,\sigma }1/2)`$.

For $`0<N<1`$, $`n_{+,}`$ and $`n_{+,}`$ initially decrease, because the degeneracy of $`ϵ^{(0)}`$ is slightly lifted by presence of pentagon-heptagon pairs making $`ϵ_{}^{(0)}<ϵ_+^{(0)}`$. Therefore $`n_{+,}=1\alpha `$ ,$`n_{+,}=1N+\alpha `$ and $`n_,=n_,=1`$. Thus, $`E_{tot}=(J_{(+,+)}W)(N\alpha )\alpha +C,`$ with $`C`$ independent of $`\alpha `$. For a given $`N`$, we define $`\alpha _{min}`$ as the value of $`\alpha `$ which minimizes $`E_{tot}`$. Thus, spin current occurs when $`\alpha _{min}=0`$. Although $`\alpha _{min}=N`$ is another solution, infinitesimal Zeeman energy favors $`\alpha _{min}=0`$. Thus the condition allowing spin down current at the first resonant peak $`(N=1)`$ is $`W|_{N=1}<J_{(+,+)}|_{N=1}`$. This condition is only satisfied when $`L>4`$; since when $`L`$ is too small, large tunneling effects mean that $`WJ`$ is positive, even without electrostatic doping. As $`L`$ increases, the value of $`N`$ over which polarization occurs also increases. This can be understood by noticing that increases in $`WJ_{(i,j)}`$ per $`N`$ become smaller as $`L`$ becomes larger, because increases in hole density in (8,0) NT also decrease.

When $`L>4`$ and $`1<N<2`$,$`(n_,,n_{+,};n_,,n_{+,})=(1,1\alpha ,1N_2+\alpha ,0)`$, where $`N_2=N1`$ and $`ϵ_,<ϵ_,`$ on account of Hund’s rule. Then

$`E_{tot}`$ $`=`$ $`(WJ_{(+,)})\left(\alpha {\displaystyle \frac{N_2}{2}}+{\displaystyle \frac{\gamma }{WJ_{(+,)}}}\right)^2+C`$ (11)
$`\gamma `$ $`=`$ $`(ϵ_{}^{(0)}ϵ_+^{(0)}+K+J_{(,)}J_{(+,)})/2.`$ (12)
When $`\gamma `$ is assumed to be positive and nearly constant, level configuration depends on the sign of $`WJ_{(+,)}(2\gamma /N_2)`$. For a negative sign, holes are only injected into spin down levels, i.e., $`\alpha _{min}=0`$. Otherwise holes are injected into both spin levels, i.e., $`\alpha _{min}=N_2/2\gamma /(WJ_{(+,)})`$. When $`L=5`$ or 6, $`WJ_{(+,)}(2\gamma /N_2)`$ increases, reversing sign as $`N_2`$ increases from 0 to 1. The increase in $`W`$ can be seen in the broadening of the peak corresponding to $`(,)`$ in Fig.5. As a result, spin up peaks overlap the second spin down peak, largely countering spin up polarization. For $`L>6`$, on the other hand, $`WJ_{(+,)}(2\gamma /N_2)`$ remains negative, and so the second spin down peak separates from the spin up peaks.

In addition to the lower limits of $`L`$ under which spin polarization can occur, upper limits also exist. As $`L`$ increases for fixed $`N`$, both thickness of the barrier $`D`$ and quantum well width $`L2D`$ increase. The former reduces $`W`$ while the latter causes a decrease in the Coulomb integral $`J_{(i,j)}`$. When $`J_{(i,j)}`$ is less than either $`W`$ or $`k_BT`$, spin polarization disappears. Because $`W`$ decreases for larger $`L`$, $`k_BT`$ becomes a more important factor in determining the upper limit of $`L`$ than it was for the lower limit. Figure 3 shows that at room temperature, the upper limit is in excess of $`39a`$.

In this letter, it has been shown that the nanotube double junction (NTDJ) functions as a nano spin filter without the use of magnetic atoms. An infinitesimal magnetic field is necessary only to determine which spin becomes unoccupied first. Conditions for spin current are: (1) Fermi level is near the band gap edges at operating gate voltage $`\stackrel{~}{V}_g`$; (2) length of the semiconducting NT region, $`L`$, is larger than a particular threshold length. In fact, there are two threshold lengths, $`L_1`$ and $`L_2`$. For $`L>L_2`$, the direction of spin can be controlled by $`\stackrel{~}{V}_g`$. For $`L_2L>L_1`$, however, spin polarization is limited to a single direction. Furthermore there is no spin current when $`L_1L`$. In the system studied, these thresholds are $`L_1=4`$ and $`L_2=6`$. Although the Coulomb interaction $`U(r)`$ of eq.(2) was chosen to be slightly smaller than estimated in Ref.  to achieve convergence of the self-consistent loop, the Author believes that the parameters used are realistic enough and results valid enough to prove that NTDJ is a useful device that combines nanotechnology and spintronics.

The author wishes to thank Prof. Masaru Tsukada for valuable advice. He also acknowledges Dr. Keiji Ohno for helpful discussions. This work was supported in part by a Grant-in-Aid for Creative Scientific Research on ”Devices on molecular and DNA levels” (No. 13GS0017) from the Japan Society for the Promotion of Science.


text:
# Implications of the KamLAND Measurement on the Lepton Flavor Mixing Matrix and the Neutrino Mass Matrix

## I Introduction

The KamLAND experiment  turns to confirm the large-mixing-angle (LMA) Mikheyev-Smirnov-Wolfenstein (MSW) solution  to the long-standing solar neutrino problem. In addition, the K2K long-baseline experiment  has unambiguously observed a reduction of $`\nu _\mu `$ flux and a distortion of the energy spectrum. These new measurements, together with the compelling SNO evidence  for the flavor conversion of solar $`\nu _e`$ neutrinos and the Super-Kamiokande evidecne  for the deficit of atmospheric $`\nu _\mu `$ neutrinos, convinces us that the hypothesis of neutrino oscillations is actually correct! We are then led to the conclusion that neutrinos are massive and lepton flavors are mixed.

The mixing of lepton flavors means a mismatch between neutrino mass eigenstates $`(\nu _1,\nu _2,\nu _3)`$ and neutrino flavor eigenstates $`(\nu _e,\nu _\mu ,\nu _\tau )`$ in the basis where the charged lepton mass matrix is diagonal:

$$\left(\begin{array}{c}\nu _e\\ \nu _\mu \\ \nu _\tau \end{array}\right)=\left(\begin{array}{ccc}V_{e1}& V_{e2}& V_{e3}\\ V_{\mu 1}& V_{\mu 2}& V_{\mu 3}\\ V_{\tau 1}& V_{\tau 2}& V_{\tau 3}\end{array}\right)\left(\begin{array}{c}\nu _1\\ \nu _2\\ \nu _3\end{array}\right).$$
(1)
The matrix elements $`|V_{e1}|`$, $`|V_{e2}|`$, $`|V_{e3}|`$ and $`|V_{\mu 3}|`$ can simply be related to the mixing factors of solar , atmospheric  and CHOOZ reactor  neutrino oscillations in the following way:

$`\mathrm{sin}^22\theta _{\mathrm{sun}}`$ $`=`$ $`4|V_{e1}|^2|V_{e2}|^2,`$ (2)
$`\mathrm{sin}^22\theta _{\mathrm{atm}}`$ $`=`$ $`4|V_{\mu 3}|^2\left(1|V_{\mu 3}|^2\right),`$ (3)
$`\mathrm{sin}^22\theta _{\mathrm{chz}}`$ $`=`$ $`4|V_{e3}|^2\left(1|V_{e3}|^2\right).`$ (4)
Taking account of the unitarity of $`V`$, one may reversely express $`|V_{e1}|`$, $`|V_{e2}|`$, $`|V_{e3}|`$, $`|V_{\mu 3}|`$ and $`|V_{\tau 3}|`$ in terms of $`\theta _{\mathrm{sun}}`$, $`\theta _{\mathrm{atm}}`$ and $`\theta _{\mathrm{chz}}`$ :

$`|V_{e1}|`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{2}}}\sqrt{\mathrm{cos}^2\theta _{\mathrm{chz}}+\sqrt{\mathrm{cos}^4\theta _{\mathrm{chz}}\mathrm{sin}^22\theta _{\mathrm{sun}}}},`$ (5)
$`|V_{e2}|`$ $`=`$ $`{\displaystyle \frac{1}{\sqrt{2}}}\sqrt{\mathrm{cos}^2\theta _{\mathrm{chz}}\sqrt{\mathrm{cos}^4\theta _{\mathrm{chz}}\mathrm{sin}^22\theta _{\mathrm{sun}}}},`$ (6)
$`|V_{e3}|`$ $`=`$ $`\mathrm{sin}\theta _{\mathrm{chz}},`$ (7)
$`|V_{\mu 3}|`$ $`=`$ $`\mathrm{sin}\theta _{\mathrm{atm}},`$ (8)
$`|V_{\tau 3}|`$ $`=`$ $`\sqrt{\mathrm{cos}^2\theta _{\mathrm{chz}}\mathrm{sin}^2\theta _{\mathrm{atm}}}.`$ (9)
Current experimental information on $`\theta _{\mathrm{sun}}`$, $`\theta _{\mathrm{atm}}`$ and $`\theta _{\mathrm{chz}}`$ allows us to get very instructive constraints on the lepton flavor mixing matrix $`V`$. One purpose of this paper is therefore to examine how accurately we can recast $`V`$ from the present KamLAND, K2K, SNO, Super-Kamiokande and CHOOZ measurements.

Another purpose of this paper is to confront two-zero textures of the neutrino mass matrix with the new KamLAND data, so as to single out the most favorable texture(s) in phenomenology. In the flavor basis chosen above, the Majorana neutrino mass matrix can be written as

$$M=V\left(\begin{array}{ccc}m_1& 0& 0\\ 0& m_2& 0\\ 0& 0& m_3\end{array}\right)V^\mathrm{T},$$
(10)
where $`m_i`$ (for $`i=1,2,3`$) are physical masses of three neutrinos. In the assumption of the LMA solution for solar neutrino oscillations, a classification of $`M`$ with two vanishing entries has been done in an analytically approximate way . Our present analysis is different from the previous ones in two important aspects: (1) we carry out a careful numerical analysis of every two-zero pattern of the neutrino mass matrix $`M`$ to pin down its complete parameter space, because simple analytical approximations are sometimes unable to reveal the whole regions of relevant parameters allowed by new experimental data; (2) we present the quantitative predictions for allowed ranges of the absolute neutrino masses, the Majorana phases of $`CP`$ violation, and the effective masses of the tritium beta decay ($`m_e`$) and neutrinoless double beta decay ($`m_{ee}`$).

The remaining part of this paper is organized as follows. With the help of current data, we derive the model-independent constraints on nine elements of the lepton flavor mixing matrix in section II. Section III is devoted to a detailed analysis of the parameter space for every two-zero texture of the neutrino mass matrix. We obtain instructive predictions for the neutrino mass spectrum and Majorana phases of $`CP`$ violation as well as $`m_e`$ and $`m_{ee}`$ in section IV. Finally a brief summary is given in section V.

## II Constraints on the lepton flavor mixing matrix

As already shown in Eq. (3), five matrix elements of $`V`$ can be determined or constrained from current experimental data. The other four matrix elements ($`|V_{\mu 1}|`$, $`|V_{\mu 2}|`$, $`|V_{\tau 1}|`$ and $`|V_{\tau 2}|`$) are entirely unrestricted, however, unless one of them or the rephasing invariant of $`CP`$ violation of $`V`$ (defined as $`J`$ ) is measured. A realistic way to get rough but useful constraints on those four unknown elements is to allow the Dirac phase of $`CP`$ violation in $`V`$ to vary between 0 and $`\pi `$ , such that one can find out the maximal and minimal magnitudes of each matrix element. To see this point more clearly, we adopt the standard parametrization $`V=UP`$ , where

$$U=\left(\begin{array}{ccc}c_xc_z& s_xc_z& s_z\\ c_xs_ys_zs_xc_ye^{i\delta }& s_xs_ys_z+c_xc_ye^{i\delta }& s_yc_z\\ c_xc_ys_z+s_xs_ye^{i\delta }& s_xc_ys_zc_xs_ye^{i\delta }& c_yc_z\end{array}\right)$$
(11)
with $`s_x\mathrm{sin}\theta _x`$, $`c_x\mathrm{cos}\theta _x`$, and so on; and

$$P=\left(\begin{array}{ccc}e^{i\rho }& 0& 0\\ 0& e^{i\sigma }& 0\\ 0& 0& 1\end{array}\right).$$
(12)
The advantage of this representation is that the neutrinoless double beta decay is associated with the Majorana phases $`\rho `$ and $`\sigma `$, while $`CP`$ violation in normal neutrino oscillations depends separately on the Dirac phase $`\delta `$. Note that three mixing angles $`(\theta _x,\theta _y,\theta _z)`$, which are all arranged to lie in the first quadrant, can be written as

$`\mathrm{tan}\theta _x`$ $`=`$ $`{\displaystyle \frac{|V_{e2}|}{|V_{e1}|}},`$ (13)
$`\mathrm{tan}\theta _y`$ $`=`$ $`{\displaystyle \frac{|V_{\mu 3}|}{|V_{\tau 3}|}},`$ (14)
$`\mathrm{sin}\theta _z`$ $`=`$ $`|V_{e3}|.`$ (15)
It is then straightforward to obtain

$`|V_{\mu 1}|`$ $`=`$ $`{\displaystyle \frac{V_{e2}V_{\tau 3}|+|V_{e1}||V_{e3}||V_{\mu 3}|e^{i\delta }|}{1|V_{e3}|^2}},`$ (16)
$`|V_{\mu 2}|`$ $`=`$ $`{\displaystyle \frac{V_{e1}V_{\tau 3}||V_{e2}||V_{e3}||V_{\mu 3}|e^{i\delta }|}{1|V_{e3}|^2}},`$ (17)
$`|V_{\tau 1}|`$ $`=`$ $`{\displaystyle \frac{V_{e2}V_{\mu 3}||V_{e1}||V_{e3}||V_{\tau 3}|e^{i\delta }|}{1|V_{e3}|^2}},`$ (18)
$`|V_{\tau 2}|`$ $`=`$ $`{\displaystyle \frac{V_{e1}V_{\mu 3}|+|V_{e2}||V_{e3}||V_{\tau 3}|e^{i\delta }|}{1|V_{e3}|^2}}.`$ (19)
Varying the Dirac phase $`\delta `$ from 0 to $`\pi `$, we are led to the most generous ranges of $`|V_{\mu 1}|`$, $`|V_{\mu 2}|`$, $`|V_{\tau 1}|`$ and $`|V_{\tau 2}|`$:

$`{\displaystyle \frac{|V_{e2}||V_{\tau 3}||V_{e1}||V_{e3}||V_{\mu 3}|}{1|V_{e3}|^2}}`$ $``$ $`|V_{\mu 1}|{\displaystyle \frac{|V_{e2}||V_{\tau 3}|+|V_{e1}||V_{e3}||V_{\mu 3}|}{1|V_{e3}|^2}},`$ (20)
$`{\displaystyle \frac{|V_{e1}||V_{\tau 3}||V_{e2}||V_{e3}||V_{\mu 3}|}{1|V_{e3}|^2}}`$ $``$ $`|V_{\mu 2}|{\displaystyle \frac{|V_{e1}||V_{\tau 3}|+|V_{e2}||V_{e3}||V_{\mu 3}|}{1|V_{e3}|^2}},`$ (21)
$`{\displaystyle \frac{|V_{e2}||V_{\mu 3}||V_{e1}||V_{e3}||V_{\tau 3}|}{1|V_{e3}|^2}}`$ $``$ $`|V_{\tau 1}|{\displaystyle \frac{|V_{e2}||V_{\mu 3}|+|V_{e1}||V_{e3}||V_{\tau 3}|}{1|V_{e3}|^2}},`$ (22)
$`{\displaystyle \frac{|V_{e1}||V_{\mu 3}||V_{e2}||V_{e3}||V_{\tau 3}|}{1|V_{e3}|^2}}`$ $``$ $`|V_{\tau 2}|{\displaystyle \frac{|V_{e1}||V_{\mu 3}|+|V_{e2}||V_{e3}||V_{\tau 3}|}{1|V_{e3}|^2}}.`$ (23)
Note that the lower and upper bounds of each matrix element turn to coincide with each other in the limit $`|V_{e3}|0`$. Because of the smallness of $`|V_{e3}|`$, the ranges obtained in Eq. (9) should be quite restrictive. Hence it makes sense to recast the lepton flavor mixing matrix even in the absence of any experimental information on $`CP`$ violation.

In view of the present experimental data from KamLAND , K2K , SNO , Super-Kamiokande  and CHOOZ , we have $`0.25\mathrm{sin}^2\theta _{\mathrm{sun}}0.40`$ , $`0.92<\mathrm{sin}^22\theta _{\mathrm{atm}}1.0`$, and $`0\mathrm{sin}^22\theta _{\mathrm{chz}}<0.1`$ at the $`90\%`$ confidence level. Namely,

$`30.0^{}`$ $``$ $`\theta _{\mathrm{sun}}\mathrm{\hspace{0.33em}39.2}^{},`$ (24)
$`36.8^{}`$ $`<`$ $`\theta _{\mathrm{atm}}<\mathrm{\hspace{0.33em}53.2}^{},`$ (25)
$`0^{}`$ $``$ $`\theta _{\mathrm{chz}}<\mathrm{\hspace{0.33em}9.2}^{}.`$ (26)
Using these inputs, we calculate the numerical ranges of $`|V_{e1}|`$, $`|V_{e2}|`$, $`|V_{e3}|`$, $`|V_{\mu 3}|`$ and $`|V_{\tau 3}|`$ from Eq. (3). Then the allowed ranges of $`|V_{\mu 1}|`$, $`|V_{\mu 2}|`$, $`|V_{\tau 1}|`$ and $`|V_{\tau 2}|`$ can be found with the help of Eq. (9). Our numerical results are summarized as

$$|V|=\left(\begin{array}{ccc}0.700.87& 0.500.69& <0.16\\ 0.200.61& 0.340.73& 0.600.80\\ 0.210.63& 0.360.74& 0.580.80\end{array}\right).$$
(27)
This result is certainly more restrictive than that obtained in Ref.  before the KamLAND measurement.

Note that the rephasing invariant of $`CP`$ violation reads as follows:

$`J`$ $`=`$ $`{\displaystyle \frac{|V_{e1}||V_{e2}||V_{e3}||V_{\mu 3}||V_{\tau 3}|}{1|V_{e3}|^2}}\mathrm{sin}\delta `$ (28)
$`=`$ $`{\displaystyle \frac{\sqrt{\mathrm{sin}^22\theta _{\mathrm{sun}}\left(\mathrm{sin}^22\theta _{\mathrm{atm}}+4\mathrm{sin}^2\theta _{\mathrm{atm}}\mathrm{sin}^2\theta _{\mathrm{chz}}\right)}}{4\mathrm{cos}^2\theta _{\mathrm{chz}}}}\mathrm{sin}\theta _{\mathrm{chz}}\mathrm{sin}\delta .`$ (29)
The term proportional to $`4\mathrm{sin}^2\theta _{\mathrm{atm}}\mathrm{sin}^2\theta _{\mathrm{chz}}`$ in $`J`$, which may correct the leading term up to $`5\%`$ (for $`\theta _{\mathrm{atm}}=45^{}`$ and $`\theta _{\mathrm{chz}}=9^{}`$), were not taken into account in Ref. . By use of Eq. (10), we find $`J0.039\mathrm{sin}\delta `$. This result implies that the magnitude of $`J`$ can maximally be 0.039, leading probably to observable $`CP`$-violating effects in long-baseline neutrino oscillations.

## III Two-zero textures of the neutrino mass matrix

The symmetric neutrino mass matrix $`M`$ totally has six independent complex entries. If two of them vanish, i.e., $`M_{ab}=M_{pq}=0`$, we obtain two constraint equations:

$`m_1U_{a1}U_{b1}e^{2i\rho }+m_2U_{a2}U_{b2}e^{2i\sigma }+m_3U_{a3}U_{b3}`$ $`=`$ $`0,`$ (30)
$`m_1U_{p1}U_{q1}e^{2i\rho }+m_2U_{p2}U_{q2}e^{2i\sigma }+m_3U_{p3}U_{q3}`$ $`=`$ $`0,`$ (31)
where $`a`$, $`b`$, $`p`$ and $`q`$ run over $`e`$, $`\mu `$ and $`\tau `$, but $`(p,q)(a,b)`$. Solving Eq. (12), we arrive at

$`{\displaystyle \frac{m_1}{m_3}}e^{2i\rho }`$ $`=`$ $`{\displaystyle \frac{U_{a3}U_{b3}U_{p2}U_{q2}U_{a2}U_{b2}U_{p3}U_{q3}}{U_{a2}U_{b2}U_{p1}U_{q1}U_{a1}U_{b1}U_{p2}U_{q2}}},`$ (32)
$`{\displaystyle \frac{m_2}{m_3}}e^{2i\sigma }`$ $`=`$ $`{\displaystyle \frac{U_{a1}U_{b1}U_{p3}U_{q3}U_{a3}U_{b3}U_{p1}U_{q1}}{U_{a2}U_{b2}U_{p1}U_{q1}U_{a1}U_{b1}U_{p2}U_{q2}}}.`$ (33)
This result implies that two neutrino mass ratios $`(m_1/m_3,m_2/m_3)`$ and two Majorana-type $`CP`$-violating phases $`(\rho ,\sigma )`$ can fully be determined in terms of three mixing angles $`(\theta _x,\theta _y,\theta _z)`$ and the Dirac-type $`CP`$-violating phase $`(\delta )`$. Thus one may examine whether a two-zero texture of $`M`$ is empirically acceptable or not by comparing its prediction for the ratio of two neutrino mass-squared differences with the result extracted from current experimental data on solar and atmospheric neutrino oscillations:

$$R_\nu \frac{\left|m_2^2m_1^2\right|}{\left|m_3^2m_2^2\right|}\frac{\mathrm{\Delta }m_{\mathrm{sun}}^2}{\mathrm{\Delta }m_{\mathrm{atm}}^2}.$$
(34)
Considering the LMA MSW solution confirmed by the KamLAND measurement, we have $`5.9\times 10^5\mathrm{eV}^2\mathrm{\Delta }m_{\mathrm{sun}}^28.8\times 10^5\mathrm{eV}^2`$  at the $`90\%`$ confidence level. In addition, we have $`1.6\times 10^3\mathrm{eV}^2\mathrm{\Delta }m_{\mathrm{atm}}^23.9\times 10^3\mathrm{eV}^2`$  at the $`90\%`$ confidence level. Thus we arrive at $`1.5\times 10^2R_\nu 5.5\times 10^2`$. The allowed ranges of three mixing angles $`\theta _x\theta _{\mathrm{sun}}`$, $`\theta _y\theta _{\mathrm{atm}}`$ and $`\theta _z\theta _{\mathrm{chz}}`$ have been given in Eq. (10). There is no experimental constraint on the $`CP`$-violating phase $`\delta `$. Hence we simply take $`\delta `$ from $`0^{}`$ to $`360^{}`$ in our numerical calculations.

There are totally fifteen distinct topologies for the structure of $`M`$ with two independent vanishing entries, as shown in Tables 1 and 2. We work out the explicit expressions of $`(m_1/m_3)e^{2i\rho }`$ and $`(m_2/m_3)e^{2i\sigma }`$ for each pattern of $`M`$ by use of Eq. (14), and list the results in the same tables . With the input values of $`\theta _x`$, $`\theta _y`$, $`\theta _z`$ and $`\delta `$ mentioned above, we calculate the ratio $`R_\nu `$ and examine whether it is in the range allowed by current data. This criterion has been used in Refs.  to pick the phenomenologically favored patterns of $`M`$ in the LMA case.

Nine of the fifteen two-zero textures of $`M`$ listed in Table 1 are found to be in accord with the LMA solution as well as the atmospheric neutrino data. They can be classified into four categories: A (with $`\mathrm{A}_1`$ and $`\mathrm{A}_2`$), B (with $`\mathrm{B}_1`$, $`\mathrm{B}_2`$, $`\mathrm{B}_3`$ and $`\mathrm{B}_4`$), C and D (with $`\mathrm{D}_1`$ and $`\mathrm{D}_2`$). The point of this classification is that the textures of $`M`$ in each category result in similar physical consequences, which are almost indistinguishable in practice. The other six patterns of $`M`$ (categories E and F) listed in Table 2 cannot coincide with current experimental data. In particular, the exact neutrino mass degeneracy ($`m_1=m_2=m_3`$) is predicted from three textures of $`M`$ belonging to category F.

Now let us focus on patterns $`\mathrm{A}_1`$, $`\mathrm{B}_1`$, $`\mathrm{C}`$ and $`\mathrm{D}_1`$ as four typical examples for numerical illustration. Our results for $`\mathrm{sin}^22\theta _{\mathrm{chz}}`$ versus $`\delta `$ and $`\theta _y`$ versus $`\theta _x`$ are shown Figs. 1 – 4. Some comments are in order.

(1) For pattern $`\mathrm{A}_1`$, arbitrary values of $`\delta `$ are allowed if $`\mathrm{sin}^22\theta _{\mathrm{chz}}`$ is large enough ($`0.014`$). The mixing angles $`\theta _x`$ and $`\theta _y`$ may take any values in the ranges allowed by current data. Therefore we conclude that pattern $`\mathrm{A}_1`$ is favored in phenomenology with little fine-tuning. A similar conclusion can be drawn for pattern $`\mathrm{A}_2`$.

(2) For pattern $`\mathrm{B}_1`$, $`\delta `$ is essentially unconstrained if $`\mathrm{sin}^22\theta _{\mathrm{chz}}`$ is extremely close to zero; and only $`\delta `$ around $`90^{}`$ or $`270^{}`$ is acceptable if $`\mathrm{sin}^22\theta _{\mathrm{chz}}`$ deviates somehow from zero. Except $`\theta _y45^{}`$, there is no further constraint on the parameter space of $`(\theta _x,\theta _y)`$. We conclude that pattern $`\mathrm{B}_1`$ with maximal $`CP`$ violation (i.e., $`\mathrm{sin}\delta \pm 1`$) is phenomenologically favored. So are patterns $`\mathrm{B}_2`$, $`\mathrm{B}_3`$ and $`\mathrm{B}_4`$.

(3) For pattern $`\mathrm{C}`$, $`\delta =90^{}`$ or $`\delta =270^{}`$ is forbidden. Furthermore, $`\theta _y=45^{}`$ is forbidden. We see that the allowed parameter space of $`(\delta ,\theta _{\mathrm{chz}})`$ and that of $`(\theta _x,\theta _y)`$ are rather large. Hence pattern C is also favored in phenomenology.

(4) For pattern $`\mathrm{D}_1`$, $`\delta `$ is restricted to be around $`0^{}`$ or $`360^{}`$. In particular, the region $`90^{}\delta 270^{}`$ is entirely excluded. $`\mathrm{sin}^22\theta _{\mathrm{chz}}>0.084`$ holds for the allowed range of $`\delta `$. Different from patterns $`\mathrm{A}_1`$, $`\mathrm{B}_1`$ and $`\mathrm{C}`$, pattern $`\mathrm{D}_1`$ requires relatively strong correlation between $`\theta _x`$ and $`\theta _y`$ (e.g., small values of $`\theta _y`$ are associated with large values of $`\theta _x`$ in the allowed parameter space). In this sense, we argue that pattern $`\mathrm{D}_1`$ is less natural in phenomenology, although it has not been ruled out by current experimental data. A similar argument can be made for pattern $`\mathrm{D}_2`$.

It is worth remarking that patterns $`\mathrm{D}_1`$ and $`\mathrm{D}_2`$ were not included into the phenomenologically allowed patterns of $`M`$ in the previous classification , where only analytical approximations were made. Our numerical analysis shows that these two patterns are marginally allowed by current data. Here we have also explored some interesting details of the parameter space for every favorable texture, which could not be seen from simple analytical approximations.

## IV Numerical predictions and further discussions

A two-zero texture of $`M`$ has a number of interesting predictions, in particular, for the absolute neutrino masses and the Majorana phases of $`CP`$ violation . With the help of Eq. (14), one may calculate the mass ratios $`m_1/m_3`$ and $`m_2/m_3`$ as well as the Majorana phases $`\rho `$ and $`\sigma `$. The absolute neutrino mass $`m_3`$ can be determined from

$$m_3=\frac{1}{\sqrt{\left|1\left({\displaystyle \frac{m_2}{m_3}}\right)^2\right|}}\sqrt{\mathrm{\Delta }m_{\mathrm{atm}}^2}.$$
(35)
Therefore a full determination of the mass spectrum of three neutrinos is actually possible. Then we may obtain definite predictions for the effective mass of the tritium beta decay,

$$m_e=m_1c_x^2c_z^2+m_2s_x^2c_z^2+m_3s_z^2;$$
(36)
and that of the neutrinoless double beta decay,

$$m_{ee}=\left|m_1c_x^2c_z^2e^{2i\rho }+m_2s_x^2c_z^2e^{2i\sigma }+m_3s_z^2\right|.$$
(37)
It is clear that the Dirac phase $`\delta `$ has no contribution to $`m_{ee}`$. Note that $`CP`$\- and $`T`$-violating asymmetries in normal neutrino oscillations are controlled by $`\delta `$ or the rephasing-invariant parameter $`J=s_xc_xs_yc_ys_zc_z^2\mathrm{sin}\delta `$. Whether $`m_e`$ and $`m_{ee}`$ can be measured remains an open question. The present experimental upper bounds are $`m_e<2.2\mathrm{eV}`$  and $`m_{ee}<0.35\mathrm{eV}`$  at the $`90\%`$ confidence level. The proposed KATRIN experiment is possible to reach the sensitivity $`m_e0.3\mathrm{eV}`$ , and a number of next-generation experiments for the neutrinoless double beta decay  is possible to probe $`m_{ee}`$ at the level of 10 meV to 50 meV.

We perform a numerical calculation of $`m_2/m_3`$ versus $`m_1/m_3`$, $`\sigma `$ versus $`\rho `$, $`m_{ee}`$ versus $`m_e`$, and $`J`$ versus $`m_3`$ for patterns $`\mathrm{A}_1`$, $`\mathrm{B}_1`$, $`\mathrm{C}`$ and $`\mathrm{D}_1`$. The results are shown in Figs. 1 – 4. Some discussions are in order.

(1) For pattern $`\mathrm{A}_1`$, $`\rho \delta /2`$ or $`\rho \delta /2180^{}`$ and $`\sigma \rho \pm 90^{}`$ hold in most cases. Two neutrino mass ratios lie in the ranges $`0.033m_1/m_30.19`$ and $`0.13m_2/m_30.28`$, and the absolute value of $`m_3`$ is in the range $`0.04\mathrm{eV}m_30.065\mathrm{eV}`$. As $`m_{ee}=0`$ is a direct consequence of texture $`\mathrm{A}_1`$, we calculate the sum of three neutrino masses $`m_i`$ instead of $`m_{ee}`$. The result is $`0.047\mathrm{eV}m_i0.093\mathrm{eV}`$, in contrast with $`0.003\mathrm{eV}m_e0.014\mathrm{eV}`$. The rephasing invariant of $`CP`$ violation $`J`$ is found to lie in the range $`0.037J0.038`$. Similar predictions are expected for pattern $`\mathrm{A}_2`$.

(2) For pattern $`\mathrm{B}_1`$, $`\rho \sigma \delta 90^{}`$ or $`\rho \sigma \delta 270^{}`$ holds in most cases. Two neutrino mass ratios $`m_1/m_3`$ may lie either in the range $`0.53m_1/m_30.99`$ or in the range $`1.01m_1/m_31.88`$, and $`m_2/m_3`$ may lie either in the range $`0.53m_2/m_30.99`$ or in the range $`1.01m_2/m_31.88`$. The value of $`m_3`$ is found to be in the range $`0.026\mathrm{eV}m_30.35\mathrm{eV}`$. Furthermore, we arrive at $`0.027\mathrm{eV}m_em_{ee}0.35\mathrm{eV}`$ as well as $`0.037J0.037`$. Similar results can be obtained for patterns $`\mathrm{B}_2`$, $`\mathrm{B}_3`$ and $`\mathrm{B}_4`$.

(3) For pattern $`\mathrm{C}`$, $`\sigma \rho `$ when $`\theta _y`$ approaches $`45^{}`$; and there is no clear correlation between $`\rho `$ and $`\sigma `$ for other values of $`\theta _y`$. Two neutrino mass ratios $`m_1/m_3`$ may lie either in the range $`0.95m_1/m_30.99`$ or in the range $`1.01m_1/m_35.4`$, and $`m_2/m_3`$ may lie either in the range $`0.95m_2/m_30.99`$ or in the range $`1.01m_2/m_35.3`$. The value of $`m_3`$ is found to lie in the range $`0.009\mathrm{eV}m_30.35\mathrm{eV}`$. It is remarkable that $`m_{ee}m_3`$ holds to a good degree of accuracy in the allowed space of those input parameters. We also obtain $`0.04\mathrm{eV}m_e0.35\mathrm{eV}`$ and $`0.037J0.037`$.

(4) For pattern $`\mathrm{D}_1`$, $`\rho \delta 90^{}`$ or $`\rho \delta 270^{}`$ and $`\sigma \rho \pm 90^{}`$ hold. Two neutrino mass ratios lie in the ranges $`7.5m_1/m_38.8`$ and $`7.35m_2/m_38.6`$, and the absolute value of $`m_3`$ is in the range $`0.005\mathrm{eV}m_30.008\mathrm{eV}`$. As for the tritium beta decay and neutrinoless double beta decay, we obtain $`0.04\mathrm{eV}m_e0.062\mathrm{eV}`$ and $`0.008\mathrm{eV}m_{ee}0.014\mathrm{eV}`$. The range of $`J`$ is found to be $`0.014J0.011`$. Similar predictions can straightforwardly be made for pattern $`\mathrm{D}_2`$.

We see that there is no hope to measure both $`m_e`$ and $`m_{ee}`$, if the neutrino mass matrix $`M`$ takes pattern $`\mathrm{A}_1`$ or $`\mathrm{A}_2`$. As for categories B and C of $`M`$, the upper limit of $`m_e`$ is close to the sensitivity of the KATRIN experiment ($`0.3\mathrm{eV}`$ ), and that of $`m_{ee}`$ is just below the current experimental bound .

## V Summary

In summary, we have discussed some implications of the KamLAND measurment on the lepton flavor mixing matrix $`V`$ and the neutrino mass matrix $`M`$. The model-independent constraints on nine elements of $`V`$ have been obtained up to a reasonable degree of accuracy. Nine two-zero textures of $`M`$ are found to be compatible with current experimental data, but two of them are only marginally allowed. Instructive consequences of these phenomenologically favored textures of $`M`$ on the absolute neutrino masses, Majorana phases of $`CP`$ violation, $`m_e`$ and $`m_{ee}`$ are numerically explored. Our results will be very useful for model building , in order to understand why neutrino masses are so tiny and why two of the lepton flavor mixing angles are so large.

Finally it is worth remarking that a specific texture of lepton mass matrices may not be preserved to all orders or at any energy scales in the unspecified interactions from which lepton masses are generated. Nevertheless, those phenomenologically favored textures at low energy scales, no matter whether they are of the two-zero form or other forms, are possible to provide enlightening hints at the underlying dynamics of lepton mass generation at high energy scales.

We would like to thank Y.F. Wang for stimulating discussions about the KamLAND measurement. One of us (Z.Z.X.) is also grateful to IPPP in University of Durham, where the paper was finalized, for its warm hospitality. This work was supported in part by National Natural Science Foundation of China.


text:
# Predicting proton-nucleus total reaction cross sections up to 300 MeV using a simple functional form

## I Introduction

Reaction cross sections from the scattering of nucleons by nuclei (stable and radioactive) are required in a number of fields of study; some being of quite current interest De01 . Those fields range from problems in basic science through many of applied nature. An example of the latter is the transmutation of long lived radioactive waste into shorter lived products using accelerator driven systems (ADS). As well, proton-nucleus ($`p`$-$`A`$) cross section values at energies to 300 MeV or more are needed not only to specify important quantities of relevance to proton radiation therapy Jo01 , but also as they are key information in assessing radiation protection for patients. As an example in basic science, total reaction cross sections are important ingredients to a number of problems in astrophysics. Also the integral observable of proton scattering from a given nucleus gives direct information on the neutron root mean square radius in nuclei Bro2000 ; a property sought in parity-violating electron scattering experiments Je00 .

It would be very utilitarian if aspects of such scattering were well approximated by a simple convenient functional form. Recently it has been shown Amos02a  that, for neutron total reaction cross sections, such a form may exist. Herein we consider that concept further to reproduce the measured total reaction cross sections from proton scattering for energies ranging from 10 MeV to 300 MeV, and from 15 nuclei ranging in mass from 9 to 238.

## II  Formalism

The total reaction cross sections for proton scattering from nuclei can be expressed in terms of partial wave scattering ($`S`$) matrices specified at energies $`Ek^2`$, by

$$S_l^\pm S_l^\pm (k)=e^{2i\delta _l^\pm (k)}=\eta _l^\pm (k)e^{2i\mathrm{}\left[\delta _l^\pm (k)\right]},$$
(1)
where $`\delta _l^\pm (k)`$ are the (complex) scattering phase shifts and $`\eta _l^\pm (k)`$ are the moduli of the $`S`$ matrices. The superscript designates $`j=l\pm 1/2`$. Total reaction cross sections then follow from

$`\sigma _R(E)`$ $`=`$ $`{\displaystyle \frac{\pi }{k^2}}{\displaystyle \underset{l=0}{\overset{\mathrm{}}{}}}\left\{(l+1)\left[1\left(\eta _l^+\right)^2\right]+l\left[1\left(\eta _l^{}\right)^2\right]\right\}`$ (2)
$`=`$ $`{\displaystyle \frac{\pi }{k^2}}{\displaystyle \underset{l=0}{\overset{\mathrm{}}{}}}\sigma _l^{(R)}(E).`$
Despite the successes with use of the $`g`$-folding approach to define $`NA`$ optical potentials Review , we note that there are discrepancies between the predictions of total reaction cross sections found thereby and actual data. Those usually are due to limitations with the structure model used to describe the ground state densities of some nuclei. However, there have been so many successes with the approach Review ; De01 ; Kara02 ; Amos02  that we believe the functional form developed here on the basis of matching the values obtained with g-folding potentials are pertinent initial guesses to start a refinement (of the parameter values) to reproduce actual measured total reaction cross sections.

As evident from the figures presented in a recent paper Amos02a , the partial total reaction cross sections, $`\sigma _l^{(R)}(E)`$, can be described by the simple functional form,

$$\sigma _l^{(R)}(E)=(2l+1)\left[1+e^{\frac{(ll_0)}{a}}\right]^1+ϵ(2l_0+1)e^{\frac{(ll_0)}{a}}\left[1+e^{\frac{(ll_0)}{a}}\right]^2$$
(3)
with $`l_0(E,A)`$, $`a(E,A)`$, and $`ϵ(E,A)`$ varying smoothly with energy and mass.

The summation giving the total reaction cross section can be limited to a value $`l_{max}`$ and the associated form tends appropriately to the known high energy limit. With increasing energy, $`l_{max}`$ becomes so large that the exponential fall of the functional form, Eq. 3, can be approximated as a straight vertical line ($`l_0=l_{max}`$). In that case, the total reaction cross section equates to the area of a triangle, and

$$\sigma _R\frac{\pi }{2k^2}l_{max}(2l_{max}+1)\frac{\pi }{k^2}l_{max}^2.$$
(4)
Then with $`l_{max}kR`$, at high energies

$$\sigma _R\pi R^2;$$
(5)
the geometric cross section as required.

## III  Results and discussions

Although on using Eqs. 2 and 3 to match values of (theoretically) calculated total reaction cross sections led Amos02a  to the three parameters $`l_0(E,A)`$, $`a(E,A)`$, and $`ϵ(E,A)`$ having smooth variations with energy, there are discrepancies between those predictions and the actual measured data. Herein we modify the method of selection of those parameter values to produce more accurate reaction cross section values, while keeping as smooth a variation with energy of those parameters as possible. Specifically, in Eq. 3, we have set the $`ϵ`$ as a constant ($`1.5`$) and so independent of energy and of mass. Further we assume that $`a(E,A)`$ varies linearly with the wave vector,

$$k=\frac{1}{\mathrm{}c}\sqrt{E^2m^2c^4},$$
(6)
and with the form

$$a(E,A)1.02k0.25.$$
(7)
Then $`l_0(E,A)`$ were adjusted to ensure that all measured total reaction cross section values are matched by using the function form, Eq. 3.

The resultant optimized values for the parameter $`l_0`$ are presented in Fig. 1. In this figure, the points represent the calculated values of $`l_0`$ at particular energies where the experimental values of total reaction cross sections are available for the considered nuclei. We have used the lower and upper uncertainty limits to get the error bars for $`l_0`$ values. We call these $`l_0`$ values with uncertainties as ”data-$`l_0`$” hereafter. The curves result on using a spline interpolation on those data points constrained to give an optimal smoothness variation of $`l_0(E,A)`$ with E. The figure displays the results for the entire energy range (to 300 MeV) in the left hand panels, while the right hand panels emphasize the low energy variations to 50 MeV since most experimental values lie in that range of energy.

In panel (a) of Fig. 1, the solid, dotted, dashed, long-dashed, and dot-dashed lines depict the values for 9Be, 12C, 16O, 19F, and 27Al respectively, as do circles, squares, diamonds, up triangles and down triangles for the data-$`l_0`$ results. The same legends apply to the 10 to 50 MeV plots shown in panel (b). In panels (c) and (d), the energy variations of $`l_0`$ values for 40Ca, 63Cu, 90Zr, 118Sn, and 140Ce are presented by solid, dotted, dashed, long-dashed and dot-dashed lines respectively. The corresponding data-$`l_0`$ points are also presented by circles, squares, diamonds, up triangles and down triangles. Finally , in panels (e) and (f), we display the variations of $`l_0`$ with energy for 159Tb, 181Ta, 197Au, 208Pb and 238U by solid, dotted, dashed, long-dashed and dot-dashed lines respectively with the corresponding data-$`l_0`$ points presented by circles, squares, diamonds, up triangles and down triangles.

Numerical values for $`l_0`$ that result from the simple functional form calculations and as displayed by the curves in Fig. 1, are listed at various energies in Table. 1. For energies 50 MeV and higher, these $`l_0`$ values monotonically increase with mass and energy as might be anticipated. At lower energies however, and while the $`l_0`$ values do still vary monotonically with energy, there is structure in their mass variation. In particular, the values of $`l_0`$ for 197Au, 208Pb and 238U in the energy range from 10 MeV to 20 MeV decrease with mass. But note that there is little data available for these nuclei at this energy range.

The mass variations given by these tabled values have extra structure as is evident from the plots given in Fig. 2.

The results are not as smooth as the counterpart variation with energy particularly at the lowest of the energy set. That may reflect the scattered (in energy) nature of the actual data from which starting values of $`l_0`$ have been taken to generate the functional forms.

In Figs. 34, and 5, we show the total reaction cross sections generated using the simple functional form and tabled values of $`l_0`$, and displayed by the solid curves, with those obtained from calculations made using $`g`$-folding optical potentials Amos02 . Dashed lines represent the predictions obtained from those microscopic optical model calculations. The experimental data Car96  are depicted by circles.

The results for scattering from 9Be, 12C, 16O, 19F, 27Al, and 40Ca are displayed in segments (a) through (f) of Fig. 3 respectively.

In segment (a), the data are well reproduced by the $`g`$-folding predictions resulting from the folding with the 9Be ground state OBDME found with (0+2)$`\mathrm{}\omega `$ spectroscopy as well as by those obtained with the simple functional form method.

Calculated p-12C reaction cross sections are compared with the experimental data in segment (b) of Fig. 3. The reaction cross sections obtained from $`g`$-folding calculations are in good agreement with the experimental data but only in the energy range above 20 MeV. On the other hand, the results obtained from the simple functional method are excellent for all energies, replicating the average trend of data well even in the lower energy range from 10 MeV. There are two data points, at 61 MeV and at 77 MeV, in disagreement with the calculated results however. But, as noted previously De01 , these data points should be discounted.

Predictions for p-16O and for p-19F scattering are compared with data in segments (c) and (d) of Fig. 3. For p-16O case, there are many data points at the energies between 20 to 40 MeV. Predictions from $`g`$-folding calculations while replicating the data well at and above 25 MeV, overestimate at lower energies. That $`g`$-folding result also underestimates the datum at 250 MeV; the sole datum above 50 MeV. In contrast, and by design, the results obtained from the simple functional method are in excellent agreement with the experimental data at all energies. For p-19F, although $`g`$-folding calculations reproduce the data well, the simple functional form method gives slightly more accurate predictions.

Total reaction cross section predictions for p-27Al and for p-40Ca are compared with the experimental data in segments (f) and (g) of Fig. 3. Again while $`g`$-folding calculations reproduce the data quite well to 200 MeV, three data points between 180 and 300 MeV are not matched. The $`g`$-folding results underestimate them noticeably. But predictions from a simple functional form approach replicate the data well at all energies. One data point at 61 MeV is exceptional in the set. With 40Ca, the folding model approach is not expected to be reliable at the energies in the range 10 to 20 MeV, as is the case with 12C, since for excitation energies of that range, both nuclei have clearly discrete spectra. That is true for most light mass nuclei but little or no total reaction cross sections have been reported for them. Indeed the reaction data from both 12C and from 40Ca show rather sharp resonance features below 20 MeV. Both the $`g`$-folding calculations and functional form calculations reproduce the rest of the 40Ca data well.

The results for scattering from 63Cu, 90Zr, 118Sn, 140Ce, 159Tb, and 181Ta are displayed in Fig. 4 in segments (a) through (f) respectively.

For 63Cu, predictions at low energies may be slightly too small and the parameter sets driven too severely by the sole datum at 30 MeV in the range 20 to 70 MeV. Also the data in the range 100 to 300 MeV are quite scattered but the simple functional form gives a good average result.

In segment (b) of Fig. 4, the predicted total reaction cross sections from p-90Zr scattering are compared with the experimental data. Results from $`g`$-folding calculations are in good agreement with the data although the data value at 30 MeV is overestimated. The results obtained from the simple functional form, and by dint of construction, are in excellent agreement with the experimental data at the few energies measured, but the shape is not optimally smooth. This we believe is the prime cause for the kink shown at $`A=90`$ in the mass variation plots in Fig. 2 and most noticeably at 200 MeV. Lack of data meant that we had to use the $`g`$-folding values to specify the functional form. That is also the case with masses 140, 159, and 181.

The p-118Sn total reaction cross section results are given in segment (c) of Fig. 4 where the two predictions again are compared with the data. Although not as good as the results found for scattering from light mass nuclei, the $`g`$-folding potential still gives reasonable shape prediction. However, the model underestimates the data by 5 to 10%. But the simple functional form model form again gives an excellent reproduction of the data at all energies except 61 MeV. This 61 MeV data point is again exceptional being much smaller than other data and the predictions as in the cases of 12C and 27Al.

Predictions for p-140Ce, p-159Tb and for p-181Ta scattering are compared with the (limited amount of) data in segments (d), (e) and (f) of Fig. 4 respectively. The $`g`$-folding calculations give good agreement with that data for p-140Ce, slightly underestimate the data for p-159Tb, and for p-181Ta, underestimate data at the energies to 20 MeV and overestimate data in the energy range 40 to 60 MeV. In all cases results predicted by the simple functional form are excellent reproductions of the experimental data.

In segments (a), (b) and (c) of Fig. 5 we compare the calculated total reaction cross sections for p-197Au, p-208Pb and p-238U scattering with the available data.

For 197Au, the $`g`$-folding optical potential calculations are in good agreement with most data; the 29 MeV datum grossly underestimated by the calculations. But that data point is also at odds with the energy trend of the other data. Save for that 29 MeV value, the simple functional form gives even better predictions of the data. The energy variation of the p-208Pb reaction cross sections is shown in segment (b) of Fig. 5 where the predictions from $`g`$-folding optical potential calculations and from the simple functional form calculations are compared with a fairly extensive set of experimental data. In making the $`g`$-folding potentials, we have used Skyrme-Hartree-Fock wave functions Bro2000  which have been shown to be more realistic Amos02  than simple oscillator model ones. Still such $`g`$-folding calculations underestimate the data up to 50 MeV. But simple functional form calculations are in excellent agreement with the experimental data, save that data values at 30, 61 and 77 MeV again are exceptional. The predictions of total reaction cross sections for p-238U scattering from the $`g`$-folding optical potential calculations and from using the simple functional form are compared with the few data in segment (c) of Fig. 5. Given the lack of data the two results are virtually identical.

As a final note regarding many of the exceptional data values so defined in the foregoing, Menet et al. Men71  argue that there may be a systematic error in the studies reported in those experiments.

## IV  conclusions

The measured reaction cross sections for 10–300 MeV proton scattering from nuclei ranging in mass from 9Be to 238U are well reproduced by calculations made using a three parameter function. The values of the parameters ($`l_0`$, $`a`$, $`ϵ`$) are set in a simple manner and when set by enough actual data, can be used to predict the total reaction cross sections at any energy for a given nucleus. The mass variations of those parameters at fixed energies also are smooth, and very much so when data exists to control the results, and we believe that the simple functional form also can be used to find reasonable estimates of the total reaction cross sections of protons from any stable nucleus in the mass range.

###### Acknowledgements.

This research was supported by a research grant from the Australian Research Council.