text:
# Tunable adsorption on carbon nanotubes

\[

## Abstract

We investigated the adsorption of a single atom, hydrogen and aluminum, on single wall carbon nanotubes from first-principles. The adsorption is exothermic, and the associated binding energy varies inversely as the radius of the zigzag tube. We found that the adsorption of a single atom and related properties can be modified continuously and reversibly by the external radial deformation. The binding energy on the high curvature site of the deformed tube increases with increasing radial deformation. The effects of curvature and radial deformation depend on the chirality of the tube.

\]

Novel mechanical, electrical and chemical properties of carbon nanotubes  have been explored actively with a motivation of finding a new technological application. A single wall carbon nanotube (SWNT) is usually described by a rolled graphene, where the hexagonal 2D lattice is mapped on a cylinder of radius $`R`$ with various helicities characterized by a set of two integers $`(n,m)`$. A SWNT can display either metallic or insulating electronic structure depending on the helicity and radius .

Recent studies  showed that the electronic properties of SWNTs can be modified by radial deformation. The energy gap of an insulating SWNT can decrease and eventually vanish at an insulator–metal transition with increasing applied radial strain. The density of states at the Fermi energy, $`𝒟(E_F)`$, of a metallized SWNT increases with the increasing radial strain. More interestingly, the radial deformation necessary to induce metallicity was found to be in the elastic range. Therefore, all strain induced changes in the electronic and also in mechanical properties are reversible.

Most noticeably, the radial strain disturbs the uniformity of charge distribution. This, in turn, may impose changes in the chemical reactivity and hence in the interaction of tube surface with foreign atoms and molecules. It is therefore anticipated that not only band gap but also chemical reactions taking place on the surface of a SWNT can be engineered through radial deformation. In this paper, we explore this feature by using the predictive power of the density functional theory and demonstrate that indeed adsorption of foreign atoms on carbon nanotubes and associated properties can be modified continuously and reversibly. Furthermore, we showed that there is a simple scaling of the adsorption energies with the radius of the SWNT. We believe that the tunable adsorption can have important implications for metal coverage and selective adsorption of foreign atoms and molecules on the carbon nanotubes, and can lead to a wide variety of technological applications, ranging from hydrogen storage to new materials .

We have investigated the electronic energy structure and charge density of bare and single atom adsorbed SWNTs with and without radial strain by performing first-principle density functional calculations . We expressed the wave functions by linear combinations of plane waves up to an energy cutoff of 500 eV within the tetragonal supercell geometry (separation between the tubes is taken to be 7.5 Å and $`c`$ is equal to the lattice constant of SWNT along its axis). With this energy cutoff and using ultra soft pseudopotentials for carbon  the total energies of nanotube-adatom systems converge within 0.5 meV/atom. The Monkhorst-Pack  special k-point scheme is used with 0.02 Å<sup>-1</sup> k-point spacing resulting 5 k-points along the tube axis. Results have been obtained within the generalized gradient approximation  (GGA) for fully relaxed geometries including all carbon and adsorbate positions and the lattice constant of the tube along the $`z`$-axis.

To investigate the effect of the radial deformation we first consider adsorption of H and a simple metal, Al, on undeformed SWNTs. In Fig. 1, we present the binding energies of H and Al adsorbed on the $`(n,0)`$ zigzag nanotubes calculated for $`n=`$7,8,9,10,12. H is adsorbed at the top site, i.e. directly above the C atoms of the tube, Al favors the hollow site i.e. above the center of the hexagon as in the graphite surface. The binding energy $`E_b`$,

$$E_b=E_T[SWNT]+E_AE_T[SWNT+A]$$
(1)
is calculated in terms of the total energy of SWNT, $`E_T[SWNT]`$, total energy of SWNT with an adsorbed atom A, $`E_T[SWNT+A]`$, and the energy of the single, free atom, $`E_A`$. Here the bare SWNT and SWNT with an adatom A are free of external stress, and all the atomic coordinates are fully relaxed. Moreover, since $`E_b`$ is calculated by using the same supercell, the spurious adatom-adatom interaction along the tube axis is substracted. The positive value of $`E_b`$ indicates that the adsorption is exothermic and hence stable. It is found that the binding energy of an adatom decreases with increasing radius (or decreasing curvature) of the tube, and eventually saturates at a value corresponding to that on graphene plane. Hence the variation of binding energies of H and Al with the radius of the zigzag tube fits to the curve given by the expression,

$$E_{b,A}(R)=E_{o,A}+\frac{C_A}{R}$$
(2)
where $`E_{o,A}`$ is the binding energy of the adatom A (H or Al) on the graphene plane. We calculated $`E_{o,H}=`$1.49 eV, and $`E_{o,Al}=`$0.02 eV. Interestingly, the fitting parameter $`C_A`$ is found to be independent of the adatom (H and Al) and is equal to $``$ 3.14 eVÅ. Currently, we are extending our calculations to see if this relation holds for other adsorbates as well. The binding energies calculated for $`(n,0)`$ SWNTs with $`n<8`$ deviate from the above simple scaling perhaps due to the fact that the singlet $`\pi ^{}`$-band, which is normally in the conduction band, falls into the band gap as a result of increased $`\sigma ^{}\pi ^{}`$ mixing at high curvature . Note that, while the band gap shows significant change with $`n`$, (for example, upon going from $`(8,0)`$ to $`(9,0)`$ $`E_g`$ changes from 0.65 eV to 0.09 eV) the binding energies vary smoothly with $`R^1`$. Increasing $`E_b`$ with decreasing $`R`$ (or with decreasing $`n`$) shows that for small $`R`$ the character of the surface deviates from that of the graphene. This finding also suggests that by creating regions of different curvature on a single SWNT by radial deformation one can attain different values of binding energies.

The radial deformation that we consider in this study is generated by applying uniaxial compressive stress on a narrow strip on the surface of the SWNT . In practice such a deformation can be realized by pressing the tube between two rigid flat surfaces. The radius is decreased in the $`y`$-direction, while it is elongated along the $`x`$-direction. As a result, the circular cross section is distorted to an elliptical one with major and minor axis, $`a`$ and $`b`$ respectively. The elliptical radial deformation can be described by the magnitude of the applied strain and it is defined as $`ϵ_{yy}=(Rb)/R`$ where $`R`$ is the radius of the undeformed nanotube (with zero strain). For different values of $`ϵ_{yy}`$ we carried out full structural optimization under the constraint that the minor axis is kept fixed at a preset value by freezing the carbon atoms at both ends of the minor axis . Then, under this constraint, the coordinates of the adatom and remaining carbon atoms and the lattice constant of the tube, $`c`$ are optimized until there remain only forces opposite to the applied strain on the fixed atoms, but all other force components are less than 0.01 eV/Å.

Figure 2a shows the variation of the binding energy $`E_b`$ of a single hydrogen atom adsorbed on the $`(8,0)`$ surface with the applied radial strain for two cases: (i) The binding energy of a single H atom adsorbed on the high curvature side of the surface (i.e. specified as the sharp-site near one of the ends of the major axis, $`x=a`$, $`y=0`$) traces the upper curve in Fig. 2a. (ii) The lower curve is for the adsorption on the low curvature side near one of the ends of the minor axis (i.e. $`x=0`$, $`y=b`$ specified as the flat-site). Here, the binding energy is defined as in Eq(1), except that $`E_T[SWNT]`$ and $`E_T[SWNT+A]`$ are calculated for radially deformed SWNTs. The binding energy of H adsorbed on the high curvature (sharp) site is increased by 0.85 eV for $`ϵ_{yy}`$=0.3. On the other hand, $`E_b`$ for the adsorption on the low curvature (flat) site first decreases with increasing $`ϵ_{yy}`$, and then saturates at an energy 0.25 eV less than that corresponds to $`ϵ_{yy}=`$0.0. The difference of binding energies of the sharp and flat sites, $`\mathrm{\Delta }E_b`$, is substantial and is equal to $``$ 1.1 eV. This value is 44 % of the binding energy of H on the undeformed SWNT. As a result of H adsorption, the $`sp^2`$ character of the bonding of the tube has changed locally and became more like to $`sp^3`$. The lengths of the C–C bonds at the close proximity of H have increased slightly.

The binding energy of Al shown in Fig. 2b exhibits a behavior similar to that of H, despite H and Al favor different sites on the $`(8,0)`$ tube: $`E_b`$ at the sharp site of the deformed SWNT increases with increasing $`ϵ_{yy}`$. For example, $`E_b`$ increases by $``$ 0.80 eV for $`ϵ_{yy}=0.3`$ which is 80% of the binding energy on the undeformed tube. For Al absorbed on the flat site, $`E_b`$ first decreases with increasing $`ϵ_{yy}`$, then gradually increases. Adsorption of Al induces local changes in the atomic and electronic structure. For example, the surface of the tube where Al is adsorbed expands.

The variation of $`E_b`$ with the radial deformation is consistent with the results illustrated in Fig. 1. In general, the higher is the curvature under deformation, the higher the binding energy. The effect of the elastic deformation is further investigated by analyzing the charge density and electronic structure. The Mulliken analysis estimates that $``$0.37 electrons is transfered from H to carbon . The C–H bond is directional and is covalent with a partial ionic component. The charge difference, $`\mathrm{\Delta }\rho =\rho [SWNT+H]\rho [SWNT]\rho [H]`$, which is calculated in terms of the undeformed (or deformed) SWNT with H (adsorbed at different sites), undeformed (deformed) clean SWNT, and single H, are shown by inset in Fig. 2a. $`\mathrm{\Delta }\rho `$ indicates no significant change with radial deformation. On the other hand, the charge transfer upon the adsorption of Al is different from H. Since Al is adsorbed on the center of a hexagon, the bond between Al and SWNT is distributed to nearest C atoms with charge accumulation between Al and those C atoms. The bond charge slightly increases at the sharp site as shown by inset in Fig. 2b. The charge transfer from Al atom is estimated to be 0.71 electrons.

