Issue 
A&A
Volume 628, August 2019



Article Number  A119  
Number of page(s)  9  
Section  Planets and planetary systems  
DOI  https://doi.org/10.1051/00046361/201833976  
Published online  20 August 2019 
Planet formation and stability in polar circumbinary discs
^{1}
Instituto de Astrofísica, Pontificia Universidad Católica de Chile,
Santiago, Chile
email: ncuello@astro.puc.cl
^{2}
Núcleo Milenio de Formación Planetaria (NPF),
Santiago, Chile
^{3}
Universidad Nacional de Córdoba, Observatorio Astronómico – IATE,
Laprida 854,
X5000BGR
Córdoba, Argentina
Received:
27
July
2018
Accepted:
25
June
2019
Context. Dynamical studies suggest that most circumbinary discs (CBDs) should be coplanar (i.e. the rotation vectors of the binary and the disc should be aligned). However, some theoretical works show that under certain conditions a CBD can become polar, which means that its rotation vector is orthogonal with respect to the binary orbital plane. Interestingly, very recent observations show that polar CBDs exist in nature (e.g. HD 98800).
Aims. We test the predictions of CBD alignment around eccentric binaries based on linear theory. In particular, we compare prograde and retrograde CBD configurations. Then, assuming planets form in these systems, we thoroughly characterise the orbital behaviour and stability of misaligned (Ptype) particles. This is done for massless and massive particles.
Methods. The evolution of the CBD alignment for various configurations was modelled through threedimensional hydrodynamical simulations. For the orbital characterisation and the analysis stability, we relied on longterm Nbody integrations and structure and chaos indicators, such as Δe and MEGNO.
Results. We confirm previous analytical predictions on CBD alignment, but find an unexpected symmetry breaking between prograde and retrograde configurations. More specifically, we observe polar alignment for a retrograde misaligned CBD that was expected to become coplanar with respect to the binary disc plane. Therefore, the likelihood of becoming polar for a highly misaligned CBD is higher than previously thought. Regarding the stability of circumbinary Ptype planets (also know as Tatooines), polar orbits are stable over a wide range of binary parameters. In particular, for binary eccentricities below 0.4 the orbits are stable for any value of the binary mass ratio. In the absence of gas, planets with masses below 10^{−5} M_{⊙} have negligible effects on the binary orbit. Finally, we suggest that mildly eccentric equalmass binaries should be searched for polar Tatooines.
Key words: planets and satellites: dynamical evolution and stability / protoplanetary disks / binaries: general / hydrodynamics / methods: numerical
© ESO 2019
1 Introduction
Planets appear as common byproducts of the star formation process within molecular clouds (Chiang & Youdin 2010). At early evolutionary stages (class 0/I), the fraction of stellar binaries and higherorder multiple systems is remarkably high: between 30 and 70% (Connelley et al. 2008; Chen et al. 2013). Consequently, protoplanetary discs can be severely affected by binary formation and stellar flybys (e.g. Bate 2018; Pfalzner 2013). Given the stochasticity of these processes, highly asymmetrical discs naturally form on timescales of the order of a few hundred kyr. Therefore, planet formation is also expected to occur in a broad variety of misaligned discs (Bate 2018; Cuello et al. 2019a). Currently, about a hundred planets have been detected in multiplestar systems (Martin 2018). There is, however, a striking tension with the more than 3700 planets detected^{1} around single stars (Batalha et al. 2013).
Planets in binary systems can either be circumbinary (Ptype, also known as Tatooines) or circumstellar (Stype). The seeming lack of Ptype planets has been ascribed to observational biases and dynamical processes. On the one hand, these systems are intrinsically difficult to observe through radial velocities (Eggenberger & Udry 2007; Wright et al. 2012) and transit methods (Martin & Triaud 2014). For transits in noncoplanar systems, the arbitrary orientation of the planet chord coupled to the fact that most binaries do not eclipse renders Ptype planet detection extremely challenging. On the other hand, recent stellartidal evolution models of shortperiod binaries show that the binary orbital period increases with time (Fleming et al. 2018). This translates into a larger region of dynamical instability around the binary, which could explain the lower frequency of Ptype planets compared to Stype planets. Based on the available data, typical planets found around binary stellar systems have a radius of the order of 10 Earth radii and orbital periods of about 160 days (Martin 2018). For instance, around a binary system with total mass equal to 1 M_{⊙}, this corresponds to an orbit with a semimajor axis of approximately 0.35 au.
Previous theoretical studies suggested that the discbinary interaction, namely the gravitational torque, should align the disc with the binary orbital plane on timescales shorter than the disc lifespan (Foucart & Lai 2013). Therefore, until recently, misaligned circumbinary discs (CBDs) were considered exotic. However, the observations of highly noncoplanar systems such as 99 Herculis (Kennedy et al. 2012), IRS 43 (Brinch et al. 2016), GG Tau (Cazzoletti et al. 2017), and HD 142527 (Avenhaus et al. 2017) suggest otherwise. For example, in the multiple system of GG Tau, Aly et al. (2018) recently explored various initial conditions to explain the disc misalignment with respect to the binary. Even more recently, Kennedy et al. (2019) reported the very first observation of a circumbinary gasrich disc in a polar configuration. These discoveries are rather puzzling and deserve to be explained.
In the field of secular dynamics, the seminal work by Ziglin (1975) presented an integrable approximation for the orbits of planets around binary systems. Later, based on numerical calculations of circumbinary polar particles, Verrier & Evans (2009) reported the coupling between the inclination and the node. Then, Farago & Laskar (2010) studied the circumbinary elliptical restricted threebody problem, giving an averaged quadrupolar Hamiltonian. In this way the authors showed that equilibrium configurations can be reached at high inclinations in the restricted and massive threebody problem. Then, relying on this analysis, it is possible to define the separatrix around the equilibrium points.
A more thorough analysis of the stability of massless particles was done by Doolin & Blundell (2011) where they consider different binary eccentricities and mass ratios between the components. This approach allows us to unambiguously identify prograde, retrograde, and polar orbits. In addition, the authors study the validity of the quadrupolar Hamiltonian approximation when the distance to the binary centre of mass is not sufficiently small. Furthermore, based on the study of hierarchical triple systems, several other studies investigated the circumbinary polar orbits for the restricted problem considering the octupole Hamiltonian (Ford et al. 2000; Naoz et al. 2013; Li et al. 2014). Interestingly, for near polar configurations, the inner binary can significantly excite the orbital eccentricity (up to about 0.3 in some cases). This excitation is weakly dependent on the distance to the binary. More recently, Naoz et al. (2017) and Zanardi et al. (2018) studied the evolution of polar orbits considering additional effects such as general relativity terms. Last but not least, the effect of the binary on the CBD can have dramatic effects on the disc alignment as described below.
Given the recent works on polar alignment and CBD dynamics (Sect. 2), it is interesting to explore further the conditions under which polar alignment of misaligned CBDs occurs (Sect. 3). Then, assuming planet formation happens in these discs, we characterise the orbits of single Ptype planets around eccentric and inclined binaries (Sect. 4). This allows us to constrain the binary eccentricity and mass ratio of systems around which we expect to observe polar circumbinary planets in the future. We discuss our findings in Sect. 5 and finally conclude in Sect. 6.
2 Summary of previous works
In the following we list the series of works that have previously addressed the topic of misaligned CBDs. We pay particular attention to nodal libration, polar CBDs, and alignment conditions.
2.1 Nodal libration, binary parameters, and stability
Doolin & Blundell (2011) numerically investigated the nodal libration mechanism in the longitude of the ascending node and in the inclination for Ptype particles around eccentric binaries. This phenomenon occurs for highly inclined orbits and for a wide range of binary mass ratios. Polar orbits were found to be remarkably stable at distances of two to three times the binary semimajor axis, which is not the case for coplanar configurations. Moreover, between the islands of stability, the authors report vertical striations of instability and hypothesise that they might be due to orbital resonances (see their Fig. 14, and Sect. 4.1). More recently, Zanazzi & Lai (2018) derived a simple analytical criterion based on the binary eccentricity and initial disc orientation that assesses whether the inclination of the CBD is expected to librate or not. It is worth noting that this criterion is symmetric with respect to the disc orientation: it does not distinguish prograde from retrograde configurations. However, as we show in Sect. 3.3, this condition is not sufficient in case of symmetry breaking.
2.2 Polar alignment of CBDs
In Aly et al. (2015) and Martin & Lubow (2017) it is demonstrated through numerical simulations how inclined gaseous CBDs around eccentric binaries can tidally evolve towards polar configuration with respect to the binary orbital plane. In this case, the dissipative torque acting on the disc is key to damping the nodal libration. This mechanism is investigated further in Zanazzi & Lai (2018) and Lubow & Martin (2018) where an analytical framework is provided to describe the polar alignment of CBDs, taking into account the gaseous viscosity of the disc.
2.3 Flybys and generation of polar alignment conditions
Stellar flybys, also known as tidal encounters, provide a natural way to tilt and twist a CBD. More specifically, a disc initially aligned with the binary orbital plane can be tilted by several tens of degrees and twisted by roughly 90° after a flybyof a massive star (Cuello et al. 2019a, Fig. D6). In particular, the smaller the pericentre distance and the more massive the perturber, the higher the resulting disc misalignment. It is worth noting that this flyby must happen in the early stages of disc formation such that enough time should remain for the disc to reach the polar configuration. In this scenario, inclined andretrograde encounters are the most efficient for disc warping (XiangGruess 2016). Alternatively, misaligned accretion onto the binary, as in Bate (2018), provides another mechanism to reach appropriate conditions for polar alignment.
2.4 CBDparameters and alignment timescales
More recently, Martin & Lubow (2018, hereafter ML18) explored the alignment of CBDs in a broad variety of configurations through hydrodynamical simulations. Specifically, they studied the dependence of the disc viscosity, temperature, size, binary mass ratio, orbital eccentricity, and inclination on the disc evolution. The authors find that polar alignment occurs for a wide range of binary–disc parameters, especially for low mass discs. Their results broadly agree with the expectation of linear theory. Interestingly, for large binary separations and/or for low mass ratios, the alignment timescale might be longer than the expected CBD lifetime. If this is the case, the CBD might remain on a misaligned configuration, neither (anti)aligned nor polar. This suggests that highly inclined discs and planets may exist around eccentric binaries (Zanazzi & Lai 2018).
3 Circumbinary disc dynamics and polar alignment
We first summarise the conditions under which we expect polar alignment. Then we present our SPH simulations, and finally describe the unexpected symmetry breaking observed.
3.1 Binary parameters and polar orbits
We consider a binary system with total mass M and individual masses M_{1} and M_{2}. We call q = M_{2}∕M_{1} the binary mass ratio. To describe the motion of a Ptype particle around the binary we use the following Jacobi orbital elements: semimajor axis a, eccentricity e, inclination i (with respect to the binary orbital plane), mean longitude λ, longitude of pericentre ϖ (alternatively argument of pericentre ω), and longitude of the ascending node Ω. The subindex B is used when referring to the binary orbit. Without loss of generality, we set initial conditions {λ_{B} = 0°, ϖ_{B} = 0°, Ω_{B} = 0°}. The angles of a given particle are measured from the direction of the binary pericentre. Finally, we set our fiducial system M_{1} = M_{2} = 0.5 M_{⊙} (i.e. q = 1), a_{B} = 0.1 au, and e_{B} = 0.5. This allows us to directly compare our results with those in Martin & Lubow (2017). Regarding the dynamics of Ptype particles, the phase space depends on Ω with two equilibrium points at {Ω = ±90°, i = 90°}, where both angles librate (see Doolin & Blundell 2011, for a detailed analysis). Interestingly, for any given couple of (Ω, i), it is possible to assert whether i librates around ± 90° or not by calculating the separatrix F (Farago & Laskar 2010): (1)
Specifically, its orbit is called polar if a particle fulfils the following condition: (2)
In Fig. 1 we show the separatrix in the plane (Ω, i) for the fiducial parameters considered. The initial conditions inside the separatrix limits correspond to the region of polar orbits for this given configuration. We recall, as described in Sect. 2.2, that the higher the value of e_{B}, the larger the region.
To sum up, in order to be polar, a Ptype particle needs to have a longitude of the ascending node close to 90° and be orbiting around an eccentric binary on a highly inclined orbit (see Fig. 1). When these conditions are met and in the absence of gas, then the particle’s inclination and the ascending node librate around 90°.
3.2 Misaligned circumbinary disc simulations
In this section, we assume that the chaotic and violent stellar cluster environment has produced an initially tilted CBD (see Sect. 2.3). We computed its resulting alignment, considering different initial disc inclinations (i_{0}) and twists (Ω_{0}). For simplicity, we only considered the set of fiducial parameters for the binary.
We ran 3D hydrodynamical simulations of CBDs by means of the PHANTOM Smoothed Particle Hydrodynamics (SPH) code (Price et al. 2018). The SPH method (see e.g. Price 2012) is well suited for misaligned CBDs simulations because there is no preferred geometry and angular momentum is conserved to the accuracy of the timestepping scheme, irrespective of the plane of the binary orbit. Both stars are represented by sink particles with accretion radii equal to 0.25 a_{B}, which interact with the gas via gravity and accretion (Bate et al. 1995). We used 10^{6} gas particles to model the disc, which is a resolution high enough to avoid artificially high numerical viscosities^{2}. The inner and outer radii of the disc were set to R_{in} = 2 a_{B} and R_{out} = 5 a_{B}, respectively.The surface density was initially set to a power law profile equal to Σ ∝ R^{−3∕2}. The disc total mass is equal to 0.001 M_{⊙}. Given this low value, we neglected the disc selfgravity in our calculations. Furthermore, we assumed that the disc is locally isothermal (i.e. that the sound speed follows c_{s} ∝ r^{−3∕4} and H∕R = 0.1 at R = R_{in}). For further details on the PHANTOM disc setup (see Price et al. 2018). Finally, we adopted a mean Shakura–Sunyaev disc viscosity α_{SS} ≈ 0.01.
Our set of PHANTOM simulations is given in Table 1. We only change the disc inclination (i_{0}) and its argument of ascendingnode (Ω_{0}). For completeness, these disc initial conditions are plotted as numbered dots in Fig. 1. Based on the linear theory in Zanazzi & Lai (2018), the disc is expected to become polar in 2 and 4, align in 1 and 6, and antialign in 3 and 5, with respect to the binary plane. From the practical point of view, we decompose the disc in 500 annuli and compute the tilt and the twist for each one of them. There is little difference between the disc warping in the inner and in the outer regions, except for r3 (see Sect. 3.3). Therefore, unless stated otherwise, we compute the tilt (i) and the twist (Ω) at r = 3 a_{B} from the PHANTOM outputs.
In Fig. 2 we plot the evolution of the tilt over time for 1000 binary orbital periods. We obtain the expected alignment for all the runs but r3. On the one hand, in r1 and r6 the disc inclination remains prograde and is slowly damped, while Ω circulates (not shown). In r5 the disc antialigns with respect to the binary (i = 180°), while Ω circulates. On the other hand, in r2, r3, and r4 the disc becomes polar (i = 90°) and Ω librates (instead of circulating). Remarkably, r2 and r4 are in phase opposition and the oscillations are quickly damped over time (compared to r1 and r6). Also, for low inclinations (r5 and r6) the disc becomes coplanar more quickly in the retrograde configuration (r5) compared to the prograde (r6) because the inner cavity, opened by resonant torques, is smaller for retrograde configurations (Miranda & Lai 2015; Nixon & Lubow 2015). Therefore, the disc responds more quickly to the binary torques. The anomalous behaviour of r3 is detailed below.
Fig. 1 Separatrix in the (i, Ω) plane constructed from Eq. (1) for e_{B} = 0.5. Initial conditions inside the separatrix limits correspond to polar orbits 2 and 4, whereas those outside correspond to prograde (i < 90°: 1, 6) and retrograde (i > 90°: 3, 5). The initial conditions of the SPH simulations of Table 1 are identified with numbers. 

