text:
# Long term stability and dynamical environment of the PSR 1257+12 planetary system

## 1. Introduction

Up to date, among about 130 extrasolar planetary systems, the one around PSR 1257+12 remains the only one which contains Earth-sized planets (Wolszczan & Frail, 1992; Wolszczan, 1994; Konacki & Wolszczan, 2003). It has been discovered with the pulsar timing technique that relies on extremely precise measurements of the times of arrival (TOA) of pulsar pulses. Such a technique in principle allows a detection of companions as small as asteroids. In the case of the PSR B1257+12 it enabled the detection of three planets A, B and C with the orbital periods of 25, 66 and 98 days and the masses in the Earth-size regime (see Table 1). Luckily, the two larger planets have the mean motions close to the 3:2 commensurability which result in observable deviations from a simple Keplerian description of the motion (Rasio et al., 1992; Malhotra et al., 1992; Peale, 1993; Wolszczan, 1994; Konacki et al., 1999). Konacki & Wolszczan (2003) applied the secular orbital theory of the PSR 1257+12 system from Konacki et al. (2000) to determine the masses and absolute inclinations of the orbits of the planets B and C. Recently, the same idea of incorporating the effects of mutual interactions between planets to remove the Doppler radial velocity (RV) signal degeneracy present in the $`N`$-Keplerian orbital model (i.e., the undeterminancy of the system inclination and the relative inclination of the orbits) has also been applied to study the strongly interacting, resonant system around Gliese 876 (Laughlin & Chambers, 2001; Rivera & Lissauer, 2001). However, the accuracy of the timing observations by far exceeds the precision of the RV measurements. In consequence, the accuracy of the PSR 1257+12 fit obtained by Konacki & Wolszczan (2003) and hence the initial condition is superior to any initial condition derived from the fits to RV observations of solar-type stars with planets. Thus, the system constitutes a particularly convenient and interesting subject for the studies of its dynamics.

There have been few and limited attempts to determine the long-term stability of the PSR 1257+12 system. Using the resonance overlap criterion and direct integrations, Rasio et al. (1992), Malhotra et al. (1992) and Malhotra (1993b) have found that the system disrupts on $`10^5`$ yr timescale if the masses exceed about 2-3 masses of Jupiter. Also Malhotra et al. (1992) and Malhotra (1993a) have estimated that if the masses were about $`2040`$ M then the system would be locked in the exact 3:2 resonances which would lead to the TOA signal very different from what is actually observed. In this paper, we investigate the stability of the PSR 1257+12 system in the Gyr-scale. We also perform a dynamical comparison of the PSR 1257+12 system to the inner Solar System (ISS) as the dynamics of the ISS is mostly driven by interactions between the telluric planets and is basically decoupled from the dynamics of the outer Solar System (Laskar, 1994, 1997). Finally, we carry out a preliminary analysis of the stability of the orbits of massless particles (e.g. dust particles or Kuiper belt type objects) to determine the zones where such particles could survive and possibly be detected.

## 2. Numerical setup

The stability analysis performed in this paper has mostly numerical character and is understood in terms of the maximal Lyapunov Characteristic Number (LCN). Hence we treat the quasi-periodic orbital motions (LCN $`0`$) as stable and the chaotic ones (LCN $`>0`$) as unstable. In order to resolve the character of the orbits efficiently, we use the so-called fast indicator MEGNO (Mean Exponential Growth factor of Nearby Orbits, Cincotta & Simó, 2000). The MEGNO is directly related to the LCN through the linear relation $`Y=at+b`$. For quasi-periodic motions $`a0`$ and $`b2`$ while for chaotic solutions $`a(\lambda /2)`$ and $`b0`$ where $`\lambda `$ is the LCN of the orbit. However, there is no general relation between the magnitude of the LCN and the macroscopic changes of the orbital elements (e.g., Murison et al., 1994; Michtchenko & Ferraz-Mello, 2001). Even if the LCN is large, the system can stay bounded for a very long time. More insight into the relation of the degree of the irregularity of motion and the macroscopic changes of the orbits one can obtain with the FFT techniques (Michtchenko & Ferraz-Mello, 2001) and its refined variant, the Modified Fourier Transform (MFT) developed by Laskar (1993). The MFT makes it possible to resolve the fundamental frequencies in a planetary system and to determine their diffusion rates (Robutel & Laskar, 2001). In our tests, we use the MFT to verify the MEGNO integrations for massless particles (Section 4) and to identify the mean motion resonances (MMRs). During the computations of MEGNO, for every planet $`k=`$A,B,C, and also for massless particles ($`k=0`$) the complex functions $`f_k(t_i)=a_k\mathrm{exp}\text{i}\lambda _k(t_i),`$ are computed at discrete times $`t_i`$, where $`t_{i+1}t_i10`$ days, over the time $`2T`$ between $`13,000P_\text{C}`$ and $`64,000P_\text{C}`$ (where $`P_\text{C}`$ is the orbital period of the outermost planet C). Here, $`a_k`$ denotes the semi-major axis of the relevant planet and $`\lambda _k`$ is its mean longitude. In general, for a quasi-periodic solution of a planetary system, the frequency corresponding to the largest amplitude, $`a_k^0`$, is one of the fundamental frequencies of motion, called the proper mean motion, $`n_k`$ (Robutel & Laskar, 2001). The 2-body MMRs can be identified through the condition $`q\nu _ip\nu _j0,ij=\text{A,B,C,0},`$ where $`p,q>0`$ are prime integers. The MFT code gives us the rates $`\nu _i/\nu _jp/q`$ and then $`p`$ and $`q`$ are found by the continued fraction algorithm. Such resonances will be labeled as $`\text{P}_1q:\text{P}_2p`$ where $`\text{P}_1`$ and $`\text{P}_2`$ denote the bodies involved in the resonance. In this work we use the code which merges the MFT and MEGNO methods. It incorporates the MFT code which has been kindly provided by David Nesvorny111www.boulder.swri.edu/~davidn/fmft/.

If a planetary system is strongly interacting, MEGNO enables us to detect chaotic behavior typically over only $`10^3`$$`10^4`$ orbital periods of the outermost planet. This feature is vital for examining large sets of the initial conditions (see e.g., Goździewski et al., 2001; Goździewski & Maciejewski, 2003). However, due to very small masses in the PSR 1257+12 system, the mutual perturbations between the planets are very week. In order to resolve the character of such orbits, much longer integration times are required. We have determined that the general purpose integrators (like the Bulirsh-Stoer-Gragg, BSG) used to compute MEGNO are not efficient enough. Therefore, we have developed a symplectic scheme of computing MEGNO based on the idea of Mikkola & Innanen (1999) (see also Goździewski (2003) and Cincotta et al. (2003)). Instead of directly solving the variational equations of the perturbed Kepler problem, one can differentiate the symplectic leap-frog scheme (Wisdom & Holman, 1991) and compute the variations (and MEGNO) using the obtained symplectic tangent map. To use such a technique, a system of canonical variables is required. In these variables, the planetary problem can be divided into integrable Keplerian motions and a perturbation which is integrable in the absence of the Keplerian part (usually it depends only on the coordinates). Our MEGNO code works internally in the Poincaré coordinates (Laskar & Robutel, 1995). As the integrator core, the second- (LR2) and third- (LR3) order schemes by Laskar & Robutel (2001) are used. These integrators are much faster than the BSG algorithm. We can efficiently perform intensive computations of MEGNO ($`10^4`$ points per stability map) for $`10^6P_\text{C}`$, keeping the fractional error of the total energy and the angular momentum below $`10^{10}`$ when the 4 d step (for the LR3 integrator) is used. The tangent map approach is particularly suitable for the pulsar system due to typically small eccentricities of the investigated configurations. The method will be described in detail elsewhere (Breiter and Gozdziewski).

However, one should be aware that even with such an extended integration time-scale, in general we can only detect the strongest chaotic behavior originating mostly from the mean motion resonances between the planets. Investigations of subtle chaos originating in the secular system would require much longer integration times, of the order of thousands of secular periods. Because these periods are in the range $`10^3`$$`10^5`$ yr (see Sect. 3) or longer, the total integration time should be counted in tens of Myr. In such a case, due to the computational limitations, direct approach applied in this paper is no longer practically possible.

In the numerical experiments, we use the initial condition given in Table 1. It represents one of the two orbital configurations obtained by Konacki & Wolszczan (2003) transformed to the classical astrocentric elements. In some of the experiments these elements have been transformed to the canonical heliocentric elements inferred from the Poincaré coordinates. In the second orbital configuration, the inclinations are such that $`i_{\text{B,C}}180^\mathrm{o}i_{\text{B,C}}`$. We assume that $`i_\text{A}`$ obeys the same rule in that solution. It follows that both these solutions are dynamically equivalent and this reflects the fact that using the secular theory as the model of TOA measurements we are still not able to determine the absolute direction of the angular momentum vector. The TOA residuals exclude opposite directions of the outermost orbits (Konacki & Wolszczan, 2003) hence the orbits are either prograde or retrograde.

## 3. Stability of the PSR 1257+12 planetary system

We computed the MEGNO signature of the PSR 1257+12 system over a few Myr ($`510^6`$ $`P_\text{C}`$). This test was repeated for many different, randomly selected initial variational vectors. The results of a few of the runs are shown in Figs 1fe. Both the regular evolution of MEGNO, $`Y(t)`$, as well as the quick convergence of its mean value, $`Y(t)`$, to 2, indicate that the system is very close to a quasi-periodic motion. Let us note, that our preliminary computations of MEGNO with the BSG integrator revealed a very slow divergence of MEGNO after $`0.510^6`$ $`P_\text{C}`$ from the theoretical value of 2. It can be explained by an accumulation of the numerical errors because this effect is absent in much longer runs employing the symplectic code. Even if such a weak instability were to be understood in terms of the chaotic behavior, it would correspond to an extremely long Lyapunov time, $`T_L=1/\lambda `$, comparable to $`10^9`$ yr (estimated through the fit of the relation $`Y=(\lambda /2)t+b`$ over $`t3.510^6`$ yr).

In Fig. 1, we present the time-evolution of the canonical heliocentric elements. Clearly, the motion of the planets B and C is tightly coupled. This is illustrated in panels b and c which show the variations of the eccentricities and inclinations. Fig. 1d demonstrates the presence of a deep, so-called secular apsidal resonance (SAR), described by the argument $`\theta =\varpi _\text{B}\varpi _\text{C}`$ (where $`\varpi =\omega +\mathrm{\Omega }`$ is the longitude of periastron; $`\mathrm{\Omega }`$ and $`\omega `$ are the longitude of ascending node and the argument of periastron of the planet, respectively), with a small semi-amplitude of the librations, $`50^\mathrm{o}`$, and the apsides on average anti-aligned. In the discovery paper, Wolszczan & Frail (1992) noticed the anti-alignment of the lines of apsides while Rasio et al. (1992) and Malhotra et al. (1992) performed first theoretical explorations of the resonance in the framework of the Laplace-Lagrange theory.

The SAR was recently found in many extrasolar systems discovered by the RV measurements. It is widely believed that the SAR, as the orbital state of multi-planetary systems involving Jupiter-size planets, is crucial for maintaining their long term stability (for an overview see, Ji et al., 2003). However in certain cases, the SAR should be understood as a typical feature of a planetary system when the secular angle $`\theta `$ oscillates about $`0^{}`$ or $`180^{}`$. In particular, this concerns the pulsar system: an inspection of Table 3 shows that the proper frequencies of the pericenter motion, $`g_\text{p}`$, are not related through any simple linear combination fulfilling the usual resonance rules. Lee & Peale (2003) and Michtchenko & Malhotra (2004) show that the orbital configurations of two planets are generically stable if the system is far from strong mean motion resonances and collision zones, i.e., when the assumptions of the averaging theorem are fulfilled. In such a case, the secular character of the system depends on the variations of the angle $`\theta `$. Three different modes of its motion are possible: circulation, oscillations around $`0^{}`$ or $`180^{}`$ and oscillations about $`0^{}`$ in the regime of large eccentricities. The first two modes are known from the classical Laplace-Lagrange linear secular theory (Murray & Dermott, 2000) while the third one was discovered by the cited authors. The third mode corresponds to the true, non-linear resonance and the secular system may be unstable in its neighborhood. It should be noted that for the two first modes there is no zero-frequency separatrix between the circulation and oscillation regimes. In this sense, such a SAR is not a true resonance.