Explanation of this remarkable and significant change of the binding energy with radial deformation is sought in the electronic energy structure and the total and partial density of states. The adsorption of H gives rise to a new state which falls in the gap at $`\mathrm{\Gamma }`$ and coincides with $`E_F`$. This state partially overlaps with the conduction band when H is adsorbed at the flat site, whereas it occurs near the valence band edge at the sharp site. A similar situation occurs with the adsorption of Al as illustrated in Fig. 3. The band gap is wide open, and a resonant state $``$ 2.5 eV below $`E_F`$ and a $`p`$-derived adsorption state which sets $`E_F`$ falls at the center of the gap for the adsorption of Al on the undistorted SWNT . In the case of adsorption on the flat site the latter state overlaps the conduction band (Fig. 3c), but it dips into the valence band at the sharp site (Fig. 3d). Accordingly, the sharp-site adsorption involves the valence band states of SWNT, while the conduction band states are involved in the flat site adsorption.

The above argument is in accord with the binding energies calculated for the adsorption of Al on the undeformed and deformed $`(6,6)`$ armchair SWNT. Note that the $`(6,6)`$ as well as all $`(n,n)`$ SWNTs are metallic with finite $`𝒟(E_F)`$ owing to the $`\pi `$ and $`\pi ^{}`$states crossing and setting $`E_F`$. Our calculations show that $`𝒟(E_F)`$ of the $`(6,6)`$ decreases slightly with increasing radial deformation, but it remains essentially finite. We found $`E_b=`$0.91 eV for the adsorption of Al on the $`(6,6)`$ SWNT. Upon the deformation of the tube with $`ϵ_{yy}=`$0.22, $`E_b`$ corresponding to sharp-, and flat-site adsorption has changed to 0.95 and 0.85 eV, respectively. Therefore, the difference between the binding energies of the sharp and flat site adsorption is only 0.1 eV, and hence is negligible as compared to the situation above discussed for the $`(8,0)`$ SWNT. In Fig. 3b the state associated with the adsortion of Al on the undeformed $`(6,6)`$ tube occurs above the Fermi level and donates its electron to the metallic $`\pi `$ and $`\pi ^{}`$bands. It appears that $`E_b`$ is almost pinned by the above mechanism.

Radial deformation induced changes in the energies of valance band edge (VBE) and conduction band edge (CBE), as well as the adsorbtion state (AE) at $`\mathrm{\Gamma }`$–point illustrated in Fig. 4 corraborate the above arguments. In the case of adsorption on the flat site, VBE is almost constant for both H and Al, while CBE and AE both shift downwards. As a result of this, more conduction band states contribute to the bond energy. By contrast, for adsorption on the sharp site, while AE behaves very similar to the previous case, CBE is now approximately constant, but VBE shifts upwards. In this situation, the contribution to the bond energy from conduction band states is decreased.

In conclusion, we showed that the chemical reactivity of a zigzag SWNT can be modified reversibly and variably by the radial deformation. The effect of deformation is significantly different for the zigzag and armchair SWNTs. It is remarkable that Al, which is not bound to the graphite surface, can be adsorbed at a high curvature site of a zigzag SWNT under radial deformation with a binding energy of 1.8 eV. This novel property may have important implications for various chemical and electronic applications, such as selective absorption and desorption of molecules and atoms, fragmentation and chemical sensors, magnetic tubes etc.

Acknowledgements This work was partially supported by the NSF under Grant No. INT97-31014 and TÜBÍTAK under Grant No. TBAG-1668(197 T 116). SC thanks Prof. S. Süzer for stimulating discussions.


format: text

text:
# References

See-saw mechanism and four light neutrino states

M.Czakon, J.Gluza and M.Zrałek
Department of Field Theory and Particle Physics

Institute of Physics, University of Silesia, Uniwersytecka 4

PL-40-007 Katowice, Poland

Abstract

A formal proof is given that in a see-saw type neutrino mass matrix with only two neutrino mass scales ($`m_Dm_R`$) and the maximal rank of $`m_{R(D)}`$ we can not get a fourth light sterile neutrino.

It has been known for long that extended gauge group models, such as $`SO(10)`$  or $`SU(2)_L\times SU(2)_R\times U(1)_{BL}`$ , naturally develop a see-saw type neutrino mass matrix. Namely ($`m_D`$ is a $`3\times n_R`$ matrix and $`m_R`$ is a $`n_R\times n_R`$ matrix)

$$m_\nu =\left(\begin{array}{cc}0& m_D\\ m_D^T& m_R\end{array}\right),$$
(1)
coupled with a large scale difference between $`m_D`$ and $`m_R`$, yields a mass spectrum containing three light Majorana neutrinos, effectively described by

$$m_{light}m_D^Tm_R^1m_D.$$
(2)
This leads to two different light $`\mathrm{\Delta }m^2`$ mass scales, enabling one to understand the solar and the atmospheric neutrino anomalies

According to the current data (now at 99 % c.l. ) a fourth sterile neutrino is not necessary to explain the Superkamiokande data. Suppose however that the LSND experiment  is confirmed and we wish also to explain its data in a language of neutrino oscillation phenomena. We then have to introduce a fourth light neutrino of sterile nature (due to the invisible width measurement at LEP). To avoid fine tunings of parameters and still have another light neutrino, one requires the theory to have additional properties, like approximate horizontal symmetry , exact parity symmetry , a discrete $`Z_5`$ symmetry , global $`S_3\times Z_2`$ symmetry , or even additional gauge ($`SU(2)_S`$) symmetry . For more examples see . Some phenomenological considerations have also appeared (see e.g. ). Interestingly enough, a see-saw type mass matrix Eq. (1) can also lead to the fourth light neutrino. This is realized by the so-called singular see-saw mechanism . The goal is achieved, by having an $`m_R`$ of rank $`n_R1`$. However, this is not enough, we still have to fine tune $`m_R`$ to the keV-MeV range. This last unwanted problem can be circumvented by building a second stage of see-saw structure. This still fits into the scheme Eq. 1, but there are in fact three scales, not two.

Here we give a formal proof that with only two scales we can not get a fourth light neutrino. The importance of this result lies in the fact that one may be tempted to believe that the larger the mass matrix the more possibilities of choosing the masses are available, and some symmetries may help getting light sterile neutrinos. The problem is defined in the following statement:

Let $`m_R`$ be a matrix of eigenvalues real positive and greater than some $`M`$, and let all of the moduli of elements of $`m_D`$ be much smaller than $`M`$, then the spectrum of $`m_\nu `$ contains the full spectrum of $`m_R`$ with corrections of the order of $`m^2/Mmax(|(m_D)_{ij}|)^2/M`$. That means that no manipulation on $`m_D`$ and/or $`m_R`$ can move a mass from the heavy $`m_R`$ matrix into the light spectrum.

The proof is a simple consequence of perturbation theory for finite matrices . There, it is shown, that if we have a matrix of the form

$$M(\beta )=M^{(0)}+\beta M^{(1)},$$
(3)
then, every non-degenerate eigenvalue of $`M^{(0)}`$ gives rise to a non-degenerate eigenvalue of $`M(\beta )`$, being an analytic function of $`\beta `$ in some surrounding of zero. Since we are interested in the heavy spectrum, the assumption of non-degeneracy is quite general. In case of degenerate eigenvalues, we still can expand the eigenvalues in series, which however will be analytic functions with branches. We limit ourselves to the former case, but the reader should understand that the theorem holds for the general case also.

We decompose $`m_\nu `$ as

$`\left(\begin{array}{cc}0& m_D\\ m_D^T& m_R\end{array}\right)`$ $`=`$ $`\left(\begin{array}{cc}0& 0\\ 0& m_R\end{array}\right)+\left(\begin{array}{cc}0& m_D\\ m_D^T& 0\end{array}\right)`$ (4)
$`=`$ $`M\left[\left(\begin{array}{cc}0& 0\\ 0& m_R/M\end{array}\right)+\beta \left(\begin{array}{cc}0& m_D/m\\ m_D^T/m& 0\end{array}\right)\right]`$
$``$ $`M(M^{(0)}+\beta M^{(1)}),`$
where $`\beta =m/M`$. The first matrix has all its elements greater than one, while the second has all elements smaller than one, both are of course dimensionless. The eigenvectors of $`M^{(0)}`$ are of the form (we chose without loss of generality a weak base for neutrinos in which $`m_R`$ is diagonal )

$$e_i=\left(\begin{array}{c}0\\ \mathrm{}\\ 1\\ \mathrm{}\\ 0\end{array}\right).$$
(5)
To find the first order correction, we expand the eigenvectors as

$$v_i=\underset{j}{}\alpha _{ji}e_j.$$
(6)
This gives us the following equation

$$(M^{(0)}+\beta M^{(1)})v_i=\lambda _iv_i,$$
(7)
which is solved into