Open with DEXTER 
Setof CBD simulations for different initial tilts (i_{0}) and twists (Ω_{0}).
Fig. 2 Tilt evolution over time, i(t), in terms of binary orbital periods. Runs: r1 (solid blue), r2 (dashed blue), r3 (solid red), r4 (dashed red), r5 (dotted red), and r6 (dotted blue). For comparison with r3, we also plot the symmetric curve of r1 with respect to i = 90° (solid green). r3 shows an anomalous behaviour because the disc breaks in two discs after ~ 150 orbits, which then evolve towards polar alignment (see Sect. 3.3). 

Open with DEXTER 
3.3 Unexpected symmetry breaking in r3
In Fig. 3 we show the CBD in r3 for three evolutionary stages: 0, 165, and 1000 binary orbital periods. Given the initial parameters of r3, the disc should antialign with respect to the binary (see Fig. 1). However, after half a tilt oscillation (~ 150 binary periods), the disc suddenly fulfils the condition for polar alignment and its tilt starts to oscillate around 90°. Additionally, its twist starts to librate instead of circulating. Both oscillations are damped, as expected for the viscous evolution of a gaseous disc (ML18). Hence, this shows that the symmetry between prograde and retrograde orbits is broken in some specific cases (e.g. r3).
In Fig. 4, we show the evolution of the ratio between the outer (a = 5 a_{B}) and the inner (a = 2 a_{B}) regions of the CBD in r3. This quantity indicates whether the disc is behaving as a solid body (~ 1) or not. By construction, this ratio is initially equal to 1. However, after a hundred binary orbits, it starts to oscillate between 0.8 and 1.2 due to the disc breaking, which is caused by the binary torque (see middle row in Fig. 3). We observe that the inner disc regions undergo a fast polar alignment, which eventually causes the entire CBD to align at i = 90°. After 400 binary orbits the entire disc has reached polar alignment and behaves as a solid body again (i.e. the tilt ratio is close to 1). This phenomenon is representative of the symmetry breaking between prograde and retrograde configurations.
The polar alignment in r3 is anomalous in the sense that it cannot be predicted by the analytical theory. More specifically, the expression of Λ in Eq. (7) in Zanazzi & Lai (2018), which shows whether the CBDs is expected to become polar or not, is symmetric for prograde and retrograde orbits (see e.g. their Fig. 2). Therefore, fulfilling the polar alignment criterion at the very beginning is not a necessary condition for polar alignment. The cavity size and further nonlinear effects play a key role in determining the final CBD alignment. Finally, as in ML18, we find that during the simulations the variations in the binary orbital parameters are only of the order of a few percent.
Polar CBDs might then be more common than previously thought. Therefore, since planet formation is expected to occur within these misaligned discs, it is relevant to characterise the planetary orbits around eccentric and inclined binaries.
4 Ptype planet orbital characterisation
We focus on the stability and evolution of single circumbinary Ptype planets (a > a_{B}) for a wide variety of binaries. To do so, we construct a regular 2D mesh of initial conditions with a pair of orbital parameters (e.g. i and e), while the remaining four parameters (e.g. a, λ, ω, and Ω) are initially set at nominal values. Then, we solve the threebody equations of motion numerically by means of a doubleprecision BulirschStoer integrator with tolerance 10^{−12}.
4.1 Dynamics of an inclined massless planet
First, we focus on the case of massless particles. A thorough study of the survival times for a regular combination of binary eccentricity 0 < e_{B} < 0.9 and mass ratio 0.2 < q < 1 in the (a, i) plane can be found in Doolin & Blundell (2011). Here, we extend their analysis by focusing on the features observed in the dynamical maps in order to better characterise inclined orbits. In particular, we use further structure indicators.
In Fig. 5, we show the dynamical maps in the (a, i) plane for Ω =0° (leftmost panels) and Ω = 90° (rightmost panels). Each orbit is integrated for 80 000 binary periods. If the particle either collides with one of the stars or escapes from the system, then it is coloured in white. In the top row, the colour scale is proportional to the Mean Exponential Growth of Nearby Orbits (MEGNO) value ⟨Y ⟩. This indicator is particularly appropriate because it identifies the chaotic orbits at a low CPUcost (Cincotta & Simó 2000). Regular orbits have ⟨Y ⟩ ~ 2 (blue), while chaotic and potentially unstable orbits have ⟨Y ⟩ > 2 (red/brown). In the bottom row, the colour scales correspond to the amplitude variation of i.
On the one hand, when Ω = 0° (first column in Fig. 5), the prograde orbits are stable only beyond 3.5 a_{B}, while retrograde orbits are stable closer to the binary system (a ~ 2.2 a_{B}). This occurs because the overlap of nearby resonances for retrograde orbits is of higher order compared to the case of prograde orbits (see Morais & Giuppone 2012). The vertical red stripes are associated with the 1 : N mean motion resonances (MMR), where eccentricity is efficiently excited. The orbits of these region have ⟨Y ⟩ ≫ 2, which indicates their chaotic nature. For prograde cases the 1 : 7 MMR delimits the stable orbits from the unstable ones at roughly a ~ 3.5 a_{B}. Alternatively, for retrograde cases this happens for the 1 : 4 MMR at a ~ 2.5 a_{B}. Interestingly, for high initial inclinations the orbits are stable for both prograde and retrograde configurations, even fairly close to 90°. It is worth noting that this does not occur for Stype configurations (see white gaps in Fig. 5 in Giuppone & Correia 2017).
On the other hand, when Ω = 90° (rightmost columns in Fig. 6), the stability maps change significantly and different stable regions appear. Specifically, there are stable orbits with high inclinations (45° ≤ i ≤ 135°) as close as retrograde orbits. Polar orbits are remarkably stable at small semimajor axis (a ~ 2.5 a_{B}). Furthermore, we analyse the impact of the binary mass ratio q on the stability of the orbits. The two rightmost columns of Fig. 5 show the stability maps for q = 0.5 and q = 0.2. As the mass ratio decreases, the MMRs become more important and severely deplete the polar regions (i ~ 90°). These vertical stripes in the e_{B} –q parameter space were first reported by Doolin & Blundell (2011). Here, we demonstrate that these peculiar patterns areindeed located at MMRs. In addition for nonequalmass binaries the libration amplitude of i (and consequently of the nodes) increases for polar orbits (see Fig. 5, bottom panel). This effect correlates well with the chaotic nature of these orbits.
To further characterise the orbits, we construct additional dynamical maps in a regular (e_{B}, i) grid. In Fig. 6, we show these maps for a fixed semimajor axis equal to a = 4.5 a_{B} and Ω =90°. This corresponds to a regular region of motion according to Fig. 5. Again, this is done for different binary mass ratios: q = 1, 0.5, and 0.2. To visually identify the region of polar orbits, we overplot a continuous red line in Fig. 6 defined by Eq. (2). The regions with the highest eccentricity excitation (Δe≳ 0.5) or with e_{B} ≥ 0.945 are unstable (coloured in white). When q = 1 (leftmost panel) the dynamical maps exhibit low eccentricity amplitudes with MEGNO values indicating regular orbits. For lower values of q, horizontal stripes of unstable and chaotic orbits appear at i = i_{N} and i = i_{c}. The values of i_{N} ~ 73°and 134° correspond to the very narrow octupole resonance wherein ϖ + ϖ_{B} − 2Ω librates about 0°, while i_{c} ~ 46°and 107°, where the frequency at ω± 180° (see e.g. Gallardo et al. 2012; Vinson & Chiang 2018, for a detailed analysis of these configurations). As a rule of thumb, the lower the mass ratio q, the larger the size of the chaotic region at i = 90°.
We also observe a diamondshaped feature centred at i = 90° and delimited by the separatrix with a vertex at e_{B} ~ 0.45. We checked that the regions where Δe ≳ 0.4 have high values of MEGNO (not shown) and are hence chaotic. In the bottom panel of Fig. 6, we show longitudinal cuts of the dynamical maps at i = 90°. The curves show the variations in eccentricity (Δe) and inclination (Δi) for different values of e_{B}. In the leftmost panel (q = 1) there is no discontinuity observed for Δe or Δi. However, for q = 0.5, we observe an abrupt change in Δe and Δi at e_{B} ~ 0.42. For q = 0.2, this occurs when 0.4 ≲ e_{B} ≲ 0.5. These curves clearly show that Δe and Δi are strongly correlated.
Fig. 3 Circumbinary disc in r3: the binary plane lies in the xyplane and the disc is initially inclined (i_{0} = 120°, Ω_{0} = 0°). At 165 binary orbits, we observe that the tilt is different between the inner and the outer regions of the disc (cf. Fig. 4). This is interpreted as the disc being broken in two discs by the binary. The inner regions experience polar alignment, causing the whole disc to become polar after roughly 400 binary orbits (cf. Fig. 2). At 1000 binary orbits, the disc is almost perfectly polar. 

