EDP Sciences
Free Access
Issue
A&A
Volume 597, January 2017
Article Number A58
Number of page(s) 15
Section Numerical methods and codes
DOI https://doi.org/10.1051/0004-6361/201629768
Published online 22 December 2016

© ESO, 2016

1. Introduction

Intermediate-mass chemically peculiar (CP) stars are A and B-type main sequence objects that exhibit unusual surface characteristics compared to normal stars of a similar spectral class. CP stars rotate slowly and show extreme chemical anomalies in their atmospheres, often corresponding to deviations from the solar chemical composition by several orders of magnitude for some heavy elements. This unusual surface chemistry is attributed to the processes of radiative levitation and gravitational settling, commonly referred to as atomic diffusion (Michaud 1970; Michaud et al. 2015), operating in stably stratified stellar atmospheres.

A subgroup of CP stars, known as magnetic Ap/Bp or CP2 stars according to the classification by Preston (1974), possesses strong, stable, globally-organised magnetic fields on their surfaces. These fields, described to a first approximation by oblique dipolar geometries, are believed to be fossil remnants of the magnetic flux generated by dynamo or acquired from interstellar medium at an earlier evolutionary phase. The kG-strength fields further stabilise stellar atmospheres and facilitate diffusive separation of chemical elements; consequently magnetic Ap/Bp stars exhibit the most extreme chemical anomalies among CP stars. Moreover, due to an anisotropic and non-axisymmetric character of the global magnetic field topologies of Ap/Bp stars, the atomic diffusion process operates differently at different parts of the stellar surface, resulting in substantial horizontal and vertical chemical abundance gradients (Michaud et al. 1981; Babel & Michaud 1991).

Chemical inhomogeneities in the atmospheres of magnetic Ap/Bp stars have profound observational manifestations in the photometric, spectroscopic and spectrophotometric behaviour of these stars. All these different types of observables show stable, strictly repeating variation with the same period as the change of the disk-integrated magnetic field characteristics. This coherent variability is therefore commonly interpreted in the framework of the so-called oblique rotator model (Stibbs 1950), postulating that the surface chemical and magnetic field distributions are static and the sole reason for the observed variability is rotational modulation of the aspect at which the surface structure is viewed by a distant observer.

Since the introduction of the oblique rotator model considerable effort was spent on deriving the geometry of chemical spots on Ap/Bp stars and understanding their relation to the surface magnetic field topology. The attempts to map the surface structure of magnetic CP stars date back to the work by Deutsch (1958, 1970), who developed a technique of fitting phase curves of magnetic, equivalent width and radial velocity measurements with surface maps expanded in spherical harmonics. A more refined star-spot mapping approach, utilising the entire information content of stellar line profiles, was introduced by Khokhlova and collaborators (Khokhlova & Riabchikova 1975; Khokhlova 1976). Their method, now commonly known as Doppler imaging (DI)1, exploits spatial resolution of the stellar surface by the rotational Doppler effect. Provided the local (intrinsic) spectral line width is substantially smaller than the stellar projected rotational velocity, any inhomogeneity on the stellar surface yields a distortion in the disk-integrated line profile displaced from the line centre according to the local rotational Doppler shift and hence the spot longitude. This distortion moves across the line profile as the star rotates. Spots located at different latitudes exhibit different temporal variability patterns. For example, a spot near the stellar rotational equator produces a signature that is visible during half of the rotation period and moves quickly across the line profile. On the other hand, spots close to the rotational pole are visible during a larger fraction of the rotational cycle; the corresponding line profile distortions are constrained to the line centre.

Initial attempts to apply DI to line profile observations of Ap/Bp stars were based on the trial-and-error approach, often leading to subjective and non-unique solutions (Megessier et al. 1979). A breakthrough contribution was made by Goncharskii et al. (1977, 1982), who recasted the DI technique as a least-squares, regularised image reconstruction inverse problem (see Kochukhov 2016, for a recent review). In this formulation a given chemical element map is specified on a discrete latitude-longitude surface grid and local abundances in each pixel of the map are iteratively adjusted by fitting model spectra to the observed spectral line profiles. The surface abundance distribution is further constrained by a regularisation function. This is a key mathematical ingredient of solving an ill-posed inverse problem. Regularisation ensures convergence to a unique solution and stability of the inversion results with respect to observational noise.

The initial and most current applications of DI to chemical spot mapping on Ap/Bp stars relied on the Tikhonov regularisation method (Tikhonov & Arsenin 1977). A smaller number of studies (e.g. Hatzes et al. 1989) used an alternative maximum entropy regularisation function. Although the external criteria imposed by the two regularisation methods are somewhat different (Tikhonov regularisation favours solutions with the smallest local gradients, maximum entropy favours solutions with the least deviation from a default value), it was demonstrated that DI maps obtained with alternative regularisation strategies are generally compatible with each other (Piskunov 1990; Strassmeier et al. 1991) provided that they are based on the observational data with sufficiently high signal-to-noise ratio (S/N) and a dense rotational phase coverage.

Early DI studies of Ap stars typically presented maps of a few chemical elements based on one or two isolated spectral lines for each element. Analytical functions were often employed to describe the local line profiles (e.g. Khokhlova et al. 1986; Hatzes 1990; Rice & Wehlau 1990). As computational resources became more affordable, DI mapping could be extended to between five and more than a dozen chemical species, using multiple lines for a given element whenever possible, and treating line blends with detailed radiative transfer calculations in realistic model atmospheres (e.g. Kochukhov et al. 2004b; Lüftinger et al. 2010a; Nesvacil et al. 2012).

