© ESO, 2013

1. Introduction

Mass loss in close-in hot-Jupiters has been a central problem since the discovery of 51 Peg b (Burrows & Lunine 1995; Guillot et al. 1996). It was first detected in HD 209458b by Vidal-Madjar et al. (2003, 2004) and in HD 189733b (Lecavelier Des Etangs et al. 2010; Bourrier et al. 2013). Several theoretical studies were devoted to this subject (e.g. Yelle 2004; Tian et al. 2005; Hubbard et al. 2007). The question whether smaller close-in super-Earths are remnants of evaporated planets has surfaced since the discovery of CoRoT-7b (Léger et al. 2009; Valencia et al. 2010) and Kepler planets (Lopez et al. 2012; Kurokawa & Kaltenegger 2013). Evidence for evaporating small-mass planets had been limited to theoretical predictions until the discovery of KIC 12557548b.

KIC 12557548b is a unique exoplanet candidate recently discovered by Rappaport et al. (2012) from Kepler long-cadence data (Borucki et al. 2011) obtained during the first two quarters. Unlike all other exoplanets, it exhibits strong variability in transit depth. For some period of time the transits even disappear. The shape of the transit is highly asymmetric with a significant brightening just before the eclipse and a sharp ingress followed by a smooth egress. The planet also has an extremely short period of 0.65356(1) days (15.6854 h). Rappaport et al. (2012) suggested that the planet has a size not larger than Mercury and is slowly disintegrating/evaporating, which creates a comet-like tail. Perez-Becker & Chiang (2013) constructed a radiative-hydrodynamic model of the atmospheric escape from such low-mass rocky planets. The hypothesis that a close-in planet can have a cometary-like tail was suggested by Schneider et al. (1998) and was revisited by Mura et al. (2011). The transit light-curve of dusty comets was first investigated by Lecavelier des Etangs et al. (1999). There is another class of objects that may look very different, but may have features analogous to this light curve. ϵ Aur is an interacting binary with the longest known orbital period, 27.1 yr. The primary star, which is the main source of light, may be either a young massive F0Ia supergiant or an evolved post-AGB star (see Guinan et al. 2012; Hoard et al. 2010 and references therein). The star is partially eclipsed by a dark dusty disk and the light curve shows a very unusual shallow mid-eclipse brightening (MEB).

Planetary transits are commonly modelled using the analytical formulae of Mandel & Agol (2002). This approach assumes a spherical shape of the objects. The JKTEBOP code can calculate and solve the transits numerically assuming the shape of bi-axial ellipsoids (Southworth 2012). The BEER algorithm (Faigler & Mazeh 2011) calculates analytical light curves including a beeming effect. The EVILMC code developed by Jackson et al. (2012) can model transit light curves assuming the Roche shape. The TAU code of Hollis et al. (2013) calculates the transmission spectrum of a spherical planet. The SHELLSPEC code of Budaj & Richards (2004) can calculate planetary light curves assuming the Roche shape, the new model of the reflection effect, and the circumstellar/planetary material. Modelling the transit light curve of KIC 012557548, however, is different and requires at least some modifications to these codes or a new approach. Shortly before the submission of this manuscript, the light curve of this planet was modelled independently by Brogi et al. (2012). The authors started with a thorough data reduction of raw Kepler long-cadence photometry from the first six quarters. They assumed a 1D model of the dust cloud in which the vertical dimension of the cloud was negligible compared with the stellar radius. They also modelled dust extinction as a free parameter and assumed analytical Henyey-Geenstein phase functions and point-source approximation for the scattered light. The particle size was estimated mainly from the overall shape of the transit.

In this paper we first revisit the light curve and orbital period using long- as well as short-cadence Kepler observations from the first 14 quadratures (Sects. 2, 3). Then we modify the code SHELLSPEC to model light curves of such objects. In Sect. 4 we calculate real opacities and phase functions of various dust grains. In Sects. 5 and 6 we construct a 3D model of the dust cloud, calculate the radiative transfer along the line of sight, and take into account a finite dimension of the source of light. We estimate the particle size from the pre-transit brightening feature, which is most sensitive to this parameter. This will enable us to verify whether the shape of the light curve agrees with the idea of a planet with the comet-like tail, put constraints on the particle size and geometry of the dust region, and it will set this interesting object, which may look like a comet, behave like an eclipsing interacting binary, but be an exoplanet, into a more general context.

2. Observations

We used the publicly available Kepler data from the first 14 quarters (Borucki et al. 2011) in the form of SAP flux¹. These are 55 657 long-cadence observations with an exposure time of about 30 min and 257 158 short-cadence observations with an exposure time of about 1 min. Each quarter has a different flux level. In our first step we scaled the flux from all quarters to the same level.

The long-cadence data were used to verify and estimate the orbital period of the planet and search for other periods in these data. We used two independent methods. The first was the Fourier method. A rough scan over the data was made to find an approximate orbital period. Consequently, we phased the data with this preliminary period and set the phase zero point such that the transits occur approximately at the phase of 0.5. Then we split the data into many chunks, each covering one orbital period, fitted a straight line to each chunk of data covering phases 0−1 and removed (divided out) the linear trend from the data. During this process, a phase interval of ⟨ 0.42,0.58 ⟩, which contains the transit, was excluded from the fitting. In this way we removed any changes on time scales longer than the orbital period from the light curve. The advantage of this particular form of detrending is that it does not introduce any additional non-linear trend to the phased light-curve data. Then we again searched for the periods and found a final value of the frequency of 1.530100(4) cyc/day, which corresponds to 0.6535521(15) days (or 15.68525 h). Transit minima can be calculated as ${\mathit{T}}_{\min }\mathrm{[}\mathit{BJD}\mathrm{-}\mathrm{2}\hspace{0.17em}\mathrm{454}\hspace{0.17em}\mathrm{833}\mathrm{]}\mathrm{=}\mathrm{495.44}\mathrm{+}\mathrm{0.6535521}\mathrm{\left(}\mathrm{15}\mathrm{\right)}\mathrm{\times }\mathit{E}\mathit{.}$(1)We estimated the error by performing Monte Carlo simulations generating and analysing about 100 artificial data-sets. This value of the period agrees very well with the value reported in the discovery paper (0.65356(1) days), however, the error of our analysis, which comprises many more data, is about seven times smaller.