Open with DEXTER 
Fig. 4 Tilt ratio between a = 5 a_{B} and a = 2 a_{B} in r3. When the disc behaves as a solid body, this ratio is ~1. The oscillations of this ratio (between 0.8 and 1.2) appear when the disc breaks in two misaligned discs. This occurs between approximately 150 and 400 binary orbits (see middle row in Fig. 3). 

Open with DEXTER 
4.2 Orbital characterisation for a massless polar planet
In the following, we exclusively focus on polar orbits (i.e. those that are initially at Ω = 90°, i = 90°). Figure 7 shows the polar orbits in the (e_{B}, q) plane. The colour scale corresponds to Δe. High values of Δe ≳ 0.3 result into chaotic behaviour. We see that some of the orbits (white dots) are unstable in less than 10^{5} binary periods. This dynamical map allows us to conclude that binaries with 0.35 ≲ e_{B} ≲ 0.5 and 0.05 ≲ q ≲ 0.4 are un likely to harbour polar planets. In other words, this corresponds to the forbidden region for polar typeP circumbinary planets.
4.3 Orbital characterisation for a massive polar planet
Here we set the eccentricity of the binary to e_{B} = 0.5 and vary themass ratio q, and the mass of the planet M_{p} for polar orbits. By doing so we aim to explore whether the stability depends on the planetary mass. In Fig. 8, we show the dynamical maps obtained for ΔΩ and for Δ(Ω − ϖ_{B}). For q > 0.7 and M_{p}< 10^{−5} M_{⊙}, the node librates (ΔΩ ~ 0°). For the same mass ratios but higher planetary masses we see that the node librates with higher amplitudes (ΔΩ > 90°). Interestingly, for massive planets M_{p} > 10^{−5} M_{⊙} (i.e. giants) the orbit precesses regardless of the value of q. This occurs because the massive planet causes the binary to start precessing. In addition, we note that the value of Δ(Ω − ϖ_{B}) does not depend on the planetary mass. Hence, in this configuration the indicator ΔΩ allows us to discern whether the binary ellipse is fixed in time or not. Figure 8 (lower panel) shows that, as expected, both the orbits of a single giant planet and of an almost massless planet remain polar around the osculating binary ellipse (i.e. independently of the binary precession).
Finally, in Fig. 9 we show the evolution of four selected orbits identified with numbers in Fig. 8. Orbit (1) corresponds to an orbit of an almost massless (10^{−15} M_{⊙}, roughly the mass of a planetary embryo) around a binary with q = 0.73. Orbit (1) is to be compared with orbit (2), which corresponds to an embryo as well, but orbiting around a binary with q = 0.67. Orbit (2) exhibits secular excitation in e and Ω, and it is stable for at least 10^{8} binary orbits. Orbit (3) corresponds to a planet of mass 10^{−3} M_{⊙} (~ Jupiter mass) orbiting around a binary with q = 0.67. We observe that the node circulates and that the secular excitations in Δe and Δi are of the same order as in orbit (2). Finally, orbit (4) corresponds to the orbit of a small moon (10^{−11} M_{⊙}) around a binary with q = 0.2. In this case the excitation in e reaches values of 0.4, which may eventually trigger instability. In fact, according to the MEGNO indicator, orbits 1, 2, and 3 are regular (⟨Y ⟩~ 2), while orbit 4 is chaotic (⟨Y ⟩ > 2). Based on the varied behaviour of these characteristic orbits, we conclude that planets with high masses (M_{p} ≳ 1 × 10^{−5} M_{⊙}) cause the binary pericentre to precess, and that the most regular orbits are obtained around binaries with q > 0.7.
In summary, the most favourable conditions to obtain stable polar orbits for circumbinary planets are reached, at moderate to high inclinations, around equalmass components and eccentric binaries for which Ω ~±90°. In addition, planets with masses above 10^{−5} M_{⊙} strongly affect the binary orbit. Assuming that such a massive planet forms in the CBD, its gravitational perturbations are likely to render the surrounding disc material unstable. The study of this effect is beyond the scope of this work and will be more thoughtfullyexplored in the future.
5 Discussion
As mentioned in Sect. 2.3, the specific initial conditions for polar alignment can either be reached after a stellar flyby (Cuello et al. 2019a) or through misaligned accretion within the molecular cloud (Bate 2018). Therefore, it is reasonable to consider polar alignment of the CBD as a possible (and potentially likely) dynamical outcome (Aly et al. 2015; Zanazzi & Lai 2018; Lubow & Martin 2018). The presence of the gaseous disc around the binary is crucial in order to viscously damp the nodal libration mechanism of the CBD (see Sect 2.2.). This is indeed a necessary condition for polar alignment, which should be distinguished from polar oscillations.
The final disc configuration, however, is highly sensitive to the binary parameters. For example, increasing the binary eccentricity e_{B} (for a fixed Ω) widens the inclination range for which polar alignment is expected. Additionally, increasing the binary mass ratio q decreases the oscillation period (cf. Eq. (16) in Lubow & Martin 2018) and makes the disc behave less rigidly. This is precisely when symmetry breaking between prograde and retrograde configurations is expected to occur. Then, for high values of q and e_{B}, the CBD is more susceptible to breaking up with at least one of the two discs eventually becoming polar.
The disc parameters are equally important as shown by ML2018: a higher disc viscosity leads to a faster damping of the nodal libration mechanism. Also, cool discs could be more easily broken by the binary torque. When the disc precession timescale is shorter than the sound crossing time, then the disc breaks at roughly 3 a_{B}. Finally, the larger the disc, the longer the timescale for the oscillations, even though these oscillations are more rapidly damped.
In this context, the very recent detection of a gasrich CBD in polar configuration by Kennedy et al. (2019) is of particular interest. In other cases, the presence of a highly inclined disc – perhaps in a polar configuration – has been proposed to explain asymmetric illumination patterns (Marino et al. 2015; Min et al. 2017; Facchini et al. 2018). In some cases, these shadows are able to trigger spiral formation in the outer disc (Montesinos et al. 2016; Montesinos & Cuello 2018), which are able to efficiently trap dust (Cuello et al. 2019b). Therefore, a shadowcasting inner disc is expected to affect the process of planet formation in the outer disc. Interestingly, the presence of shadows and spirals in an outer CBD could provide an indirect way to infer the presence of an inner CBD in a polar configuration.
Then, provided planet formation is efficient in these systems, polar and misaligned circumbinary Ptype planets should be common around mildly eccentric (e ≤ 0.4) equalmass binaries (q > 0.7). We recall that planet formation on polar orbits is more favourable when δe and δi are small since these conditions guarantee a small velocity dispersion. This kind of planet is stable at distances of a few times the binary semimajor axis (a ≳ 3a_{B}). For Ω ~ 90°, polar stability islands appear at even closer distances from the binary. These regions should be considered as sweet spots for polar Tatooine formation and evolution.
Interestingly, planets with masses above 10^{−5} M_{⊙} significantly affect the binary orbit and exhibit fast node circulation (alternatively binary pericentre precession). Observational campaigns such as BEBOP (Martin et al. 2019) provide, and will continue to provide, useful constraints on their potential existence and frequency.
Last but not least, we also expect gaps at 1 : N resonances with the binary once the gas dissipates in the remaining circumbinary debris disc. In the absence of any viscous force, these gaps naturally appear for high values of ⟨Y ⟩ in Fig. 5. This effect has profound implications for planetesimal growth and planet evolution in more evolved CBDs, which will be investigated in a future work.
Fig. 5 Dynamical maps of inclined systems with M = 1 M_{⊙} (from left to right): {q = 1, Ω = 0°}, {q = 1, Ω = 90°}, {q = 0.5, Ω = 90°}, and {q = 0.2, Ω = 90°}. The top labels correspond to the period in terms of binary periods. The regions between the horizontal red lines correspond to real polar orbits (i.e. librating around i = 90°). Longterm integrations show that particles with initial conditions in the red/brown regions do not survive more than 10^{5} periods. When Ω = 90°, stability islands appear at polar inclinations closer to the binary (between 2 and 3 a_{B}). In the bottom panel, the colour scale shows the amplitude of libration of i. The continuous red lines are defined by Eq. (2) and correspond to the analytic limits between polar and prograde or retrograde orbits. 

