On the migration-induced formation of the 9:7 mean motion resonance

TL;DR

This paper investigates how convergent planetary migration in a disc can lead to the formation of the 9:7 mean motion resonance, analyzing stability, capture regions, and the influence of planetary mass ratios on resonance dynamics.

Contribution

It provides a detailed analysis of the conditions for permanent or temporary 9:7 resonance capture, including the effects of migration parameters and planetary mass ratios, and identifies two distinct resonance modes.

Findings

Stable equilibrium surrounded by a capture region depends on migration parameters and mass ratio.

For certain parameters, equilibrium stability reverses based on mass ratio.

The 9:7 MMR has two modes separated by a separatrix, influencing post-capture evolution.

Abstract

We study formation of 9:7 mean motion resonance (MMR) as a result of convergent migration of two planets embedded in a disc. Depending on the migration parameters, initial orbits and planets' masses the system may pass through the resonance or enter it (permanently or temporarily). We illustrate that a stable equilibrium of the averaged system (a periodic orbit of the N-body model) is surrounded by the region of permanent resonance capture, whose size depends on the migration parameters and the planets mass ratio. A system located inside this region tends towards the equilibrium (the capture is permanent), while a system located outside the region evolves away from the equilibrium and leaves the resonance. We verify recent results of Delisle et al. and Xu & Lai where they show that for $m_1 \lesssim m_2$ ($m_1, m_2$ are the inner and outer planets' masses) the equilibrium is unstable when the migration is added, so the capture cannot be permanent. We show that for particular migration parameters the situation may be reversed (the equilibria are unstable for $m_1 \gtrsim m_2$). We illustrate that 9:7 MMR consists of two modes separated with a separatrix. The inner one is centred at the equilibrium and the outer one has no equilibrium in its centre. A system located outside the region of stable capture evolves from the inner into the outer mode. The evolution occurs along families of periodic orbits of the averaged system, that play a crucial role in the dynamics after the resonance capture.

Equations

$${\mathbf{f}}_i=-\frac{{\mathbf{v}}_i}{2{\tau }_{a,i}}-\frac{{\mathbf{v}}_i-{\mathbf{v}}_{c,i}}{{\tau }_{e,i}},$$

$${\tau }_{a,i}={\tau }_0{\left(\frac{r_i}{1\text{mbox}au}\right)}^{-\alpha },{\tau }_{e,i}={\tau }_{a,i}/\kappa ,$$

$${\varphi }_{5:4}^{(1)}$$

$${\varphi }_{5:4}^{(2)}$$

$${\varphi }_{4:3}^{(1)}$$

$${\varphi }_{4:3}^{(2)}$$

$${\varphi }_{9:7}^{(1)}$$

$${\varphi }_{9:7}^{(2)}$$

$${\varphi }_{9:7}^{(3)}$$

$$\overline{H}=-\frac{{\mu }_1{\beta }_1}{2a_1}-\frac{{\mu }_2{\beta }_2}{2a_2}-\frac{k^2{m}_1{m}_2}{a_2}\overline{R},$$

$$\overline{R}=\frac{1}{2\pi }{\int }_0^{2\pi }\frac{a_2}{\Vert {\mathbf{r}}_1-{\mathbf{r}}_2\Vert }dQ,Q\equiv \frac{{\lambda }_1-{\lambda }_2}{q}.$$

$${\dot{\sigma }}_i=\frac{\partial \overline{H}}{\partial I_i},{\dot{I}}_i=-\frac{\partial \overline{H}}{\partial {\sigma }_i},i=1,2.$$

$$x_i\equiv \sqrt{2I_i}\cos {\sigma }_i,y_i\equiv \sqrt{2I_i}\sin {\sigma }_i$$

$${\dot{y}}_i=\frac{\partial \overline{H}}{\partial x_i},{\dot{x}}_i=-\frac{\partial \overline{H}}{\partial y_i},i=1,2.$$

$${\dot{a}}_i^{\text{mbox}(migr)}$$

$${\dot{e}}_i^{\text{mbox}(migr)}$$

$${\dot{y}}_i$$

$${\dot{x}}_i$$

$$\dot{C}$$

$$\dot{K}$$

$${\dot{a}}^{\text{mbox}(migr)}=-\frac{a}{{\tau }_a}\left(1+\gamma \kappa e^2\right),{\dot{e}}^{\text{mbox}(migr)}=-\frac{e\kappa }{{\tau }_a}.$$

Full text

††pagerange: On the migration-induced formation of the 9:7 mean motion resonance–References††pubyear: 2017

1 Introduction

It is generally accepted that planets migrate in protoplanetary discs and that the migration plays an important role in formation of mean motion resonances (e.g. Snellgrove et al., 2001; Lee & Peale, 2002; Murray et al., 2002). It is also known that for slow and smooth migration a system of two planets evolves along families of periodic orbits (e.g., Beaugé et al., 2003, 2006; Hadjidemetriou, 2006; Hadjidemetriou & Voyatzis, 2010; Migaszewski, 2015; Voyatzis, 2016). The cited papers, though, are focused mainly on first order MMRs. Recently, Migaszewski et al. (2017) discussed the importance of the periodic orbits of the second order MMR, 9:7, in the possible formation and dynamics of a system of two super-Earths around Kepler-29.

Migaszewski (2015) showed that the system which undergoes the convergent migration not only follows the branch of periodic orbits, but also that the final system resulting from the smooth migration with efficient circularisation ends up as an exactly periodic configuration. As discussed in (Migaszewski et al., 2017), for the second order MMR, like 9:7, the final configuration is never exactly periodic, although the system may be located relatively close to the periodic system. In this work we drill the topic and try to understand why in certain situations the system passes through and in some other cases it enters the resonance. Moreover, we want to find conditions under which the system stays in the resonance either permanently or temporarily.

Our interest in the 9:7 MMR is not only motivated by the discovery of the Kepler-29 system and its resonant nature (Fabrycky et al., 2012; Jontof-Hutter et al., 2016; Migaszewski et al., 2017). It is also motivated by the fact that resonances of orders higher than one have not been studied extensively in terms of the migration. The papers devoted to the topic use various approaches and focus on different aspects of the resonance formation. For instance Xiang-Gruess & Papaloizou (2015) studied the migration of two-planet systems in wide ranges of the migration parameters and obtained several systems that ended up in higher order MMRs. The work was focused on finding constrains on the migration parameters (like the migration and circularisation rates) that lead to formation of various MMRs. Quillen (2006) and Mustill & Wyatt (2011) studied the resonance capture probability for second order MMRs as a function of masses and migration parameters. There has been also pointed out that the possibility of the capture can be determined by analysing the dynamics of the averaged system in vicinity of the stable equilibrium when the migration is added to the model (e.g., Delisle et al., 2015; Xu & Lai, 2016). It was found in the cited papers that for a few Earth mass planets and the circularisation rates comparable for both planets the capture cannot be permanent if the inner planet is less massive than the outer one.

Unlike in the cited papers, in this work we focus on particular resonance (9:7 MMR) and study its structure in details. Apart from the equilibrium of the averaged system (a periodic orbit of the unaveraged system), which was shown to play an important role in the migration-induced resonance formation, we analyse a role of periodic configurations of the averaged system in the process.

The paper is organised as follows. In section 2 we introduce a model of the system of two migrating planets. We present an example which illustrates probabilistic nature of the resonance capture. We chose three initial configurations differing only slightly one from another that evolve in three qualitatively different ways. In Section 3 we present the averaged equations of motion that will be used in most of the analysis. The study of the dynamics of the conservative model (without the migration) is presented in Section 4. We describe the structure of the phase space, i.e., equilibria, periodic orbits, chaotic regions corresponding to the separatrices. In Section 5 we study the evolution of the system of migrating planets. We start from the dynamics of the system which is already in the resonance. We explain the role of the equilibrium (in particular, we describe the region of stability around it) and the periodic orbits. Next, we focus on the process of entering into the resonance. In the last part of this section we study the importance of the planets’ mass ratio as well as a particular model of migration for the possibility of the permanent capture in the resonance. Conclusions are gathered in the last section.

