{'fraction_non_alphanumeric': 0.0348132015412967, 'fraction_numerical': 0.02395711174401072, 'mean_word_length': 4.353542600896861, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 3, 'https://': 0, 'lorem ipsum': 0, 'www.': 2, 'xml': 0}, 'pii_count': 1, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Chaotic dynamics of the planet in HD 196885 AB

Depending on the planetary orbit around the host star(s), a planet could orbit either one or both stars in a binary system as S-type or P-type, respectively. We have analysed the dynamics of the S-type planetary system in HD 196885 AB with an emphasis on a planet with a higher orbital inclination relative to the binary plane. The mean exponential growth factor of nearby orbits (MEGNO) maps are used as an indicator to determine regions of periodicity and chaos for the various choices of the planet's semimajor axis, eccentricity and inclination with respect to the previously determined observational uncertainties. We have quantitatively mapped out the chaotic and quasi-periodic regions of the system's phase space which indicate a likely regime of the planet's inclination. In addition, we inspect the resonant angle to determine whether alternation between libration and circulation occurs as a consequence of Kozai oscillations, a probable mechanism that can drive the planetary orbit to a very large inclination. Also, we demonstrate the possible higher mass limit of the planet and improve upon the current dynamical model based on our analysis.


INTRODUCTION
The study of exoplanets has made a significant leap since the first planet around a solar type star, 51 Pegasi b (Mayor & Queloz 1995), was discovered almost two decades ago. The number of confirmed exoplanets has already surpassed most peoples' expectations with a count of 1055 as of December 30, 2013 1 . With new exoplanets being added everyday into the database, the centuries old quest of discovering extrasolar planets in the Universe is closer than ever. The ultimate conquest of the contemporary exoplanet research is to find a planet within its respective habitable zone with a stable planetary orbit. Various ground based surveys and space based telescopes, such as the Kepler 2 mission and future missions such as the transiting exoplanet survey satellite (TESS) (Ricker et al. (2010)) and the James Webb 3 space telescope are set to the detection and characterization of exoplanetary systems which will deliver a promising future e-mail: ssatyal@uta.edu 1 www.exoplanet.eu 2 http://kepler.nasa.gov/ 3 http://www.jwst.nasa.gov/ in the search of terrestrial planets which could support life. Meanwhile, the search for exomoons is also an ongoing process and recent work by Kipping et al. (2013) has produced a list of candidate hosts of such moons around transiting exoplanets. Noyola et al. (2013) recently calculated the minimum required flux from the exoplanetary radio emissions in order to detect the exomoons using several ground based radio telescopes. Although an exomoon has yet to be detected, theorists have already determined relations to constrain an exomoon's habitability by considering energy flux (i.e., radiative and tidal) and orbital stability (Heller 2012;Cuntz et al. 2013).
The focus of this paper, however, is the study of the orbital stability and the chaotic dynamics of exoplanets whose orbital inclinations are significantly higher relative to the binary plane. The stability of such systems depends on various orbital parameters and the orbital long-term stability is vital for life to develop under current theories. Thus, understanding the complete dynamics of a system bears a primary importance. Various numerical tools have been developed and used in the past to address the orbital configuration leading to the stability/instability of a system. Quarles et al. (2011) used the maximum Lyapunov exponent (MLE), originally developed by Lyapunov (1907), to determine the orbital stability or instability for the CRTBP case. The Mean Exponential Growth factor of Nearby Orbits (MEGNO) (Cincotta & Simó 1999) maps have been used to study the dynamical stability of irregular satellites (Hinse et al. 2010) and extrasolar planet dynamics . Recently, Satyal et al. (2013) used MLE, MEGNO and the Hill stability time series methods (Szenkovits & Makó 2008) to study the orbital perturbation of a planet due to the stellar companion in the binary systems of γ Cephei and HD 196885.
Planets that are formed from a planetary disk are expected to have near circular orbits due to tidal circularization, yet roughly 19% of known exoplanets with well-known orbital parameters have an eccentricity greater than 0.4 4 . This estimate, however, could be biased and may not include the true population. The actual intrinsic phenomenon on how such highly eccentric orbits are formed while maintaining orbital stability is not quite fully understood, but theories exist which use a combination of phenomena (Fabrycky & Tremaine 2007) to postulate such origins. Other mechanisms have also been proposed (Zakamska & Tremaine 2004;Kozai 1962;Lidov 1962) to explain such large variations in a planet's eccentricity and inclination induced by the libration of ω. In this work, we have made use of the Lidov-Kozai mechanism to explain such highly eccentric but stable orbits of HD196885 Ab. The planet is part of a binary star system (HD 196885 AB) whose eccentricity is observationally determined to be 0.48 and the planet's inclination is unconstrained which implies a value anywhere between zero and ninety degrees with the binary plane.
The dynamics of the planet in HD 196885 have been studied in a previous work by Satyal et al. (2013). In that study, restrictions were placed on the planet's inclination, (i pl ) 25 • with the binary plane, and the dynamics of the system was not fully explored. For this paper, we have considered a full range of prograde orbits in i pl value, from 0 • to 90 • , and have made use of the chaos indicator, MEGNO, to produce maps which demonstrate regions of periodicity and chaos for a variety of initial conditions. Also, for similar initial conditions a dynamical lifetime map is produced by using the information about the planet's ejection time from the system or the collision time with the stellar host. Then, the resonant angle is analyzed for evidence of a mean motion resonance at the best-fit location of the planet and for possible alternation between libration and circulation arising due to chaos induced by Kozai oscillations. This paper is outlined as follows. In Section 2 we present the basic theory of the Kozai mechanism, our numerical approach, a discussion of MEGNO as a chaos indicator and the initial set up of the system. In Section 3 we present our results from a representative sample of singly integrated orbits and MEGNO maps followed by discussion. We conclude in Section 4 with a brief overview of our results.  (Chauvin et al. 2011). Mass of primary star, m A = 1.33 M . The planet's inclination (i pl ) is measured relative to the binary orbit, then for i pl = 0 • , the planetary orbit is coplanar with the binary orbit.

