A DOUBLE PLANETARY SYSTEM AROUND THE EVOLVED INTERMEDIATE-MASS STAR HD 4732

Precise radial velocity measurements of stars have revealed more than 500 extrasolar planets in the last 20 years. Due to the increase in velocity precision and the duration of the planet search programs, not only single planets but also many multiple-planet systems have been discovered. Around solar-type stars, including those systems with long-term radial velocity trends, ∼30%–50% of stars with giant planets are multiple-planet systems (Fischer et al. 2001; Wright et al. 2007, 2009), and systems comprising less massive planets such as super-Earths and Neptune-mass planets are more abundant (Mayor et al. 2011; Howard et al. 2012; Wittenmyer et al. 2011b).

The multiple-planet systems give us deep insight into the formation and evolution of planetary systems (e.g., Ford 2006). For example, evidence of planet–planet scattering events could be preserved in the current orbital configurations of planets (e.g., Ford et al. 2005), and the prevalence of orbital mean-motion resonances constitutes strong support for differential convergent orbital migration having occurred (e.g., Lee and Peale 2002; Kley et al. 2004). Their statistical properties may also provide us with hints about the mechanisms of these processes; the rather uniform semimajor axis distribution without the pileup of hot Jupiters and the jump at 1 AU that are seen in single planets, and systematically lower eccentricities compared to single planets (Wright et al. 2009).

Studying the long-term orbital stability of multiple-planet systems is also important. Since the best-fit orbital parameters derived from radial velocity data do not necessarily guarantee that the resulting orbits are stable over long periods of time, dynamical stability analysis can check and further constrain the system parameters, such that the solutions consistent with the data are stable (e.g., Lovis et al. 2006; Pepe et al. 2007; Niedzielski et al. 2009a; Wright et al. 2009; Johnson et al. 2011; Robertson et al. 2012b; Wittenmyer et al. 2012c).

Here we report the detection of a double giant-planet system around the evolved intermediate-mass star HD 4732 (K0 IV, M = 1.7 M_☉), one of the 300 targets of the Okayama Planet Search Program (e.g., Sato et al. 2012). This discovery was facilitated by joint precise radial velocity observations made from the Australian Astronomical Observatory (AAO). Discoveries of planets and multiple-planet systems around such intermediate-mass stars have been rapidly growing in number during the past several years. Although the number of known multiple-planet systems around intermediate-mass stars is still small, it is interesting that many of them appear to be in mean-motion resonances, such as 24 Sex (2:1; Johnson et al. 2011), HD 102272 (4:1; Niedzielski et al. 2009a), ν Oph (6:1; Quirrenbach et al. 2011; Sato et al. 2012), and BD20+2457 (3:2; Niedzielski et al. 2009b). Although detailed dynamical studies are required to determine whether these systems are truly in resonance (e.g., Robertson et al. 2012a; Wittenmyer et al. 2012b), the existence of such systems suggests that the planetary orbital migration process for intermediate-mass stars is similar to that for solar-type ones.

This paper is organized as follows. The stellar properties are presented in Section 2 and the observations are described in Section 3. Orbital solutions are presented in Section 4 and results of line shape analysis are described in Section 5. Dynamical stability analysis is presented in Section 6 and Section 7 is devoted to summary.

HD 4732 (HR 228, HIP 3834) is listed in the Hipparcos catalog (Perryman & ESA 1997) as a K0 giant with the apparent V-band magnitude V = 5.90 and the parallax π = 17.70 ± 0.99 mas, giving the distance of 56.5 ± 3.2 pc and the absolute visual magnitude M_V = +2.04 taking into account the correction of interstellar extinction A_V = 0.10 based on Arenou et al.'s (1992) table (Takeda et al. 2008).

Atmospheric parameters (effective temperature T_eff, surface gravity log g, micro-turbulent velocity v_t, and Fe abundance [Fe/H]) of all the targets for Okayama Planet Search Program were derived by Takeda et al. (2008) based on the spectroscopic approach using the equivalent widths of well-behaved Fe i and Fe ii lines of iodine-free stellar spectra. Details of the procedure and resultant parameters are presented in Takeda et al. (2002, 2008). For HD 4732, they obtained T_eff = 4959 K, log g(cgs) = 3.16, v_t = 1.12 km s⁻¹, and [Fe/H] = +0.01 as well as projected rotational velocity of vsin i_rot = 1.45 km s⁻¹.