$$\lambda _i=\lambda _i^{(0)}+\beta \frac{1}{\alpha _{ii}}\underset{j}{}\alpha _{ji}e_i^TM^{(1)}e_j,$$
(8)
where $`\lambda ^{(0)}`$ are non-degenerate nonzero eigenvalues of $`M^{(0)}`$, which are also eigenvalues of $`m_R/M`$. Obviously, only $`\alpha _{ii}`$ is $`𝒪(1)`$, and $`\alpha _{ji}`$ for $`ji`$ is $`𝒪(\beta )`$. The first order correction to $`\lambda ^{(0)}`$ is therefore

$$\lambda _i^{(1)}=\lambda _i^{(0)}+\beta e_i^TM^{(1)}e_i.$$
(9)
But this vanishes due to the nondiagonal form of $`M^{(1)}`$. Thus the first non-vanishing correction to the large eigenvalues is of the order $`M\beta ^2=m^2/M`$, which completes the proof. Remark, that the masses of the neutrinos are moduli of the eigenvalues of $`m_\nu `$. The corrections to the moduli are however of the same order.

Using simple arguments based on perturbation series, we have shown that a natural $`m_R`$ (no fine tunings and eigenvalues at the heavy scale) leads to three light neutrinos. Thus from the class of see-saw type models only a singular double see-saw mechanism can accommodate the LSND data and an additional fourth light neutrino state.

Acknowledgments

This work was supported by the Polish Committee for Scientific Research under Grant No. 2P03B04919 and 2P03B05418.


text:
# The Spectrum and Accurate Location of GRS 1758-258

## 1. Introduction

GRS 1758–258 and its sister source,  1E 1740.7–2942, were the first objects dubbed “micro-quasars”. Their X-ray spectra are typical of Galactic black hole candidates (BHCs), and they are associated with time variable cores of double-lobed radio sources, reminiscent of extra-Galactic radio sources. This morphology, seen on a parsec scale within the Milky Way, earned them their nickname.  GRS 1758–258 and  1E 1740.7–2942 are the brightest persistent sources in the Galactic bulge above $``$50 keV (Sunyaev et al., 1991). Their timing characteristics are typical of the black hole low/hard state (Main et al., 1999; Smith et al., 1997; Heindl et al., 1993; Sunyaev et al., 1991), and they consistently emit near their brightest observed levels, although they vary over times of days to years. Their X-ray emission properties are readily likened to the canonical BHC, Cyg X-1. In fact, together with Cyg X-1, they are the only known persistent, low-state BHCs, and all three sources have maximum luminosities around $`3\times 10^{37}`$ergs s-1. Radio jets have now been observed in Cyg X-1, furthering the similarity (Fender, 2000).

GRS 1758–258 and  1E 1740.7–2942 are, however, quite different from the Galactic superluminal radio sources more typically thought of as micro-quasars: GRS 1915+105 and GRO J1655-40. The X-ray emission from these objects is much brighter and more spectacularly variable. Their radio jets, too, are much brighter and are highly variable, being unresolved or absent except during exceptional ejection events which last only weeks. In contrast, the radio lobes of  GRS 1758–258 and  1E 1740.7–2942 are quite stable (Mirabel and Rodriguez, 1999).

## 2. Observations

We observed  GRS 1758–258 four times with  Chandra. Table 1 lists the observation dates and durations. Two HRC-I observations were made in 2000 September-October spanning an intermediate-to-hard spectral transition (identified with  RXTE Smith et al., 2001a, b). Another HRC-I and an ACIS/HETGS observation were made back-to-back in 2001 March following a dramatic hard-to-soft transition to a *very* low flux state. Figure 1 shows the  RXTE/PCA light curve with our  Chandra observations indicated.

## 3. Imaging and the Accurate Position

Figure 2 shows the HRC-I image from observation 400131. In this and the other HRC observations, the  GRS 1758–258 image is consistent with the HRC point spread function, indicating that the source is point-like at the sub-arcsecond level. Assuming a distance of 8.5 kpc, this says that no strong jets are present on a physical scale of a light-month, the presence of arcminute scale radio jets notwithstanding. This is perhaps not surprising, as the timescale associated with producing the parsec-sized radio lobes would be years rather than months, and in fact the radio lobes are observed to be quite stable over years(Mirabel and Rodriguez, 1999). Thus, while de facto not ruled out, it is not required that a *persistent* small-scale X-ray jet be present to produce the observed radio lobes.

Figure 3 shows the best fit source locations from our four  Chandra observations as well as the estimated 90% confidence region derived from the average of the four positions. The ACIS-S/HETGS position (400163=ch163) was based on the zeroth-order image. While the image was significantly piled up, we expect a negligible effect on the position. The  Chandra error circle is centered at (J2000)  *RA $`=`$ 18 01 12.40, Dec $`=`$ $``$25 44 36.0* and has an estimated 90% confidence radius of 0′′.6 based on  Chandra’s absolute aspect accuracy and averaging of four independent observations. The coincidence of the sub-arcsecond VLA and  Chandra error circles seals the association of  GRS 1758–258, the X-ray source, with the variable radio source (“VLA-C”, Marti et al., 1998).

## 4. Spectrum

To reduce pile-up in the ACIS/HETGS observation, we used a 1/2 chip sub-array with a 1.7 s frame time. However, the counting rates in the MEG first order spectra were still about 2 counts/s, resulting in mild ($`\stackrel{<}{}`$10%) pileup around 2 keV. For this presentation, we therefore restrict spectral fitting to the first order HEG spectra which did not suffer from pile-up. Figure 4 shows the HEG Order $`\pm 1`$ spectra rebinned to a minimum of 100 counts per bin. The data are well described by a multi-color disk-blackbody model and interstellar absorption (XSPEC: “phabs\*diskbb”). The positive deviations above 4 keV are indicative of a weak power-law component ($``$2% of the 1$``$10 keV unabsorbed flux) seen in joint fits with a contemporaneous  RXTE/PCA (3$``$15 keV) observation. Table 2 lists the best fit parameters to the HEG spectra alone.

Preliminary inspection of both the HEG and MEG spectra showed no strong emission lines. Figure 4 shows the 5$`\sigma `$ sensitivity to narrow emission lines. Because the  GRS 1758–258 spectrum is so soft, the line sensitivity is a strong function of energy, varying by nearly two orders of magnitude between 1$``$5 keV. The sensitivity was based on the best fit spectrum (see above) and the effective area of the ACIS/HETGS combination.

## 5. Discussion

During more than 5 years’ monitoring with the  RXTE prior to 2001 March, the  GRS 1758–258 hard X-ray spectrum was always dominated by a hard power law with photon index $`\mathrm{\Gamma }1.52.5`$ (Smith et al., 2001a) with occasional appearance of a weak thermal component (Mereghetti et al., 1994; Heindl and Smith, 1998; Lin et al., 2000). As shown in Figure 1,  GRS 1758–258 made an abrupt state change in 2001 March. The hard flux dropped by an order of magnitude in a few days, leaving the thermal component seen in Figure 4. Based on relative luminosity, however, the current soft state is not a *high*/soft state. Rather it is significantly less luminous than the low/hard state in this source. This can be contrasted to Cyg X-1 and the soft X-ray transients, where the *high*/soft state is more luminous. Rather, this seems to be a low-luminosity state which is fading into quiescence (Figure 1). Finally, we note that the measured column density is consistent with previous measurements (Mereghetti et al., 1994; Lin et al., 2000; Goldwurm et al., 2001)

Since strong jet ejections are generally associated with the “very high state” and transitions from the “off” to high/soft states in X-ray transients (Fender, 2001), it is perhaps not surprising that no jet emission appeared in our low/hard state observations (Sep-Oct 2000) and the recent transition observation (Mar 2001). Perhaps our best opportunity will come when (if?)  GRS 1758–258 makes a transition once again to its normal, low/hard state. We have an approved  Chandra cycle 3 proposal to monitor the morphology of  GRS 1758–258 and hope to observe a jet ejection.


text:
# Long-duration neutron star X-ray transients in quiescence: the Chandra observation of KS 1731–260

## 1. Introduction

X-ray transients spend most of their time in quiescence. The exact emission mechanisms which produce the quiescent X-rays are not yet understood. To explain the neutron star emission, several models have been put forward. So far, the model which assumes that the emission is due to the cooling of the neutron star (e.g., van Paradijs et al. 1987; Verbunt et al. 1994; Brown et al. 1998) has been most successful in reproducing the existing quiescent data and in producing testable predictions. In the most elaborate version of this model, the quiescent luminosity of a particular transient depends on the time averaged accretion rate of this source (Brown et al. 1998). Recently, this model could be tested in a new (although anticipated) way. In early 2001, the neutron star transient KS 1731–260 suddenly became quiescent again, after having actively accreted for more than a decade. A Chandra observation was performed on this source only a few months after this transition to test the Brown et al. (1998) model (Wijnands et al. 2001a). If the long duration of the outburst of KS 1731–260 is typical for this source and its quiescent time is similar to that of the other neutron star transients (of order years to a few decades), then the system was detected at a luminosity (bolometric luminosity of $`2\times 10^{33}`$ erg s-1) and a neutron star temperature ($`kT0.3`$ keV) too low compared to what would have been expected. The implications for the Brown et al. (1998) model and its validity for KS 1731–260 are discussed in detail by Wijnands et al. (2001a). If the Brown et al. (1998) model applies to KS 1731–260, its low luminosity might indicate that this system is quiescent for several hundreds of years. This would suggest that several hundreds of such systems might be present in our Galaxy and they might form a new class of X-ray transients (see Wijnands et al. 2001a).