Open with DEXTER 
Fig. 6 Top panels: dynamical maps of inclined configurations with Ω = 90° with q = 1, 0.5, 0.2 (from left to right). The colour scale corresponds to the eccentricity variation Δe of the orbit. The continuous red lines are defined by Eq. (1) and correspond to the analytic limit between polar and prograde or retrograde orbits. Bottom panels: longitudinal cuts of the dynamical maps at i = 90°. The amplitude of Δe (Δi) is shown in red (blue). The grey shaded areas correspond to chaotic orbits with ⟨Y ⟩ ≫ 2. 

Open with DEXTER 
Fig. 7 Dynamical map for polar orbits for massless planets with initial a = 4.5 a_{B} in the (e_{B}, q) plane. The colour scale corresponds to Δe. The orbits identified with white dots are unstable. Chaotic orbits identified with the MEGNO indicator are indicated with grey squares. 

Open with DEXTER 
Fig. 8 Dynamical maps of polar orbits (i = 90°, Ω = 90°) for massive planets in the (q, M_{p}) plane with e_{B} = 0.5. The colour scale in the upper and lower panels corresponds to ΔΩ and Δ(Ω − ϖ_{B}), respectively. Initially, the planet is located at the equilibrium point i = 90°, Ω = 90°. The four initial conditions identified with white circles correspond to the orbits shown in Fig. 9. Each system was integrated for 20 000 yr (~ 650 000 binary periods). 