However, some important simplifications have lingered. Looking at the bulk of published DI maps of Ap/Bp stars one can conclude that the most common approach has been to derive chemical abundance distributions from a small number, sometimes just one or two, spectral lines, using a single model atmosphere and ignoring effects of magnetic field in calculation of the model spectra. Several recent studies took advantage of the magnetic DI (MDI, Piskunov & Kochukhov 2002) extension of the indirect mapping technique to recover self-consistently magnetic field topologies and chemical spot distributions from the variation of Stokes parameter spectra (e.g. Kochukhov et al. 2015; Silvester et al. 2015; Rusomarov et al. 2016). However, only a few stars were analysed with this sophisticated MDI method; most observational constraints on the chemical element distributions are still provided by non-magnetic DI investigations. Even fewer DI analyses have addressed the problem of a feedback of chemical spots on the local atmospheric structure as part of the inversion (Lehmann et al. 2007; Kochukhov et al. 2012). Nevertheless, a series of studies (Krtička et al. 2007, 2012, 2015; Shulyak et al. 2010) has modelled photometric and spectrophotometric variation of several Ap/Bp stars based on calculations of theoretical multi-component model atmospheres with the local abundances taken from the DI analyses of these stars.

The question of reliability of abundance DI inversions was repeatedly addressed by numerical tests in which known surface maps were reconstructed from simulated observations (Hatzes et al. 1989; Rice et al. 1989; Piskunov & Wehlau 1990; Rice 1991; Wehlau & Rice 1993; Kochukhov & Piskunov 2002). These studies demonstrated that diverse, often very complicated, surface chemical over- and under-abundance patterns (circular, square and complex letter-shaped spots, rings of different sizes) are recovered satisfactorily by the same inversion codes as applied to real stellar observations. Previous simulations also provided a comprehensive assessment of the sensitivity of DI maps to the value and uncertainty of key stellar parameters (projected rotational velocity and inclination angle) and to the quality (S/N, phase coverage, resolution) of observational data. However, several sources of potentially important systematic errors inherent to non-magnetic DI inversions have not been addressed before. The goal of this paper is therefore to present a new series of DI simulations and tests using the most sophisticated and physically realistic version of the indirect abundance mapping technique currently available. Based on these calculations, we quantify the loss of information and appearance of spurious structures in DI maps due to using few and/or blended spectral lines, neglecting the moderately strong magnetic field, and using a single average model atmosphere structure. A characterisation of these error sources enables an informed judgement on the level of reliability of several dozen abundance DI studies published over the past 30 yr.

This paper is organised as follows. Section 2 reviews different versions of Doppler imaging codes applied in previous studies of Ap/Bp stars and in this paper. Section 3 describes the setup adopted for numerical DI experiments and discusses the results of the inversions using different line lists and different assumptions about the local line formation physics. In Sect. 4 we present a literature survey of DI studies, assessing in light of our simulations the reliability of chemical spot maps published for 40 CP stars. The main conclusions of our study are summarised in Sect. 5.

2. Doppler imaging (DI) and MDI inversion codes

2.1. INVERS8

This is an older non-magnetic abundance DI code (Piskunov & Rice 1993) used in many previous DI studies of Ap/Bp stars reviewed in Sect. 4. The spectrum synthesis calculations by Invers8 are based on interpolation in the precomputed local line intensity tables, which are generated by approximating abundance variation of a single chemical element by scaling its line opacity. A single model atmosphere is used and continuum intensity is considered to be independent of the local abundance. Inversions can incorporate several wavelength regions.

Despite relying on a single model atmosphere in its standard implementation, Invers8 can be easily modified to calculate the local line profile and continuum intensity tables using a grid of model atmospheres with varying abundance (Lehmann et al. 2007). This code is not used in the numerical experiments discussed here.

2.2. INVERS11 andINVERS12

These are newer non-magnetic abundance DI codes (Lueftinger et al. (Lueftinger et al. 2003); Kochukhov et al. 2004b). They are based on the same physical principles and assumptions as Invers8 but calculate the local line profiles on the fly and enable simultaneous reconstruction of several chemical maps from observations in multiple wavelength regions and from blended spectral features. Radiative transfer calculations in Invers12 are distributed over multiple processors with the Message Passing Interface (MPI) libraries. These codes are not used in the numerical experiments presented in this paper.

thumbnail Fig. 1

Surface magnetic field topology adopted for numerical tests of DI. The star is shown at five rotation phases, which are indicated above the spherical plots, and an inclination angle of i = 60°. The spherical plots show the maps of a) field modulus; b) horizontal field; c) radial field; and d) field orientation. The contours over spherical maps are plotted with a 0.5 kG step. The thick line and the vertical bar indicate positions of the rotational equator and the pole, respectively. The colour bars give the field strength in kG. The two different colours in the field orientation map correspond to the field vectors directed outwards (red) and inwards (blue).

Open with DEXTER

2.3. INVERS10

This is the main magnetic DI code currently applied to Ap/Bp stars. The physical basis and numerical methods implemented in Invers10 were described in detail by Piskunov & Kochukhov (2002); the code was thoroughly tested by Piskunov & Kochukhov (2002). The four Stokes parameter spectra computed with Invers10 agree with the calculations by other independent polarised radiative transfer codes (Wade et al. 2001) and are considered as a benchmark for stellar Stokes parameter synthesis (Carroll et al. 2008; Deen 2013).

Invers10 reconstructs simultaneously the three magnetic field vector components and an arbitrary number of chemical abundance maps from the Stokes IQUV or Stokes IV parameter spectra in multiple wavelength regions. The code can also be used for abundance mapping with Stokes I alone, using a prescribed, constant magnetic field topology or ignoring the field altogether. Invers10 uses a single model atmosphere for numerical polarised radiative transfer calculations, which are performed on the fly and parallelised with MPI. Chemical abundance variations are approximated by scaling the central opacity of relevant spectral lines. Thus, no influence of chemical spots on the continuum intensity and on the atmospheric structure is taken into account. Invers10 also includes an option to automatically adjust oscillator strengths of selected spectral lines to mitigate the line-to-line mean equivalent width scatter due to errors in the input atomic data and/or neglect of some physical effects (non-local thermodynamical equilibrium, vertical chemical stratification). All test inversions presented in this study were performed with Invers10.

2.4. INVERS13