2 A model of the system

In order to introduce the migration and circularisation of the planets into the N-body equations of motion we use a simple parametric model. For the $i$-th planet ($i=1,2$) the additional acceleration (expressed in the astrocentric reference frame) is a linear combination of the astrocentric velocity ${\boldsymbol{v}}_{i}$ and a velocity in a circular orbit at the distance of planet $i$, ${\boldsymbol{v}}_{c,i}$ (e.g., Papaloizou & Larwood, 2000; Beaugé et al., 2006; Moore et al., 2013; Voyatzis, 2016). It may be parametrised as follows:

[TABLE]

where the time-scales of migration and circularisation of planet $i$ are denoted by $\tau_{a,i}$ and $\tau_{e,i}$, respectively. We assumed (after Migaszewski et al., 2017) that

[TABLE]

where $\tau_{0}$, $\alpha$ and $\kappa$ are constant.

The model of the force given above is widely used in parametric studies of the migration that are focused on the mechanical part of the problem rather than its astrophysical aspect, although other models are used in the literature as well (e.g., Hadjidemetriou & Voyatzis, 2010). The parameters can be adopted from wide ranges of values. We will choose them to be realistic, so it would be possible, in principle, to find such a physical model of the disc that would lead to the migration governed by the parametric model of those values of the parameters.

An important limitation of the parametric model is that it assumes that if planet $i$ perturbs the mass distribution in the disc, the perturbation results in a force acting only on planet $i$ and does not affect planet $j$, $j\neq i$. For low mass planets, in a range of a few Earth masses, the assumption is relatively safe. For more massive planets it may be more risky (Podlewska-Gaca et al., 2012; Baruteau & Papaloizou, 2013). In this work we study the migration of two super-Earths, therefore the single-planet approximation may be safely applied.

2.1 Probabilistic capture

Entering into a second-order resonance is much more complex than for the first-order MMRs. Depending on initial orbits as well as the migration and circularisation rates a migrating system can pass the resonance or stay in it, permanently or temporarily (Quillen, 2006; Mustill & Wyatt, 2011). We chose an example which illustrates the complexity of the migration into 9:7 MMR, see Fig. 1. The parameters of the migration are given in the caption as are the planets’ masses and initial orbits. Each column is for one of three simulations. The top row shows the evolutions of the period ratio, $P_{2}/P_{1}$ (the outer planet’s period over the inner planet’s period), while the bottom row illustrates the behaviour of the inner planet’s eccentricity, $e_{1}$. The three simulations presented in Fig. 1 differ one from another only by the initial semi-major axes of the second planet, i.e., in the simulations shown from the left-hand to the right-hand panels, $a_{2}=0.11874,0.11875,0.11873\,\mbox{au}$. We can see that initial differences are $\sim 10^{-5}\,\mbox{au}$ for $a\sim 0.1\,\mbox{au}$. Such a small difference, though, leads to very different evolution and final configurations.

Shortly after the beginning of the simulation the system shown in the left-hand column of Fig. 1 enters into 9:7 MMR and stays there, as we show farther, permanently. The second system (middle column) enters into 9:7 MMR, stays there for a few thousands of years, leaves 9:7 MMR and enters into 5:4 MMR, in which the system stays permanently. Moreover, we observe that between $t\sim 7500\,\mbox{yr}$ and $t\sim 9500\,\mbox{yr}$, the period ratio remains around 9/7, while the eccentricity behaviour changes. The third simulation (the right-hand column of Fig. 1) shows that the system passes through 9:7 MMR, revealing only small excitation of the eccentricity and then goes into 5:4 MMR.

As already mentioned in the Introduction, the periodic orbits play a crucial role in the migration induced formation of MMRs. This fact is demonstrated in Fig. 2 which shows branches of periodic orbits in the $(P_{2}/P_{1},\log_{10}e_{1})-$plane (green dashed curves are for the branches of 5:4 and 4:3 MMRs, while the red dashed curve is for 9:7 MMR; the branches are best visible in the bottom panel). Subsequent panels from (a) to (c) show the results for the first, second and third simulations illustrated in Fig. 1. The periodic configurations are found with a help of the direct approach described in (Migaszewski et al., 2017). They were computed for initial values of the difference between the inner and the outer planets’ longitudes of pericentres $\Delta\varpi\equiv\varpi_{1}-\varpi_{2}=\pi$ and the mean anomalies for both planets $\mathcal{M}_{1}=\mathcal{M}_{2}=\pi$.

The choice of the initial values of the angles is such that the resonant angles of 9:7 MMR as well as the two first order MMRs 5:4 and 4:3 equal appropriate values. The resonant angles of 5:4 for $P_{2}/P_{1}\geq 5/4$ for which there exist periodic configurations are the following:

[TABLE]

where $\lambda_{1},\lambda_{2}$ are the mean longitudes of the inner and the outer planets, respectively. The values are reversed when the period ratio is smaller than the nominal value of 5:4 MMR. This is the case for all first order MMRs, therefore the periodic configurations of 4:3 for $P_{2}/P_{1}\leq 4/3$ occurs when

[TABLE]

It is easy to verify that for $\mathcal{M}_{1}=\mathcal{M}_{2}=\pi$ and the anti-aligned orbits the resonant angles of 5:4 and 4:3 MMR equal the values given above. The critical values of the resonant angles of 9:7 MMR read:

[TABLE]

which occurs for $\mathcal{M}_{1}=\mathcal{M}_{2}=\pi$, but also for $\mathcal{M}_{1}=\mathcal{M}_{2}=0$. In the presented illustration we chose the first case, as it covers all the resonances considered.

The evolution of the first example is shown in blue colour in the top panel of Fig. 2. The points are plotted only if both mean anomalies differ from the critical values of $\pi$ by less than half a degree. We do not check the value of $\Delta\varpi$ as shortly after the simulation starts it is already close to $\pi$. The direction of the evolution at the $(P_{2}/P_{1},\log_{10}e_{1})$-diagram is shown with black arrows. Initially the system evolves along the 4:3 MMR periodic branch. The amplitude of the oscillation around the periodic configuration decreases with time. When the system reaches the branch of 9:7 MMR, the amplitude is already very small. The system disappears from the diagram for some period of time, as there is no configurations with both $\mathcal{M}_{1},\mathcal{M}_{2}$ close to $\pi$. Next, the system appears above the 4:3/5:4 branch (the green curve), close to the 9:7 branch (the red curve) and it moves up. The zoomed plot of this fragment shows that the system tends towards the red curve, however very slowly. This system, as illustrated in the left-hand column of Fig. 1, will stay in 9:7 MMR permanently.

The second example is illustrated in Fig. 2b. Initially the system evolves in the same way as example 1 (which is expected, as the relative difference in $a_{2}$ is $\sim 10^{-4}$). Similarly to example 1, after reaching $P_{2}/P_{1}\approx 9/7$, the system goes above the green curve and evolves around the red curve. The zoomed fragment shows, though, that unlike example 1, the system tends away from the branch of periodic orbits and leaves the resonance after certain period of time. The magenta colour marks the location of the system for $t\in[7650,9500]\,\mbox{yr}$. That corresponds to the change of the behaviour of $e_{1}$ shown in Fig. 1e. For $t\gtrsim 9500\,\mbox{yr}$ the system tends towards 5:4 MMR.