Fig. 1 Comparison between the short-cadence (SC, blue) and the long-cadence (LC, red) Kepler transit light curves.

The light curve obtained with this procedure and folded with this improved period was then subjected to a running-window averaging. We used a running window with the width of 0.01 and a step of 0.001 in the phase units and calculated the averaged (phase, flux) points within each window position. The width of the window corresponds to about 0.16 h and its smoothing effect on the data is negligible in comparison with the exposure time or finite angular dimension of the star. This light curve was then used for the modelling. It is depicted in the Fig. 1 and described in more detail in Sect. 6. The standard deviation of points of this binned light curve with phases 0.9−0.1, which are not spoiled by the dusty tail, is only about 3    ×    10^-5. This illustrates that our method works very well and that there are no serious discontinuities between the individual epochs.

The phase dispersion minimization method (PDM, Stellingwerf 1978) was used next. This method is very convenient in cases with non-sinusoidal phase variability and with non-continuous sampling. The method minimizes the variance of the data with respect to the mean light curve, which is described by the parameter theta. We used the latest 100 bin version (pdm2b4.13). A period search revealed that the most significant variability has a period of about 23.0 days and an amplitude of about 1%. This period and its amplitude are not strict but consist of several similar periods that emerge in different data subsets. It might be associated with the stellar differential rotation and some spots at different latitudes. This variability makes the identification of other periods more complicated. The method confirms the orbital period found by the Fourier method.

Fig. 2 Search for a long-term period change. The minimum of the theta parameter corresponds to the period-change rate β = 0.3    ±    0.5 days/Myr, which means that there is no significant long-term period variability.

Note that this rotation period is longer than the orbital period. Consequently, the angular momentum of the planet will be transferred to the star, which would drag the planet inward. At the same time, such a close-in planet may be loosing material and one might expect all kinds of interaction that could lead to the long-term period evolution. That is why we also searched for possible long-term changes of the orbital period, P, with the PDM method. We assumed that the period changes linearly with time t according to the expression $\mathit{P}\mathrm{=}{\mathit{P}}_{\mathrm{0}}\mathrm{+}\mathit{\beta t,}$(2)where P₀ is the period obtained without the period change and β is the period change rate in d/Myr. Figure 2 displays the dependence of theta as a function of β. We fitted a parabola to this curve and obtained β = 0.3    ±    0.5 days/Myr, which means that there is no significant evidence for the long-term period change during the time span of the Kepler observations. The error was estimated by means of Monte Carlo simulations. We generated and analysed 100 artificial data sets with the same standard deviation as the original data. This constant period agrees with the expectations of Rappaport et al. (2012) and Perez-Becker & Chiang (2013).

The short-cadence Kepler observations were treated using the same method. We extracted the SAP flux, normalized it, divided it into individual epochs, fitted and subtracted the linear trend and smoothed it with the running window, which was 0.01 phase units wide with a step of 0.001. The result is plotted in Fig. 3. The comparison between the long- and short-cadence light curves is depicted in Fig. 1. The short-cadence light curve is slightly deeper and has slightly more pronounced pre-transit brightening features. An inspection of these short-cadence data also confirms the important conclusion of Rappaport et al. (2012; based on the long-cadence observations with only 2−3 points covering the transit) that the transit depth has a stochastic behaviour and can change by more than a factor of 2 on the time scale of the orbital period.

Fig. 3 Short-cadence Kepler observations (dots) phased with the orbital period. The linear trend was subtracted and data were averaged (line). See the text.

Fig. 4 Evolution of the comet-like tail. Left: there is a progressive strengthening of the absorption in the tail between Kepler quarters 3 to 6. Right: there is a progressive diminishing of the absorption in the tail between Kepler quarters 6 to 9.

3. Evolution of the tail

Apart from the known variability in the transit depth, which can easily be seen in the rough data (Rappaport et al. 2012), there may also be changes in the “cometary” tail. These changes are not that easily detected since the tail is much finer and not as deep as the core of the transit. To verify this idea we constructed an averaged light curve for the long-cadence observations during each Kepler quarter with the same method as described in the previous section. In this case, a box-car window twice as wide compared with the previous case was used. These light curves were then inspected for variability. One can indeed see in the Fig. 4 that there are apparent changes in the tail on a time scale of about 1 yr. Namely, the tail becomes more pronounced during quarters 3 to 6, then diminishes during quarters 6 to 9. This continues and the tail becomes progressively stronger during quarters 9 to 11, and again ceases during quarters 11 to 13, but this is not shown in these pictures. This indicates that the changes may be periodic with a period of about 1.3 year. Such a period indeed emerges from the period analysis. Curiously, these changes in the tail on these time scales and in these averaged light curves can be even stronger than the changes in the transit core, which happened during the quarters 3−6. The changes in the tail do not seem to correlate with the changes in the transit core.