This MDI code represents a modification of Invers10 intended primarily for mapping magnetic fields and temperature spots on active cool stars (Kochukhov et al. 2013b). Invers13 uses a grid of model atmospheres to calculate the local continuum intensity and Stokes parameter spectra according to the local temperature. It also incorporates an advanced equation of state solver suitable for treating molecular equilibrium at low temperatures (Piskunov & Valenti 2017). Thus, it can equally successfully model atomic and molecular lines. In the context of abundance mapping of Ap/Bp stars, a self-consistent treatment of the line formation and a model atmosphere structure allows one to use Invers13 for recovering abundance maps fully accounting for the influence of chemical inhomogeneities on the local atmosphere (Kochukhov et al. 2012; Oksala et al. 2015). In this case individual models in the input atmospheric grid, usually calculated with the LLmodels code (Shulyak et al. 2004), have different chemical abundance rather than different temperature. To this end, Invers13 implements the most detailed and sophisticated simultaneous magnetic and chemical abundance inversion methodology. However, its ability to model only one chemical element at a time is a drawback compared to Invers10. In this paper we used Invers13 to generate simulated observations in Sect. 3.4.5.

2.5. INVERSLSD

This MDI code, developed by Kochukhov et al. (2014), is designed for performing star spot and magnetic field inversions using least-squares deconvolved (LSD, Donati et al. 1997; Kochukhov et al. 2010) Stokes parameter profiles. The spectrum synthesis calculations by InversLSD are based on interpolation in the precomputed local LSD Stokes parameter profile tables. The latter can be generated either using simplified analytical solutions (weak-field approximation, Unno-Rachkovsky solution) of the polarised radiative transfer equation (e.g. Kochukhov & Wade 2016) or by applying the LSD line-averaging procedure to the detailed local polarised radiative transfer calculations of the entire stellar spectrum (Kochukhov et al. 2014; Rosén et al. 2015). In addition to mapping the stellar magnetic field topology, InversLSD allows one to recover a single scalar map (continuum brightness, temperature, abundance of one chemical element). This magnetic inversion code is not used in the numerical experiments presented in this paper.

3. Numerical experiments

3.1. Simulation setup

For the numerical tests of abundance DI described below we adopted the following stellar parameters: effective temperature Teff = 8000 K, surface gravity log g = 4.0, projected rotational velocity vesini = 40 km s-1, inclination angle i = 60°. The adopted vesini and i are in the middle of, a relatively wide, range optimal for DI. The chosen atmospheric parameters represent cool Ap stars to which multi-element abundance mapping has been applied in a number of recent studies (e.g. Kochukhov et al. 2004b; Lüftinger et al. 2010b; Nesvacil et al. 2012).

The corresponding stellar model atmosphere employed in spectrum synthesis was calculated with the LLmodels code (Shulyak et al. 2004) for the iron abundance log NFe/Ntot = −4.0 and solar abundance of other elements. A grid of LLmodels atmospheres with log NFe/Ntot spanning a range from −2.5 to −4.0 with a 0.25 dex step was generated for the multi-component atmosphere inversion test described in Sect. 3.4.5.

The star was assumed to have a dipolar magnetic field topology with a polar field strength of 2.5 kG and 90° inclination of the magnetic field axis relative to the rotational axis. This surface magnetic field configuration is illustrated in Fig. 1.

Two different Fe surface abundance distributions were considered in the numerical experiments. The first comprised four circular overabundance spots, with log NFe/Ntot = −2.5 in the spot centres and log NFe/Ntot = −4.0 in the background. The spots were placed at the latitudes of −30, 0, 30, and 60° and spaced equidistantly in longitude. We adopted the inner spot radius of rin = 15° and the outer radius of rout = 30°. The element abundance was linearly interpolated between rin and rout, resulting in smooth spot edges.

This multiple spot configuration is more complex than the three-spot abundance map previously considered in the numerical tests of Invers10 (Kochukhov & Piskunov 2002). In particular, it includes a low-latitude feature which is challenging to recover since it is visible close to the stellar limb in a narrow range of rotational phases. It should also be noted that this four-spot configuration is more structured than the typical chemical abundance distributions found for magnetic Ap/Bp stars: for the majority of stars single or double-spot distributions are prevalent resulting in a single or double-wave photometric variability (e.g. Bernhard et al. 2015).

The second chemical abundance distribution considered for DI experiments was represented by a ring of Fe overabundance (log NFe/Ntot = −2.5) on a uniform background (log NFe/Ntot = −4.0). Reckoning from the symmetry axis, the ring had rin = 85° and rout = 70° (yielding the width of 10° at maximum abundance and 40° at the base) and was chosen to coincide with the dipolar field equator. This relatively unusual element distribution pattern has been observed for certain elements in a few Ap stars (Rice et al. 1997; Kochukhov et al. 2004b) but is predicted to be the most common non-uniform chemical spot distribution according to the theoretical atomic diffusion calculations in the presence of a magnetic field (Alecian & Vauclair 1981; LeBlanc et al. 2009; Alecian & Stift 2010; Alecian 2015).

In one test, described in Sect. 3.4.3, we combined the two maps in a multi-element abundance mapping experiment. In that case the four-spot distribution was adopted for Fe and the ring-like abundance configuration was used for Cr.

For most of the tests described below the DI calculations were based on seven Fe i spectral lines distributed over 6 spectral intervals (the lines Fe i 5005.7 and 5006.1 Å are blended for the adopted vesini of 40 km s-1) in the 5006–5397 Å wavelength range. This group of neutral Fe lines was previously employed by Lüftinger et al. (2010b) in their DI study of HD 24712. For the experiments in Sect. 3.4.3 the spectra were calculated for a blend of Fe ii 6147.7, 6149.2 Å and Cr ii 6146.2 Å lines used by Kochukhov et al. (2004b) in the DI study of HD 83368. The main atomic parameters of the modelled transitions, including the effective Landé factors, are summarised in Table 1. This information was obtained from the vald database (Kupka et al. 1999).