## 2. Alternative models

The striking similarities between KS 1731–260 and the other systems, combined with the extra assumptions and possible adjustments which have to be made to the Brown et al. (1998) model (Wijnands et al. 2001a), could be used to argue that the cooling neutron star model might not be the correct model to explain the quiescent emission of neutron star transients. Alternative models include residual accretion (see, e.g., Campana et al. 1998; Campana & Stella 2000; Bildsten & Rutledge 2000; Narayan et al. 2001 for discussions) or the models which use a neutron star magnetic field (e.g., Campana et al. 1998; Campana & Stella 2000; Robertson & Leiter 2001111It is interesting to note that Robertson & Leiter (2001) predicted a luminosity of $`10^{33}`$ erg s-1 for KS 1731–260. Despite this success, their model was proposed to explain the power-law tail in the X-ray spectra of quiescent neutron star systems, which is not the dominate spectral component in KS 1731–260 (Wijnands et al. 2001a).). In all these models, it is expected that KS 1731–260 should be very similar to the other systems if its system parameters (i.e., the neutron star spin rate and the magnetic field properties; the binary parameters) are very similar to those of the other systems. The fact that the spin rates of KS 1731–260 and Aql X-1 (as inferred from the nearly coherent oscillations observed during type-I X-ray bursts) are very similar (524 Hz versus 549 Hz; Smith et al. 1997; Zhang et al. 1998), indicates that at least the spin rates of those systems are very similar. The orbital parameters of KS 1731–260 are not known, but the recent tentative discovery of the optical counterpart of this system (Wijnands et al. 2001b), might allow a determination of its orbital parameters. Any small differences in the details of the quiescent X-ray properties between the systems (e.g., the exact contribution of the hard power-law tail to the total 0.5–10 keV X-ray luminosity) might simply be due to small differences in the system parameters. However, the arguments why such models have problems explaining the quiescent neutron star emission in general still apply and it remains to be determined which (if any) of the models can produce accurately all of the quiescent properties (e.g., Campana et al. 1998; Bildsten & Rutledge 2000; Menou & McClintock 2001; Narayan et al. 2001).

## 3. Other quiescent long-duration neutron star transients

Besides KS 1731–260, several more systems can be identified as (possible) long-duration transients. Some of those systems have been observed in quiescence and although their results are not as clean as for KS 1731–260, it is useful to compare them with KS 1731–260 and discuss them in the context of the Brown et al. (1998) model.

### 3.1. EXO 0748–676

This source was discovered with EXOSAT in February 1985 (Parmar et al. 1986) but before that (on 22 May 1980) it was serendipitously observed with EINSTEIN with a quiescent luminosity of $`10^{34}`$ erg s-1 (Parmar et al. 1986; Garcia & Callanan 1999), at least a factor of ten larger than the other systems. Usually, this high quiescent luminosity is explained as due to a relatively high level of residual accretion during this particular observation (see, e.g., Garcia & Callanan 1999). However, a different explanation can be proposed by comparing this system with KS 1731–260. Since its discovery, EXO 0748–676 has persistently been detected at relatively high luminosities ($`>10^{36}`$ erg s-1), and, therefore, this source can be regarded as a long-duration transient which has been active now for over 15 years (although it is possible we have witnessed the birth of a new persistent source222The difference between a persistent source and a long-duration transient is somewhat arbitrary because in principle every persistent source turned on in the past and will turn off in the future.). If EXO 0748–676 is typically active and quiescent for a few decades, then from the Brown et al. (1998) model, we would expect a luminosity for this system in quiescence of $`10^{3435}`$ erg s-1, consistent with what has been observed. So, although KS 1731–260 did not behave as expected (based on the simplest version of the Brown et al.  model), EXO 0748–676 might behave as expected. High sensitivity observations of this latter system are needed when it becomes quiescent again to study its quiescent spectrum. If the luminosity is indeed due to the cooling of the neutron star, then its spectrum should resemble that of the other systems (a black-body-like spectrum), but with a higher neutron star temperature. However, the EINSTEIN quiescent data of EXO 0748–676 was fitted using a black-body (although other one-component models fitted equally well) with a $`kT`$ of $``$0.2 keV (Garcia & Callahan 1999), which is lower than expected. But, Chandra or XMM-Newton observations are needed to really constrain its quiescent spectrum.

### 3.2. MXB 1659–298