The Lidov-Kozai Mechanism
Kozai (1962) developed an analytical theory to explain the secular perturbations induced by Jupiter on the asteroids in the Solar System. Similar theory was developed by Lidov (1962) to study the evolution of orbits of artificial satellites of planets that are directly influenced by the gravitational perturbations of the Sun. This idea was later adopted to study the secular interactions between a planet in highly inclined (relative to the binary plane) exoplanetary orbit and a distant binary companion (Innanen et al. 1997;Holman et al. 1997;Takeda & Rasio 2005;Batygin et al. 2011). As a consequence of the conservation of the quantity 1 − e 2 pl cos i pl (also valid for an eccentric perturber) the time evolution of eccentricity and inclination of the planetary orbit is anti-correlated: when the eccentricity is small, then the inclination is high and vice versa. For relative orbital inclinations less than 39.2 degrees the argument of pericenter circulates and the planet's orbit is precessing. This property remains true for other values of initial eccentricity. However, for relative inclinations larger than 39.2 degrees, the planet becomes a Kozai librator (Kozai regime) with it's argument of pericentre locked and exhibiting either librations (oscillations) around 90 or 270 degrees. These libration centers are known to be stable fixed points. A Kozai librator with an initially circular orbit will undergo a large variation in eccentricity and inclination within Kozai cycles. For relative inclinations close to or around the critical value of 39.2 degrees the argument of pericentre exhibits intermittent behaviour displaying alternations between circulations and librations. In general, the coupling between eccentricity and inclination provides an effective removal mechanism. For a large initial relative inclination the amplitude of the planet's eccentricity variation increases. Following Innanen et al. (1997) the maximum extent in eccentricity for a Kozai cycle can be expressed as, where (io) pl is the initial relative inclination of the planet. Eventually, for a large enough eccentricity, the planet is either ejected from the system or collides with the primary component rendering its orbit to be unstable.
Evolution of planet's eccentricity and argument of periapsis for various choices of planet's orbital inclination with the binary plane, simulated for 1 ×10 5 years. The Savitzky-Golay smoothing function with third degree polynomial is used to smooth the high frequency data as shown in (a),(b) and (c). For the first three choices of i pl values, the ω is circulating within 0 • to 360 • . In contrast to the coplanar case, the amplitude of circulation is varying and is maximum at 90 • and 270 • (b and c). When i pl = 39.2 • (d), the ω shows circulation in the region of phase space with eccentricity amplitude maximum at 90 • and 270 • , and for even higher i pl (e) the phase space divides into two distinct regions where ω shows circulation and libration around 90 • and 270 • . This is, however, an unstable region where the amplitude of eccentricity oscillation reaches ∼1. When i pl is set to 45 • , 60 • , 70 • and 80 • (f), the ω librates around 90 • and the libration amplitude in eccentricity is found to increase with increasing i pl until the planet gets ejected from the system for i pl greater than 80 • .