The third example, which is illustrated in Fig. 2c, evolves initially similarly to the first two systems. The difference is that the third system passes through 9:7 MMR and goes directly to 5:4 MMR. The only visible result of the passage is that the amplitude of the oscillations of the system around the 4:3/5:4 MMR branch of periodic orbits is increased. The amplitude, though, is damped quickly after the 9:7 MMR encounter and next the system evolves along the 5:4 branch with almost zero amplitude. Its final position at the branch is dictated by the values of $\alpha$ and $\kappa$. For a given $\alpha$, the final eccentricity is higher when $\kappa$ is smaller (i.e., when the eccentricity damping is less efficient). Moreover, as we will show further in this work, $\kappa$ also governs the damping of the oscillations around the periodic orbits. High $\kappa$ means fast damping, so the final configuration in a first order MMR is exactly or almost exactly periodic. The damping of the oscillations around the second order MMR is more complex, which will be the subject to study in the next sections.

Figure 3 illustrates the difference between the evolution of the first and the second example systems. It shows the location of the systems at two diagrams, $(P_{2}/P_{1},e_{1})$ and $(e_{1},e_{2})$ for $t\in(1500,7500)\,\mbox{yr}$, which is the period of time counted from just after the entrance of the systems into 9:7 MMR up to the time corresponding to the transition between two regimes of motion in the evolution of example 2. The evolution of example 1 is shown in green colour, while the second example is plotted in blue. A red dashed curve shows the branch of 9:7 MMR periodic orbits. The observation we have already made is that the first system evolves towards the periodic configuration, while the second one moves away from it (see the arrows for the evolution direction). Locally both tracks can be approximated by polynomial curves (shown with black solid curves). For panel (a) we used a quadratic curve, while for panel (b) the curve is of the forth order. A black circle in each panel indicates the point towards which the system of example 1 tends. At the $(P_{2}/P_{1},e_{1})-$diagram that point corresponds to the intersection of the black and the red curves, while at the $(e_{1},e_{2})-$diagram the curves are in this point tangent one to another.

Apparently, both systems evolve along a single curve (in a vicinity of the periodic orbit). One can also notice that the initial location of the system at the black curve determines its motion along it, i.e., whether it occurs towards to or away from the black circle. It is thus straightforward that there exists such a point at the evolution track for which the system does not move either way. Such a point would indicate the critical distance from the periodic orbit below which the system stays in the resonance permanently. We will show further that the black curve corresponds to the branch of periodic orbits of the averaged system. In order to do so, we introduce in the next section the averaged model of motion.

3 The averaged system

The resonant dynamics of a two-planet system can be approximated by the averaged Hamiltonian $\overline{H}$ of two degrees of freedom with two first integrals of motion, i.e., the total angular momentum $C\equiv G_{1}+G_{2}$ and the spacing parameter $K\equiv(p+q)L_{1}+pL_{2}$ (Michtchenko & Ferraz-Mello, 2001a; Beaugé et al., 2003). The indices $p$ and $q$ indicates the (p+q):p resonance, and for 9:7 MMR $p=7,q=2$. $G_{i}$ and $L_{i}$ are the actions from the Delanuay angle-action canonical variables, i.e., $L_{i}\equiv\beta_{i}\sqrt{\mu_{i}a_{i}}$ and $G_{i}\equiv L_{i}\sqrt{1-e_{i}^{2}}$, where $\beta_{i}\equiv(1/m_{0}+1/m_{i})^{-1}$ are the reduced masses, $\mu_{i}\equiv k^{2}(m_{0}+m_{i})$ and $k$ is the Gauss’ gravitational constant. After the canonical transformation of the Delanuay variables, one can find the actions $I_{i}=L_{i}-G_{i}$ and the angles $\sigma_{i}=(1+s)\lambda_{2}-s\lambda_{1}-\varpi_{i}$, where $s\equiv p/q$ and $\lambda_{i}$ is the mean longitude of planet $i$. Therefore, the averaged Hamiltonian $\overline{H}=\overline{H}(I_{1},I_{2},\sigma_{1},\sigma_{2};C,K)$, see Michtchenko & Ferraz-Mello (2001a); Beaugé et al. (2003) for details.

The averaging of the Hamiltonian can be done numerically or analytically. For sake of self-consistency of the material, we repeat here a few informations about the averaging from (Migaszewski et al., 2017). The averaged Hamiltonian is given by a generic formula of the form of

[TABLE]

where the averaged disturbing function is given by the integral

[TABLE]

In this work we use the analytical approximation of the integral above, which is given by the sum of the secular and resonant terms, i.e., $\overline{R}=\overline{R}_{{\mbox{\scriptsize sec}}}+\overline{R}_{{\mbox{\scriptsize res}}}$, where $\overline{R}_{{\mbox{\scriptsize sec}}}$ and $\overline{R}_{{\mbox{\scriptsize res}}}$ are given in terms of power series in eccentricities and inclinations by a recipe from (Murray & Dermott, 1999). Explicit form of the expansion for 9:7 MMR up to the forth order in the eccentricities and with the inclinations assumed to be [math] is given in the appendix of (Migaszewski et al., 2017).

3.1 The equations of motion

The evolution of the Hamiltonian system (without the migration terms that will be given below) is governed by the following set of canonical equations of motion

[TABLE]

One can also introduce the non-singular equations (Michtchenko & Ferraz-Mello, 2001a), that are more suitable for small eccentricities when compared to the $(I_{i},\sigma_{i})$-set. The new variables are defined as follows:

[TABLE]

and the canonical equations of motion read

[TABLE]

The partial derivatives in the canonical equations are computed for fixed values of the integrals $C$ and $K$.

The migration can be incorporated into the equations of motion by adding extra terms to the $\dot{x}_{i}$ and $\dot{y}_{i}$ formulae as well as by adding two equations governing the time variation of $C$ and $K$. First we average out the evolution of the semi-major axes and eccentricities for each planet if the motion is governed by the N-body equations with addition of the force which mimics the migration and circularisation, Eq. 1.

If $\tau_{a}$ is given by Eq. 2 and $\tau_{e}=\tau_{a}/\kappa$, one can find approximate formulae for the averaged $\dot{a}_{i}$ and $\dot{e}_{i}$:

[TABLE]

The $\alpha$-dependence of the above equations stems from the fact that $\tau_{a}$ depends on $r$. If $\tau_{a}$ depended on $a$ instead of $r$, one would have to put $\alpha=0$ in Eq. 3.1. We will discuss the differences between the models with $\tau_{a}=\tau_{a}(r)$ and $\tau_{a}=\tau_{a}(a)$ in the last part of this work. In most of the study we use the $r$-dependent model of $\tau_{a}$.

The canonical variables $x_{1},y_{1},x_{2},y_{2}$ as well as the integrals of the Hamiltonian system $C$ and $K$ (that are not constant when formulae Eq. 3.1 are taken into account) depend functionally on $a_{1},e_{1},\sigma_{1},a_{2},e_{2},\sigma_{2}$. It is then possible to find the following equations of motion for the system of two migrating planets in a vicinity of 9:7 MMR:

[TABLE]

We note that the force which mimics the migration and circularisation, Eq. 1, does not lead to any change in $\sigma_{1}$ and $\sigma_{2}$, thus $\dot{\sigma}_{i}^{{\mbox{\scriptsize(migr)}}}=0$.

4 Dynamics of the conservative system

We start the study from an overview of the dynamics of the two-planet system without the migration. In this section we study the dynamics at energy plots. We also introduce the concept of periodic orbits of the averaged system as well as its connection with chaos.

4.1 The representative plane of initial conditions