Table 1

Spectral lines used in numerical tests of abundance DI.

The simulated observational data were produced for 20 equidistant rotational phases. Spectra were convolved with a Gaussian instrumental profile corresponding to the resolution of λ/ Δλ = 110 000 and were sampled with a step of Δλ = 0.02 Å. A normally distributed random noise component with σ = 1.9 × 10-3 was added to each spectrum in the inversions with seven Fe i lines. In the inversions with the two Fe ii lines and the blended Cr ii line we adopted σ = 10-3 to roughly offset an increase of random noise due to a smaller number of wavelength points in the simulated spectra. The adopted spectral resolution corresponds to the one provided by the HARPS-type spectrograph (Mayor et al. 2003). The assumed S/N ratio and phase coverage are realistic compared to recent observations of bright Ap stars (e.g. Silvester et al. 2012; Wade et al. 2016).

Chemical abundance maps were recovered using a 1876-element surface grid with variable number of zones in 38 latitude belts (see Fig. 5 in Piskunov & Kochukhov (2002) for an example of such a surface grid). All inversions started from a homogeneous element distribution with log Nel/Ntot = −4.0.

3.2. Choice of regularisation parameter

A considerate choice of regularisation parameter is essential for a successful Doppler imaging inversion. A too-strong regularisation will smooth out real features in a stellar surface map and will yield a poor fit to observations. Conversely, a too-weak regularisation might result in a lot of unreliable high-contrast surface details. The usual approach to finding an optimal regularisation parameter (Λ) is to ensure a balance between reaching the lowest possible chi-square (χ2) of the fit to observations on the one hand and limiting unnecessary high-contrast details in the surface maps on the other hand. Our extensive experience with different DI problems suggests that an optimal regularisation parameter is such that i) no significant improvement of the fit quality is possible to achieve by decreasing Λ and ii) the regularisation functional (i.e. regularisation function times the regularisation parameter) contribution to the total discrepancy function is a factor of between two and three below the corresponding χ2 contribution.

thumbnail Fig. 2

Illustration of the regularisation parameter determination for the inversion with two Fe ii lines. A measure of χ2 (solid black line) and the regularisation functional (solid red line) are plotted versus the Tikhonov regularisation parameter. The dotted line shows the noise level of the simulated data set. The optimal regularisation parameter is indicated with an arrow.

Open with DEXTER

Figure 2 presents a validation of this empirical prescription for determining the Tikhonov regularisation parameter. This figure shows the χ2 and regularisation functional values for a series of 12 DI inversions (corresponding to the two Fe ii line inversion test discussed in Sect. 3.4.3) starting with a large Λ and decreasing this parameter by a factor of three between consecutive inversions. The resulting dependence of χ2 on Λ is characteristic of an inverse problem regularised with the Tikhonov method. Initially, the fit to observations is poor for large Λ; it rapidly improves as Λ is let to decrease. At this stage regularisation functional is closely coupled and is comparable in magnitude to χ2. Then there is a break in the χ2 curve after which the fit quality is improving very slowly with decreasing Λ. At this point the regularisation functional decouples from the chi-square. The standard deviation of the fit reaches the nominal value of 10-3 when the regularisation functional is 2.6 times smaller than the chi-square.

This procedure of determining the Tikhonov regularisation parameter is closely related to the so-called L-curve method (Hansen 2000), which looks for the point of maximum curvature in the log Λ-log χ2 plot. In practice, some uncertainty in determining Λ is inevitable since the optimal fit quality is not known a priori and the χ2 variation with Λ might not have a well defined break. Nevertheless, the usual uncertainty in Λ is no more than a factor of between two and three; variation of the regularisation strength within these limits leads to a small and completely predictable modification of the resulting DI maps.

All DI calculations discussed in this paper were carried out with the step-wise regularisation parameter reduction, yielding 10–12 maps for different Λ values. The maps presented below correspond to the regularisation parameter chosen in the reference inversions (Sects. 3.4.1 and 3.4.2) following the procedure described above. The maps presented for other inversions correspond to the same regularisation parameter as for the reference maps.

Table 2

Numerical tests of abundance DI inversions.

3.3. Error analysis

We characterise the errors of abundance DI inversions by examining the difference between the true and reconstructed element distribution maps. Several complementary characteristics of the difference map are considered here and summarised in Table 2 for the different inversions discussed below. First, we determined the average offset between the input and reconstructed map. This offset is of no concern for DI mapping since the primary goal of this technique is to determine the relative variation of an element concentration across the stellar surface, not the absolute abundance. The average offset was subtracted from the difference map before carrying out further statistical analysis. We then determined the mean absolute abundance difference. This gives an estimate of the average inversion errors.

Of course, much higher discrepancies between the input and reconstructed maps are occasionally found locally at certain surface positions. The largest of these errors are irrelevant since they characterise a minor fraction of the stellar surface. Instead of discussing the maximum errors of individual surface elements we choose to consider discrepancies on the spatial scales comparable to the structures in adopted abundance maps. Given that the relative overabundance regions occupy from 27% (4 spots) to 34% (ring) of the stellar surface area, we adopt the 75th percentile of the unsigned difference map as a relevant measure of the maximum inversion errors.

In agreement with previous DI simulations (Berdyugina 1998; Piskunov & Kochukhov 2002), most inversions reported here showed the largest errors and coherent artefacts in the sub-equatorial, poorly visible parts of the stellar surface. This is a natural consequence of the strong latitudinal dependence of the visibility and spectral contribution of individual surface zones. To roughly assess this latitudinal trend we report the average and maximum reconstruction errors separately for latitudes −30° (essentially the entire part of the stellar surface where a reasonable DI reconstruction can be expected for the adopted i = 60°) and (the part of the surface with the best reconstruction quality).