Numerical Approach and methods
Our approach has considered the motion of a planet of mass, m pl around a star of mass, mA. When calculating the initial conditions in position and velocity, we have used the best-fit elements (a, e, i Ω, ω and M ) that are obtained from the radial velocity measurements (Chauvin et al. 2011). For any unconstrained parameters, the values are taken in a range or considered zero. The list of orbital parameters in HD 196885 AB binary system are given in Table 1. The subscripts bin and pl are used to denote the secondary and the planet respectively. The mean longitude (λ) that is used to calculate the resonant angle (Φ) is calculated from the longitude of the periapsis ( = Ω + ω) and the orbit's mean anomaly (M ).
Using the orbital integration package MERCURY (Chambers & Migliorini 1997;Chambers 1999), the builtin Radau algorithm was used to integrate the orbits of the system in astrocentric coordinates when investigating the evolution of orbital elements for a single initial condition and for producing a stability map for multiple initial conditions. MERCURY was effective in monitoring the ejection/collision of a planet due to close encounter with the secondary star and provided robust results for the purpose of determining the existence of orbital resonance. A time step of = 10 −3 year/step was considered to have high pre-cision because the error in change in total energy and total angular momentum was calculated at each time step which fell within the range of 10 −16 to 10 −13 in both cases during the total integration period of 50 Myr. The data sampling (DSP) was done per year for shorter integration periods (up to 100 kyr) and per thousand years for longer integration periods.
We used the MECHANIC 5 software (S lonina et al. 2012) optimised to N-body code to calculate the orbits of the given masses and the MEGNO maps on a multi-CPU computing environment. The MEGNO criterion is used to differentiate between the periodic and chaotic dynamics of a system. Originally developed by Cincotta & Simó (1999, 2000, MEGNO has been widely used in the study of dynamical astronomy Hinse et al. 2008). MEGNO, in general, has the parameterisation Y = α × t + β (see references above). For a quasiperiodic initial condition, we have α 0.0 and β 2.0 (or Y → 2.0) for t → ∞ asymptotically. If the orbit is chaotic, then Y → λt/2 for t → ∞. Here λ is the maximum Lyapunov exponent (MLE) of the orbit. In practice, when generating our MEGNO maps, we terminate a given numerical integration of a chaotic orbit when Y > 12.0. In contrast, quasi-periodic orbits have | Y − 2.0| 0.001. For more details on the mathematical properties of MEGNO and its relationship with Lyapunov exponents, see Hinse et al. (2010).
In our study, each map is generated with 1.5×10 5 initial conditions, for a resolution of (300 x 500) in various orbital elements and simulated up to 200 Kyr. The purple/blue color in the map denotes the quasi periodic region and the yellow is the region of chaos. The color bars in the right side of the maps indicate the amplitude of MEGNO, < Y >.

The Phase Space Evolution
To observe the system that is evolving in the phase space, the planet's eccentricity (e) is plotted versus the argument of periapsis (ω), for various choices of the planet's relative orbital inclination with the binary plane (Fig. 1). The Savitzky-Golay smoothing function is used to filter the data as we are interested in the secular dynamics of the system and the secular perturbation theory which initiates the Kozai mechanism in the case of inclined orbits. We find that, below 39 • and for our choice of i pl values, the phase space starts to evolve where ω starts to circulate between 0 • and 360 • . The eccentricity, however, oscillates with large and increasing amplitude for i pl = 20 • and 35 • while compared to the coplanar case (see Figs. 1a,b,c). When i pl begins at 39.2 • , the ω shows circulation only with eccentricity amplitude maximum at 90 • and 270 • (Fig. 1d), but the system becomes orbitally unstable due to the collision of the planet with the central body within 50 kyr. Further increment in the i pl value to 39.7 • makes the phase space divide into two distinct regions where the ω is found to circulate as well as librate (Fig. 1e). This is also a transition regime of the phase space beyond which the ω displays libration only. Holman et al. (1997) have shown a similar behaviour in the momemtum variable (1-e 2 ) plotted versus its conjugate angle (ω) for a highly inclined (i pl = 60 • ) planetary orbit.
In our case, the planet's orbit remains stable for i pl between ∼(41 • -82 • ) for up to total integration time and in this inclination regime the argument of periapsis librates only around 90 • (Fig. 1f). The libration amplitude increases with the increase in the initial i pl values, eventually entering a regime where the planet can escape from the system when i pl reaches a value greater than 82 • and as emax approaches 1.0 (Eq. 1).