Due to a very small mass of the planet A, the pulsar system may be approximated by a 3-body model with equal planetary masses, small eccentricities and inclinations. In such a case, Pauwels (1983) showed that, at least in the linear approximation, the SAR regime covers almost the entire phase space of the secular planetary system, the occurrence of a SAR is almost a certainty and the system is stable. This could explain the stability of the pulsar system. However, some of its specific dynamical features break the assumptions of the secular model. Namely, the proximity of the system to the first-order 3:2 MMR and other mean motion resonances. In order to estimate the influence of the near 3:2 MMR on the secular solution, we integrated the system over 10 Myr sampling the orbital elements every 100 yr. Subsequently, the MFT was applied to analyze the complex signals $`e_\text{p}\mathrm{exp}(\text{i}\varpi _\text{p})`$ and $`\mathrm{sin}(i_\text{p}/2)\mathrm{exp}(\text{i}\mathrm{\Omega }_\text{p})`$ derived from the heliocentric canonical elements. In this way, we can obtain an estimate of the precessional frequencies of the periastrons, $`g_\text{p}`$, and the nodal lines, $`s_\text{p}`$. These frequencies, together with the proper mean motions, $`n_\text{p}`$, compose the set of fundamental frequencies of the pulsar system. Their values and the corresponding periods are given in Table 3. The periastron frequencies, measured in arcsec/yr, are $`197.742`$, $`43.881`$ and $`13.157`$. The first one is substantially different from its analytical approximation (in terms of the Laplace-Lagrange theory). This will lead to a fast discrepancy between the analytical and numerical solutions. Clearly, the classical secular theory has to be modified in order to account for the effects of the near 3:2 MMR. Such a theory has been already developed by Malhotra et al. (1989) for the satellites of Uranus. It can also be applied to the PSR 1257+12 system however it is somewhat beyond the scope of this paper.

Using the symplectic integrator WHM (Wisdom & Holman, 1991) from the SWIFT package (Levison & Duncan, 1994) as well as the LR2 integrator, we also performed a few long term, 1 Gyr, integrations of the PSR 1257+12 system. In the first case, the time step equal to 1 d resulted in the fractional error of the integrals of the angular momentum and the total energy at the level of $`10^{10}`$. The LR2 algorithm with the 4-day integration step resulted in even higher accuracy. During the 1 Gyr interval, no signs of instability are observed. We have not detected any secular changes of the semi-major axes, eccentricities and inclinations. Their values stay within the limits shown in Fig. 1 and the SAR persists with an unchanged amplitude of the librations (Fig. 2a) over 1 Gyr.

Lest us note that from the 1 Gyr integrations it follows that the orbital elements of of the innermost planet A vary in a regular way (Fig. 1). We can also observe a time evolution of the argument $`\theta _1=\varpi _\text{A}\varpi \text{B}`$ as a semi-regular sequence of rotations alternating with irregular ”librations” about $`180^\mathrm{o}`$. This effect is preserved over the entire period of 1 Gyr (see Fig. 2c,d) and may indicate that the orbit crosses the separatrix of a resonance. This is illustrated in Figs 2c,d which show the eccentricity of the planet A (multiplied by 10,000), $`e_\text{A}`$, plotted together with the argument $`\theta _1`$. Clearly, there is a correlation between the librations of $`\theta _1`$ and small values of $`e_\text{A}`$. We also show $`e_\text{C}(\theta )`$ and $`e_\text{A}(\theta _1)`$ collected over the 1 Gyr integration (Fig. 2e,f). The first case represents a trajectory in the resonant island of the SAR whereas there is no clear sign of librations in the $`e_\text{A}(\theta _1)`$ plot. In fact, this effect is only geometrical in nature and can be explained by estimating the secular frequencies of the system using the well known Lagrange-Laplace theory. If the MMRs are absent and the disturbing function is expanded up to the first order in the masses and to the second order in the eccentricities and inclinations, the equations of the secular motion are integrable. Their solution, relative to the eccentricities and the longitudes of periastron, is given in terms of the so called eccentricity vectors, $`(h,k)=(e\mathrm{sin}\varpi ,e\mathrm{cos}\varpi )`$, by (Murray & Dermott, 2000)

$`h_j(t)`$ $`=`$ $`{\displaystyle \underset{i}{}}e_{j,i}\mathrm{sin}(g_it+\beta _i)`$
$`k_j(t)`$ $`=`$ $`{\displaystyle \underset{i}{}}e_{j,i}\mathrm{cos}(g_it+\beta _i),`$
for every planet $`j=`$A,B,C. The constant amplitudes $`e_{j,i}`$ and secular frequencies $`g_i`$ are the eigenvectors and eigenvalues of a matrix with coefficients given explicitly in terms of the masses and constant semi-major axes of the planets. The scaling factors for the eigenvectors $`e_{i,j}`$ and the phases $`\beta _i`$ are determined by the initial condition. Geometrically, the time evolution of every eccentricity vector can be described as a superposition of the eigenmodes corresponding to $`g_i`$ (Malhotra, 1993b). The parameters obtained for the PSR 1257+12 system are given in Table 2. They are consistent with the results of Malhotra (1993a) and Rasio et al. (1992) who analyzed the PSR 1257+12 system involving the two larger companions. They found that the secular evolution of the eccentricity vectors of the two outer planets is almost entirely driven by the eigenstate corresponding to $`g_1`$. Since we have $`e_{B,1}e_{C,1}`$ and the other components of $`e_{j,i}`$ are much smaller and almost equal, it follows that $`\omega _\text{B}\omega _\text{C}+180^{}`$, which corresponds to the SAR of the planets B and C. For the planet A, the eccentricity vector is a superposition of all the three eigenmodes with the leading second and third mode having comparable amplitudes. It appears that using the secular approximation, we can explain the semi-librations of $`\theta _1`$. The minima of $`e_A`$ occur when the modes corresponding to the frequencies $`g_2`$ and $`g_3`$ are in anti-phase. Because the amplitudes $`e_{A,2}`$ and $`e_{A,3}`$ have the same sign and similar magnitudes while $`e_{A,1}`$ has a much smaller magnitude, the following condition has to be satisfied to grant a minimum of $`e_A`$

$$(g_2t+\beta _2)(g_3t+\beta _3)180^{}.$$

Since $`\beta _2\beta _3180^\mathrm{o}`$ (see Table 2), we have $`(g_2g_3)t=360^\mathrm{o}`$ and the period of anti-alignment is $`P_{23}=360^{}/(g_2g_3)43,000`$ yr. Curiously, from the relation between $`g_i`$, $`g_17(g_2g_3)0.001^{}/\text{yr}`$, it follows that this period is almost commensurate with $`P_1`$ (where $`P_1=360^{}/g_1`$).

Near the moment of anti-alignment, the eccentricity of the planet A, $`e_A`$, is driven mostly by the $`g_1`$-mode that is in anti-phase with the eccentricity vector for $`e_B`$. This leads to the quasi-librations seen in Fig. 2 as they repeat with the same period as the variations of $`e_A`$, $`P_{23}`$.

### 3.1. Comparison to the inner Solar System

The results of the 1 Gyr integration, expressed in the canonical heliocentric elements, enable us to analyze the time-evolution of the so-called Angular Momentum Deficit (AMD, Laskar, 1997):

$$C=\underset{p=\text{A,B,C}}{}\frac{m_pm_{}}{m_p+m_{}}\sqrt{\mu _pa_p}(1\sqrt{1e_p^2}\mathrm{cos}i_p),$$

where $`\mu _p=G(m_p+m_{})`$, $`G`$ is the gravitational constant, $`a_p,e_p,i_p`$ is the semi-major axis, eccentricity and inclination of a planet relative to the invariant plane and $`m_{}`$ is the mass of the central body. The AMD indicates the deviation of planetary orbits from a stable circular and coplanar motion for which $`C=0`$. This quantity is preserved by the averaged equations of motion (Laskar, 1997) and its stability provides the stability of the secular system in the absence of short-period resonances. The AMD can be understood as the amplitude of the irregularity present in the averaged system. Large values of AMD lead to a chaotic motion and for a certain critical value to a crossing of the orbits and the disruption of a planetary system (Laskar, 1997, 2000; Michtchenko & Ferraz-Mello, 2001).

The PSR 1257+12 system is close to the 3:2 commensurability. Its critical argument circulates but it does not mean that the time-averaged effects of the near-resonance vanish (Malhotra et al., 1989; Malhotra, 1993a). Hence, the applicability of the AMD signature in the studies of the stability in the real, unaveraged system can be problematic. Instead, the integral obtained by averaging the quasi-resonant system should be applied (Michtchenko & Ferraz-Mello, 2001). Nevertheless, we have decided to calculate the AMD and to examine its behavior in order to compare the results with those obtained for the Solar system (which contains two planets close to 5:2 MMR, see e.g., Ito & Tanikawa, 2002). The time evolution of AMD in the pulsar system is shown in Fig. 2. The AMD stays well bounded and very regular. There is very little exchange of the AMD between the innermost planet A and the pair of the bigger companions B and C. It suggests that the motion of the planet A is decoupled from the dynamics of the planets B and C in the long-term scale. This situation is qualitatively different from the inner Solar system (ISS). The AMD of the ISS is not strictly preserved due to the perturbations of Jupiter and Saturn, nevertheless it can be considered roughly constant (Laskar, 1997; Ito & Tanikawa, 2002). Ito & Tanikawa (2002) published the results of 5 Gr integrations of the ISS which reveal rapid AMD variations of Mercury. They are much larger than the changes of AMD for the Venus-Earth pair. The variations of the AMD of Venus, Earth and Mars are also irregular and substantial. In fact, the ISS is chaotic having the Lyapunov time of about $`5`$ Myr (Laskar, 1994). Laskar (1994) found that this chaos is physically significant as it can lead to the ejection of Mercury from the inner ISS during a few Gyr. However, the source of chaos in the ISS is still not well understood (Lecar et al., 2001). In this sense, the character of the motion of the PSR 1257+12 system is quite different. The evolution of its AMD is very regular in spite of smaller distances, larger masses and thus stronger mutual interactions between the planets. Unlike Mercury’s, the AMD of the planet A is negligible when compared to the AMD of the B-C pair as it contributes only about 1/1000 of the total value. On the other hand, the orbital coupling of the pairs B-C and Venus-Earth (Ito & Tanikawa, 2002) is a similar feature of the inner Solar and PSR 1257+12 systems although it has a different nature. The long term integrations by Ito & Tanikawa (2002) revealed an anti-correlation between the changes of the orbital energies of Venus and Earth and, simultaneously, a correlation between the changes of their orbital angular momenta and the eccentricities. These effects can be explained through the influence of the Jovian planets. Hence their dynamical source is external for the coupled planets while the coupling of the B-C pair in the PSR 1257+12 system is provided by the anti-aligned SAR. Finally, the chaotic orbital evolution of the ISS may significantly depend on the weak coupling with the outer planets. An equivalent effect is obviously absent in the PSR 1257+12 system.

### 3.2. Dynamical environment of the PSR 1257+12 planetary system

In the next set of experiments, we look at the initial condition in a global manner in order to find out whether the current state of the system is robust to the changes of the formal initial condition (IC). However, thanks to a very precise determination of the initial condition such changes, if considered consistent with the TOA measurements, are very limited. For example, the formal $`1\sigma `$ error of the semi-major axes inferred from the parameter $`x^0`$ (Konacki & Wolszczan, 2003) is at the level of $`10^6`$ AU! Nevertheless, the localization of the IC in the phase space (e.g. its proximity to unstable regions) is critical to verify its character. The examination of one isolated IC does not provide a definitive answer to the question of stability.