In the further part of this work we will find it very useful a concept of the so called representative plane of initial conditions, $\Sigma-$plane, (e.g., Beaugé et al., 2003; Migaszewski et al., 2017), which can be used to construct the energy plots, that reveal the equilibria. The representative plane is a plane of the eccentricities, each point of which, $(e_{1},e_{2})$, corresponds to values of $a_{1},a_{2}$ that can be computed for given values of $C$ and $K$, or, what will be more useful when the migration is considered, for a given value of the scale-free angular momentum $c\equiv C/K$. Two remaining variables $\sigma_{1}$ and $\sigma_{2}$ are chosen such that both derivatives $\partial\overline{H}/\partial\sigma_{i}$ ($i=1,2$) vanishes for every point $(e_{1},e_{2})$ of the plane. It can be shown that it occurs for one of four combinations of the angles, i.e., $(\sigma_{1},\sigma_{2})\in\{(0,0),(0,\pi),(\pm\pi/2,\pm\pi/2),(\pm\pi/2,\mp\pi/2)\}$. Note that a simultaneous change of the signs of both angles in the two last combinations does not lead to any change in the Hamiltonian. As $\varpi= \sigma_{2}-\sigma_{1}$ and $\phi_{1}=-2\sigma_{1}$, the representative values of $(\Delta\varpi,\phi_{1})$ are $\{(0,0),(0,\pi),(\pi,0),(\pi,\pi)\}$.

A consequence of the definition given above is that a point in the plane for which $\partial\overline{H}/\partial I_{i}=0$ ($i=1,2$) is an equilibrium of the averaged Hamiltonian system (which corresponds to a periodic configuration within the N-body model of motion). Those values correspond to $\Delta\varpi=0$ or $\pi$, therefore the $\Sigma-$plane is called the symmetric representative plane. Equilibria of the averaged system, though, can also exist for $\sigma_{i}$ different from the representative values given above (e.g., Beaugé et al., 2003, for 2:1 MMR in moderate-to-high eccentricities regime). Asymmetric configurations with $\Delta\varpi$ different from [math] or $\pi$ will not be considered in this work. Such equilibria have not been found for any resonance in small eccentricities regime.

Figure 4a presents the energy levels (black solid curves) at the $\Sigma$-plane $(e_{1}\cos\Delta\varpi,e_{2}\cos\phi_{1})$. Elliptic and hyperbolic points in the plane indicate equilibria (marked with black cross/circle symbols). The energy plot has been obtained for $C$ and $K$ values corresponding to the stable equilibrium of 9:7 MMR taken from the branch of equilibria illustrated with black solid curve in Fig. 2c (almost vertical line in the centre of the plot). We chose an equilibrium with $e_{1}=0.008$ (for which $e_{2}\approx 0.0086$), as the example whose migration was discussed in Section 2 has the equilibrium eccentricities close to these values (see Fig. 3). The stable equilibrium is marked with black cross-circle symbol in the $(\pi,\pi)$-quarter of the $\Sigma$-plane. The equilibria located in the remaining three quarters are unstable.

Remaining panels of Fig. 4 presents the energy levels at different planes. They are used for better visualisation of the energy plot. Panel (b) shows the plot at the $\mathcal{S}$-plane, whose coordinates are $(e_{1}\sin\sigma_{1},e_{2}\sin\sigma_{2})$, therefore $\sigma_{i}$ can be $+\pi/2$ or $-\pi/2$. The $\mathcal{C}$-plane with $\sigma_{1},\sigma_{2}$ equal to [math] or $\pi$, $(e_{1}\cos\sigma_{1},e_{2}\cos\sigma_{2})$, is not shown here as it is irrelevant for the purpose of this work. The $\mathcal{S}$-plane possesses the same information as the two bottom quarters of the $\Sigma$-plane, while the $\mathcal{C}$-plane represents the two upper quarters. Panel (c) is a polar representation of the $\mathcal{S}$-plane, while panel (d) shows the diagram of $(P_{2}/P_{1},-e_{1}\cos\Delta\varpi)$. Thick solid curves in panel (d) show the borders of the permitted motion and is defined with $e_{2}=0$. The diagram with $e_{2}$ instead of $e_{1}$ is similar to panel (d) and was not shown.

4.2 Chaos and periodic orbits

Grey points in the panels of Fig. 4 indicate chaotic evolution of the system in terms of the spectral number (Michtchenko & Ferraz-Mello, 2001b). We integrate a given system for $2^{18}$ time-steps (which corresponds to $\sim 2500$ periods related to the resonant dynamics) and compute the Fast Fourier Transform (FFT). The spectral number is defined as a number of peaks in the power spectrum of amplitudes greater than a given limit, here $10^{-3}$ of the highest peak value. Darkest grey colour means that $SN\geq 100$, lighter grey is for systems with $SN\in[50,100)$, while the lightest denotes $SN\in[10,50)$.

Green and red curves in Fig. 4 denote branches of stable and unstable periodic orbits of the averaged system. A configuration is called periodic if there exists a period $T$ after which the system returns to its initial state, i.e., ${\boldsymbol{x}}(T)={\boldsymbol{x}}(0)$, where ${\boldsymbol{x}}=(a_{1},e_{1},\sigma_{1},a_{2},e_{2},\sigma_{2})$. Green and red points in the panels of Fig. 4 corresponds to the systems with $\delta\equiv\|{\boldsymbol{x}}(T)-{\boldsymbol{x}}(0)\|\lesssim 10^{-8}$. In order to check the stability of a given periodic configuration we change the initial eccentricities by a small value of $\pm 10^{-5}$ ($C$ and $K$ are kept unchanged) and compute $\delta$. The second step is to compute the spectral number, $SN$, for the periodic configuration. Then, for a stable periodic system both $\delta$ and $SN$ are small, while for an unstable configuration both the quantities have high values. We chose $10^{-4}$ as the limit value for $\delta$ and $100$ for $SN$.

4.3 Structure of the resonance

The structure of the resonance is best seen in Fig. 4c, i.e., at the polar version of the $\mathcal{S}$-plane. The centre of the resonance is given by the intersection of two branches of stable periodic orbits. The region of the resonance is encompassed by the separatrix (narrow regions of chaotic motions). The horizontal branch of the periodic orbits bends down at its left- and right-hand sides, passing into the branches of the unstable periodic configurations. The bifurcation takes place at the inner separatrix, which encompasses the inner part of the resonance (the one centred at the equilibrium). At both sides of the inner resonance there exists the outer part of the resonance, which is encompassed by the outer separatrix (the outer resonance is not centred at any equilibrium). The evolution of the systems in the two separated regions of the resonance is characterised by librations of both $\Delta\varpi$ and $\phi_{1}$ angles around $\pi$.

The resonance can be also investigated at the period ratio–eccentricity plane (Fig. 4d). The branch of periodic orbits which is horizontal in Fig. 4c corresponds to the track of the evolution of the system due to migration illustrated in Fig. 3 (compare Fig. 3a with Fig. 4d). The sign of the $y$-axis of Fig. 4d was changed in order to better show the correspondence. We will show in the next section that the branch of periodic orbits discussed here plays an important role in the migration induced formation of 9:7 MMR.

5 Migration

In the previous section we discussed the structure of the resonance, i.e., the equilibria, periodic orbits, division of the resonance into two parts. We will show now that after adding the migration terms to the equations of motion (Eqs. 3.1 and 3.1) the system migrates along the branches of periodic orbits (of the averaged system). We will also demonstrate that a system which enters the resonance far enough from the equilibrium evolve away from it and can leave the resonance in case of long enough migration.

5.1 Periodic orbits and equilibrium eccentricities

As we already mentioned, after entering the resonance a system evolves along certain path in the phase space (see Fig. 3a). This path is common for both types of evolution, i.e., towards and away from the resonance centre. We showed above that this path corresponds to the branch of periodic orbits (see Fig. 4d). The example from Section 2 illustrates that the system which evolves away from the equilibrium passes from one mode of the resonance to another (see Fig. 1e for the evolution of the eccentricity). We observe a change in the libration centres and amplitudes of the eccentricities. The system remains in the second mode of the resonance for limited time and leaves the resonance. Figures 4c demonstrates that the resonance indeed consists of two modes, which we called the inner (centred at the equilibrium) and the outer resonance (with no equilibrium). We will show in this section that the transition between the two modes of 9:7 MMR can be observed within the averaged model of the system as a passage from the inner to the outer MMR.

