Roche tomography of CVs - VI. Differential rotation of AE Aqr

Roche tomography of cataclysmic variables - VI. Differential rotation of AE Aqr - Not tidally locked!


We present Roche tomograms of the K4V secondary star in the cataclysmic variable AE Aqr, reconstructed from two datasets taken 9 days apart, and measure the differential rotation of the stellar surface. The tomograms show many large, cool starspots, including a large high-latitude spot and a prominent appendage down the trailing hemisphere. We find two distinct bands of spots around and latitude, and estimate a spot coverage of 15.4-17 per cent on the northern hemisphere.

Assuming a solar-like differential rotation law, the differential rotation of AE Aqr was measured using two different techniques. The first method yields an equator-pole lap time of 269 d and the second yields a lap time of 262 d. This shows the star is not fully tidally locked, as was previously assumed for CVs, but has a co-rotation latitude of . We discuss the implications that these observations have on stellar dynamo theory, as well as the impact that spot traversal across the L1 point may have on accretion rates in CVs as well as some of their other observed properties.

The entropy landscape technique was applied to determine the system parameters of AE Aqr. For the two independent datasets we find and , and , and orbital inclinations of to at optimal systemic velocities of and  kms-1.

cataclysmic variable, differential rotation, Roche tomography, magnetic activity, star spot, AE Aquarii, mass transfer variations

1 Introduction

Differential rotation of a solar-like star is thought to be a key component in stellar dynamo theory, converting and amplifying poloidal fields into toroidal fields. The magnetic activity observed on these stars is assumed to be based on a dynamo mechanism operating in the outer convection zone. The details of the generation and amplification mechanisms are currently not well understood due to the complex interactions involved and the lack of observational constraints. To test and better constrain stellar dynamo theory, we can compare models with observations of differential rotation (DR) on stars of different masses and rotation rates.

Previous measurements of DR on single stars include Donati and Cameron (1997) in their study of the rapidly rotating single K-type star AB Dor ( h). They found an equator-pole lap time of 110 d, close to solar value of 120 d. However, this displayed significant variability over several years, falling from 140 d to 70 d in the period 1988–1992, when the DR rate doubled. The Lupus post T Tauri star RX J1508.6±4423 ( h) was found by Donati et al. (2000) to have a lap time of d, and the late-type single star PZ Tel ( h) was found by Barnes et al. (2000) to have a lap time of 72–100 d. However, in the less rapidly rotating K2 dwarf star LQ Hya ( h), Kővári et al. (2005) found it to be much more rigidly rotating, with a lap time of d.

Work by Scharlemann (1982) suggested DR should be fairly suppressed in tidally locked systems such as CVs, and this was found in the tidally-locked pre-CV V471 Tau ( h) by Hussain et al. (2006). They found the surface shear consistent with solid body rotation, concluding that tidal locking may inhibit DR, but that this reduced shear does not affect the overall magnetic activity levels in this active K dwarf. Similarly, Petit et al. (2004) found weak DR on the RS CVn system HR 1099 ( h) with a lap time of d.

It is clear that different systems display a variety of DR rates, and so the study and comparison of stars of varying fundamental parameters is crucial for our understanding of the underlying dynamo mechanism. To this end, a measure of DR on cataclysmic variables would allow a test of dynamo theory in a parameter space with both rapid rotation and tidal distortion.

In addition, magnetic activity plays a critical role in the evolution of CVs, driving them to shorter orbital periods via magnetic braking (e.g. Kraft 1967; Mestel 1968; Spruit & Ritter 1983; Rappaport et al. 1983), with the transition of secondary to a fully convective state — supposedly shutting down magnetic activity — invoked to explain the period gap and dearth of CVs with 2-3 h periods. On shorter timescales, magnetic activity has been invoked to explain variations in CV orbital periods, mean brightnesses, mean outburst durations and outburst shapes (e.g. Bianchini 1990; Richman et al. 1994; Ak et al. 2001). In a study of the polar-type CV AM Her, Livio and Pringle (1994) suggested that mass-transfer variations were due to starspots traversing the L1 point. This suggestion is supported by the more recent work of Hessman et al. (2000) in their derivation of the mass transfer history of AM Her. They found that spot traversal across the L1 point is a valid mechanism for mass-transfer variation, if magnetic flux is preferentially produced around the L1 point or if there is a mechanism which forces spot groups appearing at higher latitudes to wander down towards the L1 point. Indeed, Holzwarth & Schüssler (2003) suggest tidal distortions may force starspots to form at preferred longitudes.

Clearly magnetic activity plays an important role in understanding the evolution and behaviour of CVs, and since differential rotation is a key component of such activity, a study of the secondary star in CVs allows a test of stellar dynamo theory and further insight into this class of system. In this paper we present intensity maps of the stellar surface on the secondary star in the CV AE Aqr using Roche tomography (see Watson & Dhillon (2001); Dhillon & Watson (2001); Watson et al. (2003) for details), and measure the stellar surface differential rotation rate.

2 Observations and reduction

Spectroscopic observations were carried out over 4 half-nights on 2009 August 27–28 and September 5–6. From here on we refer to nights 1 and 2 as block 1, and nights 3 and 4 as block 2. The data were acquired using the 8 m VLT, situated at Cerro Paranal in Chile. The spectroscopic observations of AE Aqr were carried out using the Ultraviolet and Visual Echelle Spectrograph (UVES). A thinned EEV CCD-44 chip with 2K x 4K pixels was used in the blue channel, and a mosaic of an EEV CCD-44 and a MIT/LL CCID-20 chip (both 2K x 4K) were used in the red channel. UVES was used in the Dichroic-1/Standard setting (390+580 nm) mode, allowing a wavelength coverage of 3259 Å–4563 Å in the blue arm and 4726 Å–6835 Å in the red arm. With a slit width of 0.9 arcsec, a spectral resolution of around 46,000 ( kms) was obtained in the blue, and a spectral resolution of around 43,000 ( kms) was obtained in the red, with the chip binning left as 1x1.

The blue and red spectra were taken quasi-simultaneously. During the first observing block, the blue and red spectra were taken using 228 s and 230 s exposures respectively, in order to minimise velocity smearing of the data due to the orbital motion of the secondary star. Each exposure corresponds to 0.65% of the orbit. Different exposure times were used in an attempt to compensate for the longer readout time of the blue CCD. The exposures were taken in blocks of 10, allowing red and blue data to remain quasi-simultaneous. For the second observing block, red and blue spectra were taken using 230 s exposures. For both observing blocks comparison ThAr lamp exposures were taken at the start and end of the night for the purpose of wavelength calibration. In the first block with this setup we obtained 118 useable spectra. 12 red spectra were lost due to a temporary shutter malfunction. In addition, spectral type templates were obtained for the purpose of determining the orbital phasing of AE Aqr via cross-correlation of absorption lines.

The seeing was around 2 arcsec at the start of night one, improving to 1.2 arcsec by mid observation, with seeing on the second night starting out around 4 arcsec, greatly improving to 1–0.9 arcsec by the middle of the night before worsening towards the end. In the second block (9 days later) we obtained 121 useable spectra in both the blue and red, giving a 100 per cent orbit coverage of AE Aqr. The seeing was 1.2–0.8 arcsec on the first night and 0.4–0.8 arcsec the second night. Table 1 gives a journal of the observations.

Object UT Date UT Start UT End Texp(s) No. spectra Comments
HD 187760 2009 Aug 27 23:15 23:19 240 1 K4V template
HD 187760 23:22 23:30 480 1
GJ 1247 23:57 00:03 360 1 K3V template
GL 653 23:40 23:44 240 1 K5V template
AE Aqr 23:49 05:11 230 57 Target spectra
GD50 2009 Aug 28 10:10 10:25 900 1 Flux standard
GJ 1192 23:15 23:22 360 1 K3V template
GL 653 23:42 23:46 240 1 K5V template
GJ 1247 23:57 00:03 360 1 K3V template
AE Aqr 2009 Aug 29 00:08 04:55 230 61 Target spectra
HR 7596 2009 Sept 05 23:31 23:31 4 1 Flux standard
HR 7596 23:34 23:34 15 1
HR 7596 23:36 23:37 30 1
AE Aqr 23:43 04:38 230 61 Target spectra
AE Aqr 2009 Sept 07 00:34 05:20 230 60 Target spectra
Table 1: Log of the VLT spectroscopic observations of AE Aqr, the spectral-type template stars and flux standards. The first column gives the object name, columns 2–4 list the UT Date and the exposure start and end times, respectively. Columns 5–6 list the exposure times and number of spectra taken for each object and the final column indicates the type of science frame taken.