The origin of this variability needs to be explored. One can speculate that it may be associated with the magnetic/coronal activity of the star or planet (see Sect. 6.6 for a discussion). One cannot also exclude that these cycles may originate from the variability in the transit core, which reflects the variability of the dense material in the close vicinity of the planet that feeds the tail. Another possibility is that the tail evolution is triggered by an interaction of the planet (or its tail) with an additional planet in the system. A TTV analysis may shed more light on the problem. Note, however, that these changes in the tail might also cause spurious changes in the orbital period due to their asymmetric nature. This indicates that the planet is not uniform and has at least two components that change on different time scales: a sort of “coma”/“inner-tail” responsible for the transit core (can change randomly on a time scale shorter than the orbital period) and the tail/“outer-tail” responsible for the egress (changes periodically on a time scale of one year).

Fig. 5 Opacities of forsterite, pyroxene (with different iron content), and iron. Blue (Fo) – forsterite, red (Py00) – pyroxene with 0% of iron (enstatite), green (Py20) – pyroxene with 20% of iron, pink (Fe) – pure iron. Solid line – absorption, dashed line – scattering. Calculated assuming a characteristic radius of the grains of 1 micron.

Fig. 6 Opacities of pyroxene with 20% of iron as a function of particle radius. Red – 0.01 micron, green – 0.1 micron, blue – 1 micron. Solid line – absorption, dashed line – scattering.

Fig. 7 Phase functions at 650 nm for different dust particle radii. Full line – phase functions assuming a point source of light. Dashed line – phase functions that take into account the finite dimension of the stellar disk. The vertical line illustrates the angular dimension of the stellar disk as seen from the planet.

4. Optical properties of the dust

According to Rappaport et al. (2012), the comet-like tail of the planet may consist mainly of dust made of pyroxene grains. Perez-Becker & Chiang (2013) argued that the planet may have lost a significant amount of mass and exposes its naked iron core. Silicates and pure iron grains are refractory species that can endure relatively high temperatures. The optical properties of iron, pyroxene, and another silicate of the olivine family, forsterite, were calculated using the Mie theory and the Mie scattering code of Kocifaj (2004) and Kocifaj et al. (2008), which is based on Bohren & Huffman (1983). By optical properties we mean opacities for scattering, absorption, and phase functions. The complex refraction index of pyroxene was taken from Dorschner et al. (1995). The complex refraction index of forsterite was taken from Jäger et al. (2003). The complex refraction index of iron was taken from two sources: Johnson & Christy (1974) for λ < 0.7 mic and Ordal et al. (1988) for λ > 0.7 mic. We assumed spherical particles of different sizes. To suppress the ripple structure that would appear in the phase function of spherical mono-dispersed particles, the poly-dispersed Deirmendjian (1964) distribution of particle sizes was assumed.

Optical properties of pyroxene (Mg_1−xFe_xSiO₃) in the optical region are quite sensitive to the amount of iron in the mineral. This is illustrated in the example of 1-micron grains in Fig. 5. Opacity (per 1 g of dust material) in the optical and near-infrared (NIR) for micron-sized particles is dominated by the scattering. Scattering is not very sensitive to the amount of iron. However, the absorption opacity in the optical and NIR is very sensitive to the amount of iron and increases with higher iron content. Its direct impact on the spectrum is thus suppressed by the strong scattering. On the other hand, stronger absorption of iron-rich pyroxene will lead to enhanced heating of the grains, which might affect its sublimation and evaporation. Opacities of forsterite (Mg₂SiO₄, another kind of silicates from the olivine family), for 1-micron grains are also shown in that picture. Its scattering properties are very similar to those of pyroxene, but it differs in absorption in the optical and NIR regions. The total (i.e. absorption plus scattering) opacity is almost grey in the optical and NIR region for 1-micron grains. Pure iron differs considerably from the silicates. It does not show the strong absorption features at 10−20 micron. For 1-micron grains its total opacity is lower then that of silicates. However, smaller iron grains have much higher true absorption in the optical region (where it overtakes the scattering) than silicates. The chemical composition of the dust in this system is not known, and therefore we carried out calculations for all these species. However, if not mentioned otherwise (in the light of the above mentioned accounts), we assumed a default model composed of the pyroxene grains with a nonzero iron content (x = 0.2) where 20% of Mg atoms were replaced by Fe.

This situation may change for particles of different sizes. Figure 6 compares the opacities of 0.01-, 0.1- and 1-micron grains of iron-contaminated pyroxene. The true absorption of smaller grains is generally slightly lower in the optical and NIR regions. The scattering opacity overcomes the absorption in the optical and NIR region for larger grains. Scattering on 1-micron grains is very grey in the optical and NIR region and governs optical properties in the 1−8 micron region. Scattering on intermediate 0.1-micron grains has a steep color dependence and governs optical properties at wavelengths shorter than 1 micron. Scattering on even smaller 0.01-micron grains is very weak and has a very strong λ^-4 wavelength dependence that is characteristic of the Rayleigh scattering regime. Absorption dominates the opacity for such small particles. Absorption also dominates scattering for wavelengths longer than 7 microns for all particle sizes smaller than 1 micron. Consequently, multi-wavelength transit observations in the optical and NIR regions might potentially constrain the particle size in this regime.

The phase functions for three different particle sizes at 6500 Å are illustrated in Fig. 7 as a function of the phase angle. The phase angle is the angle between the original and new scattered beam direction. These phase functions exhibit a strong peak near phase angle zero, which is the so-called forward-scattering peak. Larger particles and/or shorter wavelengths tend to have a stronger and narrower forward-scattering peak than smaller particles and/or longer wavelengths. The importance of this forward scattering for the eclipsing systems during the eclipse was stressed by Budaj (2011).