Open with DEXTER 
Fig. 9 Evolution of the individual polar orbits of Fig. 8: (1) blue, (2) green, (3) brown, and (4) yellow. The secular coupling between eccentricity, inclination, and node is easily identified. The periods are of roughly 5 × 10^{4} binary periods, depending on the couple of values (M_{P}, q). The evolution of giant planets (M_{P} > 10^{−4}, i.e. brown region in Fig. 9) exhibits the circulation of node. According to the MEGNO indicator, orbits 1, 2, and 3 are regular (⟨Y ⟩~ 2); instead, 4 is chaotic after ~100 000 orbital periods of the binary (⟨Y ⟩ > 2). 

Open with DEXTER 
6 Conclusions
We first studied the evolution of misaligned CBDs through hydrodynamical simulations (see Sect. 3). Then, by means of Nbody integrations,we thoroughly characterised the dynamics and stability of single Ptype planets around eccentric binaries (see Sect. 4). The main results of this study can be summarised as follows:
 1.
The viscous torque exerted by the binary on the misaligned CBD is stronger for retrograde configurations compared to progradeones. Hence, for the same relative inclination with respect to i = 90°, a retrograde CBD is more likely to become polar compared to a prograde CBD. This kind of symmetry breaking can happen in the inner disc regions under specific initial conditions (e.g. our run 3). Due to tilt and twist oscillations, the polar alignment condition of the CBD (initially not fulfilled) can be reached after a few hundred orbits. When this happens, the CBD stops behaving as a solid body (see Figs. 3 and 4) and the disc breaks. This phenomenon is not captured by the linear theory and is more likely to occur around eccentric equalmass binaries.
 2.