With the use of these atmospheric parameters together with M_V, a bolometric correction based on the Kurucz's (1993) theoretical calculation, and theoretical evolutionary tracks of Lejeune & Schaerer (2001), Takeda et al. (2008) also determined luminosity L and mass M for HD 4732 to be L = 15.5 L_☉ and M = 1.74 M_☉. The stellar radius R was estimated to be R = 5.4 R_☉ using the Stefan–Boltzmann relationship and the measured L and T_eff. Although the star is listed in the Hipparcos catalog as a giant star, the star is better regarded as a subgiant rather than giant based on the surface gravity and its position on the H-R diagram (Figure 1). The stellar properties are summarized in Table 1.

Zoom In Zoom Out Reset image size

Hipparcos observations revealed photometric stability for the star down to σ_HIP = 0.006 mag. Furthermore, the star shows no significant emission in the core of Ca ii HK lines as shown in Figure 2, which suggests that the star is chromospherically inactive. For subgiants like HD 4732, Johnson et al. (2010) reported a typical radial velocity jitter, which arises from a number of phenomena intrinsic to the star such as granulation, oscillation, and activity, of about 5 m s⁻¹. The value is comparable to that of 6 m s⁻¹ estimated for HD 4732 in Section 4.

Zoom In Zoom Out Reset image size

3.1. OAO Observations

3.2. AAT Observations

For candidate multiple-planet systems, and for planet candidates with orbital periods near one year, it is critical to obtain the most complete possible phase coverage, and it is ideal to independently confirm the signal(s) from independent observatories. A recent example is HD 38283b, a 0.34 M_Jup planet with a period P = 363.2  ± 1.6 days (Tinney et al. 2011b). That planet required 12 years of observations to confirm, and a robust detection was further aided by the fact that HD 38283 is circumpolar and so could be observed year round (although at high airmasses below the pole). Another example is HD 159868c, which has a period P = 352.3 ± 1.3 days (Wittenmyer et al. 2012c), also worryingly close to one year. Again, more than nine years of AAT data were needed for a secure detection, and the 352 day period was confirmed independently with Keck data over a four-year span (Wittenmyer et al. 2012c).

Five years of radial velocity observations of HD 4732 at OAO revealed evidence of a signal with a period near one year (∼338 days). Since this star is near the southern limit for OAO, it is only observable for five months of the year, resulting in persistent phase gaps that make the confirmation of such a period extremely difficult. AAT observations in 2010–2011 filled in the critical phase gaps and confirmed the orbital period suggested by the OAO data. Continued observations from both telescopes in 2011–2012 have also revealed a second velocity signal, first manifesting as a residual trend, and subsequently as a second, long periodicity.

To fit the two Keplerian signals evident in the HD 4732 data, we first employed a genetic algorithm. This approach has proven useful in previous work where it was necessary to fit highly uncertain orbits with long periods near the total duration of observations (Wittenmyer et al. 2012a, 2012c; Horner et al. 2012). We ran the genetic algorithm for 10,000 iterations, each consisting of 1000 individual trial fits. The best-fit set of parameters is thus the result of ∼10⁷ trial fits. The parameters of the best two-planet solution obtained by the genetic algorithm were then used as initial inputs for the GaussFit code (Jefferys et al. 1987), a nonlinear least-squares fitting routine. The GaussFit model has the ability to allow the offsets between data sets to be a free parameter. This is important because the radial velocities from OAO and AAT are not absolute radial velocities, but rather are measured relative to the iodine-free stellar template. Each data set thus has an arbitrary zero-point offset which must be accounted for in the orbit-fitting procedure (Wittenmyer et al. 2009). These offsets are 3.8 ± 3.2 m s⁻¹ (OAO) and −0.1 ± 3.7 m s⁻¹ (AAT).