However, these phase functions assume a point source of light and do not take into account the finite and non-negligible angular dimension of the stellar disk as seen from the planet, which is about 26 deg in diameter. To take this effect into account one would have to split the stellar disk into elementary surfaces and integrate over the disk. Note that this process is mathematically analogous and similar to the rotational broadening of the stellar spectra. The difference is that while in the calculations of the rotational broadening the star is split into strips with constant radial velocity, in the calculation of the scattered light, the star may be split into arcs with constant phase angle. For phase angles greater than R^∗/a, where R^∗ is the star radius and a is the semi-major axis, the arcs will resemble the strips and both integration methods will be similar. This regime also involves the pre-transit brightening (see Sect. 6). For smaller angles it will be poor. Nevertheless, this approximation will be much more realistic than the point-source approximation. Consequently, we convolved the phase functions with the “rotational” broadening function (BF, Rucinski 2002) of the star expressed in the angular units. This function is an ellipse whose width corresponds to the angular dimension of the star as seen from the planet and takes into account the linear limb-darkening of the stellar surface. We used the limb-darkening coefficients from Claret (2000; see the next section).

One has to be cautious and calculate the phase functions with a very fine step near the zero angle because of the strong forward scattering and, consequently, its convolution with the BF with a very fine step near the edge of the broadening ellipse. These phase functions, which take into account the finite dimension of the source of light, are also depicted in Fig. 7 for comparison. One can see that the finite dimension of the source of light is important for larger particles but has almost no effect for particles smaller than 0.1 micron. Broadened phase functions are re-normalized to 4π during the calculations.

5. Model and light-curve modelling

The light curves of the exoplanet were calculated using the code SHELLSPEC. This codes calculates the light curves and spectra of interacting binaries or exoplanets immersed in the 3D moving circum-stellar (planetary) environment. It solves simple radiative transfer along the line of sight and the scattered light is taken into account under the assumption that the medium is optically thin. It is written in Fortran77, but there is also a recent Fortran90 version (Šejnová et al. 2012). The code was modified so that it is possible to model dusty objects with comet-like tails. For this purpose a few new objects (structures) were introduced into the code and studied. They have the form of a cone, a ring, or an arc with a variable cross-section. Phase functions and dust opacities can be pre-calculated, loaded in the form of a table, and interpolated to a particular wavelength during the execution of the program.

In these calculations the star is assumed to be a sphere with a radius of R^∗ = 0.65  R_⊙, mass of 0.7 M_⊙, and is assumed to be subject to limb darkening. We used quadratic limb-darkening coefficients for the R filter from Claret (2000), which are based on the Kurucz models. They were interpolated from the grid assuming an effective temperature of T_eff = 4300 K, surface gravity of log g = 4.63 (cgs), metallicity of [M/H]  = 0., and micro-turbulent velocity of v_t = 2 km s^-1 (based on Rappaport et al. 2012). Given the fact that the transit is sometimes missing in the observations, the planet itself must be very small. It is its comet-like tail that gives rise to the major drop in the light curve. Consequently, the planet is modelled as a 3D object in the form of a ring (or part of the ring) with a radius of a = 2.8  R_⊙ and a non-negligible thickness (geometrical cross-section). This cross-section, C, of the ring and dust density, ρ, along the ring are allowed to change with the angle, t in [rad], in the following way: $\begin{array}{ccc}\mathit{\rho }\mathrm{\left(}\mathit{t}\mathrm{\right)}& \mathrm{=}& \mathit{\rho }\mathrm{\left(}\mathrm{0}\mathrm{\right)}\frac{\mathit{C}\mathrm{\left(}\mathrm{0}\mathrm{\right)}}{\mathit{C}}\mathrm{\left[}\mathrm{\right|}\mathit{t}\mathrm{-}\mathit{t}\mathrm{\left(}\mathrm{0}\mathrm{\right)}\mathrm{|}\mathit{/}\mathit{\pi }\mathrm{+}\mathrm{1}{\mathrm{]}}^{\mathit{A}\mathrm{1}}\mathrm{or}\\ \mathit{\rho }\mathrm{\left(}\mathit{t}\mathrm{\right)}& \mathrm{=}& \mathit{\rho }\mathrm{\left(}\mathrm{0}\mathrm{\right)}\frac{\mathit{C}\mathrm{\left(}\mathrm{0}\mathrm{\right)}}{\mathit{C}}{\mathrm{e}}^{\mathrm{|}\mathit{t}\mathrm{-}\mathit{t}\mathrm{(}\mathrm{0}\mathrm{\left)}\mathrm{\right|}\mathit{/}\mathit{\pi }\mathrm{\ast }\mathit{A}\mathrm{2}}\mathit{,}\end{array}$where ρ(0),C(0), and t(0) are the dust density, cross-section and angle at the beginning of the arc, and A1 and A2 are the density exponents to model additional phenomena, e.g. dust destruction. The ring was located in the orbital plane.

6. Results

In this section we describe the effect of various free parameters and dust properties on the light curve of this object and unravel degeneracies and correlations between the parameters. We used the above-mentioned 3D model of the dust cloud with realistic dust opacities, phase functions, finite dimension of the source of light, and place constraints (alternatives) on the physical properties of the dust cloud. No attempt was made to solve for more precise values of so many degenerate parameters using the minimization method. This approach is independent of and complementary to that of Brogi et al. (2012), who assumed an analytical 1D model of the cloud, Henyey-Geenstein phase functions, point source of the scattered light, and a few free parameters that were estimated using the minimization technique.

6.1. Inclination