Our Nbody simulations show that although “exotic” polar circumbinary Ptype planets are stable on long timescales, namely over 10^{7} binary orbital periods. Based on stability criteria, we expect to find this kind of planet around binaries with mild eccentricities (e_{B} < 0.4) and for any value of the mass ratio q (see Fig. 7).
 3.
We find that above a certain planetary mass threshold (around 10^{−5} M_{⊙} in Fig. 8), the massive polar Tatooine strongly affects the binary orbit, which then starts precessing. The impact of the binary precession and the gravitational perturbations of such a giant planet on a hypothetical disc of planetesimals will be the subject of a future work.
In conclusion, planet dynamics around misaligned and polar CBDs is dramatically different compared to similar coplanar configurations. We note that a lack of detections of polar circumbinary planets, mainly because of observational biases, does not constitute proof of nonexistence. Future radialvelocity and transit surveys around binaries will help to constrain the fraction of polar, inclined, and retrograde Ptype planets. Equalmass binaries with eccentricities ~ 0.4 should be considered as the most promising targets to discover polar Tatooines. Thus, since HD 98800BaBe has a high eccentricity (e_{B} = 0.785 ± 0.005, Kennedy et al. 2019), this system does not constitute a good candidate for harbouring such a planet.
Acknowledgements
We thank the anonymous referee and the editor for very valuable comments and suggestions that have significantly improved our work. We also thank Cristian Beaugé, Agnese Da Rocha Rolo, and Antoine Rocher for useful discussions during this project. N.C. acknowledges financial support provided by FONDECYT grant 3 170 680. The authors acknowledge support from CONICYT project Basal AFB170 002. The Geryon cluster at the Centro de AstroIngenieria UC was extensively used for the SPH calculations performed in this paper. BASAL CATA PFB06, the Anillo ACT86, FONDEQUIP AIC57, and QUIMAL 130 008 provided funding for several improvements to the Geryon cluster. Nbody computations were performed at Mulatona Cluster from CCADUNC, which is part of SNCADMinCyT, Argentina. This project has received funding from the European Union’s Horizon 2020 research and innovation programme under the Marie SkłodowskaCurie grant agreement No 823823. Figure 3 was made with SPLASH (Price 2007).
References
 Aly, H., Dehnen, W., Nixon, C., & King, A. 2015, MNRAS, 449, 65 [NASA ADS] [CrossRef] [Google Scholar]
 Aly, H., Lodato, G., & Cazzoletti, P. 2018, MNRAS, 480, 4738 [NASA ADS] [Google Scholar]
 Avenhaus, H., Quanz, S. P., Schmid, H. M., et al. 2017, AJ, 154, 33 [NASA ADS] [CrossRef] [Google Scholar]
 Batalha, N. M., Rowe, J. F., Bryson, S. T., et al. 2013, ApJS, 204, 24 [NASA ADS] [CrossRef] [Google Scholar]
 Bate, M. R. 2018, MNRAS, 475, 5618 [NASA ADS] [CrossRef] [Google Scholar]
 Bate, M. R., Bonnell, I. A., & Price, N. M. 1995, MNRAS, 277, 362 [NASA ADS] [CrossRef] [Google Scholar]
 Brinch, C., Jørgensen, J. K., Hogerheijde, M. R., Nelson, R. P., & Gressel, O. 2016, ApJ, 830, L16 [NASA ADS] [CrossRef] [Google Scholar]
 Cazzoletti, P., Ricci, L., Birnstiel, T., & Lodato, G. 2017, A&A, 599, A102 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Chen, X., Arce, H. G., Zhang, Q., et al. 2013, ApJ, 768, 110 [NASA ADS] [CrossRef] [Google Scholar]
 Chiang, E., & Youdin, A. N. 2010, Annu. Rev. Earth Planet. Sci., 38, 493 [NASA ADS] [CrossRef] [Google Scholar]
 Cincotta, P. M., & Simó, C. 2000, A&AS, 147, 205 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Connelley, M. S., Reipurth, B., & Tokunaga, A. T. 2008, AJ, 135, 2496 [NASA ADS] [CrossRef] [Google Scholar]
 Cuello, N., Dipierro, G., Mentiplay, D., et al. 2019a, MNRAS, 483, 4114 [NASA ADS] [CrossRef] [Google Scholar]
 Cuello, N., Montesinos, M., Stammler, S. M., Louvet, F., & Cuadra, J. 2019b, A&A, 622, A43 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Doolin, S., & Blundell, K. M. 2011, MNRAS, 418, 2656 [NASA ADS] [CrossRef] [Google Scholar]
 Eggenberger, A., & Udry, S. 2007, ArXiv eprints [arXiv:0705.3173] [Google Scholar]
 Facchini, S., Juhász, A., & Lodato, G. 2018, MNRAS, 473, 4459 [NASA ADS] [CrossRef] [Google Scholar]
 Farago, F., & Laskar, J. 2010, MNRAS, 401, 1189 [NASA ADS] [CrossRef] [Google Scholar]
 Fleming, D., Barnes, R., Graham, D. E., Luger, R., & Quinn, T. R. 2018, in American Astronomical Society, DDA meeting 49, 202.06 [Google Scholar]
 Ford, E. B., Kozinsky, B., & Rasio, F. A. 2000, ApJ, 535, 385 [NASA ADS] [CrossRef] [Google Scholar]
 Foucart, F., & Lai, D. 2013, ApJ, 764, 1 [NASA ADS] [CrossRef] [Google Scholar]
 Gallardo, T., Hugo, G., & Pais, P. 2012, Icarus, 220, 392 [NASA ADS] [CrossRef] [Google Scholar]
 Giuppone, C. A., & Correia, A. C. M. 2017, A&A, 605, A124 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Kennedy, G. M., Wyatt, M. C., Sibthorpe, B., et al. 2012, MNRAS, 421, 2264 [NASA ADS] [CrossRef] [Google Scholar]
 Kennedy, G. M., Matrà, L., Facchini, S., et al. 2019, Nat. Astron., 189 [Google Scholar]
 Li, D., Zhou, J.L., & Zhang, H. 2014, MNRAS, 437, 3832 [NASA ADS] [CrossRef] [Google Scholar]
 Lubow, S. H., & Martin, R. G. 2018, MNRAS, 473, 3733 [NASA ADS] [CrossRef] [Google Scholar]
 Marino, S., Perez, S., & Casassus, S. 2015, ApJ, 798, L44 [NASA ADS] [CrossRef] [Google Scholar]
 Martin, D. V. 2018, Populations of Planets in Multiple Star Systems (Berlin: Springer), 156 [Google Scholar]
 Martin, R. G., & Lubow, S. H. 2017, ApJ, 835, 2 [NASA ADS] [CrossRef] [Google Scholar]
 Martin, R. G., & Lubow, S. H. 2018, MNRAS, 479, 1297 [NASA ADS] [CrossRef] [Google Scholar]
 Martin, D. V., & Triaud, A. H. M. J. 2014, A&A, 570, A91 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Martin, D. V., Triaud, A. H. M. J., Udry, S., et al. 2019, A&A, 624, A68 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Min, M., Stolker, T., Dominik, C., & Benisty, M. 2017, A&A, 604, L10 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Miranda, R., & Lai, D. 2015, MNRAS, 452, 2396 [NASA ADS] [CrossRef] [Google Scholar]
 Montesinos, M., & Cuello, N. 2018, MNRAS, 475, L35 [NASA ADS] [CrossRef] [Google Scholar]
 Montesinos, M., Perez, S., Casassus, S., et al. 2016, ApJ, 823, L8 [NASA ADS] [CrossRef] [Google Scholar]
 Morais, M. H. M., & Giuppone, C. A. 2012, MNRAS, 424, 52 [NASA ADS] [CrossRef] [Google Scholar]
 Naoz, S., Farr, W. M., Lithwick, Y., Rasio, F. A., & Teyssandier, J. 2013, MNRAS, 431, 2155 [NASA ADS] [CrossRef] [Google Scholar]
 Naoz, S., Li, G., Zanardi, M., de Elía, G. C., & Di Sisto, R. P. 2017, AJ, 154, 18 [NASA ADS] [CrossRef] [Google Scholar]
 Nixon, C., & Lubow, S. H. 2015, MNRAS, 448, 3472 [NASA ADS] [CrossRef] [Google Scholar]
 Pfalzner, S. 2013, A&A, 549, A82 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Price, D. J. 2007, PASA, 24, 159 [NASA ADS] [CrossRef] [Google Scholar]
 Price, D. J. 2012, J. Comput. Phys., 231, 759 [NASA ADS] [CrossRef] [Google Scholar]
 Price, D. J., Cuello, N., Pinte, C., et al. 2018, MNRAS, 624 [Google Scholar]
 Verrier, P. E.,& Evans, N. W. 2009, MNRAS, 394, 1721 [NASA ADS] [CrossRef] [Google Scholar]
 Vinson, B. R., & Chiang, E. 2018, MNRAS, 474, 4855 [NASA ADS] [CrossRef] [Google Scholar]
 Wright, J. T., Marcy, G. W., Howard, A. W., et al. 2012, ApJ, 753, 5 [NASA ADS] [CrossRef] [Google Scholar]
 XiangGruess, M. 2016, MNRAS, 455, 3086 [NASA ADS] [CrossRef] [Google Scholar]
 Zanardi, M., de Elía, G. C., Di Sisto, R. P., & Naoz, S. 2018, A&A, 615, A21 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
 Zanazzi, J. J., & Lai, D. 2018, MNRAS, 473, 603 [NASA ADS] [CrossRef] [Google Scholar]
 Ziglin, S. L. 1975, Sov. Astron. Lett., 1, 194 [NASA ADS] [Google Scholar]