Similarly to the periodic orbits branch which corresponds to the path of the system which undergoes the migration, the stable equilibrium (the centre of the resonance) corresponds to the point towards which the system evolves during the migration in case of permanent capture (black circle in Fig. 3a). A system initially located in the equilibrium stays in this point during the migration only for certain values of the migration parameters (the radial dependence of the migration time-scale, $\alpha$, and the ratio between the migration and circularisation time-scales, $\kappa$). The migration time-scale is not important so long as the migration is slow compared to the time-scale for resonant libration. On the other hand, for those particular values of $\alpha$ and $\kappa$, a system which entered the resonance close enough to its centre (inside the stability zone), will move towards this point.

Figure 5a illustrates the equilibrium values of $e_{1}$ (the plot for $e_{2}$ is similar) as a function of $(\alpha,\kappa)$. Clearly, for higher values of $\kappa$ the eccentricities are lower. Additionally, for given values of the equilibrium eccentricities, when $\alpha$ is smaller, $\kappa$ also has to be smaller. The dependence of the equilibrium eccentricities on both $\alpha$ and $\kappa$ may seem unexpected, as the $\kappa$-dependence alone is usually discussed in the literature. In fact, the equilibrium eccentricities depend on the ratio between the time-scale of the period ratio variation (which corresponds to the relative migration time-scale) and the circularisation time-scale, i.e., $\kappa_{X}\equiv\tau_{X}/\tau_{e}$, where $\tau_{e}=(\tau_{e,1}+\tau_{e,2})/2$, $\tau_{X}\equiv-X/\dot{X}$ and $X\equiv P_{2}/P_{1}$. When $\alpha$ is lower, the relative migration (for the same $\tau_{1}$) is slower and $\kappa_{X}$ is higher, which leads to lower equilibrium eccentricities. When only one planet is migrating, as it is commonly studied in the literature, the distinction between $\kappa$ and $\kappa_{X}$ is not necessary, and the equilibrium eccentricities depend on $\kappa$ and not on $\alpha$.

Figure 5b presents the branch of equilibria (the red solid, almost vertical curve) together with branches of periodic orbits (only the stable part of the horizontal branch from Fig. 4c is shown) obtained for given values of the equilibrium eccentricities (or a given value of the scale-free angular momentum $c$).

As illustrated with the help of the example in Section 2, a system which enters the resonance not necessarily stays in it. Green area in Fig. 5b shows the stability zone (the region of stable capture). It means that when starting from a given point from the green area of the branch of periodic orbits (corresponding to given equilibrium eccentricities) and with the appropriate values of the migration parameters $\alpha$ and $\kappa$, the system evolves towards the equilibrium. As seen in Fig. 5a there is a continuum of combinations of $(\alpha,\kappa)$ leading to the same values of the equilibrium eccentricities. The size of the stability region, though, depends on particular values of the migration parameters. For smaller values of $\alpha$ the region is narrower. In the example shown in Fig. 5b the darker green corresponds to $\alpha=0.2$, while the lighter green denotes the stability region for $\alpha=1.8$. Naturally, if the system is located in the darker green area, it is located inside the stability zone for both values of $\alpha$.

The evolution of the period ratio and both eccentricities of three representative configurations is presented in Fig. 6. The initial orbital parameters were chosen from the branch of periodic orbits which corresponds to the equilibrium value of $e_{1}=0.008$ (see Fig. 5b). Black solid lines in each panel correspond to the evolution of the system in the equilibrium. Dashed curves are for the evolution of the system initially located in the stability region. Only epochs when $\phi_{1}=\pi$ are plotted, which means that a given parameter oscillates within the ranges given by the bottom and the upper dashed curves. As we can see, the amplitudes of oscillations decrease in time, so the system tends towards the equilibrium. The third example (shown with black points for the epochs of $\phi_{1}=\pi$ and grey area for the evolution between those epochs) corresponds to the system which is located initially outside the stability region. We observe an increase in the amplitudes of the oscillations. After $\sim 95\,$kyr the system moves to the second regime of the resonance (in which both $\phi_{1}$ and $\Delta\varpi$ keep librating around $\pi$ – not shown here) and at $t\sim 110\,$kyr the system leaves the resonance. Similar behaviour could be observed in the example simulated within the $N$-body model of motion (see Fig. 1e). Nevertheless, in the $N$-body simulation the system initially started outside the resonance, and in the example shown in Fig. 6 it starts from the branch of periodic orbits, i.e., when the system is already in the resonance. The mechanism of leaving the resonance is, though, the same in both cases.

5.2 The scale-free energy and angular momentum evolution

As we already mentioned, the structure of the equilibria, periodic orbits and the energy levels at the representative plane depend on the value of the scale-free angular momentum $c$. In order to follow the evolution of a given system one can also use another quantity, the scale-free energy $\epsilon\equiv E\,K^{2}$, where $E$ is the averaged energy given by the Hamiltonian (Eqs. 4 and 5) computed for a given point of the phase space. A system in equilibrium has a maximal value of $\epsilon$, which we denote as $\epsilon_{0}$. If a system starts away from the equilibrium but within the stability region, $\epsilon$ will grow asymptotically up to $\epsilon_{0}$.

Figure 7a illustrates the evolution of $\epsilon$ for initial systems chosen from the branch of periodic orbits for the equilibrium $e_{1}=0.008$. The systems of initial $\epsilon$ close enough to $\epsilon_{0}$ (above the green dashed line) tend towards the equilibrium. Below the green line $\epsilon$ decreases (i.e., the systems tend away from the resonance). At the level denoted by the blue dashed line the system moves from one mode of the resonance to the second one. After further decrease of $\epsilon$ the system leaves the resonance when $\epsilon$ is below the red dashed line.

Figure 7b presents the evolution of the scale-free angular momentum for the same set of initial configurations as discussed above. All the systems have the same initial value of $c=c_{0}$. Shortly after the beginning of a given simulation $c$ decreases. If it decreases below the green dashed line, $c$ will decrease further. If $c$ is above the green line after this initial decrease, the system tends towards the equilibrium ($c$ increases). Similarly to the left-hand panel, the blue and the red dashed lines denote, respectively, the transition between the two modes of the resonance and the value at which the system leaves the resonance.

5.3 The migration presented at the representative plane

In this section we illustrate the evolution of the systems with respect to the branches of periodic orbits. Figures 8 show the evolution of four systems of initial orbits chosen from the representative plane. The periodic orbits and chaotic regions are shown in the same manner as in Fig. 4c (with this difference that the chaotic orbits are marked only if $SN\geq 100$). Each system at the initial time intersects the plane in four points. Those points of intersection are denoted by circle symbols of different colours, red, green, magenta and blue, and labelled I, II, III and IV, respectively.

The points of intersection corresponding to the starting epoch are found during the integration of the equations of motion without the migration terms for the initial values of the orbital parameters. Next, we integrate the equations of the system with migration (see the parameters of the migration listed in the caption) over the time interval of $100\,$kyr. We list the osculating values of the orbital elements with the output time-step of $100\,$yr. Next, for each epoch we take the orbital elements as a starting configuration for the integration of the equations of motion without the migration terms in order to find the points of intersection of the representative plane. This way we obtain tracks of the migrating system at the plane.

Lets first describe the evolution of the example III (the magenta symbols). The system starts relatively close to the equilibrium. At first the system tends towards the horizontal branch of periodic orbits. After reaching the branch, the configuration moves towards the equilibrium. The evolution direction is shown with black arrows. Similar evolution one observes for the system II whose initial position at the plane is shown with the green points. The configuration IV (the blue symbols) also tends towards the equilibrium, however it starts very closely to the vertical branch of periodic orbits. We can then see very clearly that the evolution track is parallel to the vertical branch.