We carried out a sequence of calculations for the most straightforward assumption that the ring is in the orbital plane of the planet and that this orbital plane has an inclination of i = 90 deg with respect to the plane of sky. It turned out that the observed core of the transit is too sharp. This indicates that the transit is not entirely edge-on but may have i < 84 deg. An alternative explanation could be that the R^∗/a parameter is overestimated or that the dust-particle size in the close vicinity of the planet is larger, which could narrow the transit by filling its ingress and egress by the scattered light with sharper forward-scattering peak. There is an interesting effect that smaller inclinations or grazing transits require slightly higher dust densities and that the forward-scattering peak becomes more pronounced. Ultimately, for the non-transiting planet, the light curve shows the transit brightening only. The inclination i = 82 deg seems to reproduce the observations well and, if not mentioned otherwise, we used this as a default value in the following calculations. This value corresponds to the impact parameter b = 0.60 expressed relative to the stellar radius or b = 0.39  R_⊙ in absolute units. It agrees with the result of Brogi et al. (2012), who estimated several values of the impact parameter in the range 0.46−0.63.

6.2. Geometrical cross-section of the tail

Unfortunately, this quantity cannot be constrained very well. There is a strong degeneracy between this quantity, the dust density, and the optical depth along the line of sight. It stems from the Eqs. (3) and (4) and the best-fitting light curves tend to have the same product ρC. An increase in the cross-section by a factor of 2 results in decreasing the dust density by a factor of 2 and decreasing the optical depth through the tail by the same factor. As a default, we assumed that the tail has a cross-section of 0.1    ×    0.1  R_⊙ at the head and 0.3    ×    0.3  R_⊙ at its end. These values give reasonable optical depths and the properties can be scaled as mentioned above if necessary.

6.3. Dust density profile

The shape of the transit is highly asymmetric with a steep ingress, a short steep recovery, and a long shallow egress. This is the main argument for the idea of a comet-like tail and it suggests that the dust density decreases steeply along the tail. This density profile is described by the Eqs. (3) and (4) and is one of the most important free parameters. This is illustrated in the Fig. 8, which was calculated assuming i = 82 deg, a characteristic particle size of about 0.1 micron, and iron-contaminated pyroxene. The best results are achieved when A1,A2 are about −20 for pyroxene and −25 for iron and this was our default value when studying other effects. It means that the dust density drops by about three orders of magnitude from the head to the end of the tail, which is truncated at 60 degrees from the planet. This is a very steep function of angular distance from the planet and indicates that the dust grains are either effectively destroyed or removed from the tail. Higher exponents fit the core better, while lower exponents fit the tail better. This also suggests that the object may consist of two different components. Each density exponent requires a different density ρ(0) to fit the depth of the transit properly, and this also depends on the particle size and the chemical composition of the dust. In the Table 1 we give the dust density ρ(0) and the highest optical depth through the tail in front of the star corresponding to the particular density exponent for four types of species: pyroxene with 20% of iron, enstatite (pyroxene with 0% of Fe), forsterite, and iron. The calculations assumed an inclination of 82 degrees and the geometrical cross-section at the beginning and end of the tail of 0.01 and 0.09 ${\mathit{R}}_{\mathrm{\odot }}^{\mathrm{2}}$, respectively. The density of the Fe-contaminated pyroxene and enstatite do not differ considerably for larger particles, but enstatite requires higher densities for smaller particles. This is because the opacity at these wavelengths for larger particles is dominated by the scattering, which is not affected by the Fe content, while the opacity of smaller grains is dominated by the absorption and pyroxenes with higher Fe content have higher optical opacity. The densities of forsterite do not differ much from those of enstatite and are only slightly lower for smaller particles. Pure iron differs considerably from the above mentioned silicates. The dust density of smaller iron grains might be lower by two orders of magnitude and still would produce the same transit depth as silicates. This is because small iron grains have a very high true absorption opacity. In contrast, larger iron grains would require slightly higher dust densities than silicates. This is because in this regime the iron opacity is dominated by the scattering, which is not as high as in silicates. These dust densities and dust density profiles were obtained under the assumption that the particle size is constant along the tail. If the particle size were to change along the tail, it might result in different dust density profiles.

Fig. 8 Effect of the dust density profile on the light curve. Observations (red) are compared with calculations for 0.1-micron pyroxene grains with different density exponents (A2 = −10, −15, −20, −25).

Figure 9 illustrates the difference between the power-law behaviour (Eq. (3)) and the exponential behaviour (Eq. (4)) for A1 = A2 =  −20 and ρ(0) = 5    ×    10^-15 for an iron contaminated pyroxene with 1-micron grains. There is only a little difference and the exponential profile fits the core slightly better, while a polynomial profile fits the tail slightly better.

Fig. 9 Comparison between the exponential (exp) and power-law (power) density profile. Kepler short-cadence observations (Obs.) are polotted in red (see the text).

Fig. 10 Effect of the particle size on the light curve. Kepler short-cadence observations (red) vs. calculations for 0.01- (green), 0.1- (blue), and 1-micron grains (purple) of pyroxene. There is a strong degeneracy between this parameter and the dust density.

Fig. 11 Effect of the particle size on the light curve (a more detailed look). Long-cadence observations (red) vs. calculations for 0.01- (green), 0.1- (blue), and 1-micron grains (purple) of pyroxene. Note the significant pre-transit brightening and less significant post-transit brightening, which are caused by the strong forward scattering. It is very sensitive to the particle size (see the text).

6.4. Particle size