All Tables
All Figures
Fig. 1 Separatrix in the (i, Ω) plane constructed from Eq. (1) for e_{B} = 0.5. Initial conditions inside the separatrix limits correspond to polar orbits 2 and 4, whereas those outside correspond to prograde (i < 90°: 1, 6) and retrograde (i > 90°: 3, 5). The initial conditions of the SPH simulations of Table 1 are identified with numbers. 

Open with DEXTER  
In the text 
Fig. 2 Tilt evolution over time, i(t), in terms of binary orbital periods. Runs: r1 (solid blue), r2 (dashed blue), r3 (solid red), r4 (dashed red), r5 (dotted red), and r6 (dotted blue). For comparison with r3, we also plot the symmetric curve of r1 with respect to i = 90° (solid green). r3 shows an anomalous behaviour because the disc breaks in two discs after ~ 150 orbits, which then evolve towards polar alignment (see Sect. 3.3). 

Open with DEXTER  
In the text 
Fig. 3 Circumbinary disc in r3: the binary plane lies in the xyplane and the disc is initially inclined (i_{0} = 120°, Ω_{0} = 0°). At 165 binary orbits, we observe that the tilt is different between the inner and the outer regions of the disc (cf. Fig. 4). This is interpreted as the disc being broken in two discs by the binary. The inner regions experience polar alignment, causing the whole disc to become polar after roughly 400 binary orbits (cf. Fig. 2). At 1000 binary orbits, the disc is almost perfectly polar. 