2.1 Data reduction

The raw data were reduced using the ESO UVES pipeline. This automatically processes all calibration frames and then conducts bias subtraction and flat fielding, wavelength calibration, sky subtraction, flux calibration and optimal extraction of the target frames. The bias frames used were those closest in time to the science data. The final output consists of 1-d spectra for the blue and red arms.

It should be noted that Roche tomography cannot be performed on data which has not been corrected for slit losses. This is because the variable contribution of the secondary star to the total light of a CV forces one to use relative line fluxes during the mapping process and prohibits the usual method of normalising the spectra. We discuss our corrections for slit losses in section 5.1.

Over both observing blocks the peak signal-to-noise of the blue spectra ranged from 5–31 (typically ), and from 17–64 (typically ) in the red.

3 Ephemeris, Radial velocity curves and Continuum fitting

3.1 Ephemeris and Radial velocity curve

In order to construct an artefact-free Roche tomogram the orbital ephemeris is required. A new ephemeris for AE Aqr was determined by cross-correlation with a template star of spectral type K4V. We only considered the spectral region 6000–6500 Å for simplicity, as it contains strong absorption lines from the secondary star and reduces the probability of introducing a sloping continuum contribution from the blue accretion regions. The AE Aqr and K4V template spectra were first normalised by dividing by a constant, and the continuum was then subtracted off using a third order polynomial fit, thus preserving line strength in the spectral region of interest.

In order to perform the cross-correlation, the template spectrum was artificially broadened by an arbitrary amount to account for the rotational velocity () of the secondary star. The AE Aqr spectra were then cross-correlated with the artificially broadened template giving an initial estimate of the radial velocity of the secondary star at that orbital phase. In order to derive an improved estimate of the rotational broadening of the secondary star, the orbital motion was subtracted. These spectra were then averaged to provide one high signal-to-noise orbitally-corrected spectrum. A new value for the rotational broadening of AE Aqr was then obtained by artificially broadening the template spectrum in 0.1 kms-1 steps and optimally subtracting the broadened template from the orbitally corrected AE Aqr spectrum. A new value was then obtained by adopting the broadening that needed to be applied to the template spectrum in order to minimise the residuals in the optimal subtraction. This new broadening value was then applied to the template and the whole process repeated until the rotational broadening value no longer changed. This typically took 3 iterations. The whole procedure was carried out separately for block 1 and block 2 data.

Through the above process a cross-correlation function (CCF) was calculated for each AE Aqr spectrum. The radial velocity curve in figure 1 was then obtained by fitting a sinusoid through the CCF peaks using the equation


where is the radial velocity of the secondary at phase , the systemic velocity and the semi-amplitude of the secondary star.

From this analysis a new zero-point was obtained for the ephemeris of


with the orbital period fixed at (from Casares et al. 1996).

All data used in this work has been phased with respect to this new ephemeris. While the cross-correlation method is relatively insensitive to the use of an ill-matching template or incorrect amount of broadening, the effects of irradiation are more likely to dominate measurements, introducing systematic errors in radial velocity if not accounted for (e.g. Davey & Smith 1992). The ephemeris derived here provides a substantially improved image quality (both in terms of final reduced and artefacts in the reconstructed map) when compared against that published by Watson et al. (2006).

For the purposes of this work we have not attempted a full spectral-type and binary determination, however, for completeness we obtained a secondary star radial velocity amplitude of , a systemic velocity of and a mean from the 4 nights observations of kms-1. This assumes the K4V template star HD 187760 has a systemic velocity of , as measured by a Gaussian fit to the LSD line profile of the K4V template star (using a line-list, where lines with a central depth shallower than 10 per cent of the continuum were excluded).

Figure 1 shows the measured radial velocity for each data block as well as the fitted sinusoid. The velocities are noticeably lower around phase 0.25 due to a systemic velocity offset between blocks 1 and 2 (see section 6.2). This conventional approach to binary parameter measurements are inherently systematically biased due to surface inhomogeneities (such as irradiation and starspots) and tidal distortion, which causes the centre-of-light (COL) to no longer coincide with the centre-of-mass (COM). This inherent bias is clearly shown in the deviations from a perfect sinusoid in figure 1, and by examining the residuals after subtracting the fitted sinusoid (see figure 2), surface features can be identified. The peak around phase 0.4 and the trough around phase 0.6 correspond to a region of reduced absorption-line strength on the inner hemisphere of the secondary, and is most likely caused by irradiation from the accretion regions in addition to starspots, both moving the COL away from the COM. The large peak around phase 1 and the trough around phase 1.3 indicate a region of enhanced absorption-line strength on the outer hemisphere, again moving the COL away from the COM, and could be caused by very large starspots. The scatter in RVs between block 1 and block 2 may indicate to spot evolution or differential rotation, and is an independent indication of surface features from Roche tomography (see section 7). Due to the bias in these measurements, the parameters derived from the radial velocity curve have not been used in the subsequent analysis presented in this work.

Figure 1: Radial velocity curve of AE Aqr, which has been phase-folded for clarity. The solid line is a least-squares sinusoid fit to the RV points, assuming a circular orbit. Triangles are block 1 data, circles are block 2 data. The error bars are smaller than the points in this figure.
Figure 2: Residuals (phase-folded and repeated) after subtraction of the fitted sinusoid to the RV curve in figure 1. Red triangles are block 1 data, blue circles are block 2 data taken 9 days later. RV residuals can be traced to individual features in the Roche tomograms presented later, and differences in the residuals between blocks 1 and 2 are attributed to both spot evolution and differential rotation.

3.2 Continuum fitting

In order to generate the least-squares deconvolution line profiles required for Roche tomography, the continuum must be flattened. In order to achieve this we fit the continuum. For this and all subsequent analysis, emission lines, tellurics and regions considered too noisy were masked in the spectra (note emission lines are not required for Roche tomography since we are mapping the secondary star, not the disk). The masked spectral regions are shown in table 2.

Colour Masked region (Å) Comments
Blue Noise, emission lines
H delta emission
H gamma emission
No flux
Red lower No flux
H beta emission
He I emission
No flux
Red upper No flux,
He I emission
Na I doublet
H alpha emission
He I emission
No flux
Table 2: Spectral regions not included in continuum fitting and LSD process.

The contribution of the secondary star to the total system luminosity is constantly varying, and there is an unknown contribution to each spectrum due to variable light from the accretion regions. A master continuum fit to the data (such as that done by Collier-Cameron & Unruh 1994) is not appropriate due to the constantly changing continuum slope, caused by (for example) flaring, or from the varying aspect of the accretion region. Additionally, normalisation of the continuum would result in the photospheric absorption lines from the secondary varying in relative strength between exposures. This forces us to subtract the continuum from each spectrum.

Slit loss corrections could be undertaken by placing another star in the slit, but unfortunately this is not possible here due to the cross-dispersed nature of the echelle spectrograph. For this reason, when conducting Roche tomography it is normal to obtain photometry to monitor transparency and target brightness variations. Unfortunately for these observations we were unable to obtain simultaneous photometry, and so another means of slit correction needed to be applied. For this we employed a method of relative scaling of adjacent LSD profiles, which is discussed in detail in section 5.1.

Since we cannot construct a master continuum fit, an algorithm was developed to place spline points at regular intervals in the spectra. During this process, emission lines from the accretion disk and tellurics were masked out (see table 2 for details of excluded regions).

It is difficult to determine where the true continuum lies due to the forest of broadened and blended lines which act to lower the apparent continuum level. We therefore adopted an algorithm which tries to identify regions of line-free continuum. In order to do this we adopted the following strategy.

  • The spectra were first split into discrete wavelength windows spanning 40 Å .

  • A first order polynomial was fit to each window individually and all the data points below the fit were rejected.

  • Another first order polynomial was fit to the remaining data points in each separate window.

  • The highest flux value within of this latter fit was determined, and a constant was added to the fit so that it now passed through this flux value.

  • A spline point was then positioned on the linear fit at the central wavelength of the window in question.