This parameter is very important since it affects the absorption, scattering opacity, and emissivity, which is illustrated in Fig. 6. Unfortunately, because of this dependence there is a strong degeneracy between this parameter, the dust density, and the geometrical cross-section of the tail. This degeneracy is demonstrated in Fig. 10, which depicts a comparison between the theoretical light curves and Kepler short-cadence observations. It is possible to fit the transit with dust particles of different sizes and different dust densities. This degeneracy complicates precise estimates of the particle size and dust density from the shape of the transit. Therefore we summarized various alternatives in Table 1. It gives the edge density and density exponent for different particle radii and different grains. It can also be used to infer any possible linear combination of these parameters, if necessary. For example, 0.1-micron pyroxene particles require substantially lower dust densities to reproduce the transit than 1-micron particles, which in turn require much lower densities than 0.01-micron particles. This can be understood using Fig. 6, which shows that at 600 nm the opacity is dominated by the scattering on 0.1-micron particles followed by the scattering on 1-micron particles, and absorption on 0.01-micron grains. The situation with iron is different, as shown in Table 1.

Fig. 12 Same as Fig. 11, but for the Kepler short-cadence observations. Calculated for pyroxene grains and A2 =  −20 (see the text).

Nevertheless, the particle size does leave a few interesting imprints on the light curve. Notice the strong pre-transit brightening and a less pronounced post-transit brightening, which are quite sensitive to the particle size. It is best seen in the long-cadence observations since there are more of them with less scatter than the short-cadence data. Figure 11 compares these observations with the calculations for pyroxene. The pre-transit brightening is caused by the strong forward scattering. This forward scattering is strongly sensitive to the particle size, as demonstrated in Fig. 7. This is why the pre-transit brightening is so sensitive to the particle size and can be used to estimate the dust particles size. Larger-size particles (0.1−1 micron) or their combinations fit the transit core and the pre-transit brightening better, while smaller-size particles (0.01−0.1 micron) or their combinations fit the slow egress and post-transit brightening better. This is yet another indication that the planet is not uniform. The particle size may decrease as a function of the angular distance from the planet, which can be considered as an indication for the dust destruction in the tail. The larger particles exhibit a narrower transit than the smaller particles. This is because larger particles fill the ingress and egress by the strong forward-scattered light. This may cause another degeneracy between the particle size, orbital inclination, and density exponent since they all affect the width of the transit. Brogi et al. (2012) estimated the particle size to be about 0.1 micron in the tail using a different model and analysis.

The latest Kepler short-cadence observations are also very helpful in constraining the particle size. We compare them with the theoretical light curves for pyroxene is Fig. 12. If the calculations are to be interpreted solely in terms of the particle size, then the phases 0.3−0.4 are best reproduced by 1-micron particles. Phases at about 0.45 might be affected and screened by the streams mentioned in the next section, which are not modelled here. Phases 0.5−0.55 are best reproduced by the 1-micron particles. Phases near 0.55−0.57 are best fitted with the 0.1-micron particles, whereas even later phases towards the end of the tail are best reproduced with particles that have radii of about 0.01 micron. An alternative explanation could be that the tail has a different dust density profile or is composed of several chemical or structural components with different A1 or A2 parameters or destruction lifetimes. Figure 13 compares the short-cadence observations with the calculations for iron. A combination of 0.1−1 micron grains also fits the pre-transit brightening best. Iron has a smaller pre-transit brightening feature and requires a steeper density profile A2 =  −25 than silicates. This is because iron has a higher absorption-to-scattering opacity ratio than silicates at these wavelengths.

Finally, note that if we admit that either the particle size or chemical composition of the dust can change along the tail, it may have a strong effect on the dust density profile. For example, if the particle size of silicates were to change by one order of magnitude from 0.1 micron in the head to 0.01 micron in the tail, the dust density profile would not decline so rapidly to compensate for the more than two orders of magnitude lower opacity of smaller grains. Alternatively, according to Table 1, if the chemical composition were to change (assuming constant particle sizes of e.g. 0.01 micron) from for instance pure iron dominating in the head to pyroxene dominating in the tail, then it would also result in a considerably less steep dust density decline or vice versa.

Fig. 13 Same as Fig. 12, but for pure iron grains and A2 =  −25.

6.5. Tail morphology

The above-mentioned ring/arc morphology of the tail seems to reproduce the observations surprisingly well. The tail is about 60 degrees long. The density in the tail decreases rapidly along the tail and becomes negligible beyond 60 degrees away from the planet. Apart from the clear pre-transit brightening, there is a weak post-transit brightening in Fig. 11. This can happen if the tail is terminated by a relatively sudden drop in the density or is dispersed out of the orbital plane, reaches a vertical thickness greater than the radius of the star, and there is still a small but non-negligible fraction of particles with a size of the order of 0.1 micron. In addition to the above-mentioned tail morphology, we also tried some other geometrical models, for example, a dusty tail of the shape of a cone pointing at different angles from the planet. This geometry would produce very strong post-transit brightening that would not be compatible with the observations. Alternative dust density profiles were also considered. For example, ρ ~ e^atb with two free parameters gives encouraging results for a < 0,0 < b < 1 but we preferred to use a more simple one-parameter profile (Eqs. (3) and (4)). The theoretical light curve for 1-micron particles fits the pre-transit brightening best, but shows a slight excess compared to the observations. This is best seen in the comparison with the short-cadence observations in Fig. 12. This excess might disappear if a more realistic model of the dust morphology were considered. Hydrodynamical simulations of Bisikalo et al. (2013) and MHD simulations of Cohen et al. (2011) predicted two streams emanating from the planet. One leaves the planet via the L2 point and is deflected in the direction opposite to the orbital movement, which might result in a comet-like tail similar to the one we observed and modelled here. The other leaves the planet via the L1 point and is deflected in the direction of the orbital motion. This stream is well known in interacting binaries. It might wipe out the excess emission in the pre-transit brightening seen in the theoretical light curves. We did not model this second stream in here since we consider it too premature at this point.