#### 3.2.1 Mean motion resonances

We computed one-dimensional scans of MEGNO, $`Y`$, by changing the semi-major axis of one planet and keeping the other orbital parameters fixed at the values given in Table 1. The results of this experiment are illustrated in Fig. 3. The MEGNO scans were computed with the resolution of $`510^6`$ AU and $`10^6`$ AU for close-up scans. They reveal a large number of spikes, some of them very close to the nominal positions of the planets marked with the large filled dots. Most of these spikes can be identified as 2-body MMRs between the planets. Currently, the IC is well separated from low-order MMRs. But some unstable high-order MMRs appear close to the nominal IC. The most relevant ones seem to be C31:B21 MMR (i.e. between the planets B and C, see the bottom panel in Fig. 3) and C53:B36 MMR (shown in the scan for the planet B). Also the planet A is close to B29:A11 MMR. Yet, such resonances are extremely narrow. Their widths can be roughly estimated as less than $`510^6`$ AU what follows from their shapes and the resolution of the MEGNO scans. It is unlikely that they can affect the regular motion of the system at its nominal position. However, altering the semi-major axis of the planets by small shifts, $`510^5`$ AU, would push the system into these unstable regimes. In order to demonstrate it, we calculated MEGNO for a configuration corresponding to C31:B21 MMR (see the right column of Fig. 4). This resonance is the inclination-type MMR (according to the terminology by Peale, 1976). Clearly, its critical argument $`\sigma =31\lambda _\text{C}21\lambda \text{B}10\mathrm{\Omega }_\text{B}`$ exhibits a sequence of alternating librations and circulation giving rise to the unstable behavior. Let us note, that the chaos is formally extremely strong. The Lyapunov time, estimated by the linear fit to the MEGNO, is only about 3000 yr (see Fig. 4a).

The 2-body MMRs cannot explain all the peaks of MEGNO present in the scans. Some of them seem to be related to the 3-body resonances involving all three planets. The importance of 3-body MMRs was studied by Nesvorný & Morbidelli (1998); Morbidelli & Nesvorny (1999), Murray et al. (1998) and Murray & Holman (1999). They found that weak 3-body MMRs can strongly influence asteroidal motion and explain the short Lyapunov times which cannot be understood if only 2-body MMRs of high-order are considered. Murray & Holman (1999) proved that the chaotic behavior of the Outer Solar system is governed by the overlapping 3-body MMRs involving Jupiter, Saturn, and Uranus. It is obviously interesting to verify if some of the instabilities visible in the $`a`$-scans can be related to such 3-body resonances.

A 3-body resonance can be defined by the following condition (Nesvorný & Morbidelli, 1998):

$$i_\text{A}\dot{\lambda }_\text{A}+i_\text{B}\dot{\lambda }_\text{B}+i_\text{C}\dot{\lambda }_\text{C}0,$$

where $`\dot{\lambda }\text{p}`$ is the mean motion of the given planet. A critical argument of such resonance is a linear combination of the longitudes, arguments of pericenters and nodal longitudes:

$`\sigma _{i_\text{A}:i_\text{B}:i_\text{C}}`$ $`=`$ $`i_\text{A}\lambda _\text{A}+i_\text{B}\lambda _\text{B}+i_\text{C}\lambda _\text{C}+`$
$`p_\text{A}\varpi _\text{A}+p_\text{B}\varpi _\text{B}+p_\text{C}\varpi _\text{C}+q_\text{A}\mathrm{\Omega }_\text{A}+q_\text{B}\mathrm{\Omega }_\text{B}+q_\text{C}\mathrm{\Omega }_\text{C}`$
whose integer coefficients $`i_\text{p},p_\text{p},q_\text{p}`$ fulfill the d’Alambert rule $`_{\text{p=A,B,C}}(i_\text{p}+p_\text{p}+q_\text{p})=0`$ and the usual requirement of rotational symmetry. Due to the small ratios of the secular frequencies to the mean motions in the PSR 1257+12 system, $`10^5`$, the resonance condition may be approximated by $`i_\text{A}n_\text{A}+i_\text{B}n_\text{B}+i_\text{C}n_\text{C}0.`$ For example, one of the MEGNO peaks for $`a_\text{C}=0.46601AU`$, can be explained as the combination of the proper mean motions $`3n_\text{A}14n_\text{B}+9n_\text{C}0.01^{}/\text{yr}`$. Hence it corresponds to the 3-body MMR 3:-14:9. It is separated from the the nominal position of the planet C by only $`310^5`$ AU which is at the $`3\sigma `$ error level of $`a_\text{C}`$. The time evolution of the critical argument, $`\sigma _{3:14:9}`$, shown in the left column of Fig. 4, confirms the slow changes of this angle as well as a presence of quasi-librations alternating with circulations. This indicates a chaotic evolution. The MEGNO computed for the corresponding initial condition does not grow fast, at least over a few Myr used in this test. This mean that the chaos is moderate. Nevertheless, the long term stability of the system (likely, the chaos would affect mostly the motion of the planet A) could be confirmed only by direct numerical integrations. Our 1 Gyr integrations for the ICs of the two MMRs near the planet C do not reveal any secular changes of the elements. Likely, in order to properly investigate subtle effects of these resonances and their influence on the system, a refined model of the motion should be used. Such model should incorporate the effects of general relativity, a possible error in the pulsar mass and other factors, like the dependence on the observationally unconstrained orbital elements of the planet A.

Other low-order 3-body resonances that can be identified in the Fig. 3 are 1:-2:1, 1:-4:2, -2:10:7 (almost overlapping with B71:A27), 3:-16:12, 1:-6:5 (see the $`a_\text{A}`$-scan). Finally, also an overlapping of 2-body MMRs is possible. One such example is marked in the $`a_\text{B}`$-scan where in the neighborhood of the 3:2 MMR, a sharp peak of MEGNO at the simultaneous position of B60:A23 and C76:B51 MMRs is present. All these MMRs are very narrow since the width of 3-body MMRs is proportional to the masses in the second order and these 2-body MMRs are of high order. They are relatively distant from the nominal positions of the planets thus it is unlikely that they can affect the motion of the system.

Finally, Fig. 4 shows the effect of varying the assumed mass of the neutron star on the determination of the semi-major axis of the outermost planet and the localization of the resonances. Obviously, the sequence of MMRs does not change. A different mass of the pulsar can lead to a substantial shift of the initial semi-major axis (in our case, $`a_\text{C}`$) and the weak unstable resonances, like the C31:B21 MMR, may end up much closer to the nominal positions of the planets (for smaller than canonical mass of the parent star, $`M_{\text{psr}}=1.4M_{\mathrm{}}`$).

#### 3.2.2  Stability in the $`(a_\text{p},e_\text{p})`$\- and $`(e_\text{B},e_\text{C})`$-planes

The results of the 2-dimensional MEGNO analysis are shown in Fig. 6. In these maps we analyze the influence of the initial eccentricities of the planets on their motion. We used the initial conditions from Table 1 and changed those elements that correspond to the coordinates in the stability maps. Due to a lower resolution (400 points in the semi-major axis range), some features present in the one-dimensional scans are not so clearly visible. Nevertheless, in the $`(a_\text{B},e_\text{B})`$-plane (the middle panel), the thin strips of the MMRs can be easily identified. These maps additionally reveal the width of the strongest MMRs (like C3:B2, C13:B9 and C10:B7) and the border of the global instability. These stability maps confirm that the nominal positions of the two bigger planets, marked with filled circles, are far from strong instabilities of the motion. Also in the $`(e_\text{B},e_\text{C})`$-plane, the nominal IC lies in a wide stable region, far from zones of chaotic motion (see Fig. 7). The $`(a_\text{A},e_\text{A})`$map for the planet A shows a dense net of narrow unstable regions, some of which are close to the nominal position of the planet. They have been already identified in the 1D scans for $`a_\text{A}`$. Let us note that similar maps were computed by Ferraz-Mello & Michtchenko (2002) for the planet B with their FFT fast indicator. Their results are generally consistent with ours even though they used a different initial condition of the pulsar system.

During the MEGNO integrations we have also computed the maximum value of the SAR argument $`\theta =\varpi _\text{B}\varpi _\text{C}`$ with respect to the libration center of $`180^\mathrm{o}`$. It enables us to estimate the semi-amplitude of librations and the extent of the SAR in the space of the scanned orbital parameters. The result of such an experiment conducted in the $`(a_\text{B},e_\text{B})`$ plane are shown in Fig. 8. It uncovers an extended zone of the SAR with the smallest semi-amplitude of $`\theta `$ in the vicinity of the nominal initial condition of the PSR 1257+12 system. The near 3:2 MMR is also clearly present in this map.

#### 3.2.3 Limits on the unconstrained parameters of planet A

So far, after Konacki & Wolszczan (2003), we were using the average inclination of the orbits of the planets B and C as the orbital inclination of the planet A. However, we can try to verify whether the dynamics can provide any limits on the inclination and hence mass of the innermost planet A. To this end, we computed MEGNO for the inclination of the planet A in the range $`[30^\mathrm{o},70^\mathrm{o}]`$ and the position of its nodal line in the range $`[90^\mathrm{o}:90^\mathrm{o}]`$, with the resolution of $`1^{}`$ in both coordinates. Thus we varied both the mass (to preserve $`m\text{A}\mathrm{sin}i_\text{A}`$ determined from the TOAs) as well as the relative inclination of the planet A to the orbital planes of B and C. The result of this experiment is shown in Fig. 9. The MEGNO scan (the left panel of Fig. 9) reveals a well defined stable region. The assumed position of the planet A is in the center of this zone. To obtain this picture, an extended integration time span of $`10^6P_\text{C}`$ was used. This choice was dictated by the analysis of a few orbits of A which are substantially inclined to the orbital planes of B and C. An example corresponding to the initial $`\mathrm{\Omega }\text{A}=75^{}`$ (thus located in the unstable region) is shown in Fig. 10. For the integration time up to about 0.5 Myr ($`510^5P_\text{C}`$) MEGNO stays close to 2 but then it suddenly grows, indicating chaotic behavior. To explain its source, we computed the eccentricity, inclination $`i_\text{A}^{\text{rel}}`$ and the argument of pericenter, $`g_\text{A}`$, of the planet A relative to the invariant plane. Evolution of these elements is shown in Fig.10. Clearly, both $`e_\text{A}`$ and $`i_\text{A}^{\text{rel}}`$ exhibit long-term, large amplitude oscillations which are exactly in an anti-phase. The period of these oscillation is relatively very long, of the order of $`10^4`$ yr. Simultaneously, $`g_\text{A}`$ temporarily librates about $`90^{}`$ or $`270^{}`$. These librations indicate that the pericenter precession of A stops. Such features are typical for the well known Kozai resonance (Kozai, 1962) found for highly inclined asteroids in the Solar system. The argument of pericenter $`g_\text{A}`$ can be treated as the critical argument of this resonance. Because librations are followed by circulations of the critical argument, the orbit crosses the separatrix which explains why the motion eventually becomes chaotic.

The effects of the Kozai resonance are illustrated in the right panel of the $`e_\text{A}^{\text{max}}`$-map (the maximal eccentricity attained by the planet A during the integration time) in Fig. 9. The sharply ending zone of moderate eccentricities is narrower than the stability region. We should be aware that in the transient areas the integration time could still be to short to detect the irregular motion. Note that around the border of the stable zone, $`e_\text{A}`$ can be as large as 0.6–0.8. This resonance puts significant limits on the position of the nodal line of the planet A as well as its relative inclination to the orbital planes of B and C. Assuming that the motion of all three planets is in the same direction, the nodal line of A cannot be separated by more than about $`\pm 60^{}`$ from the nodes of B and C and the relative inclination should not exceed $`45^{}`$. Otherwise, the eccentricity of A would be excited to large values, $`>0.7`$ over only $`10^5`$ yr and for $`e_\text{A}0.8`$ the orbits of A and B would cross or close encounters would become possible.