We found that deep absorption lines could systematically lower the fit since a large number of data points would lie below the true continuum level, fooling the above procedure into placing a spline point too low. To correct for this, spline points placed in deep lines were identified and rejected by comparing it to the points on either side and to the previously included spline point. Spline points were placed on either side of masked regions, and the spline fit was allowed to vary freely in this region.

The results of this fitting process were tested and checked visually until an acceptable fit was acquired. This is obviously very subjective as it is difficult to judge where the true continuum actually lies when dealing with heavily broadened and blended lines, which effectively acts to suppress the apparent continuum level. The normalised spectra were compared visually in sequential order, and the continuum fit was found to be stable.

The impact of a poorly-fitted continuum was previously assessed by Watson et al. (2006), who found that even continuum fits which were obviously incorrect did not adversely affect the shape of the line profiles after carrying out the least squares deconvolution (LSD) process described in section 5. However, we found that the continuum was systematically fit at too high a level, leading to LSD profiles with continuum regions lying significantly below zero (as also found by Watson et al. 2006). This was solved by shifting the continuum fit to a lower level until the LSD profiles lay at zero. The line shape was not affected by this process.

4 Roche tomography

Roche tomography is analogous to Doppler imaging (e.g. Vogt & Penrod 1983) and has been successfully applied to the donor stars of CVs over the past 20 years (Rutten & Dhillon 1994, Rutten & Dhillon 1996, Watson et al. 2003, Schwope et al. 2004, Watson et al. 2006, Watson et al. 2007, Dunford et al. 2012).