The Kozai Resonance
The evolution of a planet's eccentricity and inclination was further explored for a wide range of i pl values. The conservation of the Kozai integral term (see Sec. 2.1) suggests that as the system evolves in time, the decrease/increase in planet's eccentricity is compensated by the increase/decrease in it's inclination relative to the binary plane, also shown by Holman et al. (1997). For small initial i pl values (i.e., 10 • ) we found a constant-amplitude oscillation of e pl and i pl at the average values of 0.53 and 12 • respectively (Fig. 2a) for the integration period of 50 Myr. As the initial value of i pl is increased the amplitude of oscillations of both parameters increase as well. In Fig. 2b we have shown a case for i pl = 30 • . In this double-axis plot, the y-axis on the left side indicates the time evolution of eccentricity from its best-fit value of 0.48. As the system evolves, the eccentricity varies between 0.25 and 0.55. The y-axis on right side illustrates the variation in inclination from its initial assigned value of 30 • . The i pl value oscillates between 30 • and 40 • and the variation in its amplitude from the initially assigned value behaves in an anti-correlated fashion to that of e pl . The xaxis in Figs. 2a and 2b are truncated at 20 Kyr to clearly demonstrate the oscillations of the two parameters, (e pl and i pl ).
When the initial i pl value was further increased to 39.2 • , the orbital stability is lost within 40 kyr. It is also found that there is a small instability window when i pl is set at ∼(39 • -40 • . The time series plot of e pl and i pl (Fig. 2c) shows that the planet maintains a stable orbital configuration until the e pl reaches as high as 1 and i pl closer to 80 • ; however, the planet is eventually ejected from the system. A closer look shows that when the relative inclination hits the critical angle mark (ic = 39.2 • ), the long-period oscillations between eccentricity and inclination ensue. The initial eccentricity becomes insensitive leading to forced eccentricity, which is the basis for the Kozai resonance. The other factor causing the planet to exhibit Kozai cycles is the libration of argument of periapsis (ω) around 90 • , whose significant evolution is observed in Fig. 1 for i pl greater than 39 • .
The system displays a stable configuration for initial i pl values between 40 • and 82 • . The planet gets ejected from the system only when the i pl value is increased beyond 82 • . For lower inclination values, such as, when initial i pl is set to 50 • (or, 60 • ) (Fig. 3a,b) the amplitudes of eccentricity and inclination oscillations reach as high as ∼ 0.6 and 53 • (or, ∼ 0.8 and 64 • ) respectively and the system continues to maintain the orbital stability throughout the integration period of 50 Myr but could become unstable on a longer time scale. Above the critical value of the planet's orbital inclination with the binary plane (ic = 39.2 • ), the precession of the argument of pericentre is replaced by libration around 90 • as discussed earlier. Satyal et al. (2013) encountered a similar dynamical behaviour in their study of the planet in HD 196885 AB system for the planet's orbital inclination less than 25 degrees. For higher initial i pl , as shown here, the secular perturbation causes the e pl and i pl values to oscillate with higher amplitudes, which can eventually cause the planet to eject from the system. For example, when initial i pl = 83 • the e pl is forced beyond 1 within a short integration period (15 kyr, Fig. 3c).

The MEGNO and the MMR Analysis
We have generated the MEGNO maps considering the planet's initial inclination and semi-major axis which illustrates variations in the quasi-periodic and chaotic phase space for varying i pl and a pl values (Fig. 4). We found that, for i pl less than 39 • , the best-fit value of a pl (the vertical line at 2.6 au) lies within the chaotic region. However, it's important to note that a chaotic orbit does not necessarily imply an unstable one. This is also shown above in Sec. 3 through the orbital integration using a single initial condition and through the global lifetime map of i pl vs. a pl (Fig. 5), where  the colour code is based on the planet's ejection or collision time. The lightest colour represents the stability for up to 50 Myr, comparable to MEGNO's quasi-periodic regions seen in the Fig. 4, left panel. And the dark coloured vertical bands signifies instability comparable to the MEGNO's chaotic bands. The mean motion resonance (MMR) associated with the chaotic region can be estimated by calculating its position from the perturbation theory. For the best-fit values of the planet's and binary semi-major axes, and masses of the stars, the ratio of the (p + q)/p commensurability can be calculated as (see Murray & Dermott (1999)), Where mA and mB represent masses of the primary and secondary star respectively; q is the order of the resonance and p is its degree. To search for the chaos associated with (p+q)/p MMR, we considered the resonant angle (Φ) of the form Where the coefficients follow the relations k1 = p+q, k2 = -p, and i ki = 0. The mean longitude λ is a function of the mean anomaly, M , and longitude of pericenter, . Also, the sum of k5 and k6 is even as required by symmetry in the d'Alembert rules.
For the planet at 2.6 au from the primary component, the mean motion resonance (MMR) is calculated to be near a 39:2 commensurability with the secondary component. Then, this MMR value is used to search for a librating resonant angle (Eq. 3), for k1 = 39 , k2 = -2, and various other combinations of k3, k4, k5 and k6 due to the influence of secular components. In our case, the i pl was varied with respect to the binary orbital plane and the orbits were integrated using astrocentric coordinates, considering the primary at the origin. Thus, k5 in Eq. 3 would be zero as Ω bin becomes undefined (i bin = 0 • ) and convention dictates that it be set to zero.
Initial conditions for different i pl values, but along the best-fit semi-major axis line, were picked where the measure of MEGNO is chaotic (Fig. 4, right panel) and the orbits were integrated for a short period with high data sampling frequency. We would expect the resonance angle (Φ) time series to alternate between modes of circulation and libration for initial conditions taken from the chaotic region and Φ to circulate only when the initial condition is taken from the quasi-periodic region. At 0 • , 10 • , 31 • and 35 • , Φ is found to librate and circulate as expected for chaotic behaviour (Fig. 6, a-d). This also confirms the earlier indication of a 39:2 MMR for the planet with the secondary star. One of the factors inducing the chaos in the low inclination orbit is due to this 39:2 MMR. The MMR in the case of higher inclination orbit changes its dynamical character which results in the periodic orbits. The location of the planet (vertical line at 2.6 au) is almost at the center of the MMR, thus minimizing the amplitude of the libration. The 39:2 MMR stops for i pl greater than 40 • and the planet enters a stable quasi-periodic region for i pl up to 83 • . The resonant angle  Table 1). The vertical instability stripes are seen to build up at 2.45, 2.55, 2.65, 2.75 and 2.85 au where the wider chaotic vertical regions in the MEGNO maps are observed.
is found to circulate only in these quasi-periodic regions. Φ plotted for 50 • has shown circulation in Fig. 6e. This gives a strong evidence that one of the factors responsible for the chaos below 39 • is produced by the 39:2 MMR. Also, the secular perturbation induced by the precession of ω and Ω contributes to the observed chaotic dynamics of the system. The libration and circulation of resonant angle seen in Fig.  6 is obtained only when the coefficients of pl and Ω pl are included in the Φ term.
A small region of possible quasi-periodicity appears when the initial i pl lies between 30 • and 34 • (Fig. 4 right  panel). This would seem like a periodic region but the resonance angle test (Fig. 6c) shows that the Φ is found to circulate and librate, which indicates this region as a chaotic zone. The small black and red dots in this region are indicative that many of these initial conditions have not had enough time to reveal their true nature (quasi-periodic or chaotic). What we called before a quasi-periodic region and expected the circulation of Φ has revealed itself as a region of complex dynamics with many overlapping and interacting resonances. Nevertheless, the dynamical lifetime map shows periodic orbits for up to 50 Myr.
A stripe of unstable chaos is observed when the initial i pl ranges between ∼(35 • -41 • ) (Fig. 4, right panel) and for any choice of a pl , from 2.3 au to 2.9 au. The instability at this region is also confirmed by the lifetime map but for a smaller i pl range (39 • -40 • ). The MEGNO has been successful in determining the chaotic region surrounding this instability point. The minimum critical angle (i pl 39.2 • ) required to initiate the Kozai mechanism lies within this instability region. Once the Kozai resonance phenomenon is triggered, it drives the planet into higher inclination orbits (hence higher eccentricity) causing the planet to collide with the secondary. However, beyond 40 • the planet maintains quasi-periodic orbits even though it is subject to the Kozai resonance. The observed horizontal stripe in the MEGNO map at i pl ∼ 39 • can be considered as a separatrix line which separates the phase space from circulation when i pl l 39.2 • , and libration when i pl 39.2 • The planet is found to lie in the quasi-periodic region for i pl values greater than 40 • and until it rises as high as ∼82 • above the binary plane. The system loses stability beyond this inclination. The instability was confirmed earlier through a time series plot (or the lifetime map) of e pl and i pl (Fig. 3) where e pl becomes greater than 1 and the planet gets ejected from the system. MEGNO maps have been effective at indicating the MMR locations and the chaotic regions in the system while the dynamical lifetime map has addressed whether the chaos is sufficient for an instability. The observed 39:2 MMR is when the mean best-fit value of semi-major axis is 2.6 au. When the planet's location is moved to 2.58 au and 2.68 au (but within the observational uncertainty limit), the planet exhibits oscillations due to the 20:1 and 19:1 MMRs, respectively.
The planet's stability is also analysed in the varying i pl and e pl plane (Fig. 7). The vertical line at e pl = 0.48 indicates the planet's best-fit location, and, like in the previous case (Fig. 4), the planet clearly lies in chaotic zone for e pl > 0.4 and i pl < 30 • . Driven by the Kozai resonance, a similar instability stripe appears at i pl ∼ 40 • . Periodicity continues beyond this angle, however some chaotic islands do appear at i pl ∼ 55 • and i pl ∼ 65 • . The resonance angle test for the possible chaos at these locations came out to be negative.
The discovery of the planet in HD 196885 AB system (Chauvin et al. 2011) using the radial velocity (RV) technique has only allowed the determination of the planet's minimum mass of 2.98MJ . We have performed a test to in- dicate the periodic and chaotic (or unstable) regions of the system for variations in the planet's mass (m pl ) ranging from 2MJ to 20MJ and i pl varying from 0 • to 90 • (Fig. 8). The map indicates that quasi-periodic regions exist for various combination of m pl and i pl . The system is stable for m pl mostly between 3MJ and 6MJ with a small chaotic region around 4.5MJ plus the exclusion of instability stripe at 40 • . The chaotic region continues below 30 • at the best-fit semimajor axis location of the planet. The map indicates other