Open with DEXTER  
In the text 
Fig. 4 Tilt ratio between a = 5 a_{B} and a = 2 a_{B} in r3. When the disc behaves as a solid body, this ratio is ~1. The oscillations of this ratio (between 0.8 and 1.2) appear when the disc breaks in two misaligned discs. This occurs between approximately 150 and 400 binary orbits (see middle row in Fig. 3). 

Open with DEXTER  
In the text 
Fig. 5 Dynamical maps of inclined systems with M = 1 M_{⊙} (from left to right): {q = 1, Ω = 0°}, {q = 1, Ω = 90°}, {q = 0.5, Ω = 90°}, and {q = 0.2, Ω = 90°}. The top labels correspond to the period in terms of binary periods. The regions between the horizontal red lines correspond to real polar orbits (i.e. librating around i = 90°). Longterm integrations show that particles with initial conditions in the red/brown regions do not survive more than 10^{5} periods. When Ω = 90°, stability islands appear at polar inclinations closer to the binary (between 2 and 3 a_{B}). In the bottom panel, the colour scale shows the amplitude of libration of i. The continuous red lines are defined by Eq. (2) and correspond to the analytic limits between polar and prograde or retrograde orbits. 

Open with DEXTER  
In the text 
Fig. 6 Top panels: dynamical maps of inclined configurations with Ω = 90° with q = 1, 0.5, 0.2 (from left to right). The colour scale corresponds to the eccentricity variation Δe of the orbit. The continuous red lines are defined by Eq. (1) and correspond to the analytic limit between polar and prograde or retrograde orbits. Bottom panels: longitudinal cuts of the dynamical maps at i = 90°. The amplitude of Δe (Δi) is shown in red (blue). The grey shaded areas correspond to chaotic orbits with ⟨Y ⟩ ≫ 2. 

Open with DEXTER  
In the text 
Fig. 7 Dynamical map for polar orbits for massless planets with initial a = 4.5 a_{B} in the (e_{B}, q) plane. The colour scale corresponds to Δe. The orbits identified with white dots are unstable. Chaotic orbits identified with the MEGNO indicator are indicated with grey squares. 

Open with DEXTER  
In the text 
Fig. 8 Dynamical maps of polar orbits (i = 90°, Ω = 90°) for massive planets in the (q, M_{p}) plane with e_{B} = 0.5. The colour scale in the upper and lower panels corresponds to ΔΩ and Δ(Ω − ϖ_{B}), respectively. Initially, the planet is located at the equilibrium point i = 90°, Ω = 90°. The four initial conditions identified with white circles correspond to the orbits shown in Fig. 9. Each system was integrated for 20 000 yr (~ 650 000 binary periods). 

Open with DEXTER  
In the text 
Fig. 9 Evolution of the individual polar orbits of Fig. 8: (1) blue, (2) green, (3) brown, and (4) yellow. The secular coupling between eccentricity, inclination, and node is easily identified. The periods are of roughly 5 × 10^{4} binary periods, depending on the couple of values (M_{P}, q). The evolution of giant planets (M_{P} > 10^{−4}, i.e. brown region in Fig. 9) exhibits the circulation of node. According to the MEGNO indicator, orbits 1, 2, and 3 are regular (⟨Y ⟩~ 2); instead, 4 is chaotic after ~100 000 orbital periods of the binary (⟨Y ⟩ > 2). 

Open with DEXTER  
In the text 
Current usage metrics show cumulative count of Article Views (fulltext article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.
Data correspond to usage on the plateform after 2015. The current usage metrics is available 4896 hours after online publication and is updated daily on week days.
Initial download of the metrics may take a while.