## 4. Dynamics of massless particles

While investigating the dynamics of extrasolar planetary systems one can ask whether the minor bodies can survive in the gravitational environment of the giant primary bodies. These can be small telluric planets in the systems with Jupiter-like planets. We can ask the same question about asteroids, cometary bodies or dust particles in the PSR 1257+12 system. Because the motion of the planets appears to be strictly stable, we can consider a simplified, restricted model. In this model one assumes that a probe mass moves in a gravitational tug of the primaries but does not influence their motion. In the simplest case such a model is well known as the restricted, planar circular 3-body problem (Murray & Dermott, 2000). This simplification helps us to explore the phase space of the planetary system with reasonable computational effort. In this model, the MEGNO indicator and the MFT are evaluated only for the probe mass. In the numerical experiments, we investigate the orbital stability of massless particles in a few regions of the orbital space: between the planets A and B (region I), B and C (region II) and beyond the planet C (region III).

In the first series of numerical tests, we vary the initial semi-major axis of the probe mass, $`a_0`$, and its eccentricity, $`e_\text{0}`$, while the initial inclination is constant, $`i_\text{0}=50^\mathrm{o}`$. Thus the tested orbits are slightly inclined to both orbital planes of the planets B ($`i_\text{B}=53^\mathrm{o}`$) and C ($`i_\text{C}=47^\mathrm{o}`$). The angular variables of the massless particles are set to $`\mathrm{\Omega }_\text{0}=\varpi _\text{0}=M_\text{0}=0^\mathrm{o}`$. This way we essentially follow the remarkable work of Robutel & Laskar (2001) who investigated the short-term dynamics of massless particles in the Solar System using the MFT. These authors argue that the global picture of the motion in the restricted problem does not qualitatively depend on the initial phases of the test particle. The initial plane $`(a_0,e_0)`$ is representative for the dynamics because it crosses all resonances. By changing the initial orbital phases of the massless particles, the width of the resonances may vary but their influence on the motion can still be detected.

In order to check the MEGNO signatures (e.g., to verify whether the integration time is not too short), during the integrations, we compute the diffusion rate, $`\sigma _0=1n^{(2)}/n^{(1)}`$, where $`n^{(1)}`$ and $`n^{(2)}`$ are the mean motion frequencies obtained over the intervals $`[0,T]`$ and $`[T,2T]`$, respectively, where $`T`$ (which we will call the base period from hereafter) has been set accordingly to the investigated range of the semi-major axes. The close encounters with the planets are controlled. We assume that if the distance between the probe mass and a planet is less than the Hill radius, $`r_H=[m_p/(3m_{})]^{(1/3)}a_p`$, then the particle has collided with the planet or its motion become strongly chaotic.

Within this numerical setup, in the region II we have not detected any stable motions. The region I turns out to be a much more interesting one. The results are shown in the MEGNO and $`\sigma `$ maps in Fig. 11 (left panels for the region I, middle and right panels for the region III). The integrations have been carried out for the total time, $`2T`$, of about 32,000 $`P_\text{C}`$. In both maps we can clearly recognize the planetary-crossing lines marked by large values of both indicators which correspond to strongly chaotic motions. These lines (white lines in Fig. 11) are the solutions to the equations $`a_p(1+e_p)=a_0(1e_0)`$ for $`a_0>a_p`$ and $`a_p(1e_p)=a_0(1+e_0)`$ for $`a_0<a_p`$ where $`a_p,e_p,a_0,e_0`$ denote the initial semi-major axis and the eccentricity of a planet and a test particle, respectively.

In this area and above the collision lines, the particles are scattered mostly by close-encounters and collisions with the planets A and B. Remarkably, a stable area exists under the crossing-lines for the planets A and B. It is divided by a number of MMRs with the planets. In order to verify this result, we have integrated the motion of 200 particles spread over $`[0.19,0.35]`$ AU with the same resolution in the semi-major axis, $`a_0`$, as in the 2D MEGNO test. These integrations were continued up to 15 Myr (here, we used the RMVS3 integrator from the SWIFT package). A resulting one dimensional scan of MEGNO over $`a_\text{0}`$ is shown in Fig. 12. In the regions classified with MEGNO as regular, the massless particles have mostly survived during the integration time while in the chaotic or close to chaotic areas, they have been quickly removed or collided with the planets (in this test the same criterion of 1 Hill radii for a collision event is used and in such a case we set $`\mathrm{log}\sigma _0=1`$). Some discrepancies are most likely related to the different integrators used in the experiment; the Bulirsh-Stoer integrator follows the orbits of massless particles much more accurately that the symplectic integrator. This experiment enables us also to independently estimate the proper MEGNO integration time and to ”calibrate” the scale of the diffusion rate, $`\sigma _0`$. Comparing the MEGNO and $`\sigma _0`$-maps, we conclude that for $`\mathrm{log}\sigma `$ less than about $`(7)`$$`(8)`$ the motion can be considered as regular (also compare the maps in the left column of Fig. 11). Let us note that both indicators are in an excellent accord and the diffusion rate calculated by the MFT algorithm seems to be even more sensitive to an unstable motion than the MEGNO is.

Using a similar numerical setup, we have investigated the region III. However, due to an extended range of the semi-major axes, we have divided it into two parts: the region IIIa between 0.5 AU and 0.97 AU and the region IIIb between 0.9 AU and 3.9 AU (a distance comparable to the radius of the Asteroid Belt in the Solar System). The maps of MEGNO and the diffusion rate, $`\mathrm{log}\sigma `$, in the $`(a_0,e_0)`$-plane for the region IIIa (obtained after the time $`2T32,000`$ P$`_\text{C}`$), are shown in the middle column of Fig. 11 while the maps for the region IIIb ($`2T64,000`$ P$`_\text{C}`$) are shown in the right column. For the region IIIa, in both maps, we can clearly recognize the collision line with the outermost planet while the crossing line with the planet B is in the zone of strong chaos. A number of MMRs appear as narrow vertical strips. This test shows that for a moderately low initial eccentricity, $`e_0`$, the stable zone is extended and begins just beyond the orbit of the planet C. Obviously, with a growing $`e_0`$ the border of the stable zone shifts toward a larger $`a_0`$ (at the distances of about 1 AU the zone of stability reaches $`e_00.5`$). For the region IIIb ($`a_0[0.9,3.9]`$ AU), the results are shown in the right column of Fig. 11. For efficiency reasons, in this test the step in $`e_0`$ is $`0.1`$ (the resolution in $`a_0`$ was left relatively high, at $`0.005`$ AU). The collision lines are clearly present. Otherwise, this zone is mostly regular even for very large $`e_0>0.5`$. Similarly to the region IIIa, the scans reveal some extremely narrow unstable MMRs, most of them with the two outer planets.

Finally, we examined the stability of orbits inclined to the mean orbital plane of the system. In this experiment we tested the motion of particles in the region III. The results for the region IIIa are shown in Fig. 13 (a number of additional scans, not shown here, allow us to extrapolate the results obtained for this region to the region IIIb). In the first two experiments, we scanned the $`(a_0,e_0)`$ space for the initial $`i_0`$ set to $`75^\mathrm{o}`$ and $`87^\mathrm{o}`$, respectively. These inclinations, taken relative to the mean orbital plane of the PSR 1257+12 system, are comparable to the inclinations of the Kuiper belt objects in the Solar System. Again, the collision lines and the net of MMRs clearly appear in both the MEGNO and $`\sigma _0`$ scans which are shown in the the left and the middle panels of Fig. 13. The zone of strong chaos covers the crossing zones with the planets B and C while the collision line with the planet A is much more narrow and separated by a quasi-regular area in which $`\sigma _010^6`$. Such an effect has been already observed by Robutel & Laskar (2001) — for higher inclination the close encounters with the planets are less frequent than for moderately small inclinations which explains the smaller extent of the strong chaos. Moreover, after sufficiently long time the massless bodies will be removed from above the collision lines, excluding cases when the probes are trapped within stable MMRs. In the next test illustrated in the right panel of Fig. 13, $`e_0`$ was set to 0 and the initial inclination $`i_0`$ was varied in the wide range $`[10^\mathrm{o},90^\mathrm{o}]`$. This experiments allow us to generalize the results of the previous scans. The stability of inclined particles is basically independent on the initial inclination, at least for moderately small eccentricities.

One should be aware that the results of the above analysis should be treated as preliminary ones. Our study is restricted both to the short-term dynamics and a small part of the possible volume of the parameter space. However, the timescale of the integrations is long enough to detect the most unstable regions in the phase space and to point out the regions where the particles can be long-term stable. Also in some chaotic regions the particles can still persist over a very long time but this can be verified only by direct integrations. Let us also note that in order to derive rigorous estimates of the proper mean motion frequencies one has to employ angle-action like coordinates, e.g., the Poincaré coordinates (Laskar & Robutel, 1995) or the Jacobi coordinates. Our integrations for mass-less particles were carried out in the Keplerian, astrocentric coordinates and the osculating elements analyzed by MFT are inferred from these non-canonical coordinates. Due to very small masses involved, the effects caused by the use of non-canonical coordinates are negligible.

## 5. Conclusions

In this work we carry out numerical studies of the stability of the PSR 1257+12 system using the initial condition determined by Konacki & Wolszczan (2003). The long term integrations utilizing the symplectic integrators, extended over 1 Gyr, do not reveal any secular changes in the semi-major axes, eccentricities and inclinations of the planets. Using the notion of the Angular Momentum Deficit (AMD), we do not find any substantial exchange of the angular momentum between the innermost planet and the pair of the outer, bigger planets B and C. The AMD of the planet A is negligible when compared to the AMD of the B-C pair. This is very different from the case of the inner Solar System in which the variations of the AMD of Mercury are the most significant ones. The PSR 1257+12 system has the MEGNO signature typical for a strictly regular, quasi-periodic configuration. The two outer planets are close to the 3:2 mean motion resonance and are orbitally tightly coupled. The presence of the secular apsidal resonance is quite typical for such system as demonstrated by Pauwels (1983) and recently by Lee & Peale (2003) and Michtchenko & Malhotra (2004). The semi-amplitude of the critical argument is about $`50^\mathrm{o}`$ and it persists in a wide range of the orbital initial parameters. The SAR in the PSR 1257+12 system is yet another such case among extrasolar planetary systems (Ji et al., 2003).

The neighborhood of the nominal initial condition is investigated by calculating the MEGNO signature in a few representative regimes of the semi-major axes and eccentricities of the planets. These stability maps reveal that the nominal initial condition is located in an extended stable zone, relatively far from any strong instabilities of the motion. However, numerous weak mean motion resonances can be found in close proximity to the nominal positions of the planets. These are both 2-body resonances between the planets (like the 31:21 MMR between B and C) and 3-body resonances, among which the 3:-9:14 MMR seems to be the most relevant one. Their potential influence on the motion could be investigated if the initial condition of the system was refined using a more accurate model of the dynamics, possibly including relativistic effects. However, it would most likely require much more precise TOA measurements than currently available.

These factors allows us to state that the PSR 1257+12 system is orbitally stable over the Gyr time scale. In our experiments, there are no signs of a potential instability except for a very slow divergence of the MEGNO in few of the tests. This divergence of MEGNO corresponds to the Lyapunov time of about $`1`$ Gyr. We believe that it has only a numerical character. We are aware that an alternative analytical study of the PSR 1257+12 dynamics is possible (also thanks to the accurate determination of the initial condition). Such an approach has been proposed in Malhotra et al. (1989) and Malhotra (1993a). Our numerical investigations can certainly be treated as a complement to any future analytical studies of the system.