Our model of the planet and what happens during the transit is illustrated in Fig. 14. This is a 2D image of the model that shows the intensity logarithm. One can clearly see the tail with a very steep intensity gradient that becomes much brighter close to the star due to the forward scattering and density gradient. Note also the optically thin absorption on the stellar disk.

Fig. 14 2D image (intensity logarithm, Iν, in erg   cm-2   s-1   Hz-1) of the planet with a comet-like tail during the transit. The intensity of the tail rapidly decreases because of a decrease in the dust density and forward scattering. Note also the optically thin absorption on the limb-darkened stellar disk. Calculated for 1-micron pyroxene grains.

There is one more thing that needs to be mentioned. Long-cadence observations have a relatively long exposure (30 min) that might smear any potential sharp features in the light curve. That is why the theoretical light curves were convolved with a box-car with a width corresponding to the exposure time (0.031 phase units). Short-cadence observations have exposure times (1 min) that are much shorter than the box-car smoothing applied to the observed data. Consequently, theoretical light curves used for comparison with the short-cadence data were smoothed with the box-car with a width corresponding to the box-car width applied to the observations (0.01 in phase units). This had very little effect on the result.

6.6. KIC 12557548b in a more general context

It is very interesting to note that this pre-transit brightening is analogous to the MEB observed in some long-period eclipsing binaries such as ϵ Aur (Budaj 2011) or AZ Cas (Galan et al. 2012). As an explanation of MEB in ϵ Aur it was suggested that the dusty disk eclipsing the F star has a central hole and is inclined out of the orbital plane so that the star can peek through the hole (Carroll et al. 1991). Recently, Budaj (2011) proposed that MEB might be due to the flared-disk geometry and forward scattering on dust and that one does not necessarily have to see the primary star through the hole in the disk during the eclipse. The calculation of Muthumariappan & Parthasarathy (2012) indeed indicated that the disk is flared and the hole is present, but not seen at such almost edge-on inclination. The MEB in AZ Cas may also be caused by the forward scattering during the primary eclipse. During this eclipse the hotter B star is eclipsed by the cool M super-giant. If there is dust surrounding the cool component, it would effectively scatter the light from the hidden hot star to the observer. This effect was predicted and calculated in Budaj (2011) for different geometries of the dusty component. Since very similar dust properties and physical processes can explain the light curves of all three objects, it gives stronger support to our model.

Although these are very different objects, apparently the knowledge obtained in one can be applied to the other. This analogy suggests that there might also be a tail of dust emanating from the disk of ϵ Aur, which could explain non-symmetric behaviour of gas and dust absorption during the eclipse with a stronger and longer egress absorption (Martin et al. 2013; Leadbeater et al. 2013; Sadakane et al. 2013). Moreover, the opacity of iron grains (see Fig. 5) is very “grey”, which can provide an alternative explanation for almost constant depth over a considerable wavelength range of the eclipses of ϵ Aur. On the other hand, there might be a dusty disk inside the Hill sphere of KIC 12557548b, which would feed the comet-like tail. Note, however, that the Hill radius of a planet with the masses of 0.1/0.01 Earth mass would be about 0.014/0.007 R_⊙. Consequently, the material within the Hill radius cannot be responsible for the transit (Rappaport et al. 2012) since it would produce a highest transit depth of about 0.05/0.01%.

KIC 12557548b can be viewed as the missing link between hot-Earths such as CoRoT-7b (Léger et al. 2009), Kepler-10b (Batalha et al. 2011), and comets. The discovery of the long-term comet-like tail evolution of KIC 12557548b as well as its fast transit variability leads us to invoke a possible connection with the magnetic and/or coronal activity of the star or magnetic activity of the planet. There is growing evidence of the interaction between the plasma tail of comets and our Sun, resulting in a spectacular comet-tail disconnection. Comet 2P/Encke, comet C/2001 Q4 (NEAT), and comet Lulin had their plasma tails “dismembered” by the solar wind and/or interplanetary coronal mass ejections (Vourlidas et al. 2007; Clover et al. 2010; Kuchar et al. 2008; Shi et al. 2011). MHD simulations of star-comet or star-planet interaction by Jia et al. (2009) and Cohen et al. (2011) led to the creation of variable comet-like tails of the comet/planet. There is also recent evidence for a comet disintegration, dust tail evolution and disconnection event in comet C/2011 W3 (Lovejoy), which is a member of the Kreutz-system sungrazers (Sekanina & Chodas 2012). Consequently, and in a more general context, this variability may represent a new type of star-planet interaction. Their mutual magnetic interaction may be analogous to that in the RS CVn binaries or the Jupiter-Io system and is currently intensively debated (e.g. Cuntz et al. 2000; Ip et al. 2004; Shkolnik et al. 2008; Lanza 2008; Knutson et al. 2010; Krejčová & Budaj 2012; Gurdemir et al. 2012; Miller et al. 2012; Shkolnik 2013).

7. Conclusions

The light curve of this planet candidate was reanalysed using the first 14 quarters of the Kepler data including first short-cadence observations.

The orbital period of the planet was improved. We searched for long-term period changes, but found no convincing evidence.

A periodic variability in the tail of the planet with a period of about 1.3 years was discovered.

We modelled the light curve of KIC 012557548 using the SHELLSPEC code with the following assumptions: spherical dust grains, different dust species (pyroxene, enstatite, forsterite, iron), different particle sizes, Mie absorption and scattering, finite radius of the source of light (star), and an optically thin medium. We proved that its peculiar light curve agrees with the idea of a planet with a comet-like tail. A model with a dusty ring in the orbital plane with an inclination of about 82 degrees and an exponential or power-law density profile (A1 = A2 = −20 or −25) fits the observations surprisingly well.