CONCLUSIONS
Based on our analysis the best configuration for quasiperiodicity of the planet is when the planet's orbital inclination relative to the binary plane lies between 41 • and 82 • . This range of angles gets smaller when the planet's mass is increased more than the current value of 2.98MJ . The planet can be in a 39:2 MMR with the secondary component of the binary which is responsible for the chaotic region as seen in the MEGNO maps below 39 • . The plots of the resonance angle time series for 39:2 MMR clearly demonstrates the circulation and libration of ω based on the choices of initial condition. The planet's orbital configuration below 39 • is proven chaotic, and hence lesser possibility for planet to be in that regime under the current assumptions of planetary orbits (i.e., it has been observed so it must be quasi-periodic although the Lyapunov timescale may be longer than the history of observational astronomy).
Gravitational perturbation on the planet due to the secondary companion gives rise to a highly eccentric orbit. At the best-fit value of planet's eccentricity, the system is driven towards instability when the i pl was sufficiently set to high value, as we saw going from 0 • to 39 • and from 40 • to 82 • . For the planet to be orbitally stable despite its high eccentricity its relative inclination with the binary plane is expected to remain higher. The chaotic behaviour of the planet for the smaller orbital inclinations mainly arises due to the interactions with the secondary. The stellar apastron greatly influences the accretion in the circumbinary disks and for binary periastron distance ∼10 au, up to 6 terrestrial planets could be formed around the primary with semimajor axis lesser or equal to 2.2 au (Quintana et al. 2007) . Since the periastron distance for HD 196885 AB is 12.18 au and the existing giant planet is located at 2.6 au from the primary, the possible existence of other S-type terrestrial planets could not be denied. If that's the case, the possibility that these terrestrial planets affecting the motion of the giant planet is probable. Nonetheless, HD 196885 AB system has stable orbits which could also infer that due to the perturbation from the secondary the planetesimal was forced into precession during its early formation stage, eventually driving the planet into high inclination orbits.
The planets higher mass is possibly constrained between (3 -6)MJ ; however, these choices also depend on the choice of planet's orbital inclination and it is mostly favorable when m pl is less than 6MJ and i pl less than ∼ 70 • , with some exceptions of chaotic islands (Fig. 8). The planetary mass higher than 9MJ and up to 19MJ is most likely when the planet's orbital configuration is inclined between 41 • and 67 • .