Using the MEGNO analysis, we found dynamical limits on the unconstrained elements of the planet A. Due to the destabilizing effect of the Kozai resonance, the nodal line of this planet cannot be separated by more than about $`\pm 60^{}`$ from the nodes of B and C. Otherwise, the eccentricity of A would be excited to large values permitting close encounters or collisions with the planet B. It constitutes a strong dynamical argument that the orbital plane of the planet A indeed coincides with the mean orbital plane of the system.

Finally, using the MEGNO and the MFT, we investigate the dynamics of massless particles in the PSR 1257+12 system in the framework of the restricted model. In the numerical experiments, we find a stable zone between the planets A and B extending for initially small eccentricities from 0.19 AU to 0.25 AU from the pulsar. There are no stable areas between the planets B and C. Beyond the orbit of the planet C, the stable zone begins already outside of its orbit. We find that the massless particles can move on stable orbits under the condition that their initial eccentricities and semi-major axes are located under the collision lines with the planets. The dynamics of massless particles is basically independent on the their initial inclinations. Beyond 1 AU, the motion appears to be stable except for the areas of narrow MMRs with the planets B and C. It is an encouraging result supporting the search for possible small bodies contained in a dust or Kuiper belt type disk around the PSR 1257+12 .

## 6. Acknowledgments

We are indebted to the anonymous referee for many invaluable suggestions and a detailed review which improved the manuscript. We thank Philippe Robutel for a discussion and helpful remarks. K. G. is supported by the Polish Committee for Scientific Research, Grant No. 2P03D 001 22. M. K. is a Michelson Postdoctoral Fellow.


format: text

text:
# Molecular Gas Chemistry in Starbursts and AGN

## 1. Evolution of Molecular Gas Chemistry in SB and AGN

Multi-wavelength based evidence indicates that the properties of molecular gas in SB and AGN hosts differ from that of quiescent star forming disks. The spectacular energies injected in the gas reservoirs of SB and AGN coming as strong radiation fields (UV, X-rays,…), powerful winds and jets can create a particularly harsh environment for ISM. A complete understanding of the physical/chemical evolution of molecular gas in ‘active’ galaxies requires the use of specific tracers of the relevant energetic phenomena that are at work all the way along the starburst sequence: large-scale shocks, strong UV-fields and nuclear X-ray irradiation. We have used the high-resolution/sensitivity of the IRAM 30m telescope and the PdBI to study a limited sample of prototypical SB and AGN.

### 1.1. Unveiling Shock Chemistry in SB

We have carried out a 30m multi-transition survey searching for the thermal emission of silicon monoxide (SiO) in a dozen prototypical SB galaxies including NGC 253, M 82 and NGC 1068. Studies of galactic clouds point out to SiO as a privileged tracer of shocks in the molecular gas phase. Shocks can significantly increase SiO fractional abundances to X(SiO)$``$10-8. Our survey shows that SiO emission is widespread in the circumnuclear disks (CND) of SB on scales ranging from 100 pc to 700 pc. The estimated SiO abundances vary within the sample from X(SiO)$``$10-9 in NGC 253 to 1/50 of this value in M 82. Physical parameters of SiO clouds are also very different, suggesting different scenarios for shock chemistry.

The first PdBI SiO maps obtained in NGC 253 and M 82 have given invaluable insight into the different mechanisms driving large-scale shocks in the molecular gas disks of SB (García-Burillo et al. 2000, 2001). While SiO emission is detected mainly in a 700 pc-diameter CND in NGC 253 (García-Burillo et al. 2000), the emission of SiO extends noticeably out of the galaxy plane in M 82, tracing the disk-halo interface where episodes of mass injection from the disk are building up the gaseous halo (García-Burillo et al. 2001; see Fig. 1). Large-scale shocks are driven by massive star formation and bar density waves inside the disk of NGC 253. In M 82, however, shocked molecular gas appears forming a 500 pc chimney and a giant supershell. The strikingly different distributions and average fractional abundances of SiO in NGC 253 and M 82 are suggestive of an evolutionary link between these SB; the latter pictures the M 82 starburst as a more evolved episode.

More recently, the detection of SO2,NS and NO emission in NGC 253 have confirmed that shock chemistry is at work in the nucleus of NGC 253 (Martín et al 2003, submitted).

### 1.2. Propagation of PDR Chemistry in SB

We have mapped the emission of the formyl radical (HCO) in the nucleus of M 82 (García-Burillo et al. 2002). HCO is known to be enhanced in the interfaces between the ionized and molecular gas, making it a privileged tracer of PDR. The 5” HCO map of M 82, the first obtained in an external galaxy, shows a ring-like distribution, also displayed by other molecular/ionized gas tracers in this galaxy. Most remarkably, the rings traced by HCO, CO and HII regions are nested, with the HCO ring lying in the outer edge of the molecular torus (see Fig. 2). The high overall abundance of HCO in M 82 ($``$4$`\times `$10-10) indicates that its nuclear disk has become a giant PDR of $``$650 pc size. Furthermore, the existence of a nested ring pattern with the highest HCO abundance being estimated in the outer ring ($``$0.8$`\times `$10-9), suggests that PDR chemistry is propagating in the disk.

The PdBI maps of M 82 made in SiO and HCO illustrate how two different gas chemistry scenarios can be simultaneously at play in the same galaxy. The strong UV fields of the M 82 starburst have created a giant PDR inside the disk, while the expansion of hot gas created by successive SN explosions entrains neutral gas into the halo driving shocks located in the disk-halo interface of the galaxy.

### 1.3. XDR Chemistry in AGN: the nucleus of NGC 1068

Molecular gas in the CND of AGN hosts can be exposed to strong X-ray irradiation. Contrary to UV photons, X-rays can penetrate gas column densities out to N(H2)$``$1023-24cm-2 before being attenuated. Therefore, XDR could be the dominant sources of emission for the molecular gas in the vicinity of AGN (Maloney et al. 1996). Tantalizing evidence that the chemistry of molecular gas in the CND of AGN is ‘exotic’ came from the large HCN/CO abundance ratio measured in the nucleus of the Seyfert 2 galaxy NGC 1068 (Tacconi et al. 1994), first interpreted as a signature of enhanced oxygen depletion in the NLR.

In the course of our ongoing 30m survey, we detected SiO emission in the starburst ring of NGC 1068. Most remarkably, we detected also SiO emission coming from the CND torus of NGC 1068, i.e., from a source mostly unrelated to recent star formation. We derived an SiO abundance enhanced out to $``$10-9 in the CND. Silicon chemistry in the CND of NGC 1068 is driven either by X-rays or by violent shocks near the central engine. To bring some light into the ‘obscuring torus chemistry’ we made complementary observations with the 30m telescope and PdBI for eight molecular species, purposely chosen to fully explore the predictions of XDR models for the molecular gas phase. Observations included several lines of CN, HCO, H13CO+, H12CO+, HOC+, HCN, CS and CO. A first analysis of this survey, presented by Usero et al. (2003), has shown that the bulk of the molecular gas emission in the CND of NGC 1068 can be interpreted as coming from a giant XDR.

## 2. Conclusions

The advent of highly sensitive interferometers has made possible to study the complex molecular inventory of galaxies, going beyond the classical ‘CO maps’. In particular, the information provided by molecular gas tracers of peculiar chemistry scenarios such as PDR, XDR, and Shock Chemistry is a fundamental tool to explore the physical/chemical evolution of the ISM content in SB and AGN. In the course of the combined IRAM 30m+PdBI survey we have studied a limited sample of SB and AGN. The first findings have revealed already clear study cases of large-scale shocks at work (NGC 253, M 82), propagation of PDR chemistry in SB (M 82) or XDR chemistry driven by AGN (NGC 1068).

## References

García-Burillo, S., Martín-Pintado, J., Fuente, A., & Neri, R. 2000, A&A, 355, 499

García-Burillo, S., Martín-Pintado, J., Fuente, A., & Neri, R. 2001, ApJ, 563, L27

García-Burillo, S., Martín-Pintado, J., Fuente, A., Usero, A., & Neri, R. 2002,ApJ, 575, L55

Maloney, P. R., Hollenbach, D. J., & Tielens, A. G. G. M. 1996, ApJ, 466, 561

Martín, S., Mauersberger, R., Martín-Pintado, J., & García-Burillo, S., 2003, A&A, submitted

Usero, A., García-Burillo, S., Martín-Pintado, J., Fuente, A., Rodríguez-Fernández, N., & Rodríguez-Franco, A. 2003, A&A, submitted

Tacconi, L. J., Genzel, R., Blietz, M., Cameron, M., Harris, A. I., & Madden, S. 1994, ApJ, 426, 77

Wills, K., Redman, M. P., Muxlow, W. B., & Pedlar, A. 1999, MNRAS, 309, 395


text:
# Contents

## 1 Introduction

From an analysis  using boundary CFT (Conformal Field Theory) techniques  it has been possible to deduce different properties of the open-string tachyon222 Recent progress on both open-string and closed-string tachyon dynamics include .. It was subsequently attempted to reproduce these properties from an effective field theory. The “standard form” Lagrangian333In this paper, the flat metric is $`\eta _{\mu \nu }=(1,+1,\mathrm{},+1)`$ and $`\alpha ^{}=1`$.

$$L=\frac{1}{\mathrm{cosh}\left(\frac{\alpha \stackrel{~}{T}}{2}\right)}\sqrt{1+_\mu \stackrel{~}{T}^\mu \stackrel{~}{T}}$$
(1.1)
was proposed , based on its ability to reproduce some of the properties of the tachyon known from the boundary CFT analysis. Here, $`\alpha =1`$ describes the bosonic tachyon, while $`\alpha =\sqrt{2}`$ corresponds to the superstring case.

Recently, Kutasov and Niarchos  managed to derive this action in the superstring case with an ansatz of the type $`L(T,_\mu T^\mu T)`$ using only two requirements: First, that the action should allow the known tachyon solutions. And second, that it should be consistent with certain results from boundary SFT (String Field Theory)444 A similar argument appeared already in ..

The purpose of this paper is to try to extend the derivation to the bosonic case. We give the details of this derivation in section 2. In section 3 we investigate what the properties of the resulting action are, and finally in section 4 the validity region of the action is discussed.

## 2 Deriving an Effective Action for the Tachyon

In this section, we will derive the effective action. In section 2.1 we will explain the conditions which are used to determine the Lagrangian. In section 2.2 we will then use these conditions to derive the action. The derivation is analogous to that of  for the superstring case.

### 2.1 Conditions

Let us assume that we have a Lagrangian of the form $`L=L(T,_\mu T^\mu T)`$.

We begin by discussing the EOM condition arising from requiring that the equations of motion be satisfied by the rolling tachyon solution, using the formalism of . We write

$$L=L_{even}+L_{odd}$$
(2.1)
where the even part is

$$L_{even}=\underset{n=0}{\overset{n=\mathrm{}}{}}L_{2n}$$
(2.2)
and the odd part is

$$L_{odd}=\underset{n=0}{\overset{n=\mathrm{}}{}}L_{2n+1}.$$
(2.3)
It will turn out that $`n=0`$ is sufficient to satisfy the BSFT condition (2.10).

We now introduce the parameter $`\gamma `$ to distinguish between the even ($`\gamma =0`$) and odd ($`\gamma =1`$) cases. Each term in the sums (2.2) and (2.3) goes as $`T^{2n+\gamma }`$:

$$L_{2n+\gamma }=\underset{l=0}{\overset{l=l_0}{}}a_l^{(n)}(_\mu T^\mu T)^lT^{2n+\gamma 2l}$$
(2.4)
It is known that string theory allows the “S-brane” solution

$$T(t)=T_+e^{t/\alpha }+T_{}e^{t/\alpha }$$
(2.5)
for the tachyon field $`T`$ in the spatially homogeneous case. We have here introduced another parameter, $`\beta \alpha ^2`$, which distinguishes between the bosonic ($`\beta =1`$) and the superstring ($`\beta =2`$) case.

The equations of motion which follow from (2.1), (2.2) and (2.3) are:

$$\begin{array}{cc}& \underset{l=1}{\overset{l=l_0}{}}a_l^{(n)}2l^\mu \left[(_\lambda T^\lambda T)^{l1}(_\mu T)T^{2n+\gamma 2l}\right]\hfill \\ \hfill =& \underset{l=0}{\overset{l=l_0}{}}a_l^{(n)}(2n+\gamma 2l)(_\mu T^\mu T)^lT^{2n+\gamma 2l1}\hfill \end{array}$$
(2.6)
As indicated, these should be satisfied by the solutions (2.5) for each $`n`$, as the solutions are exponentials. This renders the following recursion relation:

$$a_{l+1}^{(n)}=\frac{\beta }{2}\frac{(2l1)(2n+\gamma 2l)}{(2l+1)(l+1)}a_l^{(n)}$$
(2.7)
We see from this that for the even case, $`\gamma =0`$, we get non-zero coefficients only for $`l=0,\mathrm{},n`$, i.e. the upper summation index in (2.4) is $`l_0=n`$. For the odd terms, $`\gamma =1`$, we get non-zero coefficients for $`l=0,\mathrm{},\mathrm{}`$, i.e. $`l_0=\mathrm{}`$. Note that the lower summation index is always 0, as already indicated in (2.4).

The solution of the recursion relation is given by

$$a_l^{(n)}=\left(\frac{\beta }{2}\right)^l\frac{1}{(2l1)l!}\frac{(2n+\gamma )!!}{(2n+\gamma 2l)!!}a_0^{(n).}$$
(2.8)
In order to find the indicial coefficients $`a_0^{(n)}`$, it will be necessary to impose the boundary SFT (BSFT) condition, to which we now turn.

Using BSFT , different tachyon profiles were used in  to extract information about the tachyon Lagrangian, and rolling tachyons were studied within this framework e.g. in .

According to BSFT , the on-shell spacetime action $`S_{onshell}`$ is, both in the bosonic and the superstring  case, proportional to the disc partition function $`Z`$:

$$S_{onshell}=cZ,$$
(2.9)
where $`c`$ is a constant. The relation we will have use for is instead

$$L_{onshell}=cZ^{}(t),$$
(2.10)
where $`Z^{}(t)`$ denotes the world-sheet disc partition function with the zero mode $`t`$ unintegrated. It is consistent in the sense that (2.10) becomes (2.9) upon integrating out the $`t`$ dependence. We will take $`c=1`$, which is required for consistency with the results of  (see also ).

We will proceed by imposing this condition for the “half S-brane” case $`T_{}=0`$, i.e. $`T=T_+e^{t/\alpha }`$. $`Z^{}(t)`$ is then given by

$$Z^{}=\frac{1}{1+\frac{1}{\beta }(2\pi )^\beta T^\beta }=\underset{m=0}{\overset{m=\mathrm{}}{}}(1)^m\frac{(2\pi )^{\beta m}}{\beta ^m}T^{\beta m}.$$
(2.11)
The terms (2.4) in the on-shell Lagrangian, using (2.5), are

$$L_{2n+\gamma }^{onshell}=T^{2n+\gamma }\underset{l=0}{\overset{l=l_0}{}}a_l^{(n)}\frac{(1)^l}{\beta ^l}.$$
(2.12)
Imposing (2.10) with $`c=1`$, we see that relating equal powers of T requires

$$m=\frac{1}{\beta }(2n+\gamma ).$$
(2.13)
Apparently, if $`\beta =2`$ (superstring case) then $`\gamma =0`$, i.e. only even terms are required. This means that $`L_{super}`$ will be invariant under $`TT`$, as can be expected. If $`\beta =1`$ (bosonic case) no such restriction appears, i.e. we will need both even and odd terms.

Matching terms in (2.10) corresponding to equal powers of $`T`$ gives the relation

$$\frac{(1)^{m+1}}{\beta ^m}=\underset{l=0}{\overset{l=l_0}{}}a_l^{(n)}\frac{(1)^l}{\beta ^l},$$
(2.14)
where we have made the field redefinition $`T2\pi T`$ to get rid of the factors of $`2\pi `$.

Since $`a_l^{(n)}`$ is proportional to $`a_0^{(n)}`$ according to (2.8), this relation uniquely determines the indicial coefficient $`a_0^{(n)}`$, and that value can then be substituted back into (2.8) to give the unique solution for the expansion coefficents $`a_l^{(n)}`$.

### 2.2 Derivation

Let us begin by treating the bosonic even case, which according to section 2.1 corresponds to the parameter values $`\beta =1`$, $`\gamma =0`$, $`l_0=n`$. Performing the sum in (2.14) using (2.8) and (2.13) determines the indicial coefficent:

$$a_0^{(n)}=\frac{(2n1)!!}{n!2^n}$$
(2.15)
Substituting this back into (2.8) gives

$$a_l^{(n)}=\frac{(2n1)!!}{2^n(2l1)l!(nl)!},$$
(2.16)
which can be inserted into equations (2.1), (2.2) and (2.4) to give

$$L_{bosonic}^{even}=\underset{n=0}{\overset{n=\mathrm{}}{}}\frac{(1)^n(2n1)!!}{2^n}\underset{l=0}{\overset{l=n}{}}\frac{1}{(2l1)(nl)!l!}a^{nl}b^l,$$
(2.17)
where we have defined $`aT^2`$ and $`b_\mu T^\mu T`$. The lemma (A.1) tells us that this is equal to

$$L_{bosonic}^{even}=\frac{1}{1T^2}\sqrt{1T^2_\mu T^\mu T}.$$
(2.18)
We now continue with the bosonic odd case, which corresponds to $`\beta =1`$, $`\gamma =1`$, $`l_0=\mathrm{}`$. Performing the sum in (2.14) using (2.8) and (2.13) determines the indicial coefficent:

$$a_0^{(n)}=\frac{1}{\sqrt{\pi }}\frac{n!}{(n+\frac{1}{2})!}$$
(2.19)
Substituting this back into (2.8) gives

$$a_l^{(n)}=\frac{1}{\sqrt{\pi }}\frac{n!}{(n+\frac{1}{2}l)!}\frac{1}{(2l1)l!}.$$
(2.20)
Plugging this result into equations (2.1), (2.3) and (2.4) we get the bosonic odd Lagrangian,

$$L_{bosonic}^{odd}=\frac{1}{\sqrt{\pi }}\underset{n=0}{\overset{n=\mathrm{}}{}}n!\underset{l=0}{\overset{l=\mathrm{}}{}}\frac{1}{(n+\frac{1}{2}l)!(2l1)l!}a^{2n+1}b^l,$$
(2.21)
where we have defined $`aT`$ and $`b\frac{_\mu T^\mu T}{T^2}`$. We contend that this equals (assuming $`T=T(t)`$)

$$L_{bosonic}^{odd}=\frac{2}{\pi }\frac{1}{1T^2}\left[\sqrt{1T^2+\dot{T}^2}\mathrm{sin}^1\left(T\sqrt{1\frac{\dot{T}^2}{T^2}}\right)+\dot{T}\mathrm{sin}^1\left(\frac{\dot{T}}{T}\right)\right].$$
(2.22)
We believe that (2.21) and (2.22) are identical for the following reason: According to section 2.1, imposing the requirements that the equations of motion should be satisfied by (2.5) and that (2.10) should hold uniquely determines the Lagrangian. Consequently, since we know that (2.21) fulfills these conditions, and it can be checked that (2.22) also satisfies these conditions, they must be identical. As a check, we have confirmed that the Taylor expansion of (2.22) in $`a`$ and $`b`$ coincides with (2.21).

To get the complete bosonic Lagrangian, we need to add the even (2.18) and odd (2.22) parts:

$$\begin{array}{cc}\hfill L_{bosonic}& =L_{bosonic}^{even}+L_{bosonic}^{odd}\hfill \\ & =\frac{2}{\pi }\frac{1}{1T^2}\left[\sqrt{1T^2+\dot{T}^2}\mathrm{cos}^1\left(T\sqrt{1\frac{\dot{T}^2}{T^2}}\right)\dot{T}\mathrm{sin}^1\left(\frac{\dot{T}}{T}\right)\right]\hfill \\ & =+\frac{2}{\pi }\frac{1}{1T^2}\left[\sqrt{T^2\dot{T}^21}\mathrm{cosh}^1\left(T\sqrt{1\frac{\dot{T}^2}{T^2}}\right)+\dot{T}\mathrm{sin}^1\left(\frac{\dot{T}}{T}\right)\right].\hfill \end{array}$$
(2.23)
The first form is best suited for the field values $`0T^2\dot{T}^21`$, while the second is more appropriate when $`T^2\dot{T}^21`$ and positive $`T`$ (hence, $`T1`$). The forms are related by the relation $`\mathrm{cos}^1(1+\alpha ^2)=i\mathrm{cosh}^1(1+\alpha ^2)`$. The forms for the bosonic Lagrangian given here assume that $`T`$ only depends on $`t`$, but since $`\dot{T}`$ always appears squared in terms of the Taylor expansion, covariance is manifest.

Note that if we have a half S-brane, i.e. if either $`T_+`$ or $`T_{}`$ in (2.5) vanishes, the Lagrangian (2.23) reduces on-shell to

$$L_{bosonic}(\dot{T}=\pm T)=\frac{1}{1+T}$$
(2.24)
by virtue of (2.10). For large $`|T|`$, $`T`$ and $`\dot{T}`$ are exponentially close modulo a sign,

$$\frac{\dot{T}}{T}=\pm 1+O(e^{2t}),$$
(2.25)
using (2.5). Consequently, in the limit $`|T|\mathrm{}`$ the Lagrangian (2.23) behaves on-shell as

$$L_{bosonic}\frac{1}{T}.$$
(2.26)
This coincides with the large $`|T|`$ behaviour of (2.24), as expected since $`|\dot{T}|=|T|`$ holds exactly there.

To get (2.23) we required that (2.5) be a solution of the equations of motion. Let us mention that if we had instead imposed $`\stackrel{~}{T}=\sqrt{T_+}e^{t/2}+\sqrt{T_{}}e^{t/2}`$ as a solution for $`\stackrel{~}{T}`$ where $`\stackrel{~}{T}^2+T`$ (valid for $`|T|2\sqrt{|T_+T_{}|}`$), the resulting total bosonic Lagrangian in $`T`$, upon using the lemma (A.1), would have been

$$\stackrel{~}{L}_{bosonic}=\frac{1}{1+T}\sqrt{1+T+\frac{_\mu T^\mu T}{T}}.$$
(2.27)
This Lagrangian, which contains both even and odd parts, also has the properties (2.24) and (2.26). Its equations of motion in $`T`$ are satisfied by the half S-brane but not by the full S-brane solution (2.5), due to the large $`T`$ approximation. Making the field redefinition

$$T\mathrm{sinh}^2\left(\frac{\stackrel{~}{T}}{2}\right)$$
(2.28)
takes (2.27) to the even standard form (1.1) for $`T0`$.

Finally, let us also mention that the superstring case (non-BPS brane in type II string theory), corresponding to $`\beta =2`$, $`\gamma =0`$, $`l_0=n`$, can also be treated in an analogous way, and again using the lemma (A.1), the result is

$$L_{super}=\frac{1}{1+\frac{1}{2}T^2}\sqrt{1+\frac{1}{2}T^2+_\mu T^\mu T}.$$
(2.29)
The field redefinition

$$T\sqrt{2}\mathrm{sinh}\left(\frac{\stackrel{~}{T}}{\sqrt{2}}\right)$$
(2.30)
then transforms the Lagrangian (2.29) into the standard form (1.1).

## 3 Properties of the Bosonic Lagrangian

In the superstring case, the Lagrangian (1.1) has proven to satisfy several different consistency checks, so in this section we will investigate some of the properties of (2.23); the energy-momentum tensor in section 3.1, and the potential, closed string vacuum and lumps in section 3.2.

### 3.1 Energy-Momentum Tensor

We will now compare the behaviour of the energy-momentum tensor of the theory defined by (2.23) to that of the BCFT analysis .