For the final orbit fitting, we added radial velocity jitter in quadrature to the velocity uncertainty of each observation. This jitter arises from a number of phenomena intrinsic to the star, such as granulation, oscillations, and activity (Saar et al. 1998; Wright 2005). As HD 4732 is an evolved star, scaling relations (Kjeldsen & Bedding 1995, 2011) are not applicable. We thus determine the appropriate level of jitter empirically by performing two-Keplerian fits with varying levels of jitter. On some dates, the AAT recorded multiple consecutive exposures of HD 4732 (Table 3). This provides one way to estimate the jitter: by noting the velocity spread of the multiple exposures. On the four dates that had multiple exposures, the mean velocity spread is 2.9 m s⁻¹. This is, however, a lower limit as the timescales of many jitter sources are on the order of hours (granulation) or months (rotation). We also try a range of jitters from 4 to 10 m s⁻¹, and note the results of the fits in Table 4. The estimation of radial velocity jitter is quite uncertain, and this is reflected in Table 4, where we see that the amount of added jitter has little effect on the quality of the resulting fit. We adopt a jitter estimate of 6 m s⁻¹, as this value brings the reduced χ² of the fit close to unity. The resulting two-Keplerian planetary model is shown in Figures 3 and 4, and the parameters are given in Table 5. The uncertainties on the fitted parameters are derived in the usual way, from the covariance matrix of the least-squares fit. The phase gap for planet b has the greatest impact on the uncertainty in T₀ for that planet. As a result, the uncertainty σ_T = 18.4 days remains large even though many orbital cycles of planet b have been observed. The two-planet fit has a reduced χ² of 1.05, and the rms scatter about the fit is 7.09 m s⁻¹. Using a stellar mass of 1.74^+0.14_−0.20 M_☉ (Takeda et al. 2008), the planets have minimum masses m sin i of 2.4 ± 0.3 (inner planet) and 2.4 ± 0.4 M_Jup(outer planet).

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

In order to investigate other possible causes of observed radial velocity variations such as pulsation and rotational modulation rather than orbital motion, spectral line shape analysis was performed with the use of high-resolution stellar template spectra. Details of the analysis method are described in Sato et al. (2007, 2002).

At first, two stellar templates were extracted from five star+I₂ spectra at phases with different radial velocity level, ∼40 m s⁻¹ and ∼ −20 m s⁻¹, of the observed radial velocities based on the method by Sato et al. (2002). Cross-correlation profiles of the two templates were calculated for about 90 spectral chunks (4–5 Å width each) in which severely blended lines or broad lines were not included, and then three bisector quantities were derived for the cross-correlation profile for each chunk: the velocity span (BVS), which is the velocity difference between two flux levels of the bisector; the velocity curvature (BVC), which is the difference of the velocity span of the upper half and lower half of the bisector; and the velocity displacement (BVD), which is the average of the bisector at three different flux levels. Here we used flux levels of 25%, 50%, and 75% of the cross-correlation profile to calculate the above three bisector quantities, and obtained BVS = 3.9 ± 4.7 m s⁻¹, BVC = 2.4 ± 3.3 m s⁻¹, and BVD = −57 ± 11 m s⁻¹ for HD 4732 (Figure 5). Both of the BVS and the BVC are identical to zero and the average BVD agrees with the velocity difference between the two templates. Then the cross-correlation profiles can be considered to be symmetric and thus the observed radial velocity variations are best explained by parallel shifts of the spectral lines rather than distortion of them, which is consistent with the orbital-motion hypothesis.

Zoom In Zoom Out Reset image size

We here present dynamical studies of the HD 4732 system. Since the eccentricities of the planets are relatively small and it is thus considered that they probably have not experienced close encounters, we here assume that both the planets share the same orbital plane and are prograde. For the numerical integrations, we use a fourth-order Hermite scheme. Figure 6 shows the one-million-year evolution of the best-fit two-Keplerian model derived in Section 4, for i = 90° orbits (i.e., planet masses at minimum). The system has two secular eigenfrequencies related with eccentricity, 11'' yr⁻¹ and 28'' yr⁻¹. Accordingly, the planetary pericenter, a(1 − e), and apocenter, a(1 + e), distances oscillate with the secular timescale of 76,000 years.