The system I marked with the red symbols evolve differently. After reaching the horizontal branch it moves away from the equilibrium. That is an example of the evolution away from the resonance centre. We checked many initial systems taken from the representative plane and marked the ones that evolve this way with small blue dots. Therefore, white zone centred around the equilibrium corresponds to the region of stable capture. We can also find other white areas at the plane that are out of the resonance (in a higher eccentricity regime). They are localized around the branch of stable periodic orbits characterized by the $y$-axis value close to $\pi/4$. That means that $e_{1}\approx e_{2}$ and $\Delta\varpi=0$. We conclude that it is possible to enter the resonance by moving along this branch. The fact that the white region is structured means that the entrance into the resonance depends on the phase in which the system reaches the separatrix. On the other hand, entering the resonance along the branch of periodic orbits with $e_{1}\approx e_{2}$ and $\Delta\varpi=\pi$ (the $y$-axis value close $3\pi/4$) is not possible, when starting from the eccentricities higher than the resonant values.

We observe that after reaching the horizontal branch, the systems evolve along it (towards or away from the equilibrium). Nevertheless, an illustration of the further evolution of the system I along the branches of periodic orbits is not easy. That stems from the fact that the scale-free angular momentum changes during the migration, and because of that the branches of periodic orbits change their position at the representative plane. The shift of the branches is shown in Fig. 9. As $c$ decreases (see Fig. 7b), the equilibrium shifts towards higher values of the eccentricities. The position of the branches of periodic orbits in three representative epochs are shown with decreasing grades of green/red colours. The branches obtained for the initial epoch is marked with the darkest green/red colour. The branches for the epoch for which the system moves from the inner resonance to the outer resonance are denoted with lighter colours, while the epoch when the system leaves the resonance is illustrated with the lightest green/red colours.

The track of system I is extended here to the whole time interval (in Fig. 8 only a fragment of the evolution was shown). The system follows the horizontal branch until it reaches the bifurcation point (marked with letter A), which is a merging point of the stable and unstable periodic configurations. This is the moment of the transition between the two modes of the resonance. After the bifurcation the system moves close to the branch of periodic orbits which is the centre of the outer resonance (points B). Because the transition occurs very quickly, the system is not a periodic configuration just after that. Some amount of time is needed so the system moves closer to the branch of periodic orbits. When the system reaches the outer separatrix (points C), it leaves the resonance.

5.4 The Poincaré cross sections

The evolution of the system which is leaving the resonance can be also followed at the Poincaré cross sections. Figure 10a presents the section for the initial epoch (red points in Fig. 9). The system is relatively far from the branch of periodic orbits. The phase trajectory for the nominal system is plotted with a green curve. Other curves are phase trajectories for the same energy. All of them fulfil the condition of $\phi_{1}=\pi$ and $\dot{\phi}_{1}>0$. The fixed point towards which blue arrows point is the periodic orbits for the energy of the system. The arrows illustrate that the system will move towards the fixed point when the migration is added to the model.

Figure 10b shows the Poincaré section for the epoch shortly before the system reaches the bifurcation point (point A in Fig. 9). The section is computed for the energy and the angular momentum equal to the values of the migrating system at this epoch. The nominal system is located in the fixed stable point of the section (a green point), which means that the system is a periodic configuration. There is an unstable fixed point close to the position of the system. During the evolution, those two points get closer one to another and shortly after the epoch for which we made the section, those two points merge. Red points correspond to chaotic motions. The arrows show the direction of the evolution of this chaotic region during the migration. The chaotic area is the border of the resonance. There is one more stable fixed point (below the two mentioned above) within the resonance.

The next panel (c) in this figure corresponds to the epoch just after the bifurcation. The resonance region shrinks slightly. The two fixed points merged, and now there is only one periodic orbit for $\Delta\varpi=\pi$. The nominal system whose phase trajectory is shown with a green curve is not a periodic configuration. As we mentioned already, that stems from the fact that the bifurcation was very fast. The transition between panels (b) and (c) is the transition between the two modes of the resonance mentioned earlier in this work.

During the further migration, the nominal system tends to the periodic configuration (as pointed by the arrows). The fixed point shifts towards higher values of the eccentricities. The last panel (d) in Fig. 10 shows the Poincaré cross section shortly before the system leaves the resonance. The nominal system is already in the fixed point of the section. The chaotic region gets narrower (see arrows in panel d). After the chaotic area shrinks down and reaches the fixed point, the system leaves the resonance. In other words, for the energy of the nominal system there is no resonant region in the phase space.

5.5 The entrance into the resonance

The study presented in the previous section was devoted to the evolution of the system when it is already in the resonance. We studied on what conditions it can stay in the resonance and when it evolves away from the equilibrium and, as a consequence, leaves the resonance. As we showed, the final fate of the system depends on how close it is located to the equilibrium just after entering the resonance. One can easily guess that the initial distance from the equilibrium depends on the time-scale of migration. For slower migration the initial distance should be smaller. On the other hand, as the example shown in Fig. 1 may suggest, the dependence is not necessarily that simple. In this section we study the entrance into the resonance.

As shown in Fig. 2, a system of the initial period ratio greater than $9/7$ approaches the branch of periodic orbits of 4:3 MMR before it reaches $P_{2}/P_{1}\approx 9/7$. It is thus natural to consider the branch of equilibria of 4:3 MMR as a starting point of the migrating system which is modelled within the averaged equations of motion. A configuration that is in that equilibrium has $\Delta\varpi=\pi$ and $\phi_{4:3}^{(1)}=\pi$, while the resonant angles of 9:7 MMR rotate.

Figure 11 illustrates the entrance of the system into 9:7 MMR when starting from the branch of equilibria of 4:3 MMR. The initial period ratio is being chosen from a relatively narrow range just before $9/7$, while $\phi_{1}$ is being chosen from the whole possible range of $(0,2\pi)$. Each of the initial configurations represented by a point at the $(P_{2}/P_{1},\phi_{1})-$diagram is integrated in time over $2\,$kyr and the type of the system’s behaviour is specified. The system can pass through the resonance (then such a configuration is marked with a dark grey point), it can enter the resonance temporarily (i.e., the system moves away from the centre of the resonance; such a configuration is marked with a light grey point) or the system can enter the resonance permanently (i.e., after entering the resonance the system moves towards the equilibrium of 9:7 MMR), in such a case no grey point is plotted at the diagram.

Configurations that are not destined for being permanently captured form stripes at the diagram in its right-hand part. The strips overlap with the vector filed shown with arrows. At a given point one can compute the variation of $P_{2}/P_{1}$ and $\phi_{1}$ in a given short unit of time. The arrows in the right-hand part of the diagram are almost vertical which means that the $\phi_{1}$ variation is much faster than the migration (expressed by the period ratio variation). Closer to the centre of the diagram $\phi_{1}$ evolves slower and thus a given system stays longer in a given part of the diagram. A palette of colours code $\dot{e}_{1}$ (black means the lowest values, while yellow denotes the highest values) and shows that the eccentricity is excited the most when $\phi_{1}=\pi/2$, while is damped the most when $\phi_{1}=-\pi/2$.

For the period ratio between $\sim 1.286$ and $\sim 1.2865$ and for $\phi_{1}\approx\pi/2$ the eccentricities are accelerated the most, as $\dot{e}_{1}$ is the largest and the system stays in this region of the plane for a long time (the arrows are short). This region may be thus considered as the entry into the resonance. If a given system starts its evolution from a point which leads to that region of the diagram, the capture into 9:7 MMR is most likely.

Trajectories of two representative configurations are shown with green and red tracks. The green case corresponds to the system that happens to reach the entry into the resonance. The system is captured permanently. The red case corresponds to the system that passes through the resonance. It only differs from the green case by the initial value of $\phi_{1}$, while the initial period ratio is the same.

Because the initial $\phi_{1}$ and $P_{2}/P_{1}$ are arbitrary when the migration into the resonance is considered, one can say that the entrance into 9:7 MMR should be considered probabilistically. The probability of the stable (permanent) capture can be measured as the ratio between the area without the grey points and the total area for a given range of $P_{2}/P_{1}$ at the right-hand side of the diagram and for the whole range of $\phi_{1}$. For the case considered in Fig. 11 that probability $\approx 84.6$ per cent. Similarly the probability of the unstable (temporary) capture can be measured as the ratio between the light grey area and the total area. The probability of any capture is simply a sum of those two probabilities. In the case considered the probability $\approx 92.8$ per cent.