We begin by discussing the full S-brane case $`\rho <\tau _p`$ for which the BCFT energy-momentum tensor is given by

$$\begin{array}{cc}\hfill \rho _{BCFT}& =\tau _p\mathrm{cos}^2(\pi \stackrel{~}{\lambda })\hfill \\ \hfill (p_{BCFT})_i(t)& =\tau _p\left[\frac{1}{1+\mathrm{sin}(\pi \stackrel{~}{\lambda })e^t}+\frac{1}{1+\mathrm{sin}(\pi \stackrel{~}{\lambda })e^t}1\right]\hfill \end{array},$$
(3.1)
where $`i=1,\mathrm{},p`$. The constant $`\stackrel{~}{\lambda }`$ parameterizes the solutions, and should satisfy555 Recall that we rescaled $`T`$ in equation (2.14). It is really $`\widehat{T}=\frac{1}{2\pi }T`$ which satisfies the BSFT condition, hence $`\widehat{T_0}=\frac{1}{2\pi }T_0`$. $`\stackrel{~}{\lambda }=\widehat{T}_0+O(\widehat{T}_0^2)`$, i.e.

$$\stackrel{~}{\lambda }=\frac{1}{2\pi }T_0+O(T_0^2),$$
(3.2)
where $`T_0`$ is the turning point of the tachyon, where it has potential energy only: (2.5) becomes

$$T=T_0\mathrm{cosh}(t),$$
(3.3)
which is related by a time translation to the solution (2.5) with $`T_0^2=+4T_+T_{}`$, because for this case, $`T^2\dot{T}^2=4T_+T_{}>0`$.

The energy-momentum tensor is defined and calculated as

$$T_{\mu \nu }\frac{2}{\sqrt{\eta }}\frac{\delta S}{\delta \eta ^{\mu \nu }},$$
(3.4)
with our Lagrangian given by (2.23). The on-shell energy (3.4) is then

$$\rho T_{tt}=V(T_0),$$
(3.5)
with the potential $`V(T_0)`$ given in (3.16). In (3.1), the energy is invariant under $`\stackrel{~}{\lambda }\stackrel{~}{\lambda }`$, but (3.5) is odd for small $`T_0`$; physical tachyon energies have the wrong sign for $`T_0<0`$. For $`T_00`$, (3.5) satisfies $`0\rho 1=\tau _p`$, (using (3.17)), as expected from (3.1), and at $`T_0=\stackrel{~}{\lambda }=0`$ we get $`\rho =\rho _{BCFT}=+1`$.

Turning to the pressure, we get from (3.4) that

$$p_i(t)T_{ii}=L_{bosonic},$$
(3.6)
where $`i=1,\mathrm{},p`$ and $`L_{bosonic}`$ is given in (2.23). The numerator $`N`$ of $`p_i(t)=L_{bosonic}`$ in (2.23) on-shell is

$$N=f(T_0)Tf\left(\frac{T_0}{T}\right).$$
(3.7)
The function $`f`$ is defined by

$$f(x)\sqrt{1x^2}\mathrm{cos}^1(x)\underset{n=0}{\overset{n=\mathrm{}}{}}a_nx^n,$$
(3.8)
where we also defined the Taylor expansion coefficients $`a_n`$ for $`f(x)`$. The numerator $`N`$ vanishes for $`T=1`$, but not for $`T=1`$, so the numerator may contain a factor of $`(T1)`$. More explicitly, inserting the indicated Taylor expansion of $`f(x)`$ gives

$$N=(T1)\underset{n=0}{\overset{n=\mathrm{}}{}}a_n\frac{T_0^n}{T^{n1}}G_{n1}(T),$$
(3.9)
where $`G_n(x)\frac{x^n1}{x1}=1+x+x^2+\mathrm{}+x^{n1}`$ denotes the geometric sum with $`n`$ terms.

Consider now the well-behaved side $`T0`$. The pressure $`p_i(t)`$ is regular at $`T=+1`$, since the numerator contains a factor of $`(T1)`$ (see (3.9)). As the tachyon rolls towards the closed string vacuum $`T\mathrm{}`$ we get, using (2.26) in (3.6), that

$$p_i(t)=\frac{1}{T}0.$$
(3.10)
Thus, the pressure vanishes as it should according to (3.1), using (3.17) and making the identification $`\pi \stackrel{~}{\lambda }=\frac{T_0}{2}+O(T_0^2)`$ for sufficiently small $`T_0`$, which matches the expectation (3.2).

Consider next the singular side $`T0`$. As expected from (3.1), as the tachyon rolls down towards the singularity at $`T=1`$, the pressure blows up in a finite time: The discussion after equation (3.8) tells us that the numerator (3.7) of $`L_{bosonic}`$ in (2.23) does not vanish as $`T1`$, so $`p_i(t)=L_{bosonic}`$ has a simple pole there.

Now, consider what happens around $`T=0`$. For small $`T_0`$, equation (3.7) becomes

$$N=(1T)\left(\frac{\pi }{2}T_0\right)+O(T_0^2).$$
(3.11)
Hence, the pressure is

$$\begin{array}{cc}\hfill p_i(t)& =\left[1\frac{2}{\pi }T_0\right]\frac{1}{1+T}+O(T_0^2)\hfill \\ & =\tau _p\left(1\frac{2}{\pi }T_0\right)\left[\frac{1}{1+\frac{T_0}{2}e^t}+\frac{1}{1+\frac{T_0}{2}e^t}1\right]+O(T_0^2).\hfill \end{array}$$
(3.12)
using $`\tau _p=+1`$ from (3.17). We see that the $`t`$ dependence remains the same as in (3.1) to order $`O(T_0^2)`$, provided that we replace $`\tau _p`$ by $`\tau _p\left(1\frac{2}{\pi }T_0\right)`$, and make the identification $`\pi \stackrel{~}{\lambda }=\frac{T_0}{2}+O(T_0^2)`$, again in agreement with (3.2).

We will make one final comment before leaving the case $`\rho <\tau _p`$: Placing the bosonic tachyon at rest at the closed string vacuum, $`\stackrel{~}{\lambda }=\frac{1}{2}`$, gives vanishing pressure and energy (3.1), and the same is true for (3.5) and (3.10), because the initial conditions $`T=\mathrm{}`$, $`\dot{T}=0`$ in (2.5) corresponds to $`T_+=T_{}=\mathrm{}`$, i.e. $`T_0=\mathrm{}`$.

We now turn to the case $`\rho >\tau _p`$. For $`|\mathrm{sinh}(\pi \stackrel{~}{\lambda })|<1`$, the BCFT energy-momentum tensor is given by

$$\begin{array}{cc}\hfill \rho _{BCFT}& =\tau _p\mathrm{cosh}^2(\pi \stackrel{~}{\lambda })\hfill \\ \hfill (p_{BCFT})_i(t)& =\tau _p\left[\frac{1}{1+\mathrm{sinh}(\pi \stackrel{~}{\lambda })e^t}+\frac{1}{1\mathrm{sinh}(\pi \stackrel{~}{\lambda })e^t}1\right]\hfill \end{array},$$
(3.13)
where $`i=1,\mathrm{},p`$ and $`\stackrel{~}{\lambda }`$ parameterizes the solutions. Now, $`T^2\dot{T}^2=4T_+T_{}<0`$ on-shell using the solution (2.5), which can be time translated to

$$T=T_0\mathrm{sinh}(t),$$
(3.14)
where $`T_0^2=4T_+T_{}`$; the tachyon starts at $`T=0`$ with a velocity $`_tT=T_0`$. As will be discussed in section 3.2, the kinetic term blows up at $`T=0`$, so the energy-momentum tensor cannot be well-defined; in general, it has imaginary parts (cf. section 4). Nevertheless, a similar analysis as the one for the case $`\rho <\tau _p`$ results in a conserved energy and a pressure which blows up at $`T=1`$ and decreases as it should (3.13) for large $`T`$.

Let us now complete the discussion of the energy-momentum tensor by analysing the half S-brane $`\rho =\tau _p`$, corresponding to starting the tachyon off at rest displaced slightly from the top of the potential hill. The BCFT energy-momentum tensor for this case is666 For concreteness, we study the case $`Te^t`$, but the case $`Te^t`$ is handled similarly.

$$\begin{array}{cc}\hfill \rho & =\tau _p\hfill \\ \hfill p_i(t)& =\tau _p\frac{1}{1+\widehat{\lambda }e^t}\hfill \end{array},$$
(3.15)
where $`\widehat{\lambda }`$ is not really a parameter in this case; it can be set to unity by time translation. The energy-momentum tensor of our theory is precisely equal to the BCFT one, using $`\tau _p=1`$ from (3.17) and making the identification $`\widehat{\lambda }=T_+`$: The energy can be checked by setting $`T_0=0`$ in (3.5). The pressure is also correct, as can be seen by combining (3.6) with (2.24). This is not surprising since (2.24) is a direct consequence of having imposed the condition (2.10). The BCFT properties that the pressure is singular at $`T=1`$ and falls off exponentially for $`T\mathrm{}`$ therefore carry over automatically to our theory in this case.

### 3.2 Other properties

Let us begin by examining the potential, which follows immediately from (2.23):

$$V(T)=\frac{2}{\pi }\frac{1}{\sqrt{1T^2}}\mathrm{cos}^1(T)=\frac{2}{\pi }\frac{1}{\sqrt{T^21}}\mathrm{cosh}^1(T),$$
(3.16)
where the first form is best for $`|T|<1`$, and the second is better for $`T>1`$. The potential is shown in figure 1.

Note that at $`T=1`$ the potential explodes to $`+\mathrm{}`$, and that $`V(T\mathrm{})0`$. The energy density of the D-brane at the open-string vacuum can be defined and computed as

$$\tau _pV(0)=1.$$
(3.17)
The odd contribution makes the potential monotonically decreasing. In particular, at $`T=0`$, it is proportional to $`T+O(T^2)`$ (modulo a constant vacuum energy due to the D-brane) instead of the more traditional $`T^2+O(T^2)`$ for a maximum. To still allow for the solutions (2.5), the kinetic term blows up777 In a sense, the fact that the kinetic term blows up at $`T=0`$ is the “glue” that makes a tachyon placed at rest at $`T=0`$ stay there. and changes sign at $`T=0`$, making the initial-value problem for rolling tachyons with $`\rho >\tau _p`$ not well-defined. There cannot therefore be a field redefinition which takes (2.23) to the canonical form even locally; if there were, (2.5) would not have been restricted to $`T_+T_{}0`$. Still, the solution (2.5) tells us that (2.23) describes a tachyon of mass $`m^2=1`$.

In view of the potential (3.16) and the properties of the energy-momentum tensor discussed in section 3.1, it appears that the following regions can be identified:

* $`T=1`$ is the singularity where the pressure blows up.
* $`T=0`$ is the open-string vacuum (since (2.5) has been imposed as a solution).
* $`T=+\mathrm{}`$ is the closed-string vacuum.

An important expected property of the closed-string vacuum is that it should not admit any open-string particle excitations. As discussed in , there are basically two ways in which this property can be realized. In some field theory models  the second derivative of the tachyon potential blows up at the minimum, hence rendering infinite mass to any would-be plane-wave excitations. The second possibility, which we will investigate here analogously to what was done in , is that the absence of plane-wave solutions is due to the behaviour of the kinetic term around the minimum (i.e. $`T=\mathrm{}`$), as in .

First, note that in the large $`T`$ limit, the Lagrangian (2.23) becomes

$$\frac{\pi }{2}L=\frac{1}{T}\left[\mathrm{ln}T\left(1+\frac{1}{2}\frac{_\mu T^\mu T1}{T^2}\right)\frac{_\mu T^\mu T}{T^2}\right]+\mathrm{}.$$
(3.18)
Making the field redefinitions