It was confirmed that the light curve has a prominent pre-transit brightening. There is an indication of a less prominent brightening after the transit; both are caused by forward scattering.

The dust density in the tail is a steep decreasing function of the distance from the planet, which indicates significant destruction of the tail caused by the star. However, this dust density profile might change considerably if the particle size or chemical composition of the dust were to change along the tail.

The transit depth is a highly degenerate function of the particle size, dust density, chemical composition, and other dust properties. Various combinations of them were estimated and tabulated. There will also be a degeneracy between the inclination, particle size, and the density exponent. Simultaneous multi-wavelength transit observations in the optical and IR regions may help to put more stringent constraints on the particle size, dust density, and chemical composition of the dust.

Forward scattering and the pre(post)-transit brightening are quite sensitive to the particle size. Consequently, we estimated the particle size (radius) of the grains in the head of the tail from the pre-transit brightening to be about 0.1–1 micron. There is an indication that the particle size is larger at the head and decreases along the tail to about 0.01−0.1 micron.

We argued that there are several indications that the “planet” is not homogeneous and that it consists of several components. One component is responsible for the transit core and another for the tail. Components evolve on different time scales, may be subject to different processes, and may have different grain species, density profiles, and/or particle size.

This planet’s light curve with pre-transit brightening is analogous to light curves of some interacting binaries with mid-eclipse brightening, particularly ϵ Aur and AZ Cas, and forward scattering plays an important role in such eclipsing systems. The evolution of the planet’s tail may be connected with the magnetic/coronal activity of the star, magnetic activity of the planet, or an additional planet in the system.

Acknowledgments

I thank M. Kocifaj for his help with his code, R. Komzik and R. Stellingwerf for their help with the PDM code, the anonymous referee, J. Schneider, L. Hambalek, T. Krejcova, J. Svoreň, for their comments and discussions, and H. Mutschke for consultations about the optical data. This work was supported by the VEGA grants of the Slovak Academy of sciences Nos. 2/0094/11, 2/0038/13, by the Slovak Research and Development agency under the contract No.APVV-0158-11, and by the realization of the Project ITMS No. 26220120029, based on the supporting operational Research and development program financed from the European Regional Development Fund.

References

All Tables

All Figures

	Fig. 1 Comparison between the short-cadence (SC, blue) and the long-cadence (LC, red) Kepler transit light curves.
In the text

	Fig. 2 Search for a long-term period change. The minimum of the theta parameter corresponds to the period-change rate β = 0.3    ±    0.5 days/Myr, which means that there is no significant long-term period variability.
In the text

	Fig. 3 Short-cadence Kepler observations (dots) phased with the orbital period. The linear trend was subtracted and data were averaged (line). See the text.
In the text

	Fig. 4 Evolution of the comet-like tail. Left: there is a progressive strengthening of the absorption in the tail between Kepler quarters 3 to 6. Right: there is a progressive diminishing of the absorption in the tail between Kepler quarters 6 to 9.
In the text

	Fig. 5 Opacities of forsterite, pyroxene (with different iron content), and iron. Blue (Fo) – forsterite, red (Py00) – pyroxene with 0% of iron (enstatite), green (Py20) – pyroxene with 20% of iron, pink (Fe) – pure iron. Solid line – absorption, dashed line – scattering. Calculated assuming a characteristic radius of the grains of 1 micron.
In the text

	Fig. 6 Opacities of pyroxene with 20% of iron as a function of particle radius. Red – 0.01 micron, green – 0.1 micron, blue – 1 micron. Solid line – absorption, dashed line – scattering.
In the text

	Fig. 7 Phase functions at 650 nm for different dust particle radii. Full line – phase functions assuming a point source of light. Dashed line – phase functions that take into account the finite dimension of the stellar disk. The vertical line illustrates the angular dimension of the stellar disk as seen from the planet.
In the text

	Fig. 8 Effect of the dust density profile on the light curve. Observations (red) are compared with calculations for 0.1-micron pyroxene grains with different density exponents (A2 = −10, −15, −20, −25).
In the text

	Fig. 9 Comparison between the exponential (exp) and power-law (power) density profile. Kepler short-cadence observations (Obs.) are polotted in red (see the text).
In the text

	Fig. 10 Effect of the particle size on the light curve. Kepler short-cadence observations (red) vs. calculations for 0.01- (green), 0.1- (blue), and 1-micron grains (purple) of pyroxene. There is a strong degeneracy between this parameter and the dust density.
In the text

Fig. 11 Effect of the particle size on the light curve (a more detailed look). Long-cadence observations (red) vs. calculations for 0.01- (green), 0.1- (blue), and 1-micron grains (purple) of pyroxene. Note the significant pre-transit brightening and less significant post-transit brightening, which are caused by the strong forward scattering. It is very sensitive to the particle size (see the text).

In the text

	Fig. 12 Same as Fig. 11, but for the Kepler short-cadence observations. Calculated for pyroxene grains and A2 =  −20 (see the text).
In the text

	Fig. 13 Same as Fig. 12, but for pure iron grains and A2 =  −25.
In the text

	Fig. 14 2D image (intensity logarithm, Iν, in erg   cm-2   s-1   Hz-1) of the planet with a comet-like tail during the transit. The intensity of the tail rapidly decreases because of a decrease in the dust density and forward scattering. Note also the optically thin absorption on the limb-darkened stellar disk. Calculated for 1-micron pyroxene grains.
In the text