Zoom In Zoom Out Reset image size

For lower inclinations, since the mass increases proportional to 1/sin i, mutual perturbations become stronger. The stability map in the (a, e) plane of the outer planet for different values of i is shown in Figure 7. Here we plot the stability index D = |〈n₂〉 − 〈n₁〉| (in ° yr⁻¹) following Couetdic et al. (2010), which studied stability of HD202206 system using Laskar's frequency map analysis (Laskar 1990; Laskar et al. 2001). Following the analysis, we calculate an average of mean motion of the outer planet in 1000 Kepler periods (〈n₁〉) and subtract it from an average of mean motion of the same planet obtained in the next consecutive 1000 Kepler periods (〈n₂〉). Through comparison with the five times longer simulations, we find that the system is regular when log₁₀D ⩽ −3 (blue to navy), while the system is chaotic when log₁₀D > −1 (yellow to red). As inclination decreases, the system loses stability. There is no stable area around 1σ of the best-fitted coplanar orbits for i ⩽ 5°. Thus, in the coplanar prograde configuration, we can roughly set the upper limit on the mass of both of the two planets to be about 28 M_Jup (i > 5°), which falls in the substellar mass regime.

Zoom In Zoom Out Reset image size

In the top left panel of Figure 7, one can find the effects of mean-motion resonances as green vertical incisions. In the case of this system, based on the best-fit values of argument of pericenter, orbits closer to the mean-motion resonances are more unstable. The best-fit solution derived by the radial velocity data (white cross) sits between the 7:1 and 8:1 mean-motion resonances. From the orbital integration shown in Figure 6, we have also confirmed that their resonant angles are not librating. We thus conclude that the outer planet is not in the mean-motion resonances of the inner planet.

We report the detection of a double planet system around the evolved intermediate-mass star HD 4732 (K0 IV, M = 1.7 M_☉) by precise radial velocity measurements at OAO and AAO. The system is composed of two giant planets with minimum masses of m₂sin i = 2.4 M_J, an orbital period of 360.2 days and 2732 days, and an eccentricity of 0.13 and 0.23, respectively. The joint observations of OAO and AAO allowed us to increase the phase coverage of the nearly one-year periodicity of the inner planet and detect long-period signal of the second outer planet.

The configuration of the system is similar to those of other multiple-planet systems currently known: two giant planets in a combination of intermediate- and long-period orbits with relatively low eccentricities. HD 4732 c is the outermost planetary candidate (m_psin i < 13 M_J) ever discovered around intermediate-mass stars except for those discovered by direct imaging.⁹ The period ratio of the two planets is close to 7:1 or 8:1, but dynamical analysis showed that the system is not in the mean-motion resonance and orbits close to mean-motion resonances are unstable in this case.

The one-million-year orbital evolution of the best-fit two-Keplerian model in a coplanar prograde configuration showed no significant changes in eccentricities, but showed their slight oscillations with secular timescale of 76,000 years. Based on the dynamical stability analysis, we found that the system is dynamically unstable in the case of i ⩽ 5° for a coplanar prograde configuration. Thus, we can set the upper limit on the mass of both of the two planets to be about 28 M_J (i > 5°) in this case, which falls well into the substellar mass regime.

This research is based on data collected at OAO, which is operated by National Astronomical Observatory of Japan, and at the Australian Astronomical Observatory. We are grateful to all the staff members of OAO for their support during the observations. We thank the students of the Tokyo Institute of Technology and Kobe University for their kind help for the observations. B.S. was partly supported by MEXT's program "Promotion of Environmental Improvement for Independence of Young Researchers" under the Special Coordination Funds for Promoting Science and Technology, and by Grant-in-Aid for Young Scientists (B) 20740101 from the Japan Society for the Promotion of Science (JSPS). R.W. is supported by a UNSW Vice-Chancellor's Fellowship. M.N. is supported by Grant-in-Aid for Young Scientists (B) 21740324 and H.I. is supported by Grant-In-Aid for Scientific Research (A) 23244038 from JSPS.

This research has made use of the SIMBAD database, operated at CDS, Strasbourg, France.