$$\begin{array}{cc}\hfill T& (1+\stackrel{~}{T})e^{\stackrel{~}{T}}\hfill \\ \hfill \varphi & e^{\stackrel{~}{T}/2}\hfill \end{array}$$
(3.19)
turns (3.18) into

$$L=\frac{8}{\pi }\left(\frac{1}{2}_\mu \varphi ^\mu \varphi \frac{1}{4}\varphi ^2\right)+\mathrm{}.$$
(3.20)
In terms of $`\varphi `$, which expands around $`T=\mathrm{}`$, the Lagrangian has the right canonical form, so this is the field definition to use for analysing excitations around the closed string vacuum. The Hamiltonian

$$H=\frac{2}{\pi }\frac{1}{\sqrt{T^2+_\mu T^\mu T1}}cosh^1\left[T\sqrt{1+\frac{_\mu T^\mu T}{T^2}}\right],$$
(3.21)
corresponding to the Lagrangian (2.23), should be constant for solutions of the equations of motion, which means that

$$T^2+_\mu T^\mu T=C$$
(3.22)
for some constant $`C`$. However, the plane-wave ansatz

$$\varphi =ae^{ikx}$$
(3.23)
makes the left-hand side of (3.22) become, using (3.19),

$$\frac{1}{\varphi ^4}\left[(12\mathrm{ln}\varphi )^216k^2(1\mathrm{ln}\varphi )^2\right],$$
(3.24)
which is not constant for any $`k^2=m^2`$. We are assuming a flat closed string background, so the concept of particle excitations in itself is well-defined. We therefore conclude that open-string particle excitations are absent at the closed-string vacuum $`T=\mathrm{}`$, as promised.

Finally, let us discuss the viability of lump solutions. In the superstring case, the Dirac-Born-Infeld action has been shown to allow kink solutions  (stable D-brane configurations). In the bosonic case, there should instead exist lump solutions (unstable D-branes). However, the Lagrangian (2.23) is designed to be valid for profiles in the vicinity of rolling tachyons (cf. section 4), and we can actually show that it does not contain lumps of the type locally resembling the kink solutions of .

By analogy with the superstring case , let us begin by analysing the energy-momentum tensor, which in the bosonic time-independent case $`T=T(x)`$ becomes

$$\begin{array}{cc}\hfill T_{\alpha \beta }& =\eta _{\alpha \beta }L_{bosonic}\hfill \\ \hfill \rho T_{tt}& =L_{bosonic}\hfill \\ \hfill pT_{xx}& =\frac{2}{\pi }\frac{1}{\sqrt{1T^2T_{}^{}{}_{}{}^{2}}}\mathrm{cos}^1\left(T\sqrt{1+\frac{T_{}^{}{}_{}{}^{2}}{T^2}}\right)\hfill \\ & =\frac{2}{\pi }\frac{1}{\sqrt{T^2+T_{}^{}{}_{}{}^{2}1}}\mathrm{cosh}^1\left(T\sqrt{1+\frac{T_{}^{}{}_{}{}^{2}}{T^2}}\right)\hfill \end{array},$$
(3.25)
where $`\alpha `$,$`\beta =0,\mathrm{},p1`$, $`x^px`$. The pressure $`p`$ is given in two forms, depending on whether $`T^2+T_{}^{}{}_{}{}^{2}1`$ or $`T^2+T_{}^{}{}_{}{}^{2}1`$, respectively. Since we are looking for time-independent lump solutions, $`T=T(x)`$, we should substitute $`\dot{T}\pm iT^{}`$ in the Lagrangian (2.23). The conservation law $`^\alpha T_{\alpha \beta }=0`$ now reduces to

$$\frac{p}{x}=0.$$
(3.26)
Using (3.25), this means that to have a solution we need to have either888 This solution may also correspond to D-branes in the world-sheet theory , but we do not analyse it further here. $`T=A\mathrm{cos}(x)+B\mathrm{sin}(x)`$ (which is the Lorentz rotated version of (2.5)), or

$$T^2+T_{}^{}{}_{}{}^{2}\mathrm{}.$$
(3.27)
As in , it is clear that for example a linear profile of infinite slope would locally satsify (3.27) for $`T>1`$. However, we can never find a profile that satisfies (3.27) globally, since a lump coming in from $`T=+\mathrm{}`$ must eventually turn at some finite $`T`$. Analysing lumps therefore seems to require higher derivative corrections to the Lagrangian (2.23).

## 4 Discussion

Our main result is the effective Lagrangian (2.23) for the bosonic open-string tachyon. It is the unique covariant Lagrangian not containing second or higher order derivatives whose equations of motion are satisfied by the solutions (2.5) and which on-shell is proportional to the disc partition function (with zero modes unintegrated) on half S-branes.

In the superstring case (non-BPS brane in type II string theory), this procedure leads to (2.29) which is related by a field redefinition to the standard form (1.1). According to , it describes tachyon profiles of the type

$$T=T_+(x^\mu )e^{t/\alpha }+T_{}(x^\mu )e^{t/\alpha }$$
(4.1)
(with $`\alpha =\sqrt{2}`$) where $`|_\mu T_+|1`$ and $`|_\mu T_{}|1`$, expanded around a solution with $`T_{}=0`$ (or $`T_+=0`$). In particular, it has one free parameter $`a_0^{(n)}`$ at each order $`T^{2n}`$ (see (2.8)), precisely as the $`2n`$-point function of tachyons at order momentum squared in a spatially homogeneous background of half S-branes (i.e. $`T_{}=0`$). The bosonic Lagrangian (2.23) should then analogously be designed to describe the same type of profiles (4.1) (with $`\alpha =1`$) for all values of $`T>1`$, provided that $`|\dot{T}||T|`$. We will explain these conditions shortly.

We also noted in section 2.2 that the bosonic Lagrangian (2.27) could be derived by a similar procedure valid for $`|T|2\sqrt{|T_+T_{}|}`$ (see the discussion around equation (2.27)), and it therefore describes profiles of the type (4.1) when $`|T|2\sqrt{|T_+(x^\mu )T_{}(x^\mu )|}`$. Since it was related to the standard form (1.1) by the field redefinition (2.28), it appears that (1.1) can only be expected to describe bosonic profiles well in the vicinity of rolling tachyons for large positive $`T`$ (cf. the discussion in ).

The most important test of (2.23) was the behaviour of the physical energy-momentum tensor (section 3.1). The pressure blew up in a finite time for tachyons rolling towards the singular side of the potential, and at the other side the pressure eventually decreased exponentially, as expected. Quantitative agreement was also good for half S-branes or full S-branes with either $`T_00`$ or $`T\mathrm{}`$, i.e. close to half S-branes, which of course lies in line with the expectation (4.1).

In general, the Lagrangian (2.23), and therefore also the physical energy-momentum tensor (cf. section 3.1), acquires imaginary parts outside of the region $`T1`$ and $`|\dot{T}||T|`$. The former restriction corresponds to going past the singularity at $`T=1`$ where the pressure blows up into the unphysical region $`T<1`$. The latter restriction is related to the fact that for rolling tachyons of the type $`T=T_0\mathrm{sinh}(t)`$, the kinetic term blows up at $`T=0`$, i.e. tachyons which traverse the top of the potential hill cannot be assigned a well-defined energy. The singularity of the kinetic term can be traced back to the fact that the EOM condition in section 2.1 made the Lagrangian non-analytic in $`T`$. Indeed, the Taylor expansion in $`b=\left(\frac{\dot{T}}{T}\right)^2`$, given in (2.21), diverges by the quotient test for $`|\dot{T}|>|T|`$, leading us to restrict the validity region to $`|\dot{T}||T|`$. Analytic continuation beyond this region, while possible, introduces imaginary components999 However, note that static profiles $`T=T(x)`$ are not restricted to $`|T^{}||T|`$..

Imposing the BSFT condition (2.10) for full S-branes would rid the Lagrangian of imaginary parts. Since the freedom of the ansatz $`L(T,_\mu T^\mu T)`$ was saturated already by half S-branes, this would require adjusting the ansatz such that more freedom is provided. The right type of freedom could possibly be provided by allowing for higher derivatives. This would also change the equations of motion, possibly in such a way so as to remove the non-analyticity in $`T`$.

However, perhaps there is another ansatz than $`L=L(T,_\mu T^\mu T)`$ which contains the same amount of freedom. Suppose, for example, that we make the ansatz $`L=L(T,_\mu T^\mu T,^2T)`$. Instead of (2.4), we then get

$$L_{2n+\gamma }=\underset{l=l_0}{\overset{l=l_1}{}}\underset{m=m_0}{\overset{m=m_1}{}}a_{lm}^{(n)}(_\mu T^\mu T)^l(^2T)^mT^{2n+\gamma 2lm}$$
(4.2)
at order $`T^{2n+\gamma }`$. Requiring that the equations of motion are consistent as described in section 2.1 leads to

$$\begin{array}{cc}\hfill 0& =a_{l+1,m}^{(n)}2(m1)(l+1)(2l+1)\hfill \\ & +a_{l,m}^{(n)}\left[(2n+\gamma 2l)2l(m1)+m(2n+\gamma 2l1)(2l+1)+(2n+\gamma 2lm)\right]\hfill \\ & +a_{l1,m}^{(n)}m(2n+\gamma 2l+1)(2n+\gamma 2l)\hfill \end{array}.$$
(4.3)
The case $`m=0`$ is the case already treated in section 2.1, i.e. $`a_{l,0}^{(n)}=a_l^{(n)}`$ in (2.7) (with only one free parameter, $`a_0^{(n)}`$). But we see from (4.3) that there is another possibility, $`m=1`$, which also has only one free parameter at each order. This example therefore shows that there may be higher derivative Lagrangians which also describe profiles of the type (4.1), possibly without the restriction $`|\dot{T}||T|`$.

Acknowledgments: I am very grateful to J. Minahan for several enlightening discussions and many useful comments on the manuscript. I would also like to thank L. Freyhult, J. Gregory, F. Kristiansson, U. Lindström and V. Schomerus for conversations.

## Appendix A Appendix

We want to show that

$$L\frac{1}{1+a}\sqrt{1+a+b}=\underset{n=0}{\overset{n=\mathrm{}}{}}\frac{(1)^{n+1}(2n1)!!}{2^n}\underset{l=0}{\overset{l=n}{}}\frac{1}{l!(nl)!(2l1)}a^{nl}b^l.$$
(A.1)
This formula was basically derived as a part of , but since it is used several times in section 2.2, we repeat the derivation here. Begin by Taylor expanding the left-hand side in $`a`$ and $`b`$ and use that $`(1)!=(2)!=\mathrm{}=\mathrm{}`$ to drop terms give

$$L=\underset{n=0}{\overset{n=\mathrm{}}{}}(1)^n\underset{m=0}{\overset{\mathrm{}}{}}\frac{(2m3)!!(1)^{m+1}}{2^m}\underset{k=0}{\overset{k=p}{}}\frac{1}{k!(mk)!}a^{n+mk}b^k.$$
(A.2)
Upon substituting $`pn+m`$,

$$L=\underset{p=0}{\overset{p=\mathrm{}}{}}(1)^{p+1}\underset{k=0}{\overset{k=p}{}}\frac{1}{k!}\underset{m=0}{\overset{m=p}{}}\frac{(2m3)!!}{2^m(mk)!}a^{pk}b^k.$$
(A.3)
The last term can be rewritten using $`mm+k`$ and the identity

$$\underset{m=0}{\overset{m=pk}{}}\frac{\left[2(m+k)3\right]!!}{2^mm!}=\frac{2^{kp}(2p1)!!}{(pk)!(2k1)},$$
(A.4)
which collapses one of the sums. This leaves us with

$$L=\underset{p=0}{\overset{p=\mathrm{}}{}}\frac{(1)^{p+1}(2p1)!!}{2^p}\underset{k=0}{\overset{k=p}{}}\frac{1}{k!(pk)!(2k1)}a^{pk}b^k,$$
(A.5)
which is the same as the right-hand side of (A.1), upon replacing $`pn`$ and $`kl`$.