MXB 1659–298 was discovered in 1976 (Lewin et al. 1976) and had a clear outburst in 1978 (Lewin et al. 1978). The source was dormant until April 1999 when it was found to be in outburst again (in ’t Zand et al. 1999). Since then, the source could be detected with the RXTE all sky monitor until the writing of this paper (31 July 2001). Therefore, the source has been active for over two years and might also be considered a long-duration transient (see also Wijnands et al. 2001a), although two years is relatively short compared to the 11.5 years of KS 1731–260 and the more than 15 years of EXO 0748–676. Verbunt (2001) reported that during a 1991 ROSAT/PSPC observation, the source was not detected with an upper limit on the unabsorbed flux of $`12\times 10^{14}`$ erg s-1 cm-2 (0.5–10 keV; estimated with W3PIMMS333W3PIMMS can be found at http://heasarc.gsfc.nasa.gov/Tools/w3pimms.html assuming a black-body spectrum with $`kT0.3`$ keV and using the count rate values provided by Verbunt  but with the updated column density provided by Oosterbroek et al. ; see also the latter paper for a similar estimate but for 0.2–2.4 keV). The distance estimates to this source (as obtained from radius expansion bursts) range from 10 kpc (Muno et al. 2001) to 13 kpc (Oosterbroek et al. 2001), resulting in a 0.5–10 keV luminosity upper limit of 1–4 $`\times 10^{32}`$ erg s-1. Assuming that the upper limit on the bolometric luminosity is only a factor of a few higher (e.g., Rutledge et al. 2000), then this source had an even lower quiescent luminosity than KS 1731–260. If the time averaged mass accretion rate of the past 25 years in MXB 1659–298 is a good indication of its time averaged mass accretion rate of the past several thousands of years, then also this system is too dim in quiescence. It is also possible that the first two short outbursts are more typical for this source and the long outburst is a rare phenomenon (note that this is also possible for the other long-duration transients discussed here, including KS 1731–260). Regardless of the true averaged outburst duration of this source, high sensitivity observations are needed when it becomes dormant again to test the effect of the long period of accretion on the neutron star crust and interior.

### 3.3. 4U 2129+47

4U 2192+47 was considered to be a persistent source until in 1983 the source was not detected with EXOSAT (Pietsch et al. 1983, 1986). Although the exact beginning of its long active episode is unknown (so the duration of its active state is unclear), it is possible that 4U 2192+47 is also a long-duration X-ray transient (although its possible that we have witnessed the turn off of a persistent source). The source was detected in quiescence with the ROSAT/HRI at $`3\times 10^{33}`$ erg s-1 (0.3–2.4 keV; Garcia 1994) and the ROSAT/PSPC at $`6\times 10^{32}`$ erg s-1 (0.5–10 keV; Garcia & Callanan 1998; the bolometric luminosity would be about a factor 2 higher; Rutledge et al. 2000). The X-ray spectrum of the source in the latter data set was consistent with a black-body with a $`kT`$ of approximately 0.2 keV. If the outburst and quiescent duration so far observed for 4U 2192+47 are typical for this source, then these luminosity and black-body temperature are lower than expected from the Brown et al. (1998) model. However, during the long episode during which the source could not be detected at high levels (since September 1983; Pietsch et al. 1983) and the first detection of the source in quiescence (December 1991; Garcia 1994), the source could have cooled down. However, the exact cooling time of neutron stars is unknown and is model dependent.

### 3.4. Globular cluster source X 1732–304 (Terzan 1)

In the early eighties, Hakucho detected a bursting X-ray source in the globular cluster Terzan 1 (Makishima et al. 1981; Inoue et al. 1981). In 1985, a steady X-ray source was detected (X 1732–304) consistent with this globular cluster and it is most likely the same source as the bursting source (Skinner et al. 1987). Since then, it has persistently been detected at luminosities between a few times $`10^{35}`$ erg s-1 and $`10^{37}`$ erg s-1 (see Figure 3 of Guainazzi et al. 1999). Recently, Guainazzi et al. (1999) reported that during a BeppoSAX observation this source could only be detected at a 2–10 keV luminosity of $`1.410^{33}`$ erg s-1 (for a distance of 4.5 kpc), with a spectrum consistent with a black-body spectrum ($`kT0.34`$ keV) plus a power law (with photon index of $``$1). These spectral parameters are very similar to those of the neutron star X-ray transients in quiescence and it is likely that X 1732–304 became suddenly quiescent after having actively accreted for over 15 years (see Guainazzi et al. 1998). Thus this source is a clear example of a long-duration transient, similar to KS 1731–260.

Using the spectral parameters given by Guainazzi et al. (1998), the quiescent luminosity of X 1732–304 can be converted into a 0.5–10 keV luminosity of $`5\times 10^{33}`$ erg s-1. The bolometric luminosity will be most likely about a factor of two higher ($`10^{34}`$ erg s-1; see, e.g., Rutledge et al. 2000). This is considerably higher than the bolometric luminosity of KS 1731–260 in quiescence and might indicate that X 1732–304 behaved more as expected in quiescence due to its prolonged period of high accretion (based on the Brown et al.  model). The fact that the BeppoSAX quiescent observation of this source was taken at approximately two years after the last detection of the source above $`10^{35}`$ erg s-1 (Molkov et al. 2001), demonstrates that the exact moment when the source became dormant is unknown. The source could have been in quiescence for about two years before the BeppoSAX observation and it could be initially even more luminous and could have cooled down considerably (this depends on the exact cooling time of the neutron star crust and interior). However, due to the low angular resolution of BeppoSAX, it can also not be excluded that the detected source might not be X 1732–304 but another low-luminosity globular cluster source (see also Guainazzi et al. 1998). In that case, X 1732–304 might be less luminous than assumed in the above discussion and might be more similar to KS 1731–260. A Chandra image of Terzan 1 will reveal whether or not the BeppoSAX source is indeed X 1732–304 (its position from a ROSAT/HRI observation is known to $`5^{\prime \prime }`$; Johnston et al. 1994) and whether the luminosity if contaminated by other nearby sources. X 1732–304 might be a excellent candidate to test the cooling neutron star model proposed by Brown et al. (1998).

## 4. Conclusions

The quiescent Chandra observation of KS 1731–260 can constrain the cooling neutron star models proposed for the quiescent emission of neutron star transients. However, it is less constraining for other types of models. In those models, KS 1731–260 is expected to be very similar to the other quiescent systems if its orbital parameters and/or its neutron star parameters (i.e., magnetic field parameters, spin rate) are very similar to those of the other systems.

It is clear that KS 1731–260 is part of a group of transients which do not disappear after a few weeks to months but are active for years to decades. Several such sources can be identified for which observations were also performed in quiescence. Two of those systems, EXO 0748–676 and X 1732–304 might indeed be brighter in quiescence than the other systems, possible due to the heating of the neutron stars during the long periods of accretion. Three other systems (4U 2129+47, MXB 1659–298, and KS 1731–260) seem to have anomalously low luminosities in quiescence and might prove to be very constraining for the models dealing with the heating of the crust and the core of neutron stars in X-ray binaries. However, only for KS 1731–260 is a clear picture available about its outburst behavior (i.e., the exact mass accretion rate during outburst; the outburst duration) and was the quiescent observation performed within a few months after the source become quiescent again. For the other sources, the outburst behavior is not well known (i.e., its duration; e.g., EXO 0748–676, MXB 1659–298, 4U 2129+47), the quiescent observation were performed years after (4U 2129+47, X 1732–304) or before (MXB 1659–298, EXO 0748–676) a long-duration outburst, and source confusion might play a role (X 1732–304)

## 5. Acknowledgement

I thank Jon Miller for carefully reading this manuscript before submission. This work was supported by NASA through Chandra Postdoctoral Fellowship grant number PF9-10010 awarded by CXC, which is operated by SAO for NASA under contract NAS8-39073.

## References

Bildsten, L. & Rutledge, R. E. 2000, In: “The Neutron Star - Black Hole Connection:” (NATO ASI ELOUNDA 1999), eds. C. Kouveliotou et al. (astro-ph/0005364)

Brown, E. F. et al. 1998, ApJ, 504, L95

Campana, S. et al. 1998 A&AR, 8, 279

Campana, S. & Stella, L. 2000, ApJ, 541, 849

Garcia, M. R. 1994, ApJ, 435, 407

Garcia, M. R. & Callanan, P. J. 1999, ApJ, 118, 1390

Guainazzi, M. et al. 1999, A&A, 349, 819

Inoue, H. et al. 1981, ApJ, 250, L71

in ’t Zand, J. et al. 1999, IAUC 7138

Johnston, H. M. et al. 1995, A&A, 298, L21

Lewin, W. H. G. et al. 1976, IAUC 2994

Lewin, W. H. G. et al. 1978, IAUC 3190

Makishima, K. et al. 1981, ApJ, 247, L23

Menou, K. & McClintock, J. E. 2001, ApJ, in press (astro-ph/0010430)

Molkov, S. V. et al. 2001, Astronomy Letters, 27, 6 (see also astro-ph/0105002)

Muno, M. et al. 2001, ApJ, 553, L157

Narayan, R. et al. 2001, In the proceedings of the IX Marcel Grossman Meeting. eds. V. Gurzadyan et al. (astro-ph/0107387)

Oosterbroek, T. et al. 2001, A&A in press (astro-ph/0107388)

Parmar, A. N. et al. 1986, ApJ, 308, 199

Pietsch, W. et al. 1983, IAUC, 3887

Pietsch, W. et al. 1986, A&A, 157, 23

Robertson, S. L. & Leiter, D. J. 2001, ApJ submitted (astro-ph/0102381)

Rutledge, R. E. et al. 2000, ApJ, 529, 985

Skinner, G. K. et al. 1987, Nature, 330, 544

Smith, D. A. et al. 1997, ApJ, 479, L137

van Paradijs, J. et al. 1987, A&A, 182, 47

Verbunt, F. et al. 1994, A&A, 285, 903

Verbunt. F. 2001, A&A, 368, 137

Wijnands, R. et al. 2001a, ApJL submitted (astro-ph/0107380)

Wijnands, R. et al. 2001b, ATEL 72

Zhang, W. et al. 1998, ApJ, 495, L9


text:
# Properties of Microlensing Light Curve Anomalies Induced by Multiple Planets

## 1 Introduction

Since the first discovery of extra-solar planets around a pulsar PSR 1257+12 from the analysis of planetary perturbations on pulse arrival time by Wolszczan & Frail (1992) and Wolszczan (1994), various groups have detected planets around nearby stars using high precision radial velocity measurements, e.g., around 51 Pegasi (Mayor & Queloz 1995), 47 Ursae Majoris (Butler & Marcy1996), 70 Virginis (Marcy & Butler 1996), and 16 Cygni B (Cochran et al. 1997). Up to now, nearly 60 extra-solar planets are known (http://cfa-www.harvard.edu/planets). From close examination of the velocity residuals to the Keplerian fits to a subsample of 12 planet-bearing stars, Fischer et al. (2001) suggested that nearly 1/2 of the stars may have a second companion. Currently, six systems are identified to have multiple planets: around $`\nu `$ And (Butler et al. 1999), HD 83443 (Mayor et al. 2001), HD 168443 (Marcy et al. 2001a), G7 876 (Marcy et al. 2001b), HD 82943 and HD 74156 (http://exoplanets.org/mult.shtml).

Although majority of the best candidate extra-solar planets were discovered by the radial velocity technique, they can also be detected by using other techniques. Microlensing is one of the promising techniques, especially to search for planets located at Galactic-scale distances. Planet detections by using microlensing is possible because the event caused by a lens system composed of a planet (or planets) can exhibit noticeable anomalies in the lensing light curve when the source star passes closely to the lens caustics (Mao & Paczyński 1991). The caustic refers to the source position on which the lensing-induced amplification of a point-source event becomes infinity. The lens system with a planet (planets) form several disconnected sets of caustics. Among them, one is located very close to the primary lens, “central caustic”, and the other(s) is (are) located away from the primary lens, “planetary caustic(s)” (Griest & Safizadeh 1998). The planet-induced anomalies last only a few hours to days, and thus it is difficult to detect them from the survey-type lensing experiments (MACHO: Alcock et al. 1993; OGLE: Udalski et al. 1993; EROS: Aubourg et al. 1993; MOA: Abe et al. 1997). However, with sufficiently frequent and accurate monitoring from follow-up observations of ongoing events alerted by the survey teams, one can detect the anomalies and determine the planet parameters of the mass ratio and the instantaneous projected separation between the planet and the primary lens. Currently, several groups are actually carrying out such observations: MPS (Rhie et al. 1999), PLANET (Albrow wt al. 1998), and MOA (Abe et al. 1997).

Along with the observational efforts, planetary microlensing has also been the field of intensive theoretical studies. These include phenomenological studies of the lensing behaviors (Wambsganss 1997; Dominik 1999; Safizadeh et al. 1999), estimation of the detection probability and efficiency (Gould & Loeb 1992; Bolatto & Falco 1994; Bennett & Rhie 1996; Peale 1997, 2001; Gaudi & Sackett 2000), development of the observational strategies for optimal planet detections (Griest & Safizadeh 1998; Han & Kim 2001), and various problems and solutions in determining planet parameters (Gaudi & Gould 1997; Gaudi 1998; Griest & Safizadeh 1998). However, most of these studies were focused on lens systems composed of only a single planet, despite the fact that our solar system is composed of multiple planets. Gaudi, Naber & Sackett (1998) demonstrated that if extra-solar systems are composed of planets having orbital separations comparable to those of Jupiter and Saturn of our solar system, the joint probability of two planets having projected separations within the lensing zone of $`0.6`$ – 1.6 of the Einstein ring radius is substantial, invoking the necessity of studies about the lensing properties of multiple planetary systems.

In this paper, we investigate the properties of lensing light curve anomalies induced by multiple planets. Multiple planetary microlensing was first explored by Gaudi et al. (1998). However, since they were interested only in the anomalies occurring near the peaks of the light curves of high amplification events, for which the planet detection efficiencies are high (Griest & Safizadeh 98), they investigated anomaly patterns only in the region around central caustics. In this work, we extend the investigation in wider areas including the outer deviation regions around planetary caustics.

The paper is organized in the following way. We begin with a brief review of the basics of multiple planetary lensing in § 2. In § 3, we present the anomaly patterns of an example lens system and the resulting light curves to illustrate the effect of multiple planets. We discuss about the implications of the new findings in § 4. A brief comment and conclusion are found in § 5.

## 2 Basics of Multiple Planetary Lensing

A lens system composed of multiple planets is described by the formalism of multiple lensing with very low mass-ratio companions. If a source star is lensed by $`N`$ point-mass lenses, the locations of the resulting images are obtained by solving the mapping equation (lens equation), which is expressed in complex notations by

$$\zeta =z+\underset{j}{\overset{N}{}}\frac{m_j}{\overline{z}_j\overline{z}},$$
(1)
where $`m_j`$’s represent the mass fractions of the individual lenses such that $`_j^Nm_j=1`$, $`z_j`$’s are the positions of the lenses, $`\zeta =\xi +i\eta `$ and $`z=x+iy`$ are the positions of the source and images, and $`\overline{z}`$ denotes the complex conjugate of $`z`$ (Witt 1990). We note that all these lengths are normalized by the angular Einstein ring radius, which is related to the total mass, $`M`$, and the geometry of the lens system by

$$\theta _\mathrm{E}=\sqrt{\frac{4GM}{c^2}}\left(\frac{1}{D_{\mathrm{ol}}}\frac{1}{D_{\mathrm{os}}}\right)^{1/2},$$
(2)
where $`D_{\mathrm{ol}}`$ and $`D_{\mathrm{os}}`$ represent the distances to the lens and source from the observer, respectively. The angular separations between the microlensed images are of the order of several milli-arcsecs, which are too small to be resolved. However, the lensing event can be identified by the change of the source star flux. The amplifications of the individual images are given by the inverse of the determinant of the Jacobian of the lens equation evaluated at the image position, i.e.

$$A_i=\left(\frac{1}{|\mathrm{det}J|}\right)_{z=z_i};\mathrm{det}J=1\frac{\zeta }{\overline{z}}\frac{\overline{\zeta }}{\overline{z}}.$$
(3)
Then the total amplification is obtained by the sum of those of the individual images, i.e. $`A=_i^{N_I}A_i`$, where $`N_I`$ represents the number of images. The caustic corresponds to the source position on which the determinant becomes zero, i.e. $`|\mathrm{det}J|=0`$. The lens equation describes a mapping from the lens plane onto the source plane. Then, to find image positions $`(x,y)`$ for a given source position $`(\xi ,\eta )`$, the equation should be inverted.

For a single lens system ($`N=1`$), the lens equation can be inverted and thus is algebraically solvable, yielding two solutions (and thus the same number of images). The total amplification is expressed in a simple form of

$$A=\frac{u^2+2}{u\sqrt{u^2+4}},$$
(4)
where $`𝒖`$ is the dimensionless lens-source separation vector normalized by $`\theta _\mathrm{E}`$. The separation vector is related to the single lensing parameters by

$$𝒖=\left(\frac{tt_0}{t_\mathrm{E}}\right)\widehat{𝒙}+\beta \widehat{𝒚},$$
(5)
where $`t_\mathrm{E}`$ represents the time required for the source to transit $`\theta _\mathrm{E}`$ (Einstein time scale), $`\beta `$ is the closest lens-source separation in units of $`\theta _\mathrm{E}`$ (impact parameter), and $`t_0`$ is the time at that moment. The unit vectors $`\widehat{𝒙}`$ and $`\widehat{𝒚}`$ are parallel with and normal to the direction of the relative lens-source transverse motion.

For a multiple lens system ($`N2`$), on the other hand, the lens equation is non-linear and thus cannot be inverted. However, it can be expressed as a polynomial in $`z`$, and the image positions are obtained by numerically solving the polynomial. For a $`N`$ point-mass lens system, the lens equation is equivalent to a $`(N^2+1)`$-order polynomial in $`z`$ and there are a maximum $`N^2+1`$ and a minimum $`N+1`$ images and the number of images changes by a multiple of two as the source crosses a caustic. Due to the non-linear nature of the lens mapping for a multiple lens system, the lensing behavior of events caused by a multiple lens is not given by a simple superposition of the individual single lens events.

Unless the planet has a separation very close to $`\theta _\mathrm{E}`$, a single-planet lens system forms a single central caustic and one or two planetary caustics depending on the projected separation between the primary and the planet. Accordingly, there exist two different types of planet-induced anomalies; one affected by the planetary caustic(s) (“type I anomaly”) and the other affected by the central caustic (“type II anomaly”) (Han & Kim 2001). Compared to the frequency of type I anomalies, type II anomalies occur with a relatively low frequency due to the smaller size of the central caustic compared to that of the planetary caustic(s). However, the efficiency of detecting type II anomalies can be high because intensive monitoring is possible due to the known type of target events (i.e., high amplification events) and the predictable time of anomalies (i.e., near the peak of light curves).

As the number of lenses increases, solving the lens equation becomes very cumbersome. One commonly practiced approach to obtain amplification patterns of these lens systems is the inverse ray-shooting technique (Schneider & Weiss 1986; Kayser, Refsdal & Stabell 1986; Wambsganss 1997). In this method, a large number of light rays are shoot backwards from the observer plane through the lens plane, and then collected (binned) in the source plane. Then the amplification pattern is obtained by the ratio of the surface brightness (i.e., the number of rays per unit area) on the source plane to that on the observer plane and the light curve from a particular source trajectory corresponds to the one-dimensional cut through the constructed amplification pattern. This method has an advantage of enabling one to obtain amplification patterns for an arbitrary number of lenses, but has a disadvantage of requiring large computation time to obtain smooth amplification patterns.

In our analysis of multiple-planet lensing properties, we test lens systems composed of two Jovian-mass planets (and thus $`N=3`$). This is because if extra-solar systems have planets with orbital separations and mass ratios similar to those of our system, other planets except the two Jovian planets (Jupiter and Saturn) will have mass ratios too small to affect the amplification pattern. In addition, the probability of more than three planets being simultaneously located in the lensing zone will be very small. Besides, this allows us to directly solve the lens equation rather than using the more time-consuming ray-shooting method.

## 3 Properties of Anomalies

For the investigation of the properties of lensing light curve anomalies induced by multiple planets, we construct maps of amplification excesses. The amplification excess is defined by

$$ϵ_{\mathrm{tri}}=\frac{A_{\mathrm{tri}}A_0}{A_0},$$
(6)
where $`A_{\mathrm{tri}}`$ and $`A_0`$ represent the amplifications with and without planets. With the map, one can obtain an overview of the anomaly patterns not testing all light curves resulting from a large number of source trajectories.

In the upper panel of Figure 1, we present the constructed excess map of an example lens system. The lens system has two planets with mass ratios of $`q_1=0.003`$ and $`q_2=0.001`$, which respectively correspond to those of Jupiter- and Saturn-mass planets around a $`0.3M_{}`$ star. The individual planets have projected separations (normalized by $`\theta _\mathrm{E}`$) of $`b_1=1.3`$ and $`b_2=0.85`$ and the orientation angle between the position vectors from the primary lens to the individual planets is $`\varphi =150^{}`$. In the map, the levels of excesses are represented by contours, which are drawn at the levels of $`ϵ=10\%`$, $`5\%`$, 5%, and 10% and the regions of positive excess are distinguished by grey scales. To show the detailed excess patterns in the region around the central caustic, the map of the region is expanded and presented separately in Figure 2.

From the investigation of the map, one may notice an interesting point that the locations of the outer caustics are very similar to those of the planetary caustics of the individual planet-primary binary pairs. To show this similarity more closely, in the lower panel of Fig. 1 we present the excess map constructed by superposing the excesses of the two single-planet lens systems where the individual planet-primary pairs act as independent lens systems (binary superposition), i.e.

$$ϵ_{\mathrm{bs}}=\underset{i=1}{\overset{2}{}}\frac{A_{\mathrm{b},i}A_0}{A_0},$$
(7)
where $`A_{\mathrm{b},i}`$’s represent the binary-lensing amplifications of the individual planet-primary pairs. In Figure 3, we also present the blowup of the central region of the map. From the comparison of the excess patterns in the outer deviation regions, one finds that not only the locations of the outer caustics but also the excess patterns around them are very similar each other. This can be seen also from the comparison of the light curves of events caused by the triple lens system (solid curves in the left panels of Figure 5) to those obtained by the binary superposition (dotted curves in the same panels of the same figure), which are obtained by

$$A_{\mathrm{bs}}=A_0(1+ϵ_{\mathrm{bs}})=\underset{i=1}{\overset{2}{}}A_{\mathrm{b},i}A_0.$$
(8)
From the comparison of the maps in Fig. 2 and 3, one finds that the excess patterns in the region around the central caustic are also similar each other, although the similarities are not as exact as those in the outer deviation region.111For example, compared to the two disconnected sets of closed shape for the caustics constructed by the binary superposition, the triple lens system forms a single set of nested caustics.

To more systematically investigate how accurately the multiple-planet lensing properties can be represented as binary superposition, we compute the fractional deviation of the amplification obtained by the binary superposition from that of the exact triple lens system by

$$\delta =\frac{A_{\mathrm{tri}}A_{\mathrm{bs}}}{A_{\mathrm{tri}}}.$$
(9)
In Figure 4, we present the contour map of $`\delta `$ for the same lens system whose excess map is presented in Fig. 1. We note that in the map the contours are drawn at the levels of $`\delta =1\%`$, $`0.5\%`$, 0.5%, and 1%, which are much smaller than those of the contours in the map of $`ϵ`$. In the right panels of Fig. 5, we also present the deviation curves of $`\delta `$ for the events resulting from the source trajectories marked in Fig. 1. From the map and deviation curves, one finds that the differences occur only in small localized regions very close to the caustics. The spikes in the deviation curves are caused by the shift of the caustics. However, we note that the shift is very small ($`10^3\theta _\mathrm{E}`$), corresponding to the time shift of approximately half an hour for a typical Galactic bulge event with $`t_\mathrm{E}20`$ – 30 days. Considering the photometric precision ($`>1\%`$) and monitoring frequency ($`3`$ times/night) of the current microlensing follow-up observations, it will be difficult to notice the caustic shifts. In addition, the modified light curve can be fit by a light curve with very slightly changed lensing parameters. Therefore, the binary superposition will be a good approximation for the description of the example multiple-planet system.

Then, to what companion mass ratios is the approximation of binary superposition valid for the description of the lensing behavior of multiple-companion systems? To answer this question, we compute the fraction of the area within the Einstein ring where the binary superposition approximation deviates more than a threshold value $`\delta _{\mathrm{th}}`$ for lens systems with various combinations of mass ratios, $`f_{\mathrm{dev}}(\delta _{\mathrm{th}})`$. In Figure 6, we present the contours of constant $`f_{\mathrm{dev}}`$ in the parameter space of $`\mathrm{log}q_1`$ and $`\mathrm{log}q_2`$. We note that the presented $`f_{\mathrm{dev}}`$ are the values averaged for lens systems with planets having relative orientation angles in the range $`0\varphi 2\pi `$ and located in the lensing zone, i.e. $`0.6<b_i<1.6`$. The upper and lower panels represent the maps constructed with the threshold deviations of $`\delta _{\mathrm{th}}=0.5\%`$ and $`1\%`$, respectively. From the maps, one finds that most lens systems with heavier companion mass ratios of $`q_1<10^2`$ have $`f_{\mathrm{dev}}1\%`$, implying that the lensing behaviors of multiple-planet systems are well approximated by the binary superposition. We note, however, that if $`q_1>0.05`$, the fraction becomes considerable \[$`f_{\mathrm{dev}}(\delta _{\mathrm{th}}=0.5\%)>8\%`$ and $`f_{\mathrm{dev}}(\delta _{\mathrm{th}}=1\%)>4\%`$\], implying that the binary superposition approximation will not be appropriate to describe the lensing behavior of systems composed of very heavy planets (with masses $`>10M_\mathrm{J}`$)222For HD 168443, which was identified to have multiple planets by the radial velocity technique, one planet member has a mass of $`>17M_\mathrm{J}`$ and brown-dwarf-mass companions.

## 4 Implications

In the previous section, we show that the patterns of amplification excesses of a lens system composed of multiple planets can be well approximated by the superposition of excesses induced by the planets of the individual planet-primary binary pairs. This finding has the following implications in describing the lensing behaviors of events caused by multiple planetary systems and determining their planet parameters.

The first implication of our finding is that even if an event is caused by a multiple-planet lens system, a simple single-planet lensing model will be good enough for the description of most type I anomalies. This is because the excess pattern in the outer deviation region caused by one of the planets is hardly affected by the other planet component(s). Complication might occur if the outer deviation regions induced by the individual planets occur at a similar place. However, this case is rare since the outer excess region induced by a planet-mass companion occupies a small area.

However, the situation becomes complicated for the type II anomalies. This is because the central deviation regions of the individual planets always occur in the same region regardless of their orientations. As a result, the pattern of anomalies induced by one planet in this region can be significantly affected by the existence of other planet(s), as pointed out by Gaudi et al. (1998). Since one should apply complex multiple-planet lensing models for the description of anomalies occurred in this region, accurate determinations of the planet parameters from the observed light curves will be difficult, as pointed out also by Gaudi et al. (1998).

## 5 Conclusion

We have investigated the properties of anomalies in lensing light curves induced by multiple planets. From the analysis of the maps of excess amplification and resulting light curves, we find that the excess patterns of a multiple-planet lens system are well described by the superposition of the excesses of the single-planet systems where the individual planet-primary pairs act as independent lens systems. This finding has an important implication that a simple single-planet lensing model is good enough for the description of most type I anomalies even if an event is caused by a lens system composed of multiple planets.

Another interesting point to be mentioned is that detection of the anomalies resulting from the source trajectory passing both the outer regions of deviation induced by more than two planets (e.g., the trajectory designated by a number ‘3’ in Fig. 1) will provide a new channel of detecting multiple planets. The previously suggested method of detecting multiple planets by analyzing anomalies occurred by the source’s passage of the central caustic region (such as the trajectory ‘4’ in Fig. 1) has significant limitations in identifying the existence of multiple planets and determining their parameters due to the degeneracy of the resulting light curves from those of single-planet lensing and the complexity of multiple-planet lensing models. On the other hand, if a multiple-planet lens system is detected through the new channel, one can better determine the planet parameters because the individual anomalies are well approximated by relatively simple single-planet lensing models. We defer a more systematic analysis of this method and estimation of multiple planet detection probabilities for future work.

We thank to B. S. Gaudi for making helpful comments on multiple planetary lensing. This work was supported by a grant (2000-015-DP0449) from the Korea Research Foundation (KRF).


text:
# Analytical expressions for the spin-spin local-field factor and the spin-antisymmetric exchange-correlation kernel of a two-dimensional electron gas

## Abstract

We present an analytical expression for the static many-body local field factor $`G_{}(q)`$ of a homogeneous two-dimensional electron gas, which reproduces Diffusion Monte Carlo data and embodies the exact asymptotic behaviors at both small and large wave number $`q`$. This allows us to also provide a closed-form expression for the spin-antisymmetric exchange and correlation kernel $`K_{xc}^{}(r)`$ which represents a key input for spin-density functional studies of inhomogeneous electronic systems.

The static spin-spin response function $`\chi _\mathrm{S}(q)`$ of a paramagnetic electron gas (EG) can be written in terms of the Lindhard function $`\chi _0(q)`$ by means of the spin-antisymmetric many-body local field $`G_{}(q)`$ through the relationship

$$\chi _\mathrm{S}(q)=\mu _B^2\frac{\chi _0(q)}{1+v_qG_{}(q)\chi _0(q)},$$
(1)
where $`\mu _B`$ is the Bohr magneton. Thus $`G_{}(q)`$ is a fundamental quantity for the determination of many properties of a general electron system. It is a key input, together with the charge-charge local-field factor $`G_+(q)`$, in the spin-density functional theory (SDFT) of the inhomogeneous electron gas and in studies of quasiparticle properties (such as the effective mass and the effective Landè g-factor) in the electronic Fermi liquid.

For what concerns SDFT calculations, a common approximation to the unknown exchange-correlation energy functional $`E_{\mathrm{xc}}[n,\zeta ]`$ of electron density $`n`$ and spin-polarization $`\zeta `$ appeals to its second functional derivatives,

$$K_{\mathrm{xc}}^{\sigma \sigma ^{}}(\overline{n},\overline{\zeta };|𝐫𝐫^{}|)\frac{\delta ^2E_{\mathrm{xc}}[n,\zeta ]}{\delta n_\sigma (𝐫)\delta n_\sigma ^{}(𝐫^{})}|_{n=\overline{n},\zeta =\overline{\zeta }},$$
(2)
where $`\overline{n}=\overline{n}_{}+\overline{n}_{}`$ is the average local density of the EG and $`\overline{\zeta }=(\overline{n}_{}\overline{n}_{})/\overline{n}`$ is the average local spin-polarization. The local field factors $`G_{\sigma \sigma ^{}}(q)`$ and the matrix of the exchange-correlation kernels are simply related in Fourier transform by

$$\stackrel{~}{K}_{\mathrm{xc}}^{\sigma \sigma ^{}}(q)d^dre^{i𝐪𝐫}K_{\mathrm{xc}}^{\sigma \sigma ^{}}(r)=v_qG_{\sigma \sigma ^{}}(q),$$
(3)
where $`d`$ is the dimensionality of the system and $`v_q`$ is the Fourier transform of the Coulomb potential $`e^2/r`$. In the paramagnetic state the following relation holds between the charge-charge and spin-spin local-field factors $`G_\pm (q)`$ and $`G_{\sigma \sigma ^{}}(q)`$:

$$G_\pm (q)=\frac{G_{}(q)\pm G_{}(q)}{2}.$$
(4)
Thus, from Eq. (3) and (4), one obtains that the spin-symmetric and the spin-antisymmetric exchange and correlation kernels $`K_{\mathrm{xc}}^\pm (r)`$ are related to $`G_\pm (q)`$ by

$$\stackrel{~}{K}_{\mathrm{xc}}^\pm (q)\frac{\stackrel{~}{K}_{\mathrm{xc}}^{}(q)\pm \stackrel{~}{K}_{\mathrm{xc}}^{}(q)}{2}d^dre^{i𝐪𝐫}K_{\mathrm{xc}}^\pm (r)=v_qG_\pm (q).$$
(5)
In what follows we shall only consider $`K_{\mathrm{xc}}^{}(r)`$ in the case of two spatial dimensions, with $`d=2`$ and $`v_q=2\pi e^2/q`$. An analytical expression for $`K_{\mathrm{xc}}^+(r)`$ has already been given in Ref. \[ \].

A number of exact asymptotic properties of $`G_{}(q)`$ in two dimensions are readily proven. In particular,

$$\underset{q0}{lim}G_{}(q)=A_{}\frac{q}{k_F}$$
(6)
with

$$A_{}=\frac{1}{r_s\sqrt{2}}\left(1\frac{\chi _0}{\chi }\right),$$
(7)
where $`k_F=\sqrt{2\pi n}=\sqrt{2}/r_sa_B`$ is the Fermi wave number, $`r_s=\sqrt{\pi na_B^2}`$ is the usual EG density parameter with $`a_B`$ the Bohr radius, $`\chi _0=\mu _B^2m/\pi \mathrm{}^2`$ is the (Pauli) spin susceptibility of the ideal Fermi gas while $`\chi `$ is the spin susceptibility of the interacting system. By making use of the thermodynamic definition of $`\chi `$ we can write

$$\frac{\chi _0}{\chi }=1\frac{\sqrt{2}}{\pi }r_s+\frac{r_s^2}{2}\frac{^2ϵ_c}{\zeta ^2}|_{\zeta =0}$$
(8)
where $`ϵ_c`$ is the correlation energy per particle in the spin-polarized fluid. Unfortunately the correlation energy at finite spin polarization is poorly known. By extrapolating the Diffusion Monte Carlo (DMC) data for $`G_{}(q)`$ to $`q=0`$ and by observing from Eq. (8) that $`A_{}1/\pi `$ for $`r_s0`$, we obtain the following parametrization of $`A_{}`$ in the range $`0r_s10`$:

$$A_{}=\frac{1}{\pi +1.4954r_s+0.3193\sqrt{r_s}}.$$
(9)
The asymptotic behaviour of $`G_{}(q)`$ at large $`q`$ is also known exactly,

$$G_{}(q)C_{}\frac{q}{k_F}+B_{}$$
(10)
for $`q\mathrm{}`$, where $`C_{}`$ is proportional to the difference in kinetic energy between the interacting and the ideal gas,

$$C_{}=\frac{tt_0}{2\pi ne^2}k_F=\frac{r_s}{2\sqrt{2}}\frac{d}{dr_s}\left(r_sϵ_c(r_s)\right).$$
(11)
The correlation energy per particle $`ϵ_c(r_s)`$ in the paramagnetic fluid is available from Monte Carlo data. Finally $`B_{}=g(0)`$, $`g(0)`$ being the value of the pair-correlation function at the origin.

In this work we fit the values of $`G_{}(q)`$ originally obtained by DMC in Ref. \[\] in such a way as to obtain analytical expressions for both $`\stackrel{~}{K}_{\mathrm{xc}}^{}(q)`$ and $`K_{\mathrm{xc}}^{}(r)`$. The functional form previously employed for the charge-charge local-field factor $`G_+(q)`$ is also suitable for $`G_{}(q)`$:

$$G_{}(\overline{q})=A_{}\overline{q}\left[\frac{\kappa _{}(r_s)}{\sqrt{1+\left(A_{}\kappa _{}(r_s)\overline{q}/B_{}\right)^2}}+\left(1\kappa _{}(r_s)\right)e^{\overline{q}^2/4}\right]+C_{}\overline{q}\left(1e^{\beta _{}\overline{q}^2}\right)+\mathrm{P}_{}(\overline{q})e^{\alpha _{}\overline{q}^2},$$
(12)
where $`\overline{q}=q/k_F`$, $`\kappa _{}(r_s)=\sqrt{1+0.0082r_s^2}`$ and $`\mathrm{P}_{}(\overline{q})`$ is the polynomial $`\mathrm{P}_{}(\overline{q})=h_2\overline{q}^2+h_4\overline{q}^4+h_6\overline{q}^6+h_8\overline{q}^8`$. The free parameters $`\alpha _{}`$, $`\beta _{}`$ and $`h_2,\mathrm{},h_8`$ are fitted so as to minimize the differences from the DMC numerical results. As for $`G_+(q)`$, it proves useful to have a continuous parametrization of these parameters at least in the range $`0r_s10`$. We propose the following:

$`\alpha _{}(r_s)=\mathrm{exp}\left({\displaystyle \frac{0.2231+81.2115(r_s/10)}{1+54.6665(r_s/10)+50.7534(r_s/10)^2}}\right),`$ (13)
$`\beta _{}(r_s)=0.80890.4025(r_s/10)^30.0941(r_s/10)^{1/2},`$ (14)
$`h_2(r_s)=\mathrm{exp}\left({\displaystyle \frac{12.6262+20.9673(r_s/10)}{1+12.4002(r_s/10)}}\right),`$ (15)
$`h_4(r_s)=0.0531\left(1e^{2.154r_s^2}\right)0.4984(r_s/10)^{3/2}+0.4021(r_s/10)^2,`$ (16)
$`h_6(r_s)=0.0076\left(1e^{r_s}\right)+0.0977(r_s/10)^{3/2}0.0726(r_s/10)^2,`$ (17)
$`h_8(r_s)=0.0027(r_s/10).`$ (18)
In Figure 1 we compare the fit given by Eq. (12) with the DMC data for $`r_s=1,2,5`$ and $`10`$. In Figure 2 we show the local field factor $`G_{}(q)`$ as from Eq. (12) for various values of $`r_s`$.

We turn next to the evaluation of $`K_{\mathrm{xc}}^{}(r)`$. From Eq. (5) and (12), the expression of the spin-antisymmetric exchange-correlation kernel in real space (in Ryd) is readily obtained as

$$K_{\mathrm{xc}}^{}(r)=N_1\frac{\delta ^{(2)}(𝐫)}{k_F^2}+N_2\frac{\mathrm{exp}\left(B_{}k_Fr/(A_{}\kappa _{}(r_s))\right)}{k_Fr}+N_3e^{(k_Fr)^2}+N_4e^{(k_Fr)^2/4\beta _{}}+\underset{n=1}{\overset{4}{}}N_{5,2n}F_{2n}(\alpha _{},k_Fr),$$
(19)
where $`N_1=4\pi \sqrt{2}C_{}/r_s`$, $`N_2=2\sqrt{2}B_{}/r_s`$, $`N_3=4\sqrt{2}A_{}(1\kappa _{}(r_s))/r_s`$, $`N_4=\sqrt{2}C_{}/r_s\beta _{}`$ and $`N_{5,n}=2^{\frac{3}{2}}h_n/r_s`$. The function $`F_n(\alpha ,x)`$ is given by

$$F_n(\alpha ,x)=_0^{\mathrm{}}𝑑yy^n\mathrm{J}_0(xy)e^{\alpha y^2}=\frac{1}{2}\alpha ^{(1+n)/2}\mathrm{\Gamma }(\frac{1+n}{2})_1\mathrm{F}_1(\frac{1+n}{2};1;\frac{x^2}{4\alpha })$$
(20)
where $`\mathrm{\Gamma }(z)`$ is Euler’s Gamma function and $`{}_{1}{}^{}\mathrm{F}_{1}^{}(a;b;z)`$ is Kummer’s function. In practice the function $`F_n(\alpha ,x)`$ can be obtained via the recursive relation

$$F_{n+2}(\alpha ,x)=\frac{dF_n(\alpha ,x)}{d\alpha },F_2(\alpha ,x)=\frac{\sqrt{\pi }}{16\alpha ^{5/2}}\left[(4\alpha x^2)\mathrm{I}_0\left(\frac{x^2}{8\alpha }\right)+x^2\mathrm{I}_1\left(\frac{x^2}{8\alpha }\right)\right]e^{x^2/8\alpha }$$
(21)
where $`\mathrm{I}_n(z)`$ is the modified Bessel function of order $`n`$. It is also useful to recall that $`d\mathrm{I}_0(z)/dz=\mathrm{I}_1(z)`$ and that $`d\mathrm{I}_1(z)/dz=\mathrm{I}_0(z)\mathrm{I}_1(z)/z`$.

In Figure 3 we show the matrix elements $`K_{\mathrm{xc}}^{\sigma \sigma ^{}}(r)`$ of the exchange-correlation kernel for various values of $`r_s`$. These kernels describe the local structure (Pauli-Coulomb hole) of the paramagnetic electron fluid around an electron of a given spin. $`K_{\mathrm{xc}}^{}(r)`$, being the sum of $`K_{\mathrm{xc}}^+(r)`$ and $`K_{\mathrm{xc}}^{}(r)`$, has the same behavior of these two components, namely no structure at both small and high $`r_s`$. The situation is different for $`K_{\mathrm{xc}}^{}(r)`$, which is instead given by the difference of $`K_{\mathrm{xc}}^+(r)`$ and $`K_{\mathrm{xc}}^{}(r)`$. Although these two components are structureless at small $`r_s`$, the difference between them around $`k_Fr1`$ is large and this is reflected in the large structure shown in Figure 3.

In conclusion, we have presented an analytic parametrization of the local field factor entering the spin response of the two-dimensional electron gas in the paramagnetic state, incorporating the known asymptotic behaviors and giving an accurate representation of the available quantum Monte Carlo data. We have obtained from it and from our earlier results on the dielectric response analytic expressions of the exchange-correlation kernels for spin-density-functional calculations on inhomogeneous two-dimensional electronic systems.

###### Acknowledgements.

This work was partially supported by MURST through the PRIN 1999 program. We are grateful to Dr. S. Moroni for fruitful discussions and for providing us with the results of the Diffusion Monte Carlo study. Three of us (B. D., M. P. and M. P. T.) wish to thank the Condensed Matter Group of the Abdus Salam International Center for Theoretical Physics in Trieste for their hospitality during the final stages of this work.