Naturally the probabilities depend on the time scale of migration, $\tau_{0}$, and the $\alpha$ and $\kappa$ values (or equivalently on $\alpha$ and $e_{1}^{{\mbox{\scriptsize(eq)}}}$). Figure 12 shows the capture probability as a function of $\tau_{1}$ (the top label) and $\tau_{X}\equiv-X/\dot{X}$, where $X\equiv P_{2}/P_{1}$ (the bottom label). The green curve corresponds to the stable capture probability, while the red dashed curve denotes the probability of any capture (stable or unstable). Predictably, the capture probability increases for slower migration. In the next section we try to verify the results shown in Fig. 12 with the help of the N-body simulations.

5.6 N-body verification of the capture probability

In general, the capture probability depends not only on the migration parameters and planets’ masses but also on the initial orbits. Moreover, the probability depends on the distribution of the initial orbits, both its character and ranges for particular orbital elements. In the previous analysis, we assumed that the systems initially start from the equilibrium of 4:3 MMR and that the initial period ratio is close to $9/7$. As the equilibria in the averaged model correspond to the periodic orbits in the full N-body model, a natural choice for the N-body verification would be that the systems start from the branch of periodic orbits of 4:3 MMR. The periodic configurations require $\Delta\varpi=\pi$ and particular combination of the mean anomalies. We mentioned already that $\mathcal{M}_{1}=\pi,\mathcal{M}_{2}=\pi$ as well as $\mathcal{M}_{1}=\pi,\mathcal{M}_{2}=0$ both lead to appropriate values of the resonant angles of 4:3 MMR, when $P_{2}/P_{1}<4/3$. We can choose the initial phase freely as far as we keep the eccentricities at the values corresponding to the periodic orbits.

Unfortunately, when starting from the exactly periodic system, the capture does not occur for $\tau_{1}=10\,$kyr (the averaged model gives almost certain capture for this value of $\tau_{1}$). It is necessary to increase $\tau_{1}$ by a factor of $3$ in order to make the capture possible. Direct comparison of the N-body results with the results presented in Fig. 11 cannot be done, as a given value of $\phi_{1}$ can be achieved for different combinations of $(\mathcal{M}_{1},\mathcal{M}_{2})$. Moreover, a choice of $(e_{1},e_{2})$ other than the one for a periodic orbit would be arbitrary. Instead, we performed simulations in which the orbital elements are being chosen randomly from given ranges, i.e., $e_{1},e_{2}\in[0,0.002]$, $\mathcal{M}_{1},\mathcal{M}_{2}\in[0,2\pi]$, $P_{2}/P_{1}\in[1.288,1.29]$ – for the first set of simulations and $P_{2}/P_{1}\in[1.29,1.292]$ – for the second set. Initial $a_{1}=0.1\,\mbox{au}$, $\Delta\varpi=\pi$ for both sets. For each given $\tau_{0}$ we integrate $3000$ initial systems over time $t=\tau_{0}$ ($P_{2}/P_{1}$ approaches 9/7 in $t\sim 0.1\,\tau_{0}$ when starting from the initial values given above) and note final $P_{2}/P_{1}$. If the period ratio is below $1.275$ at the end of the simulation we count the system as non-resonant. Otherwise, the system is considered as captured into resonance. Therefore, we compute the probability of any capture (stable or unstable). The results are presented in Fig. 12. The black points connected with lines correspond to the first range of the period ratio, while the grey symbols indicate the capture probability obtained for the second range of $P_{2}/P_{1}$.

The probabilities obtained on the ground of the N-body model agree with the results of the averaged model (black and grey curves should be compared with the dashed red curve), although the differences for fast migration ($\tau_{1}\sim 4\,$kyr) are noticeable. The systems that start farther from 9:7 MMR enter the resonance less likely. The explanation of this fact is the following. As it was shown in Fig. 2, periodic orbits act as attractors in the phase space when the migration is included in the model. A system which starts farther from a given resonance approach the branch of periodic configurations closer when compared with a system whose initial $P_{2}/P_{1}$ is lower, i.e., closer to the nominal resonant value. As we mentioned already, the closer a given system is from the branch of periodic orbits, the less likely is the capture.

5.7 Dependence of the planets mass ratio

Recently Xu & Lai (2016) presented a study of the migration induced formation of second order MMRs. Although they were focused on 3:1, 5:3 and 7:5 MMRs, their results (especially the ones for 5:3 and 7:5) should be qualitatively transferable into 9:7 MMR. They used an analytic model of the system (similar to the averaged model used in this work, but of a lower order in the eccentricities) and found that a system initially in the equilibrium (from the branch of equilibria which exists for $\Delta\varpi=\pi$ and $\phi_{1}=\pi$) of the second order resonance not necessarily stays in the resonance when the migration is added. They show that when the time-scales of circularisation for both planets are comparable $\tau_{{\mbox{\scriptsize e,1}}}\approx\tau_{{\mbox{\scriptsize e,2}}}$, the systems can remain resonant only for $m_{1}\gtrsim m_{2}$. Similar criterion was found by Delisle et al. (2015). They found that the equilibrium is stable when $\tau_{{\mbox{\scriptsize e,1}}}/\tau_{{\mbox{\scriptsize e,2}}}>(e_{1}^{{\mbox{\scriptsize(eq)}}}/e_{2}^{{\mbox{\scriptsize(eq)}}})^{2}$. Because $e_{1}^{{\mbox{\scriptsize(eq)}}}/e_{2}^{{\mbox{\scriptsize(eq)}}}\approx m_{2}/m_{1}$ for 9:7 MMR in the low eccentricity limit, both criteria are equivalent. Below we try to verify those stability limits.

Figure 13 illustrates the dependence of the stability of equilibria when the migration is added to the model on the mass ratio and on $\alpha$. For each point $(e_{1},e_{2})$ from a given branch of equilibria (black/green and red/white points) which was found for a given value of $m_{1}/m_{2}$ and $\alpha$ we find such $\kappa$ for which those $e_{1}$ and $e_{2}$ are the equilibrium eccentricities. Next, we integrate the equations of motion and verify if the system evolves towards or away from the equilibrium. If the system moves towards the centre of the resonance, we mark the point as black/green, otherwise the point is red/white.

Although our parametrisation is different from the ones used in (Delisle et al., 2015; Xu & Lai, 2016), we can easily find that $q_{{\mbox{\scriptsize e}}}\equiv\tau_{{\mbox{\scriptsize e,2}}}/\tau_{{\mbox{\scriptsize e,1}}}=(7/9)^{\alpha}$. For a given value of the ratio $q_{{\mbox{\scriptsize e}}}$ one can find the critical value of $q\equiv m_{1}/m_{2}$ below which the system cannot stay in the resonance permanently. For all the $\alpha$ values considered here this critical value is between $\sim 0.8$ and $\sim 0.9$. The dashed and solid red lines plotted in each panel of Fig. 13 corresponds to the branch of equilibria for this critical value of the mass ratio computed on the ground of the criteria from (Delisle et al., 2015; Xu & Lai, 2016).

The results of Delisle et al. (2015); Xu & Lai (2016) are confirmed in our study for $\alpha=0.2$ (Fig. 13a), although we can see that for $e\gtrsim 0.01$ the equilibria are stable against the migration even if for lower eccentricities they were unstable. For $\alpha=1$ (Fig. 13b) the border which stems from (Delisle et al., 2015; Xu & Lai, 2016) does not overlap with our results. For slightly higher $\alpha=1.2$ (Fig. 13c) the equilibria are stable within the whole range of the mass ratio. For $\alpha=1.8$ (Fig. 13d) the situation reverses when compared to panels (a) and (b). The equilibria are unstable against the migration for high mass ratio, while for low $m_{1}/m_{2}$ they are stable. This last case shows that the permanent capture of the planets with $m_{1}/m_{2}$ as small as $0.1$ is possible for high values of $\alpha$.