In addition to the numerical error assessment reported in Table 2, we present a comparison between the true and recovered abundance maps in a series of spherical plots (Figs. 311). Each such figure is composed of three rows of panels showing the input abundance map, reconstructed map and the corresponding difference map after subtracting the mean offset. The true and reconstructed maps are rendered using the same linear colour scale, without applying any cuts to the maximum and minimum values. The difference map is shown using a different colour table scaled to the same abundance range – this allows one to judge the significance of error map structures in comparison to the input and recovered abundance distributions. Patterns in the difference maps are further highlighted with contours plotted every 0.3 dex.

thumbnail Fig. 3

Results of the Fe abundance inversion for the reference test case (4 circular spots, 7 Fe i lines, including magnetic field). The star is shown at five rotation phases, which are indicated above the spherical plots, and an inclination angle of i = 60°. The spherical plots show a) the true abundance map; b) reconstructed abundance map; c) difference map. The contours over the difference map are plotted at ±0.3,0.6,0.9 dex. The thick line and the vertical bar indicate positions of the rotational equator and the pole, respectively. The colour bars give the local abundance in log NFe/Ntot units.

Open with DEXTER

thumbnail Fig. 4

As Fig. 3 but for the reconstruction of the overabundance ring.

Open with DEXTER

thumbnail Fig. 5

As Fig. 3 but using two Fe ii lines.

Open with DEXTER

thumbnail Fig. 6

As Fig. 3 but for the simultaneous reconstruction of Cr (panels a)c)) and Fe (panels d)f)) abundance maps using two Fe ii lines and one blended line of Cr ii.

Open with DEXTER

thumbnail Fig. 7

As Fig. 3 but ignoring magnetic field. A constant offset is subtracted from the difference map shown in panel c).

Open with DEXTER

thumbnail Fig. 8

As Fig. 4 but ignoring magnetic field. A constant offset is subtracted from the difference map shown in panel c).

Open with DEXTER

thumbnail Fig. 9

Spurious chemical abundance distribution resulting from chemical mapping ignoring magnetic field. The true abundance map is a homogeneous distribution with log NFe/Ntot = −4.0. Reconstructed Fe abundance map is shown in the top panel adopting the same abundance range as in Fig. 7. The bottom panel shows the same distribution after subtracting a constant offset.

Open with DEXTER

thumbnail Fig. 10

As Fig. 3 but ignoring effect of local abundance variation on the model atmosphere structure.

Open with DEXTER

thumbnail Fig. 11

As Fig. 4 but ignoring effect of local abundance variation on the model atmosphere structure.

Open with DEXTER

3.4. Inversion results

3.4.1. Multiple circular spots

We started by assessing accuracy of the DI reconstruction of the multiple circular spot abundance map under the optimal conditions, taking into account the stellar magnetic field. The true and recovered maps for this test are shown in Fig. 3. As is evident from Fig. 3 and Table 2, very accurate inversion results are achieved. There is no significant mean offset between the true and recovered maps. The average reconstruction errors are 0.07–0.09 dex while the maximum errors are 0.11–0.12 dex. The largest errors are found for the region occupied by the sub-equatorial overabundance spot. This structure is inferred to have a somewhat larger extension and a lower contrast than in the input map. A careful inspection of the difference map shows a qualitatively similar but much less pronounced smoothing of higher latitude spots. These reconstruction artefacts are, however, entirely negligible compared to the 1.5 dex contrast of the input non-uniform Fe distribution.

3.4.2. Overabundance ring

In the second test we assessed the ability of our DI code to recover a ring-like overabundance structure. Similar to the previous test, magnetic field is fully accounted for in spectrum synthesis and no systematic errors are introduced in the DI reconstruction. The results of this experiment are illustrated in Fig. 4. Again, we find a successful recovery of the true overabundance structure, with average reconstruction errors of 0.06 dex, maximum errors of 0.08 dex, and no systematic offset between the input and recovered maps. As before, the largest discrepancy between the maps is found for the lowest visible latitude belts (−30°), which contribute negligibly to the disk-integrated spectra and, therefore, are poorly constrained by the inversion.

3.4.3. Inversions based on a few spectral lines

Among many chemical elements studied by the DI technique in Ap/Bp stars only a few species (Si, Cr, Fe) present a large number of useful diagnostic lines in the spectra of moderately and fast-rotating (vesini = 30–100 km s-1) stars. The majority of elements can be studied using only one or two lines because other lines are unsuitable due to blending or absent altogether in the optical wavelength range typically covered by spectroscopic observations. In this context, we are interested to ascertain the usefulness of a one or two-line DI inversion compared to using five or more lines with diverse parameters. To this end, we tested reconstruction of the four-spot Fe abundance distribution using the Fe ii 6147.7 and 6149.2 Å lines (originating from the same multiplet and having the same strength) instead of the 7 Fe i lines (having different strengths and excitation potentials) employed above. The purpose of this experiment is not to analyse a degradation of the DI images with increased observational noise. Therefore, we have increased the S/N of simulated observational data to compensate for the reduction of the number of wavelength points.

The results of the inversion using the two Fe ii lines are presented in Fig. 5. One can see, based on this figure and information in Table 2, that the reconstruction errors slightly increase compared to the inversion documented in Sect. 3.4.1. In particular, the Fe concentration is underestimated at the centre of the −30° latitude spot. The spot at + 30° is shifted to lower latitudes by 5–10°. Nevertheless, the average errors are still only 0.10–0.11 dex and the maximum errors are 0.13–0.15 dex, which is small compared to the 1.5 dex abundance contrast in the input map.

We have carried out a further test of multi-element DI mapping using a single blended line by adding the Cr ii 6146.2 Å line to the list of the two Fe ii lines. The same ring-like element distribution as studied in Sect. 3.4.2 was adopted for Cr. This experiment represents a challenging test of the simultaneous recovery of two abundance maps since for the adopted vesini the Cr ii line is an unresolved, weaker component of the blend with Fe ii 6147.7 Å.

