Roche tomography of cataclysmic variables - VII. The long-term magnetic activity of AE Aqr
C.A. Hill, C.A. Watson, D. Steeghs, V.S. Dhillon, T. Shahbaz
MMon. Not. R. Astron. Soc. , 1–19 (2015) Printed 1 October 2018 (MN L A TEX style file v2.2)
Roche tomography of cataclysmic variables - VII.The long-term magnetic activity of AE Aqr
C.A. Hill , ∗ , C.A. Watson , D. Steeghs , V.S. Dhillon , and T. Shahbaz , Astrophysics Research Centre, Queen’s University Belfast, Belfast, BT7 1NN, Northern Ireland, UK IRAP, Observatoire Midi-Pyrénées, University of Toulouse, 14 avenue Edouard Belin, 31400, Toulouse, France Department of Physics, University of Warwick, Coventry, CV4 7AL, UK Department of Physics & Astronomy, University of Sheffield, Sheffield, S3 7RH, UK Instituto de Astrofísica de Canarias (IAC), E-38200 La Laguna, Tenerife, Spain Departamento de Astrofísica, Universidad de La Laguna (ULL), E-38206 La Laguna, Tenerife, Spain
ABSTRACT
We present a long-term study of the secondary star in the cataclysmic variable AE Aqr,using Roche tomography to indirectly image starspots on the stellar surface spanning8 years of observations. The 7 maps show an abundance of spot features at both highand low latitudes. We find that all maps have at least one large high-latitude spot re-gion, and we discuss its complex evolution between maps, as well as its compatibilitywith current dynamo theories. Furthermore, we see the apparent growth in fractionalspot coverage, f s , around ◦ latitude over the duration of observations, with a per-sistently high f s near latitudes of ◦ . These bands of spots may form as part of amagnetic activity cycle, with magnetic flux tubes emerging at different latitudes, simi-lar to the ‘butterfly’ diagram for the Sun. We discuss the nature of flux tube emergencein close binaries, as well as the activity of AE Aqr in the context of other stars. Key words: stars: novae, cataclysmic variables – stars: starspots – stars: activity –stars: magnetic field – stars: individual: AE Aqr – dynamo
Understanding the behaviour of stellar magnetic activity,and the nature of the underlying dynamo mechanism, aresome of the most pressing challenges in solar and stellarphysics. It is well known that the Sun displays an 11 yrsun spot cycle. Since the first detection of cyclic magneticbehaviour in solar-like stars (e.g. Wilson 1978), there hasbeen great interest in determining which parameters, suchas binarity, spin, or convective zone depth (and hence stellartype), are pivotal to both the duration and amplitude ofmagnetic activity cycles.In a survey of stellar activity on 111 lower main-sequence stars, Baliunas et al. (1995) used chromosphericCa ii HK measurements as a proxy of the surface magneticfields. They found that, of the stars with solar-like activ-ity cycles, the measured activity cycle periods P cyc rangedfrom 2.5 yr to the 25 yr maximum baseline of observa-tions. They also found that G0–K5V type stars show changesin rotation and chromospheric activity on an evolutionarytimescale, with stars of intermediate age showing moder-ate levels of activity and occasional smooth cycles, whereas ∗ E-mail: [email protected] young rapidly-rotating stars exhibit high average levels ofactivity and rarely display a smooth, cyclic variation.In other work, Saar & Brandenburg (1999) used alarge and varied stellar sample (including evolved stars andcataclysmic variable secondaries) to explore the relation-ships between the length of the activity cycle P cyc andthe stellar rotation period P rot . They parameterized the re-lationships using the ratio of cycle and rotation frequen-cies ω cyc / Ω (= P rot /P cyc ) , as well as the inverse Rossbynumber Ro − ( ≡ τ c Ω , where τ c is the convective turnovertimescale). They found that stars with ages >0.1 Gyr layon two nearly parallel branches, separated by a factor of ∼ in ω c / Ω , with both branches exhibiting increasing ω c / Ω with increasing Ro − . Furthermore, they found that, if thesecondary stars in close binaries can be used as proxies foryoung, rapidly rotating single stars, the cycles of these starspopulate a third ‘superactive’ branch, that shows the oppo-site trend of decreasing ω cyc / Ω with increasing Ro − .Elsewhere, Radick et al. (1998) found that the lumi-nosity variation of young stars was anti-correlated withtheir chromospheric emission, in the sense that young starsare fainter near their activity maxima. This suggests thatthe long-term variability of young stars is spot-dominated,whereas older stars are faculae-dominated (Lockwood et al.2007). Such behaviour has been observed in a number of © 2015 RAS a r X i v : . [ a s t r o - ph . S R ] M a r C.A. Hill et al. young single stars (e.g. Berdyugina et al. 2002; Messina &Guinan 2002; Järvinen et al. 2005) as well as in binary sys-tems (e.g. Henry et al. 1995). Indeed, magnetic activity cy-cles have been found in several systems using photometrictechniques, with some systems appearing to show preferredlongitudes for spot activity. The increase and correspond-ing decrease of spot activity on opposite stellar longitudeshas been interpreted as a so-called ‘flip-flop’ magnetic ac-tivity cycle (e.g. Berdyugina & Tuominen 1998; Berdyug-ina & Järvinen 2005). By tracking the number and posi-tion of spots on the stellar surface using Doppler imaging,such activity was also found on the RS CVn star, II Peg, byBerdyugina et al. (1999).While the magnetic activity of single stars and detachedbinaries is reasonably well studied, studies of magnetic ac-tivity cycles on interacting binaries are critically lacking.Cataclysmic variables (CVs) are semi-detached binaries con-sisting of a (typically) lower main-sequence star transferringmass to a white dwarf (WD) primary via Roche-lobe over-flow. These systems, with both rapid rotation and tidal dis-tortion, provide a unique parameter regime to allow criticaltests of stellar dynamo theories. In addition, CVs form thefoundation of our understanding of a wide range of accretiondriven phenomena, and in turn, the secondary stars are keyto our understanding of the origin, evolution and behaviourof this class of interacting binary. The secondary star regu-lates the mass transfer history and is intimately tied in withthe orbital angular momentum transport that determinesthe evolutionary timescales of the various accretion stages.In particular, magnetic braking is thought to drain angularmomentum from the system, sustaining the mass transferthat causes CVs to evolve to shorter orbital periods. Thishas been a standard ingredient of compact binary evolutiontheory for several decades.Furthermore, magnetic activity cycles in secondarystars have been invoked to explain the variations in or-bital periods in interacting binaries caused by the Applegate(1992) mechanism. This causes angular momentum changeswithin the secondary star throughout the activity cycle tobe transmitted to the orbital motion, resulting in cyclicalorbital period variations. In addition, an increase in thenumber of magnetic flux tubes on the secondary star duringa stellar maximum is thought to cause the star to expand(Richman et al. 1994) and to result in enhanced mass trans-fer – giving rise to an increased mass transfer rate throughthe disc and a corresponding increase in the system luminos-ity. Additional mass transfer also reduces the time requiredto build up sufficient material in the disc to trigger an out-burst, resulting in shorter time intervals between consecu-tive outbursts (e.g. Bianchini 1990). On shorter timescales,starspots are thought to quench mass transfer from the sec-ondary star as they pass the mass-losing ‘nozzle’, resultingin the low-states observed in many CVs (see Livio & Pringle1994; King & Cannizzo 1998; Hessman et al. 2000). Previ-ous surface maps of the secondary stars in the CVs BV Cen(Watson et al. 2007) and AE Aqr (Watson et al. 2006; Hillet al. 2014) show a dramatic increase in spot coverage on theside of the star facing the WD. This suggests that magneticflux tubes are forced to emerge at preferred longitudes, aspredicted by Holzwarth & Schüssler (2003b), and is possiblyrelated to the impact of tidal forces from the nearby compactobject. If these particular spot distributions are confirmed to be long-lasting features, they would require explanation bystellar dynamo theory (e.g. Sokoloff & Piskunov 2002; Mosset al. 2002), and would provide evidence for the impact oftidal forces on magnetic flux emergence. In addition, sincethe number of star spots should change dramatically overthe course of an activity cycle, the density of spots aroundthe mass transfer nozzle may also vary. This would providean explanation for the extended high and low periods seenin polar type CVs such as AM Her (Hessman et al. 2000).Thus, the magnetic activity of CV secondary stars iscrucial to the long and short term behaviour of these sys-tems. Furthermore, it is clear that comparisons of the long-term magnetic activity across a range of stellar types, indifferent systems, are crucial to understanding the nature ofthe stellar dynamo, how it evolves, and what system param-eters are most important in its operation. In light of this,we present a study of the long-term magnetic activity of thesecondary star in the CV, AE Aqr ( P rot = 9 . h) by usingRoche tomography to map the number, size, distributionand variability of starspots on the surface. This is the firsttime a CV secondary has been tracked with this type of acampaign, and given that CVs with both rapid-rotation andtidal distortion provide unique test-beds for dynamo the-ories, we can better understand what parameters are mostimportant to the behaviour of the underlying dynamo mech-anism, and the duration and amplitude of magnetic activitycycles Simultaneous spectroscopic and photometric data of AE Aqrwere taken in 2001, 2004, 2005 and 2006 (hereafter D01, D04,D05, D06), with spectroscopic data only in 2008 and 2009(hereafter D08, D09a, D09b), where D9a and D09b weretaken 9 d apart. As D01 and D09a & D09b have previouslybeen published in Watson et al. (2006) and Hill et al. (2014),respectively, we refer the reader to these works for detailsof the reduction methods for both the spectroscopic andphotometric data. Logs of the observations are shown inTables 1 and 2.
For D04, D05 and D06, spectroscopic observations were car-ried out using the dual-beam Magellan Inamori KyoceraEchelle spectrograph (MIKE, Bernstein et al. 2003) on the6.5 m Magellan Clay telescope, situated at the Las Cam-panas Observatory in Chile. The standard set-up was used,allowing a wavelength coverage of 3330–5070 Å in the bluearm and 4460–7270 Å in the red arm, with significant wave-length overlap between adjacent orders. A slit width of 0.7arcsec was used, providing a spectral resolution of around38,100 ( ∼ . kms − ) and 31,500 ( ∼ . kms − ) in theblue and red channels, respectively. A Gaussian fit to sev-eral arc lamp lines gave a mean instrumental resolution of ∼ kms − , which was adopted for use in Roche tomogra-phy in Section 4. Exposure times of 250 s (0.7 per cent ofthe orbital period) were used in order to minimize veloc-ity smearing of the data due to the orbital motion of the © 2015 RAS, MNRAS , 1–19 he long-term magnetic activity of AE Aqr secondary star. ThAr lamp exposures were taken every 10exposures for the purpose of wavelength calibration.The data were reduced using the MIKE pipeline writtenin python by Kelson (2015). This automatically conductsbias subtraction, flat-fielding, blaze correction and wave-length calibration. The final output provides 1-D spectrasplit into orders, for both the blue and red arms. After re-duction, it was found that each extracted order was not fullyblaze-corrected, and so we applied an additional correctionusing a flux standard star. After flux calibration, the or-ders in the blue and red arms, respectively, were combinedinto continuous spectra by taking a variance-weighted meanacross the spectral range. The blue spectra were then scaledto match the red spectra by optimally-subtracting the over-lapping spectral regions (where the blue spectra was scaledand subtracted from the red spectra, with the optimal scal-ing factor being that which minimizes the residuals). Finally,a variance-weighted mean was made by combining the spec-tra from both arms, creating a single spectrum for each ex-posure. For D08, spectroscopic observations were carried out usingthe Ultraviolet and Visual Echelle Spectrograph (UVES,Dekker et al. 2000) on the 8.2-m UT2 of the VLT, sit-uated on Cerro Paranal in Chile. UVES was used in theDichroic-1/Standard setting (390+580 nm) mode, allowinga wavelength coverage of − Å in the blue arm and − Å in the red arm. A slit width of 0.9 arcsecwas used, providing a spectral resolution of around 46,000( ∼ . kms − ) and 43,000 ( ∼ kms − ) in the blue andred channels, respectively – an instrumental resolution of7 kms − was adopted for use in Roche tomography in Sec-tion 4. Exposure times of 230 s (0.65 per cent of the orbitalperiod) were used, with ThAr lamp exposures taken at thestart and end of the night. The data were taken from the Eu-ropean Southern Observatory (ESO) data products archiveafter being reduced automatically using version 5.1.5 of theESO/UVES pipeline. The final output consisted of 1-D spec-tra for both the red and blue arms. Simultaneous photometry was carried out for D04, D05 andD06 using a Harris V-band filter on the Carnegie Institu-tion’s Henrietta Swope 1-m telescope, situated at the LasCampanas Observatory in Chile. The data were reduced us-ing standard techniques. The master bias frame showed noramp or large scale structure across the CCD, and so thebias level of each frame was removed by subtracting the me-dian value of pixels in the overscan region. Pixel-to-pixelvariations were corrected by dividing the target frames by amaster flat-field taken at twilight. Optimal photometry wasperformed using the package photom (Eaton et al. 2009),where three suitable comparison stars were identified usingthe catalogue of Henden & Honeycutt (1995) to performdifferential photometry. The light curves of D04, D05 andD06 are shown in Figure 1. Flaring activity is clearly evi-dent over both slow and rapid timescales, with amplitudesof up to ∼ . mag. This most likely stems from accretionvariability rather than the secondary star. V m agn i t ude V m agn i t ude V m agn i t ude Figure 1.
The light curves of AE Aqr for D04 (top panel), D05(middle) and D06 (bottom). The points are phase folded for clar-ity, and the typical uncertainties (not shown) are given in Table 2.Rapid and frequent flaring is apparent in all plots, and is due toaccretion variability.
The analysis carried out in this section was completed forthe sole purpose of determining a revised ephemeris in or-der to improve the quality of the Roche tomograms in Sec-tion 7. New ephemerides were determined from the radial- © 2015 RAS, MNRAS000
The analysis carried out in this section was completed forthe sole purpose of determining a revised ephemeris in or-der to improve the quality of the Roche tomograms in Sec-tion 7. New ephemerides were determined from the radial- © 2015 RAS, MNRAS000 , 1–19
C.A. Hill et al.
Table 1.
A log of the spectroscopic observations of AE Aqr. Columns 1-3 list the UT date, the start, and end times of observations,respectively. Column 4 lists the instrument and telescope used. Columns 5-8 show the exposure time, the number of spectra taken, thepeak signal-to-noise ratio around the central wavelength of each spectrum (with the typical value in parentheses), and the phase coverageachieved. Column 9 gives the abbreviation used throughout the text to refer to that specific data set. ut date ut start ut end Instrument T exp (s) No. spectra SNR Phase coverage Abbreviation2001 Aug 09 21:01 04:22 UES+WHT 200 88 −
44 0 . − . D012001 Aug 10 20:49 04:37 UES+WHT 200 95 −
44 0 . − . −
129 ( ∼
96) 0 . − . D042004 Jul 10 03:08 08:57 MIKE+Magellan 250 64 −
116 ( ∼
83) 0 . − . −
121 ( ∼
81) 0 . − . −
99 ( ∼
62) 0 . − . D052005 Aug 06 02:49 07:29 MIKE+Magellan 250 57 −
131 ( ∼
98) 0 . − . −
100 ( ∼
95) 0 . − . D062008 Aug 06 00:01 05:02 UVES+VLT 230 65 −
185 ( ∼ . − . D082008 Aug 07 00:08 05:13 UVES+VLT 230 66 −
176 ( ∼ . − . −
150 ( ∼ . − . D09a2009 Aug 28 00:08 04:55 UVES+VLT 230 61 −
147 ( ∼ . − . −
125 ( ∼ . − . D09b2009 Sept 06 00:34 05:20 UVES+VLT 230 60 −
158 ( ∼ . − . Table 2.
A log of the photometric observations taken of AE Aqr taken with the Henrietta Swope 1-m telescope. Columns 1-3 list thedate, the start and end times of observations, respectively. Columns 4-6 give the exposure time, the number of exposures taken, and thetypical uncertainty in the measured magnitude. Column 7 gives the abbreviation used throughout the text to refer to that specific dataset. ut date ut start ut end T exp (s) No. exp. σ mag Abbreviation2004 Jul 09 02:42 09:22 10 416 ∼ . D042004 Jul 10 02:54 09:02 10 384 "2004 Jul 11 03:14 08:50 10 348 "2005 Aug 05 01:55 07:33 20 290 ∼ . D052005 Aug 06 02:56 07:32 15 272 "2006 Jul 04 03:27 07:33 15 228 ∼ . D06Data set Star γ (kms − )D04 HD 214759 . ± . D05, D06 HD 24916 − . ± . D08 HD 187760 − . ± . Table 3.
Spectral-type templates used to calculate newephemerides for AE Aqr. Columns 1-3 list the data for which thetemplate star was used, the star’s designation, and its systemicvelocity. velocity curves independently for each AE Aqr data set, bycross-correlation with a spectral-type template star, follow-ing Watson et al. (2006) and Hill et al. (2014). The details ofthe template star used for each data set are shown in Table 3,where the systemic velocity was measured by a Gaussian fitto the least-squares deconvolution (LSD, see Section 5) lineprofile for each star (using a line-list where lines with a cen-tral depth shallower than 10 per cent of the continuum wereexcluded).For this, we restricted ourselves to the spectral regionslying between − Å and − Å, as thesecontain strong absorption lines from the secondary star, andreduce the probability of introducing a continuum slopefrom the blue primary. Both the AE Aqr and K4V tem-plate spectra were normalized by dividing by a constant,and the continuum was fit using a third-order polynomial,and subtracted, thus preserving the line strength. The tem-plate spectrum was then artificially broadened (initially by 100 kms − ) to account for the rotational velocity ( v sin i )of the secondary, multiplied by a constant, and subtractedfrom an averaged high-signal-to-noise orbitally-corrected AEAqr spectrum. These latter three steps were repeated, ar-tificially broadening the template spectrum in 0.1 kms − steps until the scatter in the residual spectrum was mini-mized. This typically took two to three iterations. Throughthe above process, a cross-correlation function (CCF) wascalculated for each AE Aqr spectrum, and the peak of theCCF was found using a parabolic fit. A radial velocity (RV)curve was then derived by fitting a sinusoid through theCCF peaks, obtaining new zero-point ephemerides for eachdata set (shown in Table 4), with the orbital period fixed at P orb = 0 . d (from Casares et al. 1996). All subse-quent analysis of each data set have been phased with re-spect to these new ephemerides. Separate ephemerides werecalculated for each data set as the RV curves are affectedby systematics, and so combining all data to calculate a sin-gle global ephemeris may not be optimal. Furthermore, thescatter in the O–C values of all data is relatively small, witha standard deviation of 0.12 per cent of the orbital period.The RV measurements obtained from the cross-correlation method described above are relatively insensi-tive to the use of a poorly-matched template, or an incor-rect amount of template broadening. Surface features suchas irradiation or star spots, as well as the tidal distortionof the secondary are more likely to introduce systematic er-rors in RV measurements, if not properly accounted for (e.g. © 2015 RAS, MNRAS , 1–19 he long-term magnetic activity of AE Aqr Data T (HJD) γ (kms − ) K (kms − ) v sin i (kms − )D01 . ± . − . ± . . ± . -D04 . ± . − . ± .
024 169 . ± . . ± . − . ± .
04 169 . ± . . ± . − . ± .
18 165 . ± . . ± . − . ± .
018 166 . ± . . ± . − . ± .
017 167 . ± . Table 4.
The new ephemerides for each data set of AE Aqr based on the radial velocity analysis described in Section 3. Columns 1-4list the data set, the ephemeris and associated statistical uncertainty, as calculated from the fit to the RV curves, the systemic velocity,and the radial velocity semi-amplitude. Only the ephemerides are adopted for the Roche tomography analysis.
Table 5.
Spectral regions excluded from analysis.Masked region (Å) Comments < Noisy, He i & H emission − He ii emission − H β emission − He i emission & Na i doublet − Tellurics − H α emission − He i emission − Tellurics − He i emission − Tellurics − Tellurics
Davey & Smith 1992). In addition, no detailed attempt wasmade to determine the best-fitting spectral-type or binaryparameter determination in this analysis, however, for com-pleteness we include the systemic velocity ( γ ), the radialvelocity semi-amplitude ( K ) and the projected rotationalvelocity ( v sin i ) in Table 4.Figure 3 shows the measured radial velocities for eachdata set, the fitted sinusoid, and the residuals after subtract-ing the sinusoid. The inherent systematic biases are clearlyevident as deviations from a perfect sinusoid, and as such,the binary parameters derived from this RV analysis havenot been used in the subsequent analysis presented in thiswork.The small variation in γ , K and v sin i between datasets (see Table 4) may be due to instrumental offsets be-tween observations, but the spread in values is most likelydominated by the systematic biasing of RV measurementsdue to surface features such as irradiation and starspots.Such features may alter the slope of the RV curve, changing K , and due to their non-uniform distribution, surface fea-tures may cause a shift in the measured γ . We note that the σ spread in ephemerides is ∼ s, suggesting the periodwas stable over the duration of all observations. Roche tomography is a technique analogous to Dopplerimaging (e.g. Vogt & Penrod 1983), and is specifically de-signed to indirectly image the secondary stars in close bi-naries such as CVs (Rutten & Dhillon 1994, 1996; Schwopeet al. 2004; Watson et al. 2003, 2006, 2007; Dunford et al.2012; Hill et al. 2014), pre-CVs (e.g. Parsons et al. 2015) andX-ray binaries (e.g. Shahbaz et al. 2014). The technique as-sumes that the secondary is locked in synchronous rotation with a circularized orbit, and that the star is Roche-lobefilling. We refer the reader to the references above and thetechnical reviews of Roche tomography by Watson & Dhillon(2001) and Dhillon & Watson (2001) for a detailed descrip-tion of the axioms and methodology.
Least squares deconvolution (LSD) was applied to all spec-tra in the same manner as in Watson et al. (2006), producingmean line profiles with a substantially increased signal-to-noise ratio (SNR). LSD requires that the spectral continuumbe flattened. However, the contribution to each spectrumfrom the accretion regions is unknown, and a constantlychanging continuum slope due to, for example, flaring (seeFigure 1) or the varying aspect of the accretion regions,means a master continuum fit to the data cannot be used.In addition, as the contribution of the secondary star to thetotal light of the system is constantly varying, normalizingthe continuum by division would result in the photosphericabsorption lines from the secondary star varying in relativestrength from one exposure to the next. Hence, we are forcedto subtract the continuum from each spectrum. This wasachieved by fitting a spline to the data. As the spectral typeof AE Aqr has been determined to lie in the range K3-K5V(Crawford & Kraft 1956; Chincarini & Walker 1981; Tanziet al. 1981; Bruch 1991), we generated a stellar line list fora K4V type star ( T eff = 4750 K and log g = 4 . , the clos-est approximation available) using the Vienna Atomic LineDatabase (VALD, see Kupka et al. 2000), adopting a de-tection limit of 0.2. The normalized line depths were scaledby a fit to the continuum of a K4V template star so eachline’s relative depth was correct for use with the continuumsubtracted spectra.Emission lines and telluric lines were masked in thespectra and line list – the excluded spectral regions are de-tailed in Table 5. This meant that over the 4600–7700 Åspectral range for D04, D05 and D06, 2354 lines were avail-able over which to carry out LSD. Similarly, 1558 lines wereavailable for D08 in the spectral range 4780–6810 Å. Aftercarrying out LSD, a small continuum slope was present inthe LSD profiles. This was removed by masking out the linecentre and subtracting a second-order polynomial which wasfit to the continuum. Details for D01 and D09a & D09b maybe found in Watson & Dhillon (2001) and Hill et al. (2014),respectively.The variable light contribution of the secondary meanswe cannot normalize the data in the usual way. Instead, weare forced to use relative line fluxes, requiring the spectra © 2015 RAS, MNRAS , 1–19 C.A. Hill et al. R ad i a l v e l o c i t y ( k m s − ) R ad i a l v e l o c i t y ( k m s − ) R ad i a l v e l o c i t y ( k m s − ) R ad i a l v e l o c i t y ( k m s − ) Figure 2.
The radial velocity curves of AE Aqr for D04 (top left), D05 (bottom left), D06 (top right) and D08 (bottom right). The pointsare phase folded for clarity using the ephemerides in Table 4, and a least-squares sinusoid fit to the RV points (assuming a circular orbit)is shown as a solid line. The lower panels of each plot show the residuals after subtracting the fitted sinusoid, as well as the statisticaluncertainies of the measured RVs. to be slit-loss corrected. For D04, D05 and D06 we used si-multaneous photometry to monitor transparency and targetbrightness variations (see Section 2.2). We corrected for slitlosses by dividing each LSD profile by the ratio of the fluxin the spectrum (after integrating the spectrum over thephotometric filter response function) to the correspondingphotometric flux. The value of photometric flux used wasthe mean over the duration of the spectroscopic exposure.As we were unable to obtain simultaneous photometry forD08, we used the fits from Roche tomography to iterativelyscale the LSD line profiles in the same manner as carried outfor D09a & D09b in Hill et al. (2014). The resulting scaledprofiles were visually inspected and found to be consistent.The final LSD profiles, the computed fits, and the resid-uals after subtracting the fits from the LSD profiles, aretrailed for each data set in Figures 3 to 6, where the or-bital motion has been removed. Starspots and surface fea-tures are clearly visible as emission bumps moving throughthe profiles from negative to positive velocities as AE Aqrrotates. The variation in v sin i due to the tidal-distortionis also clearly apparent. Furthermore, the residuals of D04, D05 and D08 (see Figures 3, 4 and 6) show narrow emissionfeatures that are seen to lie outside the stellar absorptionprofile, and appear to move in anti-phase with respect tothe secondary. Similar features were also seen in the trailsof D09a & D09b (see Figure 3 of Hill et al. 2014), and asthey are visible at all phases, they may be due to accretionmaterial in the system. However, as this emission is weak,we did not carry out any further analysis. The residuals alsoshow the relatively poor fit to the wings of the LSD profiles,resulting from adopting a spherical limb darkening law fora non-spherical object. We note that the SNR of D08 is sig-nificantly higher than for any other data set (see Table 1),resulting in a relatively better fit with Roche tomography.This means that visually, the fits to the wings of the LSDprofiles appear of lower quality, but the absolute level of theresiduals is indeed lower than that of the other data sets. © 2015 RAS, MNRAS , 1–19 he long-term magnetic activity of AE Aqr Figure 3.
Trailed LSD profiles of AE Aqr for D04. The orbital motion has been removed assuming the binary parameters found inSection 6.2, allowing individual starspot tracks across the profiles and the variation in v sin i to be more clearly observed. Panels show(from left to right) the observed LSD profiles, the computed fits to the data using Roche tomography, and the residuals (increased bya factor of 10). Starspots and surface features appear bright in these panels, where a grey-scale of 1 corresponds to the maximum linedepth in the reconstructed profiles. The system parameters (systemic velocity γ , orbital inclina-tion i , primary star mass M and secondary star mass M ) ofAE Aqr were determined using the standard methodology ofRoche tomography (e.g. Watson et al. 2006, 2007). Adopt-ing incorrect system parameters when carrying out Rochetomography reconstructions results in spurious artefacts inthe final image. These artefacts are well characterised (seeWatson & Dhillon 2001), and always increase the amount ofstructure (information content) of the map, decreasing themap entropy. We can constrain the binary parameters bycarrying out map reconstructions for many pairs of compo-nent masses, fitting to the same χ . This can be visualized asan entropy landscape, with an example shown in Figure 7,where the optimal masses are the pair that produce the mapof maximum entropy (least information content). Entropylandscapes are then repeated for different values of i and γ , with the optimal set of parameters those which producethe map containing least structure (the map of maximumentropy).The optimal system parameters are unique to each dataset, as systematic effects may result in different optimal pa-rameters between data sets. Hence, we do not adopt themean values across all data sets for our analysis, as to do somay increase the number of artefacts reconstructed in themaps. This is further discussed in Section 7. Table 6.
Limb-darkening coefficients.Coefficient D04 D05 D06 D08 a a -0.759 -0.768 -0.764 -0.764 a a -0.415 -0.411 -0.412 -0.412 Following Hill et al. (2014), we adopted the four-parameternon-linear limb darkening model of Claret (2000). The stel-lar parameters closest to that of a K4V star were adopted,which for the PHOENIX model atmosphere were log g = 4 . and T eff = 4800 K. The adopted coefficients for each dataset are shown in Table 6, where different (but very simi-lar) values were used for each data set due to the differentcentral wavelengths of the spectra. The treatment of limbdarkening for D01 and D09a & D09b may be found in thepreviously published work of Watson et al. (2006) and Hillet al. (2014), respectively.
All data sets were fit independently. For each, we con-structed a series of entropy landscapes for a range of orbitalinclinations i and systemic velocities γ . For given values of © 2015 RAS, MNRAS000
All data sets were fit independently. For each, we con-structed a series of entropy landscapes for a range of orbitalinclinations i and systemic velocities γ . For given values of © 2015 RAS, MNRAS000 , 1–19 C.A. Hill et al.
Figure 4.
The same as Figure 3 but for D05.
Figure 5.
The same as Figure 3 but for D06. © 2015 RAS, MNRAS , 1–19 he long-term magnetic activity of AE Aqr Figure 6.
The same as Figure 3 but for D08. The fits to the wings of the LSD profiles appear relatively poor, but due to the higherSNR of the data, the absolute level of the residuals is actually lower (see text for discussion).
Figure 7.
The entropy landscape for AE Aqr using D05, as-suming the parameters given in Table 7. Dark regions indicatemasses for which no acceptable solution could be found. The crossmarks the point of maximum entropy, corresponding to compo-nent masses of M = 0 .
86 M (cid:12) and M = 0 .
56 M (cid:12) . i and γ we selected the pair of masses that produced themap of maximum entropy. The results of this analysis arepresented here. Figure 8 shows the map entropy (after adopting optimalvalues of M , M & i ) as a function of systemic velocity,for each data set. Crosses mark the peak of the ‘entropyparabola’, giving the optimal values of γ , as listed in Ta-ble 7. The measured values of γ are consistent with thatof previous work, falling within the uncertainties given byboth Welsh et al. (1995) and Casares et al. (1996). The sig-nificant difference between the γ found here and that foundby Echevarría et al. (2008) stems from the uncertainties inthe absolute radial velocities of the template stars used inthe latter authors’ analysis. The spread in γ as measured byRoche tomography may be explained by instrumental off-sets between different instruments, as well as for the sameinstrument over different observation periods. Values of γ obtained by using entropy landscapes was found to be in-dependent of the assumed inclination, as previously foundby Watson et al. (2003, 2006) and Hill et al. (2014). In ad-dition, the values obtained using the radial velocity curvesare similar, although these will be biased, as discussed insection 3. Figure 9 shows the maximum entropy obtained as a functionof inclination, for each data set, assuming the systemic ve-locities derived in Section 6.2.1. Crosses mark the maximumentropy obtained for a given data set, and the correspondinginclinations are listed in Table 7. All values of i are consis- © 2015 RAS, MNRAS000
56 M (cid:12) . i and γ we selected the pair of masses that produced themap of maximum entropy. The results of this analysis arepresented here. Figure 8 shows the map entropy (after adopting optimalvalues of M , M & i ) as a function of systemic velocity,for each data set. Crosses mark the peak of the ‘entropyparabola’, giving the optimal values of γ , as listed in Ta-ble 7. The measured values of γ are consistent with thatof previous work, falling within the uncertainties given byboth Welsh et al. (1995) and Casares et al. (1996). The sig-nificant difference between the γ found here and that foundby Echevarría et al. (2008) stems from the uncertainties inthe absolute radial velocities of the template stars used inthe latter authors’ analysis. The spread in γ as measured byRoche tomography may be explained by instrumental off-sets between different instruments, as well as for the sameinstrument over different observation periods. Values of γ obtained by using entropy landscapes was found to be in-dependent of the assumed inclination, as previously foundby Watson et al. (2003, 2006) and Hill et al. (2014). In ad-dition, the values obtained using the radial velocity curvesare similar, although these will be biased, as discussed insection 3. Figure 9 shows the maximum entropy obtained as a functionof inclination, for each data set, assuming the systemic ve-locities derived in Section 6.2.1. Crosses mark the maximumentropy obtained for a given data set, and the correspondinginclinations are listed in Table 7. All values of i are consis- © 2015 RAS, MNRAS000 , 1–19 C.A. Hill et al.
Figure 8.
The points show the maximum entropy obtained ineach data set as a function of systemic velocity for AE Aqr. Theoptimal inclination and masses were adopted for each data set, asfound in Sections 6.2.2 & 6.2.3, respectively. The points are offsetin the ordinate for clarity. Crosses mark the optimal value of γ ,and solid lines are shown only as a visual aid. tent with previously published work by Welsh et al. (1995)and Casares et al. (1996), although all values lie below thatfound by Echevarría et al. (2008). Furthermore, all values of i lie between the previously determined inclinations of D01( i = 66 ◦ , Watson et al. 2006) and D09a & D09b ( i = ◦ ,Hill et al. 2014).The consistency of i = ◦ across the four data setsis a reassuring result, as inclination is the worst constrainedparameter when using Roche tomography. We do not havea clear explanation for the discrepancy between the inclina-tions found here and that of D01 and D09a & D09b. The component masses of AE Aqr were determined usingentropy landscapes, with an example shown for D05 in Fig-ure 7. For each data set we assumed i and γ , as derived inSections 6.2.1 & 6.2.2. Our derived masses shown in Table 7are consistent (within the uncertainties) of those found byEchevarría et al. (2008), Welsh et al. (1995), and Casareset al. (1996), once the masses have been adjusted to ac-count for the change in inclination. The differences betweenthe masses determined here, and those found in previousstudies of AE Aqr using Roche tomography, are simply dueto the use of slightly different inclination values. Indeed the Figure 9.
The points show the maximum entropy obtained foreach data set as a function of inclination, assuming the optimalvalues of γ and masses, as found in Sections 6.2.1 & 6.2.3, re-spectively. Crosses mark the optimal value of i , and a solid lineis shown as a visual aid. mass ratios q are in excellent agreement with previous work,and are typically 6 per cent larger than those found byEchevarría et al. (2008), Welsh et al. (1995), and Casareset al. (1996). We note that the masses determined in thiswork are (in principle) the most reliable, as we correct forthe systematic effects of surface features that may bias themeasured RVs used to determine the system parameters inother work.The target reduced χ to which our data were fit waschosen as the point where the entropy of the reconstructedmaps dramatically decreased when fits to a lower χ wereperformed. An increase in small scale structure contributesto a dramatic decrease in entropy, indicative of mappingnoise in the Roche tomograms. Fitting to a higher reduced- χ caused fewer features to be mapped, and thus the systemparameters were less well defined as more map pixels wereassigned the default map value. Figure 10 shows how themap entropy depended on the aim χ , where the adoptedvalue is circled. The absolute value of χ is not a good re-flection of the quality of fit, as a value above or below 1indicates our error bars were systematically under or overestimated.Assigning uncertainties to any of the derived systemparameters ( i, γ, M , M ) is not trivial. As previously dis-cussed in Watson & Dhillon (2001) and Watson et al. (2006),it would require using a Monte Carlo style technique com- © 2015 RAS, MNRAS , 1–19 he long-term magnetic activity of AE Aqr D08 E n t r op y ( a r b r i t r a y un i t s ) Reduced χ Figure 10.
The reconstructed-map entropy as a function of re-duced χ for the Roche tomograms of each data set (marked intop left of panels). The system parameters derived in Section 6.2were adopted for the fits. The selected aim χ is circled in eachplot, and is taken as the point where there is a dramatic decreasein map entropy, and a corresponding increase in small scale fea-tures in the reconstruction. The final reduced χ for D04, D05,D06 and D08 are 0.95, 1.2, 1.3 & 0.26, respectively, where a χ < indicates that we overestimated the size of our propagated uncer-tainties. See Section 6.2.3 for further discussion. bined with bootstrap resampling to generate synthetic datasets drawn from the same parent population as the observeddata. Then, the same analysis carried out in this work wouldneed to be applied to the hundreds of bootstrapped datasets, requiring an unfeasible amount of computation. Hencewe do not assign strict uncertainties to our derived systemparameters. Roche tomograms of the secondary star in AE Aqr wereconstructed for each data set using the system parametersderived in Section 6.2. The corresponding fits to the dataare shown in Figures 3 to 6, and the Roche tomograms areshown in Figures 12 to 15. For ease of comparison, the previ-ously published Roche tomograms of D01 and D09a & D09bare shown in Figures 11, 16 & 17, respectively. These werepreviously analysed by Watson et al. (2006) and Hill et al.(2014), and we highlight the relevant features here. In theanalysis which follows, the map coordinates are defined suchthat ◦ longitude is the centre of the back of the star, with in-creasing longitude towards the leading hemisphere, and withthe L point at ◦ . We note that, due to the inclination ofthe system combined with limitations in the technique, weonly consider features mapped in the Northern hemisphereto be reliable. Hence, the Southern hemisphere is excludedfrom our analysis. A Figure 11.
The Roche tomogram of AE Aqr using D01. Darkgrey-scales indicate regions of reduced absorption-line strengththat is due to either the presence of starspots, the impact of irradi-ation, or gravity darkening. The absolute grey-scales are relativeand are not necessarily comparable between maps. The orbitalphase is indicated above each panel. Roche tomograms are shownwithout limb darkening for clarity.
We have adopted the optimal system parameters de-termined for each data set rather than the mean valuesacross all data sets. This is due to the fact that system-atic effects between data sets may result in different opti-mal system parameters being determined. Hence, adoptingthe mean values may lead to an increase in the number ofartefacts reconstructed in the maps. Nevertheless, we haveassessed the impact of adopting the mean system parame-ters by additionally carrying out the analysis in this sectionfor Roche tomograms reconstructed using the mean valuesof i = 57 . ◦ , M = 0 . M (cid:12) and M = 0 . M (cid:12) (wherethe component masses were calculated from the mean of theconstant M (1 , sin i ). We find that, for each data set, thespot features reconstructed using the mean parameters arenot significantly different to those reconstructed using theoptimal parameters. Furthermore, there are no significantdifferences in the fractional spot coverage as a function oflongitude and latitude (see Section 7.3) on the maps recon-structed using the mean and the optimal parameters. Thisshows the robustness of the surface features reconstructedagainst incorrect inclination and component masses. How-ever, to obtain the same map entropy using the mean pa-rameters, as compared to the optimal parameters, we wererequired to fit the data to a higher aim χ . Thus, we adoptedthe optimal parameters for the spot analysis as they providea better fit to the data. Spot features, both large and small, are clearly prevalent inall tomograms. Common to all maps are the dark regionsaround the L point, primarily due to the time-averagedeffects of irradiation (ionising weak metal lines preventingphoton absorption), as well as gravity darkening. As bothof these effects appear dark in the maps, we are unable to © 2015 RAS, MNRAS000
We have adopted the optimal system parameters de-termined for each data set rather than the mean valuesacross all data sets. This is due to the fact that system-atic effects between data sets may result in different opti-mal system parameters being determined. Hence, adoptingthe mean values may lead to an increase in the number ofartefacts reconstructed in the maps. Nevertheless, we haveassessed the impact of adopting the mean system parame-ters by additionally carrying out the analysis in this sectionfor Roche tomograms reconstructed using the mean valuesof i = 57 . ◦ , M = 0 . M (cid:12) and M = 0 . M (cid:12) (wherethe component masses were calculated from the mean of theconstant M (1 , sin i ). We find that, for each data set, thespot features reconstructed using the mean parameters arenot significantly different to those reconstructed using theoptimal parameters. Furthermore, there are no significantdifferences in the fractional spot coverage as a function oflongitude and latitude (see Section 7.3) on the maps recon-structed using the mean and the optimal parameters. Thisshows the robustness of the surface features reconstructedagainst incorrect inclination and component masses. How-ever, to obtain the same map entropy using the mean pa-rameters, as compared to the optimal parameters, we wererequired to fit the data to a higher aim χ . Thus, we adoptedthe optimal parameters for the spot analysis as they providea better fit to the data. Spot features, both large and small, are clearly prevalent inall tomograms. Common to all maps are the dark regionsaround the L point, primarily due to the time-averagedeffects of irradiation (ionising weak metal lines preventingphoton absorption), as well as gravity darkening. As bothof these effects appear dark in the maps, we are unable to © 2015 RAS, MNRAS000 , 1–19 C.A. Hill et al.
Table 7.
System parameters. Columns 1-6 list the data set or paper from which the parameters were taken, the systemic velocity asmeasured by Roche tomography, the inclination, the mass of the primary star, the mass of the secondary star and the mass ratio. Thesignificantly higher component masses found in D09a & D09b are due to the lower inclination found in that study.Author / Data γ (kms − ) i (degrees) M (M (cid:12) ) M (M (cid:12) ) q = M / M D01 −
66 0.74 0.50 0.68D04 − .
60 0.84 0.55 0.65D05 − .
59 0.86 0.56 0.65D06 − .
59 0.87 0.56 0.64D08 − .
57 0.94 0.64 0.68D09a − . ± .
50 1.20 0.81 0.68D09b − . ± .
51 1.17 0.78 0.67Echevarría et al. (2008) -63 ± . ± .
05 0 . ± . − . ± . ± . ± .
16 0 . ± . − ± . ± . . ± .
23 0 . ± . BC Figure 12.
The same as Figure 11 but for D04. clearly distinguish between the two. Likewise, we are unableto disentangle any starspots that may inhabit the affectedregion. Nevertheless, the impact of gravity darkening on thefractional spot coverage, f s , was assessed, and is discussedin Section 7.3.The map of D01 (see Figure 11) clearly shows a sin-gle large spot (labelled ‘A’) extending 60–80 ◦ latitude and260–320 ◦ longitude, with a second prominent spot centredat 50 ◦ latitude and ◦ longitude. Also clearly apparentis a spot extending from ∼ ◦ latitude down to the L point, becoming indistinguishable in the irradiated region.The map of D04 (see Figure 12) shows two separate largespots, with the largest (labelled ‘B’) extending 40–75 ◦ inlatitude and 160–220 ◦ in longitude, and the second largest(labelled ‘C’) covering 60–70 ◦ latitude and 310–340 ◦ longi-tude. The map of D05 (see Figure 13) has a single large spot(labelled ‘D’), extending 55–75 ◦ in latitude and 140–210 ◦ inlongitude. The poor phase coverage of AE Aqr in D06 re-sulted in the reconstructed map having a relatively feature-less leading hemisphere, with spots on the rest of the starbecoming smeared out over the image, resulting in a lowercontrast (see Figure 14). Despite this, at least one, possi-bly two large spots (labelled ‘E’) are evident above ∼ ◦ latitude, centred on ∼ ◦ longitude, with another promi-nent feature extending from ◦ down to the L point. The D Figure 13.
The same as Figure 11 but for D05. map of D08 (see Figure 15) exhibits a single large spot (la-belled ‘F’) extending 65–85 ◦ and 340–060 ◦ in latitude andlongitude, respectively. Finally, the maps of D09a and D09b(see Figures 16 and 17) show one large spot (labelled ‘G’)spanning 65–90 ◦ latitude and 340–050 ◦ longitude.The latitude of the largest spot in each map (labelledA, B, D, E, F, G) remains fairly constant over the 8 yearsbetween the first and last observations, ranging between 60–80 ◦ . However, the longitude of the dominant spot is notfixed. In D01, the dominant spot (A) lies ∼ ◦ longi-tude, whereas for D04, D05 and D06, the dominant spot(labelled B, D and E, respectively) lies ∼ ◦ longitude,although the position of spot E in D06 is less certain due tosmeared features. In contrast, the largest spot in D08 (la-belled F) and D09a & D09b (labelled G) lies at ∼ ◦ longi-tude. Such monolithic spots are prevalent in many Dopplerimaging studies of rapidly-rotating solar-type stars such asLQ Hya (Donati 1999), and in CVs such as BV Cen (Watsonet al. 2007). The high-latitude spots imaged here, and theirpossible evolution, are further discussed in Section 8.Starspots are also prevalent at low to mid latitudes in allmaps, and in order to make a more quantitive assessment oftheir properties and the underlying dynamo mechanism, wemust consider their size and distribution across the stellarsurface. © 2015 RAS, MNRAS , 1–19 he long-term magnetic activity of AE Aqr E Figure 14.
The same as Figure 11 but for D06. F Figure 15.
The same as Figure 11 but for D08. G Figure 16.
The same as Figure 11 but for D09a. G Figure 17.
The same as Figure 11 but for D09b.
To determine the spot coverage in the Roche tomograms,it was first necessary to define the pixel intensity of theimmaculate photosphere as well as that of a totally spot-ted pixel. We do not adopt a two-temperature model whenfitting with Roche tomography, where a spot filling factoris predetermined (e.g. Collier-Cameron & Unruh 1994), assecondary stars in CVs are expected to exhibit large tem-perature differences due to irradiation by the primary. Ourmethod of determining a totally spotted pixel was to simplyselect the lowest pixel intensity at the centre of the largestspot feature. The adopted value of a totally spotted pixel foreach map is shown as a dotted line on the left side of the his-tograms of pixel intensities in Figure 18, where the brightestpixel is assigned an intensity of 100 and other pixels scaledlinearly relative to this. Pixels with a lower intensity arepresent in the Roche tomograms, but these are confined tothe region around the L point, and as discussed above, arenot likely to be due to a spot feature. Pixels on the Southernhemisphere are not included in the histograms or any of theanalysis below for two reasons. Due to the inclination of thebinary, a large portion of the surface is not visible, and so asubstantial number of pixels on the Southern hemisphere areassigned the default map value. Furthermore, as RVs cannotconstrain whether a feature is located in the Northern orSouthern hemisphere, features may be mirrored about theequator, reducing their contrast as they are smeared over alarger area.The intensity of the immaculate photosphere was moredifficult to define due to the growth of bright pixels – anartefact known to affect maps that are not thresholded (e.g.Hatzes & Vogt 1992). To assess the extent of bright pixelgrowth in our Roche tomograms, we carried out reconstruc-tions of simulated maps with a random spot distribution,adding in varying levels of noise. We found that the brightestreconstructed pixel was between 4–15 per cent higher in in-tensity than that of the original map, and that up to 12% ofpixels were classed as ‘bright’ (typically ∼ per cent). Datawith a higher SNR increased the number of bright pixels,however, the most dramatic increase was found when the © 2015 RAS, MNRAS , 1–19 C.A. Hill et al.
Figure 18.
Histograms of the pixel intensities in the Roche to-mograms of AE Aqr for each data set, where the pixel density isthe fraction of the total number of pixels. Pixels on the Southernhemisphere (latitude < ◦ ) were not included (see Section 7.2 fordetails). The brightest pixel in each map was assigned an intensityof 100 and all other pixel intensities were scaled linearly againstthis. The definition of the pixel intensity representing the immac-ulate photosphere is shown as a dashed line on the right side ofthe histogram (a dotted line for D09b), and that representing atotally spotted pixel is shown as a dashed line on the left side (seetext for details). data were fit to progressively lower reduced χ . Indeed, themaps with a large number of bright pixels were clearly over-fit, exhibiting substantial reconstructed noise. Given thatour Roche tomograms of AE Aqr were fit to an aim re-duced χ that limited the amount of reconstructed noise (seeFigure 10), we assumed a percentage growth of bright pixelsof 3% of the total number of pixels in the map. This some-what conservative estimate, given the results of our simula-tions, means we likely underestimate f s . Hence, we definedthe immaculate photosphere as the lowest pixel intensitythat includes 97% of all pixels, and is shown in Figure 18 asa dotted line on the right side of the histograms.The histograms of pixel intensities in Figure 18 showbroad peaks, with long tails towards lower intensities. Anidealised histogram would have a significantly bimodal dis-tribution of pixel intensities, where spotted pixels and thoserepresenting the immaculate photosphere would be moreclearly separated. The large number of intermediate pixel intensities found here may be explained by a lack of contrastin the maps, stemming from both a population of unresolvedspots, as well as spots that have been smeared in latitude,increasing their areal coverage. Figures 19 and 20 show the fractional spot coverage f s asa function of longitude and latitude for each map, respec-tively. We calculate f s using Equation 1, where I is the pixelintensity, I p is the intensity of the immaculate photosphere,and I s is the intensity of a totally spotted pixel. f s = max (cid:20) , min (cid:20) , I p − II p − I s (cid:21)(cid:21) (1)The most prominent feature common to all plots is thelarge value of f s around ◦ longitude in Figure 19, as wellthe increased coverage below ◦ latitude in Figure 20. Theseregions include the features around the L point which aredominated by the effects of irradiation and gravity darken-ing (as discussed above). The impact of these phenomenaon the maps was assessed by simulating blank maps usingthe corresponding system parameters and limb darkeningcoefficients specific to each data set. Additionally, a gravitydarkening coefficient of β = 0 . was adopted, as this is rep-resentative of that measured for late-type secondary starsin close binaries (Djurašević et al. 2003, 2006). The simu-lated maps were used to create synthetic line profiles withthe same orbital phases, exposure times and instrumentalresolutions of the original data. These synthetic line profileswere then reconstructed in the same manner as the originaldata, and f s was calculated using the same definitions of atotally-spotted pixel and that representing the immaculatephotosphere, as determined above. The value of f s of thesesynthetic maps was then subtracted from that of the originalmaps, effectively removing the systematic effects of inclina-tion, phase sampling, and incorrect limb and gravity dark-ening in the spot coverage of the original maps. However, wecannot distinguish between a spot and the effects of gravitydarkening or irradiation, as all three appear as dark regionsin the maps. Hence, when we apply the correction describedabove, spots that are located in the regions most affected byirradiation and gravity darkening are not preserved (as theregions are made brighter, regardless of spots being present;see discussion in Section 7.1), and thus we may underesti-mate f s in those regions (see Figures 19 and 20). The fractional spot coverage f s as a function of longitudevaries significantly for each map (see Figure 19). Even aftersubtraction of the simulated maps there still remains a sig-nificant spot coverage around ◦ longitude for D04, D05,D06 and D09a & D09b. The map of D01 shows an increasein f s around ◦ longitude, with the maps of D04 and D06showing a significant increase in f s between 280–360 ◦ lon-gitude, and the maps of D08 and D09a & D09b showing alarger f s over a broader range of 280–060 ◦ longitude. Theincrease in f s between 0–120 ◦ longitude for D06 is not realand results from poor phase coverage, leading to the pixels © 2015 RAS, MNRAS , 1–19 he long-term magnetic activity of AE Aqr Figure 19.
The fractional spot coverage f s as a function of lon-gitude for the Northern hemisphere of AE Aqr for each data set,where f s is normalized by the number of pixels within a ◦ lon-gitude bin. The solid line shows f s for the original map, and thedotted line shows f s after subtracting f s for the simulated map(removing the effects of irradiation and gravity darkening, seeSection 7.3 for details). The shaded region for D06 indicates thatfeatures in this region are not reliable due to a lack of phasecoverage. in this region being assigned a value similar to that of thedefault map (the mean of all map pixels).The distributions of f s in longitude suggests the ex-istence of two longitude regions, separated by ∼ ◦ ,with significantly higher f s in the maps of D04, D06 andD09a & D09b. However, no such distributions are presentin D01 and D05, with the high f s around ◦ in D08 be-coming much lower once the effects of irradiation and grav-ity darkening have been removed. Such ‘active longitudes’are observed in single stars such as LQ Hya, AB Dor andEK Dra (Berdyugina 2005) as well as in RS CVn bina-ries (e.g. Rodonò et al. 2000). In shorter period systems( P rot < d), active regions are preferentially located atquadrature longitudes (e.g. Olah et al. 1994; Heckert et al.1998), whereas longer period systems show no such prefer-ence. Indeed the RS CVn binaries HD 106225 (Strassmeieret al. 1994) and II Peg (Rodonò et al. 2000) show a mi-gration of active longitudes with respect to the companionstar, ascribed to differential rotation. However, the lack ofclearly defined active longitudes in the maps of D01, D05and D08 suggest that such active regions are not permanentfeatures in AE Aqr, or are at least not fixed with respectto the companion. We discuss possible explanations of theobserved stellar activity in Section 8. The fractional spot coverage f s as a function of latitude isshown in Figure 20 for all maps. Latitudes above ◦ showhigh f s due to the large high-latitude spots. The distribu-tion of f s above ◦ latitude is similar in D01, D08, andD09a & D09b, with maximum f s at latitudes > ◦ . In com-parison, f s peaks around ◦ latitude for D04 and D06, withthe flatter distribution of D05 exhibiting a peak that is con-sistent with D04 and D06.At mid-latitudes of ∼ ◦ we see the apparent growthof f s over the 8 years of observations, with a peak becomingmore pronounced in the plots of D08 and D09a & D09b.This may be indicative of increasing spot coverage at thislatitude, however, it could also be explained by a relative de-crease in the spot coverage at surrounding latitudes. Givenour lack of an absolute spot filling factor, the total spot cov-erages between maps (given as a percentage in the top-leftof each panel in Figure 20) are not necessarily comparableon an absolute scale, rather, they indicate the relative spotcoverage for that particular map. However, given the mea-sured spot coverage increases between D05 and D09b, wecan be confident the increase of f s at ◦ latitude is indeeddue to an increase in spot coverage localised to this latitude.A common feature to all maps, excluding that of D01, isthe apparent lack of spots around ◦ latitude. This is mostobvious for D04, D05 and D08, as f s is lowest at this latitude,and becomes more pronounced for D06 and D09a & D09bwhen the systematic effects of irradiation and gravity dark-ening are subtracted. The lower quality of data in D01 meansfeatures at lower latitudes are more smeared out comparedwith other maps in our sample, and so while this featuremay exist, we cannot resolve it to the same degree as in theother maps.At lower latitudes, we see a clear increase in f s around ◦ for all maps. As previously discussed, the effects of ir-radiation and gravity darkening are significant at these lat-itudes, which makes it difficult to determine how f s variesbetween maps. However, even after these effects have beensubtracted, we still see a persistently high-level of f s around ◦ latitude, suggesting spots at this latitude are common.The reliability of the reconstruction of small scale fea-tures in the maps was tested with reconstructions of simu-lated data sets. Test maps were created using the same pa-rameters as those for D04, D05 and D08, with a large polarspot, bands of spots at ◦ and ◦ latitude, and gravity andlimb darkening. The bands of spots each contained at least17 individual spots, with sizes ranging between 5–10 ◦ in lat-itude and longitude, separated by 10–20 ◦ longitude. Trailedspectra were created using the same phases as the originaldata, and representative noise was added. Maps were recon-structed from the synthetic data in the same manner as thatcarried out for the original data. The resulting maps of allthree simulated data sets show that the polar spot and thespot band at ◦ are clearly recovered, with the spot band at ◦ becoming moderately smeared in latitude, reducing thecontrast. However, the latter spots are still clearly distinctfrom the reconstructed noise, and so we are confident thatwe can reliably reconstruct features of this size at these lati-tudes. Furthermore, the latitudinal spread in f s due to thesesmeared spots has been taken into account in the analysispresented here. © 2015 RAS, MNRAS , 1–19 C.A. Hill et al.
Figure 20.
The fractional spot coverage f s as a function of latitude for the Northern hemisphere of AE Aqr for each data set, normalizedby the surface area at that latitude. The solid line shows f s for the original map, and the dotted line shows f s after subtracting f s for thesimulated map (removing the effects of irradiation and gravity darkening, see Section 7.3 for details). The total spot coverage for eachmap is given as a percentage in the top left of each panel, where the value in parenthesis is the total spot coverage after subtracting thesimulated map. Several high-latitude spots are seen in the maps of AE Aqr,and are labelled A-G in Figures 11-17. There are severalpossible scenarios that may explain the behaviour of thedominant high-latitude spot; The largest spot in each mapmay be the same feature, with only its position changingbetween observations. Alternatively, the largest spot maynot be the same feature, but instead is evolving significantly,disappearing and appearing elsewhere between maps.In the first scenario, if the largest spot in each map (la-belled A, B, D, E, F, G) is the same feature, and is ableto move position, then stable spots may have a lifetime of ∼ yr. However, if the largest spot must remain in a fixedposition, then spots B, D and E in Figures 12–14 implythat a large, stable spot can live for ∼ yr. The chang-ing position of the dominant spot may be explained by thepresence of differential rotation (DR) on the surface of thestar, as measured by Hill et al. (2014) using the maps ofD09a & D09b. By using the measured shear-rate, we havedetermined the longitude-shift for latitudes between 60–80 ◦ over the intervals between observations. By comparing theshift in longitude due to DR with the observed change inlongitude of the spot in question, we find that spot A inD01 may have shifted to the position of spot B in D04, andspot E in D06 may have shifted to the position of spot Fin D08. Only the movement of the spots in these two casesmay be explained by DR, with all other observed shifts be-ing incompatible with this mechanism. It is possible thatthe spot C in D04 grew in area and rotated due to DR to the same position as spot D in D05, however, the spot B inD04 would then have to disappear completely in D05, whichseems unlikely. Furthermore, the fact that spot F in D08 isat the same location as spot G in D09a & D09b suggests DRis unlikely to be the cause of the observed shifts. Indeed DRmay primarily affect smaller spots that have magnetic fluxtubes anchored closer to the surface, whereas larger spotsare less affected as their flux tubes may be anchored deeperin the stellar interior.The second scenario requires the evolution of these dom-inant spots over a relatively short time period. Similar be-haviour, observed in other systems, has been explained bya phenomenon known as the ‘flip-flop’ activity cycle. In thisscenario, active longitudes are present on the star. The ‘flip-flop’ cycle occurs when the active longitude with the highestlevel of activity (i.e. the most spotted) switches to the oppo-site longitude, with cycles taking years to decades to occur(see Berdyugina 2005 for a summary). The disappearance ofspot C in D04, as compared to D05, may imply the switch-ing of dominant longitudes between observations of AE Aqr.However, the lack of at least two clearly defined active lon-gitudes in D01, D05 and D08 suggests that this type of cycleis not present in AE Aqr. In any case, a robust detection ofa ‘flip-flop’ cycle in AE Aqr would require a much shortergap between observations than obtained for our sample inorder to clearly track the emergence of spots at preferredlongitudes. © 2015 RAS, MNRAS , 1–19 he long-term magnetic activity of AE Aqr To understand why surface features are distributed as theyare, and why they evolve in the observed manner, we cancompare the Roche tomograms to numerical simulationsof emerging magnetic flux tubes in close binary systems.Holzwarth & Schüssler (2003a) carried out such simulations,assuming that starspots are formed by erupting flux tubesthat originate from the bottom of the stellar convective zone(by analogy with the Sun). Magnetic fields, believed to beamplified in the rotational shear layer (tacholine) near thebase of the convection zone, are stored in the form of toroidalflux tubes in the convective overshoot layer (Schussler et al.1994). By studying the equilibrium and linear stability prop-erties of these flux tubes, Holzwarth & Schüssler (2003a,b)examined whether the influence of the companion star wasable to trigger rising flux loops at preferred longitudes, sincethe presence of the companion breaks the rotational symme-try of the star.The authors find that while the magnitudes of tidal ef-fects are rather small, they nevertheless lead to the forma-tion of clusters of flux tube eruptions at preferred longitudeson opposite sides of the star, a phenomenon resulting fromthe resonating action of tidal effects on rising flux tubes.Pertinently, the authors establish that the longitude distri-bution of spot clusters on the surface depends on the initialmagnetic field strength and latitude of the flux tubes in theovershoot region, implying there is no preferred longitudein a globally fixed direction. Moreover, flux tubes that areperturbed at different latitudes in the convective overshotregion, show a wide latitudinal range of emergence on thestellar surface, with considerable asymmetries appearing ashighly peaked f s distributions or broad preferred longitudes.In a binary with P orb = 2 d, the authors find it takes severalmonths to years for a flux tube, perturbed from the convec-tive overshoot region, to emerge on the stellar surface. Overthis time, the Coriolis force acts on the internal gas flowand causes the poleward deflection of the tube (Schuessler& Solanki 1992), with the largest deflection for flux tubesstarting at lower latitudes, and those starting > ◦ latitudeshowing essentially no deflection. Such simulations are con-sistent with the large high-latitude spots found on AE Aqr.Clearly the behaviour of flux tubes in a binary system iscomplex. However, the results of the models by Holzwarth& Schüssler (2003b) may explain the varying size, distri-bution and evolution of the starspots imaged in AE Aqr;Namely, the low value of f s around ◦ latitude in the mapsof AE Aqr is consistent with the fact that some latitudesare avoided by the simulated erupting flux tubes. Further-more, the variable peak of f s at high latitudes in AE Aqrmay be a redistribution of magnetic energy, changing thefield strength of perturbed flux tubes and causing them toemerge at different latitudes, as well as shifting in longitude.In addition, flux tube eruption at high latitudes, due to flowinstabilities, leads to spots emerging over a broad longituderegion – similar to the large high-latitude spots observed inAE Aqr. It is unclear what the dominant dynamo mechanism is inAE Aqr. There is no clear evidence that we see a ‘flip-flop’ cycle, especially given the lack of active longitudes in threemaps of our sample. Indeed the most prominent evidence ofan activity cycle is the increase in f s at ◦ latitude overthe course of the 8 years of observations, combined with thepersistently high f s around ◦ latitude. The growth in spotcoverage around ◦ latitude may be indicative of an emerg-ing band of spots, forming part of an activity cycle similar tothat seen on the Sun. Furthermore, the increase in f s around ◦ latitude may be a second band of spots that form partof a previous cycle. In the case of the Sun, the latitude ofemergence of flux tubes gradually moves towards the equatorover the course of an activity cycle, taking ∼ years, withlittle overlap between consecutive cycles of flux tube emer-gence. However, simulations by Işık et al. (2011) show thatstronger dynamo excitation may cause a larger overlap be-tween consecutive cycles. Given that the high spot coveragein AE Aqr suggests a strong dynamo excitation, the presenceof two prominent bands of spots may be indicative of suchan overlap between cycles. However, if such a cycle were toexist, we would expect to see the higher-latitude peak movetowards lower latitudes over the course of our observations.Given that we do not clearly see this, any solar-like activitycycle must take place over a timescale longer than 8 years.Saar & Brandenburg (1999), in their study of a large sam-ple of stars in a range of systems (including single stars andbinaries), found several correlations between the duration ofthe magnetic activity cycle and the rotation period. Perti-nently, using the fit to all stars in their sample, we estimatethat AE Aqr would have a magnetic activity cycle lasting ∼ years. Furthermore, we estimate a longer cycle periodof ∼ years by using the fit to stars defined as ‘superactive’by Saar & Brandenburg (2001). If the correlations found forother systems are also true for AE Aqr, then we may haveobserved less than half of an activity cycle.This is the first time the number, size and distributionof starspots has been tracked in a CV secondary. While anyspecific interpretation of the long term behaviour of the im-aged starspots is somewhat challenging, the presence andevolution of two distinct bands of spots may indicate anongoing magnetic activity cycle in the secondary star inAE Aqr. Hence, it is crucial we continue our study of itsmagnetic activity. Future maps would allow us to track theevolution of the large high-latitude spots to determine iftheir long term behaviour is periodic. Moreover, by trackingthe evolution of the spot bands we may determine if theyform part of a periodic activity cycle, and if so, the length ofsuch a cycle could be measured, providing a unique insightinto the behaviour of the stellar dynamo in an interactingbinary. In addition, shorter, more intensive campaigns wouldallow the position of specific spot features to be tracked. Thiswould allow further measurements of differential rotation aswell as determining if meridional flows are present. The im-pact of tidal forces on magnetic flux tube emergence andpossible quenching of mass transfer could then be assessed.
10 CONCLUSIONS
We have imaged starspots on the secondary star in AE Aqrfor 7 epochs, spread over 8 years. This is the first time suchas study has been carried out for a secondary star in a CVand, in some cases, the number, size and distribution of spots © 2015 RAS, MNRAS , 1–19 C.A. Hill et al. varies significantly between maps. In particular, the chang-ing positions of the large high-latitude spots cannot be ex-plained by differential rotation, nor by the ‘flip-flop’ activitycycle. At lower latitudes, we see the emergence of a band ofspots around ◦ latitude, as well as a persistently high spotcoverage around ◦ latitude. These bands may form partof an activity cycle similar to that seen in the Sun, wheremagnetic flux tubes emerge at progressively lower latitudesthroughout a cycle. Furthermore, the complex distributionand behaviour of spots may be attributed to the impact ofthe companion WD on flux tube dynamics. ACKNOWLEDGMENTS
We thank Tom Marsh for the use of his molly software pack-age in this work, and VALD for the stellar line-lists used.We thank the staff at Carnegie Observatories for their assis-tance with the MIKE pipeline and for access to the HenriettaSwope Telescope. C.A.H. acknowledges the Queen’s Univer-sity Belfast Department of Education and Learning PhDscholarship, C.A.W. acknowledges support by STFC grantST/L000709/1, and D.S. acknowledges support by STFCgrant ST/L000733/1. This research has made use of NASA’sAstrophysics Data System and the Ureka software packageprovided by Space Telescope Science Institute and GeminiObservatory.
REFERENCES
Applegate, J. H., 1992, ApJ, 385, 621Baliunas, S. L., et al., 1995, ApJ, 438, 269Berdyugina, S. V., 2005, Living Reviews in Solar Physics,2, 8Berdyugina, S. V., Järvinen, S. P., 2005, Astron. Nachr.,326, 283Berdyugina, S. V., Tuominen, I., 1998, A&A, 336, L25Berdyugina, S. V., Berdyugin, A. V., Ilyin, I., Tuominen,I., 1999, Astron. Astrophys., 350, 626Berdyugina, S. V., Pelt, J., Tuominen, I., 2002, A&A, 394,505Bernstein, R., Shectman, S. A., Gunnels, S. M., Mochnacki,S., Athey, A. E., 2003, in Iye, M., Moorwood, A. F. M.,eds., Instrument Design and Performance for Opti-cal/Infrared Ground-based Telescopes, vol. 4841 of
So-ciety of Photo-Optical Instrumentation Engineers (SPIE)Conference Series , p. 1694Bianchini, A., 1990, AJ, 99, 1941Bruch, A., 1991, A&A, 251, 59Casares, J., Mouchet, M., Martinez-Pais, I. G., Harlaftis,E. T., 1996, MNRAS, 282, 182Chincarini, G., Walker, M. F., 1981, A&A, 104, 24Claret, A., 2000, VizieR Online Data Catalog, 336, 31081Collier-Cameron, A., Unruh, Y. C., 1994, MNRAS, 269,814Crawford, J. A., Kraft, R. P., 1956, ApJ, 123, 44Davey, S., Smith, R. C., 1992, MNRAS, 257, 476Dekker, H., D’Odorico, S., Kaufer, A., Delabre, B., Kot-zlowski, H., 2000, in Iye, M., Moorwood, A. F., eds., Opti-cal and IR Telescope Instrumentation and Detectors, vol. 4008 of
Society of Photo-Optical Instrumentation Engi-neers (SPIE) Conference Series , p. 534Dhillon, V. S., Watson, C. A., 2001, in Boffin, H. M. J.,Steeghs, D., Cuypers, J., eds., Astrotomography, IndirectImaging Methods in Observational Astronomy, vol. 573 of
Lecture Notes in Physics, Berlin Springer Verlag , p. 94Djurašević, G., Rovithis-Livaniou, H., Rovithis, P., Geor-giades, N., Erkapić, S., Pavlović, R., 2003, A&A, 402, 667Djurašević, G., Rovithis-Livaniou, H., Rovithis, P., Geor-giades, N., Erkapić, S., Pavlović, R., 2006, A&A, 445, 291Donati, J.-F., 1999, MNRAS, 302, 457Dunford, A., Watson, C. A., Smith, R. C., 2012, MNRAS,422, 3444Eaton, N., Draper, P. W., Allan, A., 2009, Starlink UserNote, 45Echevarría, J., Smith, R. C., Costero, R., Zharikov, S.,Michel, R., 2008, MNRAS, 387, 1563Hatzes, A. P., Vogt, S. S., 1992, MNRAS, 258, 387Heckert, P. A., Maloney, G. V., Stewart, M. C., Ordway,J. I., Hickman, A., Zeilik, M., 1998, AJ, 115, 1145Henden, A. A., Honeycutt, R. K., 1995, PASP, 107, 324Henry, G. W., Eaton, J. A., Hamer, J., Hall, D. S., 1995,Astrophys. J. Suppl. Ser., 97, 513Hessman, F. V., Gänsicke, B. T., Mattei, J. A., 2000, A&A,361, 952Hill, C. A., Watson, C. A., Shahbaz, T., Steeghs, D.,Dhillon, V. S., 2014, MNRAS, 444, 192Holzwarth, V., Schüssler, M., 2003a, A&A, 405, 291Holzwarth, V., Schüssler, M., 2003b, A&A, 405, 303Işık, E., Schmitt, D., Schüssler, M., 2011, A&A, 528, A135Järvinen, S. P., Berdyugina, S. V., Tuominen, I., Cutispoto,G., Bos, M., 2005, Astron. Astrophys., 432, 657Kelson, D., 2015, Website, http://code.obs.carnegiescience.edu/mike/index.html
King, A. R., Cannizzo, J. K., 1998, ApJ, 499, 348Kupka, F. G., Ryabchikova, T. A., Piskunov, N. E., Stem-pels, H. C., Weiss, W. W., 2000, Baltic Astronomy, 9, 590Livio, M., Pringle, J. E., 1994, ApJ, 427, 956Lockwood, G. W., Skiff, B. A., Henry, G. W., Henry, S.,Radick, R. R., Baliunas, S. L., Donahue, R. A., Soon, W.,2007, ApJS, 171, 260Messina, S., Guinan, E. F., 2002, Astron. Astrophys., 393,225Moss, D., Piskunov, N., Sokoloff, D., 2002, A&A, 396, 885Olah, K., Budding, E., Kim, H.-L., Etzel, P. B., 1994, A&A,291, 110Parsons, S. G., et al., 2015, A high precision study of thepost common envelope binary qs vir, in prep.Radick, R. R., Lockwood, G. W., Skiff, B. A., Baliunas,S. L., 1998, ApJS, 118, 239Richman, H. R., Applegate, J. H., Patterson, J., 1994, Pub-lications of the Astronomical Society of the Pacific, 106,1075Rodonò, M., Messina, S., Lanza, A. F., Cutispoto, G., Teri-aca, L., 2000, A&A, 358, 624Rutten, R. G. M., Dhillon, V. S., 1994, A&A, 288, 773Rutten, R. G. M., Dhillon, V. S., 1996, in A. Evans &J. H. Wood, ed., IAU Colloq. 158: Cataclysmic Variablesand Related Objects, vol. 208 of
Astrophysics and SpaceScience Library , p. 21Saar, S., Brandenburg, A., 2001, in Mathys, G., Solanki,S. K., Wickramasinghe, D. T., eds., Magnetic Fields © 2015 RAS, MNRAS , 1–19 he long-term magnetic activity of AE Aqr Across the Hertzsprung-Russell Diagram, vol. 248 of
As-tronomical Society of the Pacific Conference Series , p. 231Saar, S. H., Brandenburg, A., 1999, ApJ, 524, 295Schuessler, M., Solanki, S. K., 1992, A&A, 264, L13Schussler, M., Caligari, P., Ferriz-Mas, A., Moreno-Insertis,F., 1994, A&A, 281, L69Schwope, A. D., Staude, A., Vogel, J., Schwarz, R., 2004,Astronomische Nachrichten, 325, 197Shahbaz, T., Watson, C. A., Dhillon, V. S., 2014, MNRASSokoloff, D., Piskunov, N., 2002, MNRAS, 334, 925Strassmeier, K. G., Handler, G., Paunzen, E., Rauth, M.,1994, A&A, 281, 855Tanzi, E. G., Chincarini, G., Tarenghi, M., 1981, PASP,93, 68Vogt, S. S., Penrod, G. D., 1983, Publications of the As-tronomical Society of the Pacific, 95, 565Watson, C. A., Dhillon, V. S., 2001, MNRAS, 326, 67Watson, C. A., Dhillon, V. S., Rutten, R. G. M., Schwope,A. D., 2003, MNRAS, 341, 129Watson, C. A., Dhillon, V. S., Shahbaz, T., 2006, MNRAS,368, 637Watson, C. A., Steeghs, D., Shahbaz, T., Dhillon, V. S.,2007, MNRAS, 382, 1105Welsh, W. F., Horne, K., Gomer, R., 1995, MNRAS, 275,649Wilson, O. C., 1978, ApJ, 226, 379 © 2015 RAS, MNRAS000