5.8 Dependence on the migration model

The differences discussed above stem from different models of the migration used in this work and in the cited papers. Here we use the force given by Eq. 1 together with the assumption of $\tau_{a}=\tau_{a}(r)$. That leads to Eq. 3.1 for the $a$ and $e$ migration induced evolution (in the mean sense). Delisle et al. (2015) and Xu & Lai (2016) use slightly different formulae:

[TABLE]

The force given by Eq. 1 and $\tau_{a}=\mbox{const}$ (or $\tau_{a}=\tau_{a}(a)$) leads to similar expressions, apart from the numeric factor $\gamma$. Instead of $2$, the averaging over the mean anomaly gives $5/8$. In other words, if the force which mimics the migration is a linear combination of the ${\boldsymbol{v}}$ and ${\boldsymbol{v}}_{c}$ vectors (e.g., Papaloizou & Larwood, 2000; Beaugé et al., 2006; Moore et al., 2013; Voyatzis, 2016), then $\gamma=5/8$. When $\tau_{a}$ is a function of $r$, the averaging gives $\gamma=\frac{5}{8}(1-\frac{4}{5}\alpha)$. As recently shown (Goldreich & Schlichting, 2014; Deck & Batygin, 2015) a numeric value of this factor may be important for the stability of the equilibrium. In the literature $\gamma=0$ is used as well (e.g., Lee & Peale, 2002; Kley et al., 2004; Libert & Tsiganis, 2009; Bodman & Quillen, 2014). Below we study the stability of the equilibrium for different mass ratios and different models of the migration.

Models Ia and Ib are characterized by $\gamma=0$ (Eq. 11), where $\tau_{a}=\mbox{const}$ for model Ia and $\tau_{a}=\tau_{a}(a)$ for model Ib. The factor $\gamma=5/8$ for models IIa and IIb (IIa means that $\tau_{a}=\mbox{const}$, while in model IIb $\tau_{a}=\tau_{a}(a)$). The last model studied (III) has $\dot{a}$ and $\dot{e}$ given by Eq. 3.1.

Figure 14 shows the differences in the stability zone around the equilibrium resulting from the five models. The blue curves in Fig. 14 denote the borders between the dotted-blue and white regions inside the resonance as presented in Fig. 8. For the mass ratio $m_{1}/m_{2}=0.5$ (panel a) the stability zone disappears within models IIa,b. Models Ia,b and model III give similar results, all of them indicate that the equilibrium is stable and the zone spans significant part of the resonant region. When $\alpha=1$, model III has $\gamma=1/8$ which is only slightly more than [math] (as for models Ia,b). For equal masses (panel b) all models give the stable equilibrium, although the sizes of the stability zone differs from one model to another. Still the zone spans significant part of the resonance. In the last case (panel c;for $m_{1}/m_{2}=2$), models IIa,b leads to the stability zone as wide as the whole resonance. The other three models with smaller values of $\gamma$ gives narrower stability zones, although again for all the models the equilibrium is stable when the migration is considered.

Figure 13 shows that for small values of $\alpha=0.2$ the equilibrium is stable if $m_{1}/m_{2}\gtrsim 1$, while for higher $\alpha=1.8$ the situation reverses. For an intermediate value of $\alpha=1.2$ the equilibrium is stable in the whole range of the mass ratio. Because $\gamma=\gamma(\alpha)$ in the model with $\tau_{a}=\tau_{a}(r)$, this transition can be understood as the dependence of the stability on $\gamma$. Figure 15 illustrates the dependence. For the mass ratio $\gtrsim 0.9$ the equilibrium is stable for $\gamma$ greater the the critical value, while for $m_{1}/m_{2}\lesssim 0.9$ the stability requires $\gamma$ to be smaller than the critical value. Although the critical value depends on the mass ratio, it is close to [math] almost in the whole range, apart from $m_{1}/m_{2}\approx 0.9$. Negative values of $\gamma$, though, require an inverse profile of the surface density of the disc, when classic Type I migration is considered (Tanaka et al., 2002), which can be achieved at the disc inner border. An interesting observation is that for $\gamma=0$ the equilibrium is stable in the whole range of the mass ratio.

6 Conclusions

We have studied the formation of 9:7 mean motion resonance between two migrating planets of a few Earth masses. We used a simple parametric model of the migration and circularisation incorporated into the N-body as well as the averaged equations of motion. We confirm and explain the probabilistic nature of the entrance into the resonance (e.g., Quillen, 2006; Mustill & Wyatt, 2011). We show that the capture can be permanent or temporary, which is in agreement with previous studies (e.g., Delisle et al., 2015; Xu & Lai, 2016).

We demonstrate that the equilibrium for $\Delta\varpi=\pi$ and $\phi_{1}=\pi$ of the averaged conservative model (which corresponds to the periodic configuration within the N-body model of motion) plays a crucial role in the capture. The equilibrium is surrounded by the stability zone, the existence and the size of which depend on the mass ratio as well as on a particular model of migration. When $m_{1}<m_{2}$ and the disc surface density decreases with the distance from the star, the stability zone does not exist and the capture cannot be permanent. When the disc profile is reversed, the stability zone exists for $m_{1}<m_{2}$ and does not exist when $m_{1}>m_{2}$. For particular disc profiles the equilibrium may be stable in the whole range of the mass ratio. Generally, the stability zone spans only some part of the resonance, a region centred at the equilibrium. If the system, which just entered the resonance, is located inside the stability zone, it will move towards the equilibrium and the capture is permanent. Otherwise the system leaves the resonance despite the fact that the stability zone exists.

We show that 9:7 MMR consists of two modes that are separated one from another with the separatrix. We called them the inner mode and the outer mode. The inner mode is centred at the equilibrium of the averaged model, while the outer mode is centred at a branch of periodic configurations. We illustrate that a system in resonance evolves along branches of periodic orbits of the averaged system. Shortly after entering the resonance the system achieves the branch which goes through the equilibrium (the one with $e_{1}/e_{2}\approx\mbox{const}$, i.e., the horizontal branch at the polar $\mathcal{S}$-plane) and evolves towards or away from the equilibrium. In the latter case, when reaching the separatrix the system moves from the inner mode of the resonance to the outer one. The system evolves further along another branch of periodic orbits, until it passes through the outer separatrix and leaves the resonance. We illustrate the process of leaving the resonance at the representative plane as well as at the Poincaré cross sections.

We study the probability of permanent and temporary capture into the resonance and show that the system needs to enter the resonance in particular phase (a particular value of the resonant angle $\phi_{1}$). Therefore, even small changes in the initial orbits may lead to different final configurations. The phase-sensitivity of the process of entering the resonance has its consequence in the probabilistic nature of the capture. We illustrate that the system which starts the migration with $P_{2}/P_{1}>9/7$ tends towards the branch of periodic orbits of 4:3 MMR. The further from 9:7 MMR the system starts, the closer to the branch of 4:3 MMR it is when $P_{2}/P_{1}$ approaches $9/7$. We show that the proximity to the branch of 4:3 MMR periodic orbits makes the capture into 9:7 MMR more difficult. Therefore, the capture probability depends strongly on the initial orbits.

The results presented in this paper illustrate the complexity of the process of the migration induced formation of 9:7 MMR (and possibly other second order resonances). Further studies based on simplified parametric as well as hydrodynamical models are needed in order to understand the formation of such compact resonant configurations.

acknowledgements

CM thanks Alexander Mustill for valuable comments that helped to improve the manuscript. This work was supported by Polish National Science Centre MAESTRO grant DEC-2012/06/A/ST9/00276.