Figure 6 shows the Cr and Fe abundance maps reconstructed simultaneously by Invers10. The quality of the Fe map reconstruction (average error 0.11–0.13 dex, maximum error 0.15–0.16 dex) is comparable to the inversion using the Fe ii lines alone. At the same time, the errors of the Cr ring-like map recovery (average error 0.08, maximum error 0.11) are almost the same as those achieved in Sect. 3.4.2 for the inversion based on seven Fe i lines.

Based on these tests, we conclude that there is no appreciable gain, apart from decreasing observational noise, in using many spectral lines for chemical spot DI of Ap/Bp stars. Abundance mapping with a single line, even a blended one, is not significantly inferior to an inversion based on five or more lines. This is true, provided that different blend contributors are identified correctly and modelled accurately.

3.4.4. Ignoring a moderately strong magnetic field

As mentioned above, a significant fraction of previous abundance DI studies of Ap/Bp stars have ignored the magnetic field in chemical spot mapping. This potentially leads to systematic errors and artefacts in the reconstructed maps. The goal of this section is to quantify these effects and compare possible errors with typical horizontal abundance contrasts.

The influence of the magnetic field on Stokes I line profiles and therefore its possible impact on chemical abundance DI strongly depends on the magnitude of the field and on the sensitivity of a given spectral line to the Zeeman effect. A line with a small or null effective Landé can be reliably studied with non-magnetic spectral modelling methods even in strong-field targets. On the other hand, typical spectral lines with effective Landé factors 0.9–1.5, such as the group of 7 Fe i lines studied here, exhibit a noticeable magnetic local profile modification for the field strength of 2 kG. However, these profile shape distortions are not observable for moderately or rapidly rotating Ap/Bp stars typically targeted by DI. Instead, the primary manifestation of the Zeeman effect is a modification of the local equivalent widths of spectral lines through the magnetic intensification mechanism. This effect changes across the stellar surface according to the local field strength and orientation. This magnetic line strengthening might be misinterpreted as chemical inhomogeneities by a non-magnetic DI code.

It is essential to adopt a representative magnetic field strength and geometry for an assessment of systematic errors of an inversion ignoring magnetic field. To this end, the 2.5 kG polar field strength of the dipolar topology assumed in our experiments corresponds to the average field strength of a volume limited sample of nearby Ap/Bp stars (Power et al. 2008). On the other hand, the average dipolar field strength for a sample of non-magnetic abundance DI targets discussed below (see Sect. 4) is 2.0 kG.

To study the effects of neglecting magnetic field in chemical abundance mapping we have repeated the inversions discussed in Sects. 3.4.1 and 3.4.2 using the same simulated data but assuming zero field instead of the correct 2.5 kG dipolar field topology. Furthermore, since the magnetic intensification of individual spectral lines leads to a line-to-line scatter of the mean equivalent width, which in real DI applications is difficult to separate from other systematic errors (for instance, erroneous oscillator strengths), we have enabled an automatic oscillator strength correction for all Fe i lines from the list in Table 1 except Fe i 5005.7 and 5006.1 Å.

The results of abundance inversions ignoring the magnetic field are illustrated in Fig. 7 for the multiple circular spot distribution and in Fig. 8 for the chemical overabundance ring. In both cases we detect a systematic offset on the order of 0.3 dex between the true and reconstructed maps. This offset is subtracted before further analysis of the difference map. Not surprisingly, Figs. 7 and 8 show that inversion errors increase when the magnetic field is ignored in abundance mapping. The contrast of reconstructed abundance maps is somewhat larger than for the case of correct magnetic field treatment. The maps also show a systematic latitude dependence of reconstruction errors: the relative abundance is underestimated for higher latitudes but overestimated for lower latitudes. The fit to simulated observations is also not as good as before, with the standard deviation at convergence larger by some 20% in the ring-structure inversion compared to the test in Sect. 3.4.2. These issues are reflected in the average errors of 0.19–0.22 dex and maximum errors of 0.28–0.33 dex. These errors are about a factor of two larger than those found for the optimal inversions in Sects. 3.4.1 and 3.4.2 but still not overwhelming compared to the 1.5 dex input abundance contrast and certainly do not preclude a correct identification of the main chemical overabundance features at the stellar surface.

We have carried out an additional inversion to assess a spurious component of the abundance DI maps resulting from neglecting magnetic field. In this test simulated observations were generated for a uniform Fe distribution with log NFe/Ntot = −4.0 and the same 2.5 kG dipolar field as before. Then an abundance inversion was carried out assuming zero magnetic field. In this case all structures in the resulting chemical spot map, shown in Fig. 9, are spurious. Quantitative assessment indicates an overestimation of the mean abundance by 0.34 dex, the average errors of 0.13 dex, and the maximum errors of 0.18–0.19 dex. These errors are smaller than in the inversions with inhomogeneous surface abundance distributions. The structures in the difference map in Fig. 9 are also far less conspicuous than in Figs. 7c and 8c. This suggests that the primary effect of an unaccounted magnetic field on abundance DI is to distort a real chemical abundance distribution rather than to produce a completely spurious abundance structure. Thus, the non-magnetic abundance inversion errors include two contributions: spurious abundance structures mimicking the effect of local magnetic line intensification (~0.1 dex for the 2.5 kG dipolar field considered here) and abundance map distortions due to an incorrect recovery of chemical inhomogeneities present on the stellar surface. The latter contribution scales with the surface element distribution contrast, reaching ~0.2 dex for the 1.5 dex spot to photosphere abundance contrast considered in our experiments.

3.4.5. Effect of local atmospheric structure

Apart from the Si and Fe abundance inversions discussed by Lehmann et al. (2007) for HR 7224 and by Kochukhov et al. (2012) for α2 CVn, respectively, and mapping of He spots in a few hot CP stars (Yakunin et al. 2015; Oksala et al. 2015), no systematic effort was made by previous DI studies to achieve a self-consistency between chemical spot inversions and lateral variation of the atmospheric structure and continuum brightness. Since the metal overabundance spots are brighter in the optical than the surrounding atmosphere (e.g. Lüftinger et al. 2010a; Krtička et al. 2015), an abundance DI based on a single, average stellar model atmosphere may overestimate the contrast of surface abundance distribution and lead to other artefacts. However, the studies by Lehmann et al. (2007) and Kochukhov et al. (2012) did not reveal any major effect of this problem on DI maps.