Roche tomography is specifically designed to map Roche-lobe-filling secondary stars in interacting binaries such as CVs and X-ray binaries (e.g. Shahbaz et al. (2014), and assumes the secondary is locked in synchronous rotation and has a circularised orbit. Rather than repeat a detailed description of the methodology and axioms of Roche tomography here, we refer the reader to the references above and the technical reviews of Roche tomography by Watson & Dhillon (2001) and Dhillon & Watson (2001).

Relevant points of note with respect to this work are that we select a unique map by employing the maximum-entropy MEMSYS algorithm developed by Skilling & Bryan (1984). We use a moving uniform default map, where each element is set to the mean value of the reconstructed map. We do not adopt a two-temperature or filling-factor model (e.g. Collier-Cameron & Unruh 1994), since donor stars are expected to exhibit large temperature differences due to irradiation by the primary. The lack of a two-temperature model means our Roche tomograms may be prone to the growth of bright pixels which make the tomograms quantitatively more difficult to analyse.

5 Least squares deconvolution

Spot features appear in line profiles as an emission bump (actually a lack of absorption), and are typically a few percent of the line depth. This means that very high signal-to-noise data are required, which cannot directly be achieved for a single line with AE Aqr due to its faintness and the requirements for short exposures to avoid orbital smearing. Thus, the technique of Least Squares Deconvolution (LSD) was employed which effectively stacks the ’s of photospheric absorption lines observable in a single echelle spectrum to produce a ‘mean’ profile of greatly increased signal-to-noise. Theoretically, the multiplex gain in signal-to-noise is the square root of the number of lines observed. The technique is well documented and was first applied by Donati et al. (1997) and has since been used in many Doppler imaging studies (e.g. Jeffers et al. 2002; Barnes et al. 2004; Marsden et al. 2005) as well as in the mapping of starspots on CVs (e.g. Watson et al. 2007; Dunford et al. 2012). For more details on LSD, see these references as well as the review by Cameron (2001).

The technique assumes all photospheric lines have the same local line profile shape, and that starspots affect all of the rotationally broadened line profiles in the same way, so the morphology of the characteristic starspot bump is identical. The positions and relative strengths of observed lines in each echelle spectrum must be known. For this work, a line list generated by the Vienna Atomic Line Database (VALD) was used (see Kupka & Ryabchikova 1999; Kupka et al. 2000). The spectral type of AE Aqr has been determined to lie in the range K3–K5 V (Crawford & Kraft 1956; Chincarini & Walker 1981; Tanzi et al. 1981; Bruch 1991), and so a line-list for a stellar atmosphere with and (the closest approximation available to a K4V spectral type) was obtained and used in the LSD process.

Initially the LSD line profiles exhibited a slope in the continuum. After trials it was found that including more lines in the deconvolution acted to reduce the slope in the continuum region of the deconvolved line profiles. A detection limit of 0.2 (of the normalised line depth, below which all lines with a smaller central depth were excluded) was adopted, giving 3202 lines over which to carry out LSD. Including more lines resulted in minimal improvement but significantly increased computation time. Since the line-list obtained from VALD contains normalised line-depths, where as Roche tomography uses continuum subtracted spectra, each line-depth was scaled by a fit to the continuum of a K4V template star, meaning each line’s relative depth was correct for use with continuum subtracted data.

The exact choice of line-list is unlikely to affect the results presented here as Barnes (1999) found that the LSD process was insensitive to the use of an incorrect line-list (see Watson et al. 2006 for more detail).

After this process, there was still evidence of a continuum slope in the LSD profiles. This was removed by masking out the centre of the line and subtracting a second-order polynomial which was fit to the continuum. This smooth function did not significantly alter the profiles’ shape.

5.1 Scaling of profile depths

As previously mentioned, Roche tomography cannot be performed on data which has not been corrected for slit losses. Normally, when conducting Roche tomography, simultaneous photometry is taken to monitor transparency and target brightness variations, and due to its absence here, the LSD line profiles were scaled relative to each other using an adjacent-profile scaling approach. For this we assumed that the line profile does not change significantly from profile to profile, apart from small changes due to rotational effects. We also assume that any large changes are due to slit losses and transparency variations, rather than changing spot features.

For the line profile scaling we employed an optimal subtraction method, where one profile is scaled and subtracted from its neighbouring profile, the optimal scaling factor being that which minimises the residuals. This scaling factor was found between a benchmark profile (chosen for its near-Gaussian shape) and its neighbour. The factor was applied to the neighbouring LSD profile and the process repeated for the entire time series of LSD profiles, using the newly scaled profile as the comparison for the next. It was found that applying a slight Gaussian smoothing to the line profile allowed for better scaling as this reduced the effect of noise and of a changing line profile shape. The reader should be reminded that Roche tomography uses relative fluxes, and so only the shape, not the absolute depth of the line profiles matter.

The above scaling technique was extensively tested in simulations. For this we constructed a synthetic dataset which resembled the final AE Aqr tomogram. We took each line profile of our scaled dataset (from above) and scaled them by varying amounts to simulate slit losses and transparency variations. We then applied the adjacent-profile scaling method and compared our calculated scaling factors to those that were injected at the start of the simulation.

The initial adjacent-profile scaling was typically correct to within 4 per cent for both the simulated block 1 and block 2 datasets, with a particular worst case scenario, discussed below, that leads to a 16 per cent error. This occurs around phase 0.5 and is due to rapidly changing line depths caused by the irradiation around the L1 point. These poorly scaled profiles showed up in the Roche tomograms as prominent artefacts in the the form of two large dark spots on the southern hemisphere. The systematic errors in this method were corrected as follows

  • A map was reconstructed from the data.

  • Random scaling factors were applied to the computed fits to the line profiles.

  • The adjacent-profile scaling method above was applied to these line profiles (starting from a benchmark profile).

  • The scaling factors between the computed profiles and those from the previous step were found.

  • These scaling factors were divided by the random scaling factor assigned to the benchmark profile, giving the percentage to which they were correct.

  • Then, the scaling factors between the observed line profiles and the reconstructed line profiles were found.

  • These latter scaling factors were divided by the percentage above, giving the correction factors that need to be applied after using the adjacent-profile scaling method.

The systematic corrections were found to be typically 5 per cent, and 13 per cent in the worst case. The resulting profiles were used in the reconstruction process in Roche tomography. The computed line profiles were compared with the input profiles, and we found profiles around phase 0.5 were poorly fitted. We then scaled the input profiles to match the computed profiles. We note that, even including the poorly scaled profiles, the spot features in the northern hemisphere are not significantly altered.

The final LSD profiles are trailed in figure 3, showing the emission bumps due to starspots moving from blue (negative velocities) to red (positive velocities) through the profile as AE Aqr rotates, in addition to the secondary stars orbital motion and variations in the projected equatorial velocity, .

Figure 3: Trailed spectra of AE Aqr (top panel is block 1 data, bottom panel is block 2 data). Panels show (from left to right) the observed LSD data, computed data from the Roche tomography reconstruction and the residuals (increased by a factor of 10). Starspot features in these panels appear bright. A grey-scale wedge is also shown, where a value of 1 corresponds to the maximum line depth in the reconstructed profiles. The orbital motion has been removed assuming the binary parameters found in section 6, which allows the individual starspot tracks across the profiles and the variation in to be more clearly observed.

The residuals in figure 3 show narrow emission features in the form of sinusoids across the entire orbitally-corrected trail for both observation blocks. These appear as vertical lines when the orbital motion has not been removed, and the most prominent of these features can be visually tracked through sequential line profiles. These emission bumps are also seen lying outside the stellar absorption profile, off the stellar limb, and therefore cannot lie on the stellar surface. We suggest that the most likely cause for these features is the presence of circumstellar material in the form of prominences. A similar feature was also observed in the LSD line profiles of BV Cen, with the conclusion that this was due to a slingshot prominence (see Watson et al. 2007 for further discussion). The emission features seen here, however, are too faint to carry out any further conclusive analysis.

6 System parameters

Adoption of incorrect system parameters such as systemic velocity, component masses and orbital inclination when carrying out a Roche tomography reconstruction results in spurious artefacts in the final map (see Watson & Dhillon 2001 for details). Artefacts normally appear as bright and dark streaks, increasing the amount of structure mapped in the final image, which in turn leads to a decrease in the entropy regularisation statistic.

We can constrain the binary parameters by carrying out reconstructions for many pairs of component masses (iterating to the same ), creating an ‘entropy landscape’ (see figures 9 & 10) where each grid tile corresponds to the map entropy value for the reconstruction for that pair of component masses. Entropy landscapes can also be repeated for different values for the orbital inclination and the systemic velocity . One then adopts the set of parameters () that produce the map containing the least structure (the map of maximum entropy).

6.1 Limb Darkening

After trial reconstructions, it was found that the square-root limb darkening law as used in previous work (see Watson et al. 2006) resulted in a poor fit to the wings of the observed LSD profiles. The mismatch to the wings is evident in this dataset due to the superior signal-to-noise compared to previous Roche tomography datasets.

For this work we adopted a four-parameter non-linear model (see Claret 2000), given by


where ( is the angle between the line of sight and the emergent flux), and is the monochromatic specific intensity at the centre of the stellar disk.

In order to calculate the correct limb darkening coefficients we determined the effective central wavelength of the UVES data using,


where is the line depth at wavelength , and the error in the data at .

Since there are effectively three wavelength ranges in the UVES setup (the blue CCD and the two halves of the red CCD, red-upper and red-lower) we calculated for all the exposures separately, for each CCD. We then averaged these for each CCD, and averaged the resulting three values, giving the final central wavelengths as  Å  and  Å  for block 1 and block 2, respectively.

Once was found, the limb darkening coefficients were then determined by linearly interpolating between the tabulated wavelengths given in Claret (2000). We adopted stellar parameters closest to that of a K4V star, which for the PHOENIX model coefficients were and . We also trialled limb darkening coefficients from the ATLAS model, but found that maps reconstructed using the PHOENIX coefficients resulted in reconstructions to higher entropy values (by around 10%). We therefore adopted the PHOENIX coefficients for the rest of the work, which are listed in table 3 for both the separate observing blocks.

Coefficient Block 1 Block 2
a1 0.7253 0.7254
a2 -0.7380 -0.7348
a3 1.3314 1.3288
a4 -0.4264 -0.4281
Table 3: Limb-darkening coefficients used in equation (3).

6.2 Systemic velocity, inclination and masses

The binary parameters were found independently for each observation block. A sequence of entropy landscapes were constructed for a range of orbital inclinations i and systemic velocities . For each and we chose the map of maximum entropy in the corresponding entropy landscape. We present the results of the entropy landscape analyses below.

Systemic velocity

Figures 4 & 5 show the systemic velocities yielded; kms-1 and kms-1, where the error bars are estimated from the spread of velocities which produce maps of similar entropy (the “entropy plateau” seen in figures 4 & 5). These measurements are consistent with that of previous work (see table 4). The difference between the two blocks of 1.8 kms-1 can be explained by an instrumental offset, since offsets of ms-1 are frequently reported for wavelength stabilised high-precision echelle spectra such as SOPHIE (see Simpson et al. 2011).

The value obtained using Roche Tomography is similar to that obtained from the radial velocity curve of , although the RV curve measurement will be biased (see section 3.1). The value obtained for the systemic velocity from the entropy landscapes is independent of the assumed inclination, as has been found previously by Watson et al. (2003, 2006).

Author Systemic velocity (kms-1) Inclination (degrees) (M) (M) Mass ratio
Block 1 (this work) 50 1.20 0.81 0.68
Block 2 (this work) 51 1.17 0.78 0.67
Echevarría et al. (2008) 70 0.63 0.37 0.60
Watson et al. (2006) 66 0.74 0.50 0.68
Casares et al. (1996) 0.63
Welsh et al. (1995) 0.64
Table 4: System parameters as found by the respective authors.
Figure 4: Points show the maximum entropy value obtained in each entropy landscape for block 1 data as a function of systemic velocity, assuming an orbital inclination of . The solid line shows the general trend, with a maximum of kms-1.
Figure 5: As figure 4, but for block 2, assuming an orbital inclination of . The solid curve shows the general trend, with a maximum of kms-1.


Figures 6 & 7 show the maximum entropy value obtained as a function of inclination from the entropy landscape analyses, assuming systemic velocities of kms-1 and kms-1 as derived above for block 1 and block 2, respectively. We obtain inclinations of for block 1 and for block 2.

The inclination measurements found here are lower than previous estimates, but still in agreement with those found by Welsh et al. (1995) and Casares et al. (1996) (see table 4).

Work by Echevarría et al. (2008) suggests a higher inclination up to but, perhaps more importantly, a previous Roche tomography study of AE Aqr by Watson et al. (2006) yielded a higher inclination of , despite using the same methods and approaches as this work.

It should be noted, however, that the inclination is the worst constrained parameter when using Roche tomography, but it is reassuring that the two independent data sets agree to within on the inclination of to , and that this lies below the limit inferred from the lack of eclipses (Chanan et al. 1976). In addition, the reconstructed maps (see section 7) at both and and at the previously determined value of show very similar features, and the calculated differential rotation rate (see section 8) is similar. This shows that the surface features mapped by Roche tomography, as well as the inferred differential rotation rate calculated later, are robust even when the orbital inclination is varied over the full span of previously reported values.

Figure 6: Points show the maximum entropy value obtained in each entropy landscape as a function of inclination for block 1 data, assuming a systemic velocity of kms-1. The solid curve shows a sixth order polynomial fit through these points, as a guild only to highlight the maximum of .
Figure 7: The same as figure 6, but for block 2 data, assuming a systemic velocity of kms-1. The solid curve shows a sixth order polynomial fit through these points, highlighting the maximum of .


The entropy landscapes for AE Aqr are shown in figures 9 & 10 for block 1 and block 2, respectively. For block 1 we assume and kms-1, and from this we derive a secondary mass of 0.81 M and a primary mass of 1.20 M. For block 2 we assume and kms-1, and from this we derive a secondary mass of 0.78 M and a primary mass of 1.17 M. These masses are consistent within this work, but are much larger when compared with previous estimates (see table 4), due to the inclination difference. It should be noted, however, that the mass ratios derived here are in good agreement with previous work.

In all reconstructions we fit to a reduced , which indicates that our propagated error bars are systematically overestimated for the scaled profiles, but this does not affect the reconstructions. This aim was chosen as the entropy of the reconstructed maps dramatically decreases when fits to lower than this are performed (see figure 8). This is due to a marked increase in the presence of small-scale structure due to the mapping of noise features in the Roche tomogram. On the other hand, fitting to a higher reduced causes the reconstructed maps to have less distinct structure, and the inclination and masses are not constrained as more pixels are assigned the default map value.

Figure 8: Reconstructed-map entropy as a function of reduced , for both observation blocks, using the system parameters derived in section 6.2. Note the sharp decrease in entropy below a reduced of 0.4.

Assigning errors to the derived system parameters obtained from the entropy landscapes is not straight forward. As discussed in Watson et al. (2006) and Watson & Dhillon (2001), it would require using a Monte-Carlo style technique combined with a bootstrap resampling method (Efron 1979; Efron & Tibshirani 1993) to generate synthetic datasets drawn from the same parent population as the observed dataset. Then the same analysis carried out in this work would need to be applied to the hundreds of bootstrapped datasets, requiring several months of computation. This is unfeasible and so we do not assign strict error bars to our derived system parameters.

Figure 9: The entropy landscape for AE Aqr using data from the first observation block, assuming an orbital inclination of and a systemic velocity of kms-1. Dark regions indicate masses for which no acceptable solution could be found. The cross marks the point of maximum entropy, corresponding to component masses of and .
Figure 10: The same as figure 9, but this time showing the entropy landscape for AE Aqr using data from the second observation block, assuming an orbital inclination of and a systemic velocity of kms-1. The point of maximum entropy corresponds to component masses of and .

7 Surface maps

7.1 Global features

Using the system parameters derived in section 6 we have constructed Roche tomograms of the secondary star in AE Aqr for both observation block 1 and block 2 of our dataset (see figures 11 & 12). The corresponding fit to the datasets obtained on both nights are displayed in figure 3.

Figure 11: The Roche tomogram of AE Aqr using block 1 data. Dark grey scales indicate regions of reduced absorption line strength that is due to either the presence of starspots or the impact of irradiation. The orbital phase is indicated above each panel. Roche tomograms are shown without limb darkening for clarity.
Figure 12: The same as figure 11 but displaying the Roche tomogram of AE Aqr reconstructed using the data obtained on observing block 2 (9 days after observing block 1).

Dark spots features are evident in both tomograms. Common to both maps is the large, high-latitude starspot located above latitude. Similar high-latitude spots are commonly found in Doppler imaging studies of rapidly-rotating single stars such as AB Dor (Hussain et al. 2007), LQ Hya (Donati 1999) and PZ Tel (Barnes et al. 2000), and have also been found in Roche tomograms of other CVs including RU Peg (Dunford et al. 2012), BV Cen (Watson et al. 2007) and a previous study of AE Aqr (Watson et al. 2006). An explanation for these high-latitude spots is given by Schuessler & Solanki (1992) who suggest that strong Coriolis forces in these rapidly rotating stars drive the magnetic flux tubes towards the polar regions. Alternatively, spots may migrate polewards after forming at lower latitudes, as has been found for the RS CVn binary HR1009 (Vogt et al. 1999).

We note that in all previous Roche tomograms of CV secondaries (RU Peg, BV Cen and AE Aqr) it appeared that the dominant high latitude spot feature was located towards the trailing hemisphere — raising the prospect that these large spots form at preferential longitudes. This does not appear to be the case for the high-latitude spot in either of the Roche tomograms of AE Aqr presented here.

There are two more prominent spot features in the Roche tomograms, particularly evident in figure 11. The first is the large ‘appendage’ extending from near the high-latitude spot at a latitude of down to on the trailing hemisphere (labelled ‘A’ in figure 11). We note that this spot feature shows a distinct morphological change between the two observing runs. The second feature is the large spot near the rear of the star (labelled ‘B’ in figure 11).

Finally, a ‘chain’ of starspots is evident that lead down from latitude towards the L1 point, and is offset towards the leading hemisphere by around in longitude. A similar chain of spots down to the L1 point was seen on previous Roche tomograms of BV Cen and AE Aqr (Watson et al. 2007, 2006), and a higher spot coverage was seen on the hemisphere facing the white dwarf in the Doppler image of the pre-CV V471 Tau (Hussain et al. 2006). This ‘chain’ of spots suggests a mechanism which forces magnetic flux tubes to preferentially arise at these locations. This may be due to tidal forces, which are thought to be able to force spots to form at preferred longitudes (Holzwarth & Schüssler 2003), and may also be due to the tidal enhancement of the dynamo action itself (Moss et al. 2002). In addition, surface flows may drag emergent magnetic flux tubes towards the L1 point (discussed later in section 9.2).

In order to make a more quantitative estimate of the spot parameters on AE Aqr, we have examined the pixel intensity distribution in the Roche tomograms. We discarded all pixels in the southern hemisphere as, since this hemisphere is least visible, pixels are mainly assigned the default intensity, which in this case is the average intensity of the map. Histograms of the pixel values are shown in figure 13 where we have assigned the brightest pixel in the map an intensity of 100 and linearly scaled the other pixel intensities relative to this.

Unlike the previous Roche tomogram of AE Aqr by Watson et al. (2006), we do not see a clear bimodal distribution in pixel intensities from which we can confidently distinguish between immaculate photosphere and spotted photosphere. Instead, the histograms of pixel intensities show broad peaks, with long tails towards high and low pixel intensities. Pixels with a scaled intensity of around 50-64 may be explained by the smearing of spots across latitudes, increasing areal coverage and reducing contrast, and by unresolved spots which cause a decrease in pixel intensity without being clearly reconstructed in the tomograms. Above the peak in the histograms we see a tail of high intensity pixels. It is unlikely that these represent the immaculate photosphere, as the growth of bright pixels in maps which are not ‘threshholded’ is a known artefact of Doppler imaging techniques (e.g. Hatzes & Vogt 1992 – also see section 4). Instead, we have defined an intensity of 81 and above (the upper end of the peak of pixel intensities) to represent the immaculate photosphere.

We have defined a spot by taking the pixel intensity at the centre of the large, high-latitude spot to represent the intensity of a 100 per cent spotted region (in the same way as Watson et al. 2007). Pixels with a lower intensity are present in the Roche tomogram, but these are confined to the irradiated region near the L1 point. These are significantly lower in intensity, confirming our interpretation that this dark feature is not a spot.

Figure 13: Histograms of the pixel intensities from the Roche tomograms of block 1 and block 2 datasets of AE Aqr after pixels on the southern hemisphere (latitude ) were discarded. The brightest pixel in each map is assigned an intensity of 100 and other pixel intensities are scaled linearly against this. For the purposes of estimating the global spot properties of AE Aqr, the immaculate photosphere has been defined as pixels with intensities of 81 or greater, and the spotted photosphere those regions where pixel intensities are less than 81.

If we assume a spot simply alters the continuum level and not line depths (a secondary influence on the LSD profiles), we can assume a blackbody scaling. This gives a temperature contrast between the photosphere and spot of K for block 1 and K for block 2. This is similar to K found for BV Cen (see Watson et al. 2007) and fits well with models by Frasca et al. (2005) who found a temperature differences of K, and Biazzo et al. (2006) who found K for stars in various locations in the HR diagram. The small disagreement in the spot temperatures calculated for the two maps (and also the overall spot filling factors determined later) are most likely due to slight differences in the quality of the maps achieved. This is largely due to the technical fault which affected the block 1 data, which means that overall, small scale features are less readily reconstructed.

Each pixel in our Roche tomograms was given a spot-filling factor between 0 (immaculate) to 1 (totally spotted) depending on its intensity between our predefined immaculate and spotted photosphere intensities. After removal of the region near the L1 point (which is irradiated and would cause us to overestimate the total spot coverage), we estimate that 15.4 per cent and 17 per cent of the northern hemisphere of AE Aqr is spotted, for block 1 and block 2 datasets, respectively. These values should be taken with caution since it depends on our classification of a spot, and is likely to be a lower limit due to the presence of unresolved spots, but they match a previous estimate of 18 per cent spot coverage by Watson et al. (2006), and is similar to the 22 per cent spot filling factor found by Webb et al. (2002) in a TiO study of the CV SS Cyg.

7.2 Spot coverage as a function of longitude and latitude

In the analysis which follows, the map coordinates are defined such that longitude is the centre of the back of the star, with increasing longitude towards the leading hemisphere, and with the L1 point at .

When plotted as a function of longitude (see figure 14) we find a large fractional spot coverage around longitude for both maps. The largest fractional coverage is between longitude. This longitude range contains the irradiated region around the L1 point, but even without including this, the trail of spots down to the L1 point have a large fractional spot coverage.

The two prominent features mentioned earlier, a large appendage labelled ‘A’ and a prominent spot labelled ‘’B’, are readily visible at longitudes of and , respectively, in figure 14. These become longitudinally blended in block 2 data due to differential rotation between maps.

Figure 14: Spot coverage on the northern hemisphere as a function of longitude for block 1 and block 2, normalised by the number of pixels within a 10 degree longitude bin. The region enclosed within longitude and latitude has been masked off for the L1 point.

As well as preferred longitudes, spots appear to form at preferred latitudes, as shown in figure 15. There are two distinct bands around and latitude with a larger fractional spot coverage, in addition to the large spot coverage at latitudes above . In the Roche tomogram of AE Aqr by Watson et al. (2006) from observations in 2001, the authors reported evidence of increasing spot coverage towards lower latitudes, but their observations were not of sufficiently high signal-to-noise to determine whether these spots formed in a distinct latitude band. This work confirms that spots do indeed form at lower latitudes in distinct bands on AE Aqr. One key difference, however, between the 2001 map and the maps presented here is that we find a weak band of spots at a latitude centred around 43 degrees, whereas the 2001 observations show a distinct paucity of spots at the same latitude. This may be indicative of a solar-like magnetic activity cycle in operation in AE Aqr, where the latitude of spot formation may change in a manner that mimics the butterfly-diagram for the Sun. We tentatively speculate that the two latitude bands may be explained by a cross over between activity cycles, with the higher band just emerging and the lower band diminishing, as seen on the sun. The band is similar to that found on BV Cen (Watson et al. 2007).

The larger spots mentioned above are reconstructed independently on each tomogram, and are very similar in morphology (with some evolution). This similarity shows they are not artefacts due to noise or over fitting of the data.

Figure 15: Spot coverage on the northern hemisphere as a function of latitude for block 1 and block 2 data, normalised by the surface area at that latitude. Since the grid in our Roche tomograms is not aligned along strips of constant latitude, some interpolation between grid elements is required in order to produce these plots, resulting in their slightly noisy nature. The region around the L1 point has been excluded.

8 Differential rotation

The Roche tomograms produced in section 7 were then used to measure the differential rotation rate on the secondary star. The maps are separated by 9 days, over which we expect minimal spot evolution, apart from changing morphology due to differential rotation.

We assumed AE Aqr follows a solar-like differential rotation law of the form


where is the rotation rate at latitude , is the rotation rate at the equator, and is the difference between the rotation rates at the pole and the equator.

Two different methods were used to measure the differential rotation rate, and for consistency both maps were reconstructed at an inclination of , with M and M.

For the first method, we cross-correlated strips of constant latitude for each map from block 1 and block 2, taking the peak of the cross-correlation function (CCF) as the shift for that latitude. The tiling of the Roche lobe as implemented in Roche tomography, however, does not yield surface elements on strips of constant latitude, as this would lead to a discontinuity at the inner-Lagrangian point. For this reason we cannot incorporate a differential law into the reconstruction itself, as has been done by (for example) Hussain et al. (2006) and Petit et al. (2002). Instead, we interpolated the maps using the ‘griddata’ module from the SciPy library (Oliphant 2007) to produce an equirectangular-projection with resolution, which was then used for the cross-correlations, as these interpolated maps now contain constant latitude strips.

In order to estimate the error on our derived differential rotation measurement, we adopted a Monte Carlo approach to assessing the impact of noise on the Roche tomography reconstructions. The usual method of ‘jiggling’ the observed data in accordance with the error bars cannot be implemented in Roche tomography. This is because this process adds noise to the data, which means that the synthesized data sets are not being drawn from the same parent population as the observed data set and leads to maps dominated by noise. To surmount this issue, we took the reconstructed line profiles for each observation block and varied the flux of each data point about its value by the error bar in the observed data multiplied by a number output by a Gaussian random-number generator with zero mean and unit variance. In addition, to account for systematic effects in the observed data that were not reconstructed, we calculated a moving average for the residuals (thus removing unfitted noise from the observed data but keeping the large scale variation). This was added to the newly created line profiles, which were then visually checked to make sure the noise and shape of the line profiles were representative of the observed data. These new line profiles were then reconstructed to a similar map entropy as the observed data (the same reduced could not be used due to the variation introduced in the data). This whole process was repeated 100 times for each observation block and the resulting maps were cross-correlated as before, resulting in a distribution of shifts for each latitude. The uncertainty in the shifts found for the original maps was taken to be of this distribution.

The CCFs for each latitude may be found in an additional figure online. It was found that below latitude no significant shift was found, most likely due to the broad CCF peaks stemming from a lack of well defined features on the map. Note, the region around the L1 point up to latitude was ignored in the CCF due to the large, unchanging irradiation feature. Additionally, no shift was found above latitude, where the broad CCF peaks stem from the large polar-region spot appearing unchanged between block 1 and 2. For these reasons, we only included latitudes between in our analysis. A minimised- fit of equation 5 was made to the CCF peaks, iteratively adding an offset to account for the co-rotation latitude. This gave rad d-1 (corresponding to a equator-pole lap time of 269 days) and a co-rotation latitude of . The limits of this fit gave rad d-1 and rad d-1, with co-rotation latitudes of and , respectively (see figure 16).

Figure 16: The peaks of the cross-correlation of constant-latitude strips in the interpolated maps of block 1 and 2. The dashed line shows a minimised- best-fit of a solar-like differential rotation law, with rad d-1 and a co-rotation latitude of . The dotted lines show fits, with a minimum and maximum of 0.0176 rad d-1 and 0.0291 rad d-1, respectively.

In addition, to check the results from the CCF analysis, we implemented a second approach to measure the shear rate in the tomograms. This consisted of applying different levels of differential rotation to the 1st map (reconstructed from block 1 data and assuming a solar-like differential law shown in equation 5), and subtracting this ‘sheared map’ from the second map and analysing the residuals. In this process, the level of applied shear that results in the least residuals is deemed to best represent the true differential rotation rate. As before, since this requires strips of constant latitude, we used the map interpolated using the ‘griddata’ module as described earlier. The residuals found after map subtraction were weighted by (where is the latitude) since the map projection increases the weighting of higher latitudes if not corrected for. The residuals were then squared and summed and are shown in figure 17 as a percentage of the summed pixel intensities in the interpolated map of block 1 data. The least-residuals fit yielded rad d-1 (corresponding to a equator-pole lap time of 262 days) and a co-rotation latitude of , in good agreement with the CCF analysis.

Figure 17: The residuals after applying a solar-like DR law to the interpolated block 1 data and subtracting the interpolated block 2 data (see text for full description). The minimum gives rad d-1 and a co-rotation latitude of .

To more clearly show the movement of surface features, segments of the constant latitude strips are shown in figure 18, showing the surface intensity across the longitude range for various regions of the maps.

This analysis was also carried out after varying the masses by per cent, yielding equator-pole lap times of 278-354 d by changing , and 336-345 d by changing . However, these are not representative of the uncertainty on the shear rate since incorrect masses cause large artefacts to dominate the maps, thus diluting the contribution made by real features. To provide a more realistic error estimate, the systemic velocity was varied by kms-1, yielding a lap time of 267-296 d. The analysis was also carried out at an inclination of using M and M, as found by Watson et al. (2006). This yielded equator-pole lap times of 322 d and 316 d and co-rotation latitudes of and for the first and second methods, respectively. These are similar to the results found above, and show the robustness of the shear measurement against an incorrect inclination and systemic velocity. We therefore estimate the shear rate uncertainty to be rad d-1 and the corresponding lap time to be d.

Figure 18: The map intensity as a function of longitude. The solid line is block 1 data and the dashed line is block 2 data. The left column shows the unshifted data and the right column shows the same data after applying the shift found from the fit to the CCF peaks (see figure 16). The panels from top to bottom are: a bright region at latitude; the top of appendage ’A’ at latitude; a prominent spot ’B’ at latitude; the bottom of appendage ’A’ at . Note, the second and fourth panels down on the right are visually not shifted enough, reflecting the fit shown in figure 16.

9 Discussion

The measurement of differential rotation (DR) in a CV donor is the first of its kind and shows that the secondary star is (technically) not tidally locked. Theoretically, DR is predicted to be weak in tidally distorted systems (Scharlemann 1982) and the result found in this work confirms this observationally.

When compared to other systems we find the DR rate determined in this work is much lower than that for single stars, although LQ Hya is very similar, despite its rotation being four times slower. The pre-CV V471 Tau was found by Hussain et al. (2006) to have a drastically different DR rate (essentially that of a solid body), despite having a similar rotation period and mass, suggesting that the filling of the star’s Roche lobe may impact the DR rate. Finally, HR 1099 is significantly different in every parameter. The lack of systems with similar parameters with which to compare means no significant conclusions can be drawn, although we can say DR rates vary wildly in the systems observed so far.

Author Object System type Stellar type (hrs) (M) Lap time (days)
This work AE Aqr CV secondary K4V 9.88 0.78-0.8 269
(cross-correlation method)
This work AE Aqr CV secondary K4V 9.88 0.78-0.8 262
(map subtraction method)
Hussain et al. (2006) V471 Tau Tidally-locked pre-CV K2V 12.5 3900
Petit et al. (2004) HR 1099 Tidally locked RS CVn K1 68.1 480
Cameron (2001) AB Dor Single star K0V 12.36 0.76 70-140
(in a distant trinary system)
Donati et al. (2000) RX J Single star G2 7.4
Barnes et al. (2000) PZ Tel Single star G9IV/V 22.68 1.1
Kővári et al. (2005) LQ Hya Single star K2 38.4 0.8 280

(3) Guinan & Ribas (2001); (4) Fekel (1983)

Table 5: Comparison of system types and differential rotation rates, as found by the respective authors.

9.1 Mass transfer variations

The presence of differential rotation on AE Aqr has significant implications for theories invoking starspots to explain rapid variations in accretion rates. The passage of spots across the mass-transfer nozzle may cease mass-transfer due to the depression of the photosphere by strong magnetic fields in a starspot, decreasing thermal pressure and density (Livio & Pringle 1994). Hessman et al. (2000), in their study of mass-transfer variations in AM Her, concluded that about half the surface of the star near the L1 point needs to be covered with spots to produce the observed variations. They note this scenario requires a mechanism which preferentially produces magnetic flux around the L1 point and/or one which forces spot groups, which appear at higher latitudes, to wander down towards the L1 point. The trail of spots leading down to the L1 point, which is evident in several previous Roche tomograms as well as this work, lends observational support for this scenario (for further discussion see section 9.2).

King & Cannizzo (1998) note that the most rapid transitions to low states occur over the course of a day. Using our measurement of differential rotation, we can estimate how long a spot would take to traverse the mass-transfer nozzle and cut off mass-transfer. We estimate the mass-transfer nozzle radius of AE Aqr to be kms-1 (see section 4 of Hessman et al. 2000). Taking the average radius of the star as km, combined with the equatorial longitude-shift over the 9 day observation gap, we estimate the equator rotates at a speed of ms-1. This means it takes days for a spot to fully cover the mass-transfer nozzle, assuming the spotted region is large enough, and that the surface flow is constant across the L1 region. This is similar to the length of time taken for the large dips in brightness observed in AM Her, which has a much shorter orbital period and is of a different spectral type, suggesting that this is a viable mechanism.

In addition, the highly magnetised region around the mass-transfer nozzle may cause the formation of magnetised blobs of material in the accretion stream, causing the flicker seen in CV lightcurves (Bruch 1992).

9.2 Surface flows: Impact on starspot movement

To understand why surface features are located where they are, and to determine a mechanism which preferentially produces magnetic flux around the L1 point, we can compare the Roche tomograms to theoretical predictions and numerical simulations of surface mass flow on Roche lobe filling secondaries.

Oka et al. (2002) define the ‘surface’ as the region where appreciable mass flow occurs, and is not necessarily the photosphere. Lubow & Shu (1975) state that mass flow is directed around contours of equal pressure, proceeding from high to low pressure, with the highest pressure at the poles. They predict that surface flow may be astrostrophic (i.e. parallel to isobars) on each equipotential surface, but that matter may be transferred from the entire surface layer to the equator if there is an outward flux at the bottom of the surface layer.

In their numerical simulations, Oka et al. (2002) find gas ascends in higher density regions, and has a circulating flow in lower density regions. They find a high pressure circulatory flow around the north pole (denoted the H-eddy), where gas seems to drift gradually downward, increasing its velocity with decreasing latitude. Near the L1 region they find a strong low pressure, due to outflow through the L1 point, and a circulatory flow around the L1 point. In a counter-cockwise rotating secondary (as viewed from the north pole), this L1 eddy rotates in a coulter-clockwise direction. The L1 eddy is fixed to the L1 region, but as this forces the H-eddy towards leading hemisphere, the system is not steady. To stabilise the eddy system, a third circulation is needed in the form of a L2 eddy on the back of the star, which circulates counter-clockwise and acts to convect the H-eddy towards the trailing hemisphere direction. These three eddy features are present in simulations covering a mass-ratio range (Oka et al. 2004).

These simulations may provide an explanation for the trail of spots leading from the pole to the L1 point, as seen in Roche tomograms of AE Aqr and BV Cen (Watson et al. 2006, 2007; this work). These may be caused by the combined effect of the H-eddy and the L1 eddy, bringing magnetic flux tubes which emerge near polar regions down the leading hemisphere side of the L1 region. In addition, a combination of the H-eddy and L2 eddy may cause magnetic flux tubes to be dragged down the trailing hemisphere, creating the large appendage we see in the AE Aqr tomograms in this work, and additionally causing the morphological change we see between the two blocks. This would lead to a skewed measure of differential rotation, since the surface would be circulating rather than truly rotating longitudinally. However, this explanation is lacking as the simulations by Oka et al. do not account for magnetic effects, so we do not know how magnetic flux tubes are influenced by this “surface” flow.

9.3 DR rate variations and future opportunities

Substantial changes in differential rotation on magnetically active binaries may cause quasi-cyclic orbital period changes (Applegate 1992). A stellar activity cycle may alter the viscous transport of angular momentum sufficiently to produce substantial changes in DR, altering the gravitational quadropole moment of the star sufficiently to produce orbital period changes (e.g. Donati 1999). Collier Cameron & Donati (2002) suggest that the large modulation of surface DR in AB Dor would alter the stars oblatness such that if it were a close binary system it would be expected to produce observable long-term orbital period changes. This effect could be tested in AE Aqr with future DR measurements, testing if the Applegate mechanism exists in a Roche lobe filling star. In addition, variation in DR could alter the mass transfer rate (see section 9.1) such that it would affect the evolution of the binary system as a whole. This is the first measure of DR for a CV and it would be interesting to compare this to another with different system parameters. A good candidate is the secondary in BV Cen, as it is Roche lobe filling, has a longer orbital period ( hrs) and is more massive (). We intend to carry out this study in the near future.

10 Conclusions

We have unambiguously imaged starspots on two Roche tomograms of AE Aqr, separated by 9 days. We find a spot coverage of per cent. A comparison with a previous Roche tomogram of AE Aqr by Watson et al. (2006) shows many similarities, confirming the highly spotted nature of the secondary. The spot distributions indicate preferred latitudes for spot formation, with a ‘chain’ of spots from the pole to the L1 point and the two distinct latitude bands of spots around and latitude.

Using these two Roche tomograms, we have measured the surface differential rotation using two different techniques — cross-correlation of constant-latitude strips, and subtraction of maps after injecting a solar-like differential law into the first tomogram. The first yields an equator-pole lap time of 269 days and the second yields a lap time of 262 days, with a co-rotation latitude of . This shows that the star is not tidally locked (as was previously assumed for CVs), and is in stark contrast with the near solid body rotation of the tidally locked pre-CV V471 Tau, found by Hussain et al. (2006). The discovery of differential rotation has large implications for stellar dynamo theory, and the cessation of mass transfer due to spot traversal of the L1 point. The highly magnetised region around the L1 point may create magnetised blobs of material, explaining the ‘blobby’ accretion observed in CVs.


We thank Tom Marsh for the use of his molly software package in this work, and VALD for the stellar line-lists used. C.A.H. acknowledges the Queen’s University Belfast Department of Education and Learning PhD scholarship, C.A.W. acknowledges support by STFC grant ST/L000709/1. T.S. acknowledges support by the Spanish Ministry of Economy and Competitiveness (MINECO) under the grant (project reference AYA2010-18080). This research has made use of NASA’s Astrophysics Data System and the Ureka software package provided by Space Telescope Science Institute and Gemini Observatory. Finally, we would like to thank the referee, Gaitee Hussain, for their comments that helped to substantially improve the paper.


  1. pagerange: Roche tomography of cataclysmic variables - VI. Differential rotation of AE Aqr - Not tidally locked!References
  2. pubyear: 2013


  1. Ak, T., Ozkan, M. T., Mattei, J. A., 2001, A&A, 369, 882
  2. Applegate, J. H., 1992, ApJ, 385, 621
  3. Barnes, J., 1999, PhD thesis, University of St. Andrews
  4. Barnes, J. R., Collier Cameron, A., James, D. J., Donati, J.-F., 2000, MNRAS, 314, 162
  5. Barnes, J. R., Lister, T. A., Hilditch, R. W., Collier Cameron, A., 2004, MNRAS, 348, 1321
  6. Bianchini, A., 1990, AJ, 99, 1941
  7. Biazzo, K., Frasca, A., Catalano, S., Marilli, E., Henry, G. W., Tǎs, G., 2006, Memorie della Societa Astronomica Italiana Supplementi, 9, 220
  8. Bruch, A., 1991, A&A, 251, 59
  9. Bruch, A., 1992, A&A, 266, 237
  10. Cameron, A. C., 2001, in H. M. J. Boffin, D. Steeghs, & J. Cuypers, ed., Astrotomography, Indirect Imaging Methods in Observational Astronomy, vol. 573 of Lecture Notes in Physics, Berlin Springer Verlag, p. 183
  11. Casares, J., Mouchet, M., Martinez-Pais, I. G., Harlaftis, E. T., 1996, MNRAS, 282, 182
  12. Chanan, G. A., Middleditch, J., Nelson, J. E., 1976, ApJ, 208, 512
  13. Chincarini, G., Walker, M. F., 1981, A&A, 104, 24
  14. Claret, A., 2000, VizieR Online Data Catalog, 336, 31081
  15. Collier Cameron, A., Donati, J.-F., 2002, MNRAS, 329, L23
  16. Collier-Cameron, A., Unruh, Y. C., 1994, MNRAS, 269, 814
  17. Crawford, J. A., Kraft, R. P., 1956, ApJ, 123, 44
  18. Davey, S., Smith, R. C., 1992, MNRAS, 257, 476
  19. Dhillon, V. S., Watson, C. A., 2001, in H. M. J. Boffin, D. Steeghs, & J. Cuypers, ed., Astrotomography, Indirect Imaging Methods in Observational Astronomy, vol. 573 of Lecture Notes in Physics, Berlin Springer Verlag, p. 94
  20. Donati, J.-F., 1999, MNRAS, 302, 457
  21. Donati, J.-F., Collier Cameron, A., 1997, MNRAS, 291, 1
  22. Donati, J.-F., Semel, M., Carter, B. D., Rees, D. E., Collier Cameron, A., 1997, MNRAS, 291, 658
  23. Donati, J.-F., Mengel, M., Carter, B. D., Marsden, S., Collier Cameron, A., Wichmann, R., 2000, MNRAS, 316, 699
  24. Dunford, A., Watson, C. A., Smith, R. C., 2012, MNRAS, 422, 3444
  25. Echevarría, J., Smith, R. C., Costero, R., Zharikov, S., Michel, R., 2008, MNRAS, 387, 1563
  26. Efron, B., 1979, The Annals of Statistics, 7, 1
  27. Efron, B., Tibshirani, R., 1993, An introduction to the bootstrap, Monographs on statistics and applied probability, Chapman & Hall
  28. Fekel, Jr., F. C., 1983, ApJ, 268, 274
  29. Frasca, A., Biazzo, K., Catalano, S., Marilli, E., Messina, S., Rodonò, M., 2005, A&A, 432, 647
  30. Guinan, E. F., Ribas, I., 2001, ApJ, 546, L43
  31. Hatzes, A. P., Vogt, S. S., 1992, MNRAS, 258, 387
  32. Hessman, F. V., Gänsicke, B. T., Mattei, J. A., 2000, A&A, 361, 952
  33. Holzwarth, V., Schüssler, M., 2003, A&A, 405, 303
  34. Hussain, G. A. J., Allende Prieto, C., Saar, S. H., Still, M., 2006, MNRAS, 367, 1699
  35. Hussain, G. A. J., et al., 2007, MNRAS, 377, 1488
  36. Jeffers, S. V., Barnes, J. R., Collier Cameron, A., 2002, MNRAS, 331, 666
  37. Kővári, Z., Weber, M., Strassmeier, K. G., 2005, in Favata, F., Hussain, G. A. J., Battrick, B., eds., 13th Cambridge Workshop on Cool Stars, Stellar Systems and the Sun, vol. 560 of ESA Special Publication, p. 731
  38. King, A. R., Cannizzo, J. K., 1998, ApJ, 499, 348
  39. Kraft, R. P., 1967, ApJ, 150, 551
  40. Kupka, F., Ryabchikova, T. A., 1999, Publications de l’Observatoire Astronomique de Beograd, 65, 223
  41. Kupka, F. G., Ryabchikova, T. A., Piskunov, N. E., Stempels, H. C., Weiss, W. W., 2000, Baltic Astronomy, 9, 590
  42. Livio, M., Pringle, J. E., 1994, ApJ, 427, 956
  43. Lubow, S. H., Shu, F. H., 1975, ApJ, 198, 383
  44. Marsden, S. C., Waite, I. A., Carter, B. D., Donati, J.-F., 2005, MNRAS, 359, 711
  45. Mestel, L., 1968, MNRAS, 138, 359
  46. Moss, D., Piskunov, N., Sokoloff, D., 2002, A&A, 396, 885
  47. Oka, K., Nagae, T., Matsuda, T., Fujiwara, H., Boffin, H. M. J., 2002, A&A, 394, 115
  48. Oka, K., Matsuda, T., Hachisu, I., Boffin, H. M. J., 2004, A&A, 419, 277
  49. Oliphant, T. E., 2007, Computing in Science & Engineering, 9, 10
  50. Petit, P., Donati, J.-F., Collier Cameron, A., 2002, MNRAS, 334, 374
  51. Petit, P., et al., 2004, in A. Maeder & P. Eenens, ed., Stellar Rotation, vol. 215 of IAU Symposium, p. 294
  52. Rappaport, S., Verbunt, F., Joss, P. C., 1983, ApJ, 275, 713
  53. Richman, H. R., Applegate, J. H., Patterson, J., 1994, Publications of the Astronomical Society of the Pacific, 106, 1075
  54. Rutten, R. G. M., Dhillon, V. S., 1994, A&A, 288, 773
  55. Rutten, R. G. M., Dhillon, V. S., 1996, in A. Evans & J. H. Wood, ed., IAU Colloq. 158: Cataclysmic Variables and Related Objects, vol. 208 of Astrophysics and Space Science Library, p. 21
  56. Scharlemann, E. T., 1982, ApJ, 253, 298
  57. Schuessler, M., Solanki, S. K., 1992, A&A, 264, L13
  58. Schwope, A. D., Staude, A., Vogel, J., Schwarz, R., 2004, Astronomische Nachrichten, 325, 197
  59. Shahbaz, T., Watson, C. A., Dhillon, V. S., 2014, MNRAS
  60. Simpson, E. K., et al., 2011, MNRAS, 414, 3023
  61. Skilling, J., Bryan, R. K., 1984, MNRAS, 211, 111
  62. Spruit, H. C., Ritter, H., 1983, A&A, 124, 267
  63. Tanzi, E. G., Chincarini, G., Tarenghi, M., 1981, PASP, 93, 68
  64. Vogt, S. S., Penrod, G. D., 1983, Publications of the Astronomical Society of the Pacific, 95, 565
  65. Vogt, S. S., Hatzes, A. P., Misch, A. A., Kürster, M., 1999, The ApJSupplement Series, 121, 547
  66. Watson, C. A., Dhillon, V. S., 2001, MNRAS, 326, 67
  67. Watson, C. A., Dhillon, V. S., Rutten, R. G. M., Schwope, A. D., 2003, MNRAS, 341, 129
  68. Watson, C. A., Dhillon, V. S., Shahbaz, T., 2006, MNRAS, 368, 637
  69. Watson, C. A., Steeghs, D., Shahbaz, T., Dhillon, V. S., 2007, MNRAS, 382, 1105
  70. Webb, N. A., Naylor, T., Jeffries, R. D., 2002, The ApJ, 568, L45
  71. Welsh, W. F., Horne, K., Gomer, R., 1995, MNRAS, 275, 649
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description