In this section we further explore and characterise the systematic DI errors and artefacts related to neglecting an impact of chemical abundance spots on the local atmospheric structure. To this end, we have generated simulated observational data with the same characteristics and for the same magnetic field and abundance geometries as discussed above using forward calculations with the Invers13 code and a grid of LLmodels atmospheres covering the full range of Fe abundances present in the input maps. The local line and continuum spectra were thus calculated fully self-consistently with the adopted local Fe abundance. We then analysed these simulated observations with the Invers10 code, using a single model atmosphere calculated for log NFe/Ntot = −4.0. The resulting comparison between the true and recovered abundance maps is presented in Figs. 10 and 11, for the multiple spots and overabundance ring surface distribution scenarios respectively.

According to Table 2, the average abundance reconstruction errors stemming from an inconsistent treatment of chemical spots and local atmospheric structure are about 0.2 dex while the maximum errors reach 0.3 dex. On the other hand, no significant mean abundance offset is detected in the reconstructed maps. As is evident from Figs. 10 and 11, the spot maps obtained with a single model atmosphere exhibit a somewhat higher abundance contrast, with spurious underabundance zones found close to or in between the true overabundance features. In addition, one can find a couple of small areas in the difference map for the multiple circular spot test case (Fig. 10c) where the Fe abundance appears to be significantly overestimated. These artefacts correspond to a latitudinal smearing by 10–20° of real overabundance spots.

In general, results of the experiments presented in this section lead to the conclusion that the atmospheric structure variation due to the 1.5 dex Fe abundance contrast adopted in our tests has a marginal impact on the classical abundance DI. On the other hand, abundance anomalies of elements other than Fe are less important for the atmospheric structure of A- and B-type stars (Khan & Shulyak 2007) and hence their impact on DI maps is even smaller than discussed here. These considerations suggest that a self-consistent treatment of metal overabundance spots and local atmospheric structure is not a critical factor determining reliability of CP star abundance DI maps.

4. Discussion

Table 3

Doppler imaging studies of chemically peculiar stars.

In the preceding sections we have examined an impact of several systematic error sources, such as neglecting the magnetic field, on the Doppler imaging reconstruction of chemical spot maps in early-type stars. What are the implications of these results for the interpretation of published chemical abundance distributions? To address this question, we have compiled a summary of all chemical abundance DI studies which were based on reasonably modern, high-quality observations (i.e. high S/N spectra recorded with solid-state detectors) published between 1989 and 2016. This summary, presented in Table 3, lists 36 magnetic Ap/Bp stars investigated by 39 individual DI and MDI studies. A few DI publications superseded by later studies of the same chemical elements in the same stars were omitted.

For each star, Table 3 provides an estimate of the peak surface magnetic field strength, Bmax. This parameter was taken from the MDI magnetic field maps whenever available or adopted according to the dipolar field strength reported in the literature. For a few stars the existing magnetic measurements are insufficient to constrain even the simplest dipolar field model. In these cases we adopted 3 | ⟨ Bz ⟩ | max as an approximation of Bmax.

For completeness, Table 3 also includes information on four HgMn stars investigated with abundance DI. These stars are qualitatively different in their spectral appearance and characteristics of non-uniform element distributions (low-contrast and evolving) from magnetic Ap/Bp stars. Although extremely weak fields of the type recently found in Vega and Sirius (Petit et al. 2010, 2011) cannot be excluded for HgMn stars, these objects certainly lack the organised global fields typical of classical magnetic CP stars (Makaganiuk et al. 2011b; Kochukhov et al. 2013a).

thumbnail Fig. 12

Distribution of peak surface magnetic field strengths for DI (filled histogram bins) and MDI (hatched histogram bins) targets.

Open with DEXTER

The distribution of peak surface magnetic field strength for 36 magnetic Ap/Bp stars with available chemical spot maps is shown in Fig. 12. The MDI targets span a wide range of Bmax from ~1 to >30 kG. On the other hand, the stars studied with non-magnetic DI techniques are concentrated at lower field strengths. Their average Bmax is equal to 2.0 kG and exceeds 3.0 kG for only three out of 25 stars. HD 3980 with its 6.9 kG dipolar field (Nesvacil et al. 2012) is a notable outlier. However, the field strength and topology of this star were constrained using only five mean longitudinal magnetic field measurements and thus are by no means definitive. The presence of a strong field in HD 3980 needs to be confirmed with high-resolution spectropolarimetric observations.

Given the results of our numerical tests in Sect. 3.4.4, which demonstrated only minor effects of ignoring a 2.5 kG dipolar field, there is no evidence that a systematic neglect of the magnetic field is detrimental for the quality of chemical maps for the majority of DI targets with fields of a few kG or less. In fact, only the non-magnetic DI analysis of HD 3980 by Nesvacil et al. (2012) may potentially suffer from non-negligible artefacts for some elements provided that the strong magnetic field of this star is real. One can ascertain, though, that the majority of chemical abundance maps of HD 3980 cannot be dominated by spurious effects of an unaccounted Zeeman line broadening and intensification since the derived element distributions are diverse and do not follow a clear common pattern. In any case, quantifying the probable inversion errors requires a separate element-by-element assessment that takes into account the magnetic field response of specific spectral lines employed for the DI mapping of this star.

5. Summary and conclusions

In this study we presented a series of numerical experiments addressing the question of reliability of chemical abundance maps derived for early-type magnetic stars with the Doppler imaging technique. In particular, we focus on assessing the systematic errors of DI reconstructions coming from common simplifying approximations, such as neglect of a magnetic field in radiative transfer and ignoring horizontal variation of the stellar atmosphere structure. These issues were investigated for a set of atmospheric parameters corresponding to a cool Ap star. At the same time, none of the studied systematic effects is expected to have a significant temperature dependence, suggesting that our conclusions are also valid for hotter Ap/Bp stars.

The main results of our investigation can be summarised as follows.

  • Our DI code achieves an average accuracy of0.06–0.09 dex and maximum errors of0.12 dex for the chemical map reconstructionwhich includes magnetic field effects in spectrum synthesis and isbased on about half a dozen of individual spectral lines.

  • DI reconstruction is equally successful for a complex element distribution comprised of multiple circular element overabundance spots and for a narrow overabundance ring located at the magnetic equator of a dipolar field.

  • A reduction of the line list to one or two spectral features, including a blended one, leads to only a marginal increase of the average inversion errors to about 0.1 dex and maximum errors to ~0.15 dex. Therefore, apart from counteracting random observational noise, abundance DI inversions do not gain significantly from modelling a large number of spectral lines.

  • Ignoring a moderately strong, 2.5 kG dipolar magnetic field in abundance DI introduces a mean offset of about 0.3 dex in the recovered chemical maps. The average relative reconstruction errors increase to ~0.2 dex while the maximum relative errors reach ~0.3 dex. These errors correspond to distortions of real surface abundance inhomogeneities. The spurious component of abundance maps (i.e. Zeeman broadening and intensification misinterpreted as abundance variations) does not exceed ~0.15 dex on average.

  • Ignoring local atmospheric structure variations in the areas of element overabundance leads to average reconstruction errors of ~0.2 dex and maximum errors of ~0.3 dex. These numbers correspond to the effect of overabundance spots of the element (iron) providing the most important contribution to the atmospheric opacity. The largest artefacts resulting from using a single average model atmosphere for DI inversions correspond to the distortion of real chemical inhomogeneities; no spurious spots are produced.

  • Reviewing results of several dozen non-magnetic DI studies in the light of our findings, we conclude that the absolute majority of these investigations have targeted stars with sufficiently weak fields and therefore could not have been adversely affected by the neglect of these fields in abundance mapping. In general, the systematic DI reconstruction errors (up to 0.2–0.3 dex) inferred in this paper are small or entirely negligible when compared to the horizontal element abundance contrasts of 2–5 dex typically derived with DI for Ap/Bp stars.


1

The term “Doppler imaging” was suggested later by Vogt & Penrod (1983).

Acknowledgments

The author thanks Drs. T. Lüftinger and J. Silvester for comments on this paper. This research is supported by the Knut and Alice Wallenberg Foundation, the Swedish Research Council, and the Swedish National Space Board.

References

All Tables

Table 1

Spectral lines used in numerical tests of abundance DI.

Table 2

Numerical tests of abundance DI inversions.

Table 3

Doppler imaging studies of chemically peculiar stars.

All Figures

thumbnail Fig. 1

Surface magnetic field topology adopted for numerical tests of DI. The star is shown at five rotation phases, which are indicated above the spherical plots, and an inclination angle of i = 60°. The spherical plots show the maps of a) field modulus; b) horizontal field; c) radial field; and d) field orientation. The contours over spherical maps are plotted with a 0.5 kG step. The thick line and the vertical bar indicate positions of the rotational equator and the pole, respectively. The colour bars give the field strength in kG. The two different colours in the field orientation map correspond to the field vectors directed outwards (red) and inwards (blue).

Open with DEXTER
In the text
thumbnail Fig. 2

Illustration of the regularisation parameter determination for the inversion with two Fe ii lines. A measure of χ2 (solid black line) and the regularisation functional (solid red line) are plotted versus the Tikhonov regularisation parameter. The dotted line shows the noise level of the simulated data set. The optimal regularisation parameter is indicated with an arrow.

Open with DEXTER
In the text
thumbnail Fig. 3

Results of the Fe abundance inversion for the reference test case (4 circular spots, 7 Fe i lines, including magnetic field). The star is shown at five rotation phases, which are indicated above the spherical plots, and an inclination angle of i = 60°. The spherical plots show a) the true abundance map; b) reconstructed abundance map; c) difference map. The contours over the difference map are plotted at ±0.3,0.6,0.9 dex. The thick line and the vertical bar indicate positions of the rotational equator and the pole, respectively. The colour bars give the local abundance in log NFe/Ntot units.

Open with DEXTER
In the text
thumbnail Fig. 4

As Fig. 3 but for the reconstruction of the overabundance ring.

Open with DEXTER
In the text
thumbnail Fig. 5

As Fig. 3 but using two Fe ii lines.

Open with DEXTER
In the text
thumbnail Fig. 6

As Fig. 3 but for the simultaneous reconstruction of Cr (panels a)c)) and Fe (panels d)f)) abundance maps using two Fe ii lines and one blended line of Cr ii.

Open with DEXTER
In the text
thumbnail Fig. 7

As Fig. 3 but ignoring magnetic field. A constant offset is subtracted from the difference map shown in panel c).

Open with DEXTER
In the text
thumbnail Fig. 8

As Fig. 4 but ignoring magnetic field. A constant offset is subtracted from the difference map shown in panel c).

Open with DEXTER
In the text
thumbnail Fig. 9

Spurious chemical abundance distribution resulting from chemical mapping ignoring magnetic field. The true abundance map is a homogeneous distribution with log NFe/Ntot = −4.0. Reconstructed Fe abundance map is shown in the top panel adopting the same abundance range as in Fig. 7. The bottom panel shows the same distribution after subtracting a constant offset.

Open with DEXTER
In the text
thumbnail Fig. 10

As Fig. 3 but ignoring effect of local abundance variation on the model atmosphere structure.

Open with DEXTER
In the text
thumbnail Fig. 11

As Fig. 4 but ignoring effect of local abundance variation on the model atmosphere structure.

Open with DEXTER
In the text
thumbnail Fig. 12

Distribution of peak surface magnetic field strengths for DI (filled histogram bins) and MDI (hatched histogram bins) targets.

Open with DEXTER
In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.