© ESO 2020

1. Introduction

Magnetic fields play a central role in the surface activity of cool stars. It is now well established that magnetism is responsible for such phenomena as dark spots, flares, coronal mass ejections, enhanced chromospheric, and X-ray emission. Magnetic fields directly affect stellar evolution by altering the mass loss and governing redistribution of angular momentum between different parts of stellar interiors. Planets orbiting cool stars are influenced by the stellar magnetic activity in many different ways. This makes understanding stellar magnetism essential for studying evolution, atmospheres, and the habitability of terrestrial exoplanets.

Despite the availability of a massive body of circumstantial observations of stellar magnetic activity, direct detections and measurements of magnetic fields on stellar surfaces are still very challenging. This type of research relies on exploiting the signatures of the Zeeman effect in stellar spectra, which requires high-resolution, high signal-to-noise ratio spectroscopic and spectropolarimetric observational data. Two complementary approaches to studying the Zeeman effect are commonly used to infer the presence of a magnetic field and derive its characteristics. The first method relies on line polarisation measurements with high-resolution spectropolarimetry. The second technique extracts information on the magnetic broadening and splitting of spectral lines from the usual intensity spectra.

Analyses of weak circular polarisation signals in spectral lines, which are often enhanced with a multi-line technique, have been very successful in studying cool-star magnetic fields (Donati et al. 1997). The polarimetric method has yielded magnetic field detections for hundreds of stars using snapshot circular polarisation observations (e.g. Fossati et al. 2013; Marsden et al. 2014; Moutou et al. 2017). It has also enabled the reconstruction of detailed magnetic field maps for dozens of objects (e.g. Donati et al. 2003; Petit et al. 2008; Morin et al. 2008; Rosén et al. 2016; Folsom et al. 2016) with the Zeeman Doppler imaging (ZDI, Kochukhov 2016) inversion technique applied to time-resolved spectropolarimetry. The success of polarimetric diagnostic methods stems from an unambiguous nature of the magnetic field detection and a relative simplicity of the theoretical modelling of weak stellar circular polarisation signatures, which can be carried out relying on a very basic line formation treatment (Folsom et al. 2018).

However, it is also understood that polarimetry is only able to recover a small fraction of the magnetic field energy, vastly underestimating the true strength of magnetic structures on the surfaces of cool stars (Vidotto 2016; Kochukhov et al. 2017; Lehmann et al. 2019). This polarimetric bias is caused by the topological complexity of a typical cool-star magnetic field geometry, which comprises many unresolved magnetic features with opposite field polarities. Circular polarisation signals corresponding to these regions have opposite signs and mostly cancel out in any disc-integrated polarimetric observable. Thus, polarimetry is only sensitive to a large-scale magnetic field component, particularly for slowly rotating stars, which exhibit no significant rotational Doppler broadening of their line profiles. At the moment, it is unclear how this large-scale field component is related to fields at smaller spatial scales where the bulk of the magnetic energy is concentrated.

On the other hand, the Zeeman splitting and broadening observed in stellar intensity spectra is proportional to the absolute value of magnetic field strength and thus includes contributions from magnetic structures at all spatial scales. In this way the Zeeman broadening diagnostic provides an unbiased estimate of the total surface magnetic field strength. The relative intensities of Zeeman components are only weakly sensitive to the field orientation, making it impossible to infer the vector field maps from intensity spectra considering typical field strengths of ∼1 kG encountered in cool stars. Moreover, a challenging aspect of this type of magnetic field analysis is that magnetic broadening has to be separated from many other broadening mechanisms affecting spectral lines. This makes the field detections from intensity spectra more ambiguous compared to the polarimetric method and often requires observational data of exceptional quality. Additionally, the Zeeman response of the intensity profiles of spectral lines is more complex and diverse than the circular polarisation in the same lines, impeding application of multi-line techniques and requiring the use of sophisticated polarised radiative transfer codes.

These problems have been overcome, with varying degree of success, by a number of pioneering studies which inferred the presence of magnetic fields in different types of cool active stars (Robinson 1980; Saar et al. 1986; Basri & Marcy 1988; Valenti et al. 1995). The Zeeman broadening method was particularly successful in application to T Tauri stars (Johns-Krull et al. 1999; Johns-Krull 2007; Yang & Johns-Krull 2011; Lavail et al. 2017, 2019) and active M dwarfs (Johns-Krull & Valenti 1996; Shulyak et al. 2014, 2017, 2019; Kochukhov & Lavail 2017; Kochukhov & Shulyak 2019), which typically have rather strong fields and narrow lines. At the same time, relatively little progress has been made for Sun-like G dwarfs (see review by Reiners 2012). The majority of Zeeman broadening field detections and measurements for these stars come from historic publications by S. Saar and collaborators (Saar 1987, 1996, 2001; Saar & Linsky 1986a,b; Saar & Baliunas 1992), with a few measurements contributed by other studies (Basri & Marcy 1988; Rüedi et al. 1997; Anderson et al. 2010). To summarise, despite numerous recent polarimetric investigations of global magnetic fields of Sun-like stars (e.g. See et al. 2019, and references therein), the properties of their overall magnetic fields, including typical surface field strengths, their rotational modulation and cyclic variation, relationship to large-scale fields and different indirect magnetic activity indicators, remain largely unexplored.

This unsatisfactory situation is largely due to the absence of an efficient Zeeman broadening diagnostic technique that can be applied to moderate quality high-resolution optical spectra of solar-type stars. In this paper we develop such a technique and present its application to a sample of G dwarf stars with different activity levels. The rest of this paper is organised as follows. Section 2 details various methodological aspects of our study, including a review of different manifestations of the Zeeman effect in stellar intensity spectra, description of our line profile modelling codes, motivation of the choice of key diagnostic lines, target selection and discussion of observational data. This is followed in Sect. 3 by the presentation of magnetic field measurement results for each target star. We discuss our results in Sect. 4, where we establish correlations between magnetic field characteristics, stellar parameters and magnetic activity proxies. Finally, Sect. 5 summarises main conclusions of our investigation.

2. Methods

2.1. Zeeman broadening and intensification of spectral lines

The presence of a magnetic field in stellar atmosphere leads to splitting of each spectral line into three groups of differently polarised Zeeman components. The linearly polarised π components are distributed symmetrically around the line centre. Elliptically polarised σ components form two groups, one shifted bluewards, and another redwards of the line centre. The consequence of this Zeeman effect on spectral lines observed in the usual intensity spectrum is twofold. First, lines are broadened (if the magnetic splitting is less than the non-magnetic line width) or split (if the field strength is large enough) due to a wavelength separation of the π and σ components. The magnitude of magnetic broadening increases linearly with the field strength and can be expressed in km s⁻¹ units as

$$\begin{array}{c}\hfill \mathrm{\Delta }{v}_B=1.4\times {10}^{-4}{g}_{\mathrm{eff}}\lambda B\end{array}$$(1)

for the field strength in kG and wavelength in Å. The effective Landé factor, g_eff, expresses the relative span of Zeeman splitting for lines with different magnetic sensitivity. Considering that the most magnetically sensitive lines one can find in stellar spectra have g_eff  ≈  3, we get Δv_B ≲ 2 km s⁻¹ for a 1 kG field and a line at λ = 5000 Å. This magnetic broadening is comparable to the intrinsic line width, dominated by the ∼2−3 km s⁻¹ turbulent broadening, and is smaller than the instrumental broadening (3−6 km s⁻¹ for the resolving power of R = λ/Δλ= 0.5−1 × 10⁵) of most of the actively used night-time spectrographs. Furthermore, any significant rotational broadening, with v_e sin i exceeding a few km s⁻¹, effectively renders Zeeman broadening unobservable in cool stars. All these factors limit practical applications of the Zeeman broadening diagnostic to active, very slowly rotating stars observed at R ≳ 10⁵ with high signal-to-noise ratio spectra (Robinson 1980; Marcy 1984; Saar et al. 1986; Basri & Marcy 1988; Johns-Krull & Valenti 1996; Rüedi et al. 1997; Anderson et al. 2010). These restrictions can be partly relaxed with the help of observations at near-infrared wavelengths (Saar & Linsky 1985; Valenti et al. 1995; Johns-Krull 2007; Yang & Johns-Krull 2011; Lavail et al. 2017, 2019; Flores et al. 2019) thanks to a faster increase of Zeeman splitting with wavelength (grows as λ²) compared to other broadening mechanisms (grow as λ). For instance, Δv_B ≈ 10 km s⁻¹ for a 1 kG field and a g_eff = 3 line at λ = 2.3 μm, enabling magnetic measurements of moderately fast rotators using lower quality data.

Another, less commonly discussed, consequence of Zeeman effect is the overall strengthening of absorption lines. This Zeeman intensification effect occurs due to a separation of Zeeman components and the resulting desaturation of strong spectral lines. This effect is well known from studies of much stronger magnetic fields encountered in the early-type chemically peculiar stars (Babcock 1949; Hensberge & De Loore 1974; Mathys 1990; Takeda 1991; Mathys & Lanz 1992; Kupka et al. 1996; Stift & Leone 2003; Kochukhov et al. 2004, 2013). It was first studied in cool stars by Basri et al. (1992) and Basri & Marcy (1994). They demonstrated that Zeeman intensification is a complex function of both magnetic field parameters (field intensity and orientation) and spectral line characteristics (line strength and Zeeman splitting pattern). Unlike Zeeman broadening, magnetic intensification is most effective for strong spectral lines with a large number of uniformly spaced Zeeman components but not necessarily for lines with the largest Landé factors. A major advantage of Zeeman intensification analysis over Zeeman broadening is that the former does not require observational data of exceptional quality and can be applied to rapid rotators provided that diagnostic lines remain free from blends. On the other hand, quantitative interpretation of spectral line strengths in terms of Zeeman intensification requires detailed modelling of the magnetic desaturation process with a realistic polarised radiative transfer code and is more sensitive to errors of atomic parameters, particularly transition probabilities, which determine relative line intensities in the absence of a magnetic field.

A series of recent studies of magnetic fields in active M dwarf stars using the Ti I multiplet at λ 9647−9788 Å implemented a combined Zeeman broadening and intensification approach (Shulyak et al. 2017, 2019; Kochukhov & Lavail 2017; Kochukhov & Shulyak 2019). These investigations demonstrated that a detection and analysis of ≳2 kG fields in stars rotating as fast as v_e sin i = 30−40 km s⁻¹ is within reach. However, this methodology cannot be directly applied to Sun-like stars since this particular group of Ti I lines becomes too weak at T_eff ≳ 4500 K. Here we aim to develop an equivalent magnetic field measurement procedure for hotter stars based on a different set of diagnostic lines.

2.2. Magnetic spectrum synthesis

We used the polarised spectrum synthesis code SYNMAST (Kochukhov 2007; Kochukhov et al. 2010) to model magnetic field effects on absorption lines in the spectra of active stars. This code solves the polarised radiative transfer equation with an efficient numerical algorithm (de la Cruz Rodríguez & Piskunov 2013) using realistic stellar model atmospheres. It treats ionisation of atomic species and dissociation of molecules with an up-to-date equation of state package shared with the SME code (Piskunov & Valenti 2017)¹. Theoretical spectra in four Stokes parameters are computed by SYNMAST for a given limb angle and a depth-independent magnetic field vector. Information on atomic and molecular line parameters is obtained from the VALD database (Ryabchikova et al. 2015)². In addition to the usual set of line parameters required for spectrum synthesis (central wavelength, excitation potential, oscillator strength and damping constants), VALD supplies Landé factors and J quantum numbers of the upper and lower atomic levels necessary for calculation of Zeeman splitting patterns.

A magnetic broadening and intensification analysis of cool active stars does not require a detailed geometrical model of the surface magnetic field distribution. We follow previous Zeeman broadening studies (e.g. Valenti et al. 1995; Johns-Krull 2007; Yang & Johns-Krull 2011; Shulyak et al. 2017; Kochukhov & Lavail 2017; Lavail et al. 2019) in assuming that the field is uniform and oriented normally with respect to the stellar surface. The assumption of a radial field orientation is justified by the analogy with solar flux tubes (e.g. Valenti et al. 1995; Johns-Krull et al. 1999, 2004) and represents an intermediate case, in terms of the impact on line profiles, between (clearly unrealistic in the stellar case) extremes of the magnetic field vectors strictly parallel and strictly perpendicular to the observer’s line of sight. In this way, a uniform radial field geometry provides a good mixture of field lines with a range of orientations to the line of sight (Yang & Johns-Krull 2011).

Using this simple field geometry model, the local intensity spectra are computed with SYNMAST at seven limb angles and then convolved with appropriate kernels to take into account the radial-tangential macroturbulence, rotational Doppler broadening, and instrumental broadening. The final theoretical stellar flux spectrum is obtained by adding these local calculations with the weights corresponding to relative areas of the seven annular regions and normalising by the continuum flux integrated over the visible stellar hemisphere in the same way. Further details on this disc integration procedure can be found in Valenti & Piskunov (1996).

We treat the field strength distribution with a standard two-component model adopted by many previous Zeeman broadening studies of active FGK stars (Basri & Marcy 1988; Rüedi et al. 1997; Saar 2001; Anderson et al. 2010; Kochukhov & Lavail 2017). This approach is inspired by the notion that small-scale magnetic fields are concentrated in distinct surface elements, referred to as flux tubes in solar physics (e.g. Stenflo 1973, 1994; Solanki & Stenflo 1984; Solanki et al. 2006). In this case, the total observed stellar spectrum is given by the weighted superposition of magnetic and non-magnetic contributions

$$\begin{array}{c}\hfill S(\lambda )=(1-f)·S_0(\lambda )+f·S(\lambda ,B),\end{array}$$(2)

where S₀(λ) and S(λ, B) is the non-magnetic, continuum-normalised stellar flux spectrum and the spectrum calculated with the field strength B, respectively, and f is the filling factor of magnetic regions. This is undoubtedly a highly simplified approximation of the actual continuous distribution of magnetic field strengths. Nevertheless, it has proven to be a successful practical approach to the problems of diagnosing small-scale solar magnetic fields (Stenflo 1994 and references therein) and inferring mean magnetic field strength from the stellar intensity spectra. Information on the magnetic filling factors determined within this framework (provided that f can be reliably separated from B) is important for understanding a range of processes taking place in active cool stars (e.g. Montesinos & Jordan 1993; Cranmer & Saar 2011; Cranmer 2017; See et al. 2019).

The solar small-scale strong-field regions, such as flux tubes, are known to have a distinctly different thermodynamic structure. However, there are no reliable analytical models or numerical simulations that quantify this difference for active stars with much stronger mean magnetic fields than those observed in the quiet Sun. In the absence of suitable grids of one-dimensional models of magnetised cool star atmospheres we resort to employing the same normal model atmosphere for computing both S₀(λ) and S(λ, B). This approach is also obligatory to enable a meaningful comparison with the results of previous studies, most of which did not take a difference between the structures of magnetic and non-magnetic atmospheres into account. Nevertheless, expecting that this difference is primarily reflected in unequal temperatures of the regions with different magnetic fields, we explore below the impact of adopting different normal model atmospheres for calculation of S₀(λ) and S(λ, B). Since our approach relies on the differential magnetic intensification of atomic lines from the same multiplet, the main effect of a temperature difference is unequal continuum brightness of the spectra corresponding to magnetic and non-magnetic regions. This difference is subsumed, to first order, by the magnetic filling factor f. This parameter thus represents continuum-intensity-weighted fraction of the stellar surface covered by a magnetic field.

2.3. Zeeman-sensitive lines in the optical solar spectrum

Several neutral Fe lines (Fe I 5250, 6173, 6302, 8468 Å) are commonly used for the analysis of polarisation and Zeeman broadening in the atmosphere of the Sun and solar-type stars. These lines are distinguished by large values of effective Landé factors and thus exhibit the largest profile shape modification when a magnetic field is present. However, as discussed in Sect. 2.1, these lines are not necessarily the most useful diagnostics when Zeeman intensification is considered. To assess the latter in a comprehensive and systematic manner, we have carried out magnetic spectrum synthesis calculations with SYNMAST for the entire optical spectrum (400 nm to 1 μm) covered by modern echelle spectropolarimeters. These calculations were based on the solar model atmosphere from the MARCS grid (Gustafsson et al. 2008), employed a line list retrieved from VALD and adopted the solar chemical abundances (Grevesse et al. 2007), microturblent velocity v_mic = 0.85 km s⁻¹, macrotubulent velocity v_mac = 3 km s⁻¹, projected rotational velocity v_e sin i = 5 km s⁻¹, and instrumental resolution R = λ/Δλ = 10⁵. Two theoretical calculations were produced: one without a magnetic field and another one with a B = 1 kG radial field covering the entire stellar surface. This choice of v_e sin i and B corresponds to a moderately active star rotating significantly faster than the Sun. This parameter combination is in the regime where a Zeeman broadening analysis is already quite challenging since the signatures of magnetic line broadening are largely washed out by the stellar rotation.

The resulting magnetic and non-magnetic synthetic spectra and their difference are shown in Fig. 1. The spikes in the difference plot identify spectral features exhibiting the strongest intensification for the 1 kG magnetic field considered in this calculation. Setting an arbitrary threshold of 5% for the change of the central line depth, we compiled the list of eight unblended lines potentially useful for a Zeeman intensification analysis. These lines, all belonging to Fe I, are listed in Table 1. The well-known magnetic diagnostic lines Fe I 5250 and 8468 Å are among the lines with the largest magnetic intensification. However, by far the strongest Zeeman response is found for the Fe I 5497.5 Å line. It has a moderately large, though not exceptional, effective Landé factor of 2.25. Interestingly, the nearby Fe I 5506.8 Å (g_eff = 2.00) is also present in Table 1 and another line from the same multiplet, Fe I 5501.5 Å (g_eff = 1.87), exhibits a weaker but still significant Zeeman intensification response.

Fig. 1. Response of the continuum normalised spectrum in the 4000−10 000 Å wavelength range to a 1 kG magnetic field covering the entire stellar surface. Lower panel: synthetic spectra calculated for the solar atmospheric parameters and abundances with (light curve) and without (dark curve) magnetic field. Upper panel: difference between these two spectra. The lines showing the strongest Zeeman response are identified according to the numbering adopted in Table 1.

Stenflo et al. (1984) discussed the solar Stokes I and V spectra of the Fe I 5497.5−5506.8 Å lines in a strong plage, noting that these lines are substantially more polarised than Fe I 5250 Å. Rosén & Kochukhov (2012) used these three Fe I lines for numerical tests of ZDI, finding that the Stokes I Zeeman intensification signal in the 5497.5 Å line aids reconstruction of stellar magnetic field geometries provided that the intensity and circular polarisation spectra are modelled self-consistently. Morgenthaler et al. (2012) correlated the widths of the 5497.5 and 5506.8 Å lines with chromospheric emission indicators for the active Sun-like star ξ Boo A (HD 131156A). Apart from these few studies, to the best of our knowledge, the Fe I 5497.5 and 5506.8 Å lines have not been systematically utilised for either solar or stellar magnetic field diagnostic.

2.4. Determination of magnetic field parameters

The magnetic intensification of the three Fe I lines 5497.5, 5501.6, and 5506.8 Å represents a promising tool for measuring surface magnetic field strength in Sun-like stars, provided one can find a suitable reference spectral feature with weak or no magnetic field sensitivity. According to Nave et al. (1994), these three lines belong to the Fe I multiplet 88, also known as multiplet 15 in the older tables by Moore (1959), formed by transitions between the a ⁵F and z ⁵F^o energy levels in a neutral iron atom. This multiplet includes a handful of other strong unblended Fe I lines, one of which, 5434.5 Å, has a very small effective Landé factor and is therefore essentially insensitive to a magnetic field. Thus, a suitable method of extracting information on stellar magnetic fields from the Fe I 5497.5−5506.8 Å lines is to compare them with Fe I 5434.5 Å. The latter line can be employed to establish the stellar projected rotational velocity and Fe abundance. Then, the three magnetically sensitive lines can be used to determine magnetic field parameters. This two-step procedure comprises the new magnetic diagnostic method advanced in this paper.

Parameters of the four Fe I lines are summarised in Table 2. The oscillator strengths listed in this table come from high-precision laboratory measurements (Fuhr et al. 1988; O’Brian et al. 1991). Other sources of atomic data may provide oscillator strengths with a different overall scale, yet the relative strengths of the four lines are going to be identical since all these transitions come from the same multiplet. This alleviates the problem of oscillator strength uncertainties which complicated many previous attempts to measure stellar magnetic fields with Zeeman intensification and broadening (Basri et al. 1992; Valenti et al. 1995; Rüedi et al. 1997; Anderson et al. 2010).

Figure 2 illustrates response of the four Fe I diagnostic lines to a uniform radial magnetic field increasing in strength from 0 to 4 kG. These calculations were carried out with SYNMAST for the same set of stellar parameters as used in Sect. 2.3 and treating each of the four lines in isolation, ignoring possible blends. Two sets of calculations are shown: one without any broadening applied to the theoretical spectra and another one using representative values of v_mac, v_e sin i and R adopted in Sect. 2.3. The Zeeman splitting patterns are schematically shown in Fig. 2 for the 4 kG field. It is evident that the 5497.5−5506.8 Å lines are strongly influenced by the magnetic field compared to Fe I 5434.5 Å. The largest magnetic broadening and intensification effect is shown by Fe I 5497.5 Å and is due to its wide Zeeman splitting pattern composed of five groups of overlapping π and σ components.

Similarity of the oscillator strengths, wavelengths and excitation potentials of the four Fe I lines studied here translates into their similar formation physics in a non-magnetic cool-star atmosphere. These lines exhibit little differential response to variation of thermodynamic structure, allowing one to disentangle magnetic intensification from other effects. A series of spectrum synthesis calculations documented by Fig. 3 compares the relative change of equivalent widths of the 5497.5–5506.8 Å lines with respect to the 5434.5 Å line caused by variation of magnetic field and stellar parameters. The reference model parameters adopted for these calculations were T_eff  =  5750 K, log g  =  4.5, v_mic  =  0.85 km s⁻¹, and B  =  0 kG. The relative normalised equivalent width change was computed for each magnetically sensitive line by dividing its equivalent width change ΔW by the initial equivalent width W₀ and subtracting the same ratio for the Fe I 5434.5 Å line. The purpose of these line formation calculations was to explore to what extent inevitable uncertainties of stellar parameters can interfere with the determination of magnetic field strength using our method. The three upper panels in Fig. 3 suggest that an uncertainty of ∼50 K for T_eff, ∼0.05 dex for log g, and ∼0.1 km s⁻¹ for v_mic, typical of modern spectroscopic analyses of solar-type stars (Valenti & Fischer 2005), have little impact (ΔW/W₀ − ΔW/W₀(5434.5)≲0.01) on the magnetic measurements unless B is much smaller than ∼200 G. The studied lines show no significant differential response to a moderate variation of Fe abundance. One can also note that the differential equivalent width change due to temperature variation (top panel in Fig. 3) is small and has a positive slope. This means that large cool spots, known to be present at the surfaces of the most active stars in our sample, should yield a relative equivalent width change opposite to that of the magnetic intensification effect.

Fig. 2. Theoretical profiles of the four Fe I spectral lines studied in this paper. Calculations for the solar parameters and magnetic field strengths ranging from 0 to 4 kG are shown for representative rotational, macroturbulent and instrumental broadening (solid lines) and without broadening (dotted lines). The bar plots above line profiles schematically show the Zeeman splitting patterns for a field strength of 4 kG.

Fig. 3. Equivalent width change of the Fe I lines λ 5497.5, 5501.5, and 5506.8 Å relative to the Fe Iλ 5435.5 Å line in response to the variation (from top to bottom) of the effective temperature, surface gravity, microturbulent velocity, Fe abundance, and magnetic field strength.

Among the three nuisance spectroscopic parameters considered here, the microturbulent velocity is characterised by the largest relative error. The choice of v_mic can influence the strengths of the 5497.5−5506.8 Å lines relative to the 5434.5 Å line because the latter is about 40% stronger than any of the former and is more sensitive to v_mic. Throughout this paper we follow Valenti & Fischer (2005) and Brewer et al. (2016) in using the same fixed microturbulent velocity, v_mic = 0.85 km s⁻¹, for all G-type targets. This is a reasonable assumption for weakly or moderately active stars within a narrow parameter range around the solar values. One can suspect that this assumption does not hold for stars significantly more active than the Sun due to a modification of their convective turbulent spectrum by a magnetic field. However, such targets also exhibit a larger Zeeman intensification due to stronger fields, making the v_mic uncertainty less of a concern.

The bottom panel in Fig. 3 demonstrates that the magnetic intensification curves of the three Fe I lines follow theoretically expected ∝B² dependence (Landi Degl’Innocenti & Landolfi 2004) up to 200−300 G and then behave linearly with B. It may be more problematic to recognise the presence of the field using these Fe I lines in the quadratic regime, corresponding to the hypothetical situation of ⟨B ⟩≲300 G and f ≈ 1. However, we did not encounter such situations in our analysis.

The increase in the equivalent width of the Fe I 5497.5 Å line is steeper than for the other two magnetically sensitive lines. This difference of the equivalent width responses, coupled with the Zeeman broadening of the 5497.5 Å line detectable in slower rotators, enables disentangling, to some extent, the field strength B from the magnetic filling factor f. Nevertheless, as we show below, the filling factor and the field strength are still partially degenerate in our approach, yielding considerably larger individual errors of B and f compared to the uncertainty of their product, the mean field strength ⟨B⟩ = B ⋅ f.

Drawing from the forward theoretical spectrum synthesis calculations described above we proceed to the analysis of the four Fe I lines in the spectrum of the Sun. For this purpose we consider the HARPS solar flux spectrum calibrated in wavelength with the help of a laser frequency comb (Molaro et al. 2013). The atomic and molecular data were extracted from VALD, now taking into account all known absorption features around the four Fe I lines studied here.

A comparison between the observed solar spectrum and the best fit model calculation is shown in Fig. 4. For this spectrum synthesis we adopted v_mic = 0.85 km s⁻¹, v_e sin i = 1.63 km s⁻¹ (Valenti & Fischer 2005) and determined v_mac = 2.83 km s⁻¹ and log(N_Fe/N_tot) = − 4.58 from the 5434.5 Å line. Then, the fit to the three magnetically sensitive Fe I lines was optimised by changing both B and f. This was accomplished with a straightforward grid search approach, varying f from 0 to 1 with a step of 0.01 and B from 0 to 7 kG with a 0.1 kG step and looking for a B, f pair yielding the lowest chi-square. The chi-square probability statistics was employed to establish 68.3% confidence limits for B, f, and their product B ⋅ f. We found $B=2.6_{-1.2}^{+3.2}$ kG, $f=0.{07}_{-0.03}^{+0.07}$, and $⟨B⟩=0.{18}_{-0.05}^{+0.11}$ kG. The 99.7% confidence limits for ⟨B⟩ are 0.12 to 0.32 kG. This mean field strength is comparable to the average field strength of 130−220 G reported by some Zeeman and Hanle studies of the quiet Sun magnetism (Trujillo Bueno et al. 2004; Shchukina & Trujillo Bueno 2011; Danilovic et al. 2010; Orozco Suárez & Bellot Rubio 2012). The solution with f = 1 (one-component model) increases the chi-square by a factor of 1.9 relative to the best fitting two-component model and is well outside the 99.7% confidence range.

Fig. 4. Solar flux spectrum obtained with the HARPS spectrograph (symbols connected by thin solid lines) compared to best fitting synthetic profiles (solid lines) for the four diagnostic Fe I lines studied in this paper.

As part of the analysis of the solar spectrum we adjusted oscillator strengths of several weak lines adjacent to the four Fe I diagnostic features. Only one of these lines, Y II 5497.4 Å, for which the oscillator strength had to be reduced by 0.12 dex, directly affects the wing of one of the studied Fe I lines in the solar spectrum. Some other blending features contribute to the wings of these Fe I lines in the spectra of faster rotators. As can be seen from Fig. 4, we have not succeeded in reproducing the shallow diffuse absorption in the far red wing of Fe I 5501.6 Å. Several weak C₂ lines contribute to the solar spectrum at those wavelengths. It is possible that the VALD C₂ line data are incomplete or inaccurate, explaining the discrepancy between the observed solar spectrum and the model calculation. This disagreement is unimportant for the Sun due to its low v_e sin i. But this issue becomes progressively more problematic as v_e sin i increases. Consequently, we systematically excluded the far red wing of the Fe I 5501.6 Å line from the set of wavelength intervals employed for chi-square calculations.

Here we also present a detailed account of the magnetic parameter inference for the active star HD 129333 (EK Dra). For this series of calculations we considered the average spectrum of this star corresponding to the epoch 2009.02. This spectrum was obtained by co-adding four individual observations obtained over nine nights (see Sect. 2.6 for further details on the observational data employed in this study). The Fe I 5434.5 Å line profile was used to determine v_e sin i = 17.0 km s⁻¹ and log(N_Fe/N_tot) = − 4.46 adopting v_mac = 3.0 km s⁻¹. This macroturbulent broadening parameter follows from the v_mac(T_eff, log g) calibration by Doyle et al. (2014) with the solar value replaced by $v_{\text{mac}}^{\odot }=2.83$ km s⁻¹ determined above. The macroturbulent velocity was calculated in this manner for all G-type stars in this study. However, the exact choice of the macrotubulent velocity value is unimportant in the context of our analysis. For narrow-line stars considered here a change of v_mac can be compensated by modification of v_e sin i, resulting in the same chi-square of the fit to observations and identical magnetic field parameters.

After constraining the Fe abundance and v_e sin i of HD 129333 with the 5434.5 Å line we proceeded to determination of B and f using the same grid search technique as was employed for modelling the solar spectrum. The resulting chi-square surface is shown in Fig. 5. The best fitting parameters are B = 3.5 kG, f = 0.40, and B ⋅ f = 1.40 kG. The corresponding 68.3% confidence limits are 2.7−4.4 kG, 0.32−0.54, and 1.33−1.47 kG respectively. As expected, there is a significant anti-correlation between B and f along the B ⋅ f = const line. The error determination procedure adopted here fully takes this anti-correlation into account, allowing to derive realistic constraints on B and f. In this particular case, the best one-component model (f = 1, B = ⟨B⟩ = 1.6 kG) is clearly excluded as it yields a chi-square increase by a factor of 1.8. We emphasise that the 68.3% confidence limits discussed in this section and elsewhere in the paper should not be treated as one-σ error bars of the normal distribution. For example, for the 2009.02 observation of HD 129333 considered here, the 95.5% confidence limits (corresponding to the second contour in Fig. 5) are 2.6−4.5 kG, 0.31−0.56, and 1.33−1.48 kG for B, f, and ⟨B⟩ respectively.

Fig. 5. Illustration of the determination of magnetic field strength B and filling factor f for HD 129333 (2009.02 epoch). The greyscale corresponds to the χ2 of the fit to the three magnetically sensitive Fe I lines, as quantified by the side bar. The best fitting parameters, B = 3.5 kG and f = 0.4, are indicated with the cross. Solid lines show 68.3, 95.5, and 99.7% confidence limits. The dashed line corresponds to the B ⋅ f = const curve.

The two-step grid search procedure described above, with the initial determination of v_e sin i and the Fe abundance using the 5434.5 Å line followed by the measurement of B, f, and ⟨B⟩ from magnetically sensitive lines, is largely equivalent to a general chi-square optimisation using all four lines simultaneously. However, our approach is computationally faster, more straightforward, reproducible and less prone to degeneracies thanks to a clear separation of the information content of lines with different magnetic sensitivity. Anyway, we have verified that application of a general least-squares fitting algorithm to the observation of HD 129333 discussed above yields the same set of magnetic field parameters, v_e sin i, and Fe abundance as was obtained with our two-step grid search procedure.

The observed spectrum of HD 129333 and the best fitting magnetic model spectrum are displayed in Fig. 6. This plot also shows the non-magnetic theoretical calculation for the same set of stellar parameters. The magnetic intensification of the 5497.5–5506.8 Å lines is readily apparent. The equivalent width and the residual central depth of the Fe I 5497.5 Å line increase by 28% and 9%, respectively. Such an effect can be easily detected for this very active star even using moderate quality spectra. Figure 6 also demonstrates that the rotational variability of the Fe I lines in HD 129333 induced by cool spots is much smaller than the Zeeman intensification signature.

Fig. 6. Magnetic field measurement for HD 129333 (2009.02 epoch). The average observed profiles of the four Fe I lines are shown with black symbols connected by thin solid lines. The grey curve in the background corresponds to ± twice the standard deviation for each pixel of the observed spectra. The thick red solid line shows theoretical spectrum for the magnetic model yielding the lowest χ2 according to Fig. 5. The dashed blue line shows the corresponding non-magnetic spectrum.

We conclude the assessment of the new magnetic field measurement methodology with investigation of the sensitivity of the analysis results to variation of stellar parameters. For this purpose we use the same observation of HD 129333 as was discussed above. Determination of v_e sin i, log(N_Fe/N_tot) from the 5434.5 Å line, followed by measurement of the magnetic field parameters from the 5497.5−5506.8 Å lines, was repeated varying T_eff by ±100 K, log g by ±0.1 dex, and v_mic by ±0.2 km s⁻¹. These uncertainties are about a factor of two larger than the typical errors of T_eff, log g, and v_mic reported by spectroscopic studies of Sun-like stars and, in this respect, represent a conservative estimate of possible systematic errors. The outcome of this error analysis is summarised in Table 3. We found that B changes by up to 0.3 kG, f by up to 0.03, and ⟨B⟩ by up to 0.11 kG in response to the variation of stellar parameters. In all cases these changes are compatible with the formal error bars of the magnetic field parameters obtained using the reference set of T_eff, log g, and v_mic. This suggests that our magnetic field measurements are not strongly affected by the uncertainties of stellar parameters adopted from the literature.

Finally, we study the impact of neglecting the multi-component nature of active star atmospheres on our magnetic field analysis. The work by Järvinen et al. (2018) demonstrated that HD 129333 has several cool spots with 200−1000 K temperature contrast occupying 14% of the stellar surface. We therefore repeated determination of the magnetic field parameters for the 2009.02 observation of HD 129333 assuming that 20% of the star is 500 K cooler than the rest of the surface. In this test we assumed that both hot and cool atmospheric components have the same distribution of small-scale magnetic field strengths, that is both have the same B and f. Table 3 shows that adopting this multi-component model has a very minor impact. Both B and f remain within 68% confidence limits of the reference determination whereas ⟨B⟩ is altered by just 0.02 kG.

A systematic temperature difference between magnetic and non-magnetic regions causes a more significant modification of our spectrum fitting results. To test this effect, we repeated analysis of the 2009.02 spectrum of HD 129333 assuming that (a) magnetic regions are 100 K cooler relative to the mean stellar effective temperature T_eff = 5845 K and non-magnetic regions are hotter by the same amount (T_mag = 5745 K, T₀ = 5945 K) and (b) the temperature difference is reversed (T_mag = 5945 K, T₀ = 5745 K). In both of these situations the 5434.5 Å line is affected by the choice of the filling factor f, so our two-step procedure was iterated until convergence was achieved for f and Fe abundance. The results of these tests are reported in the last two rows of Table 3. A 200 K temperature contrast modifies B, f, and ⟨B⟩ systematically by the amount comparable to formal error bars. This change is explained primarily by the difference of continuum brightness of the two spectral contributions, as predicted in Sect. 2.2. If magnetic regions are cooler, their contribution to the total spectrum is diminished (case a) and a larger ⟨B⟩ is required to fit the observations. The opposite is happening in the case (b). In fact, one can analytically predict ⟨B⟩ correction factors of 1.08 and 0.91 for the cases (a) and (b), respectively, considering the continuum brightness at λ = 5500 Å of the spectra corresponding to T_eff = 5745 and 5945 K. The actual change of the mean field strength according to Table 3 is a factor of 1.07 and 0.92 for the scenarios with cooler and hotter magnetic regions, respectively.

2.5. Target selection and stellar parameters

The targets for our study were selected according to the following criteria (i) fundamental stellar parameters close to solar values, (ii) information on the global magnetic field is available from previous ZDI studies, and (iii) high-quality optical spectra are available. The first of these constraints is motivated by our general goal of expanding the number of early-G dwarfs with reliable magnetic field measurements and certain limitations of our magnetic diagnostic method (see Sect. 2.4), which relies on solar calibration of the line list and turbulent velocities. Starting from the summaries of ZDI studies published by Vidotto et al. (2014) and See et al. (2019), we identified 14 dwarfs in the spectral type range from G0 to G7 which satisfy these criteria. This list includes 6 stars (HD 1835, HD 20630, HD 39587, HD 72905, HD 129333, HD 206860) from the well-studied “Sun in Time” reference sample (Ribas et al. 2005; Güdel 2007; Rosén et al. 2016; Fichtinger et al. 2017; Pognan et al. 2018) composed of solar twins at different evolutionary stages. To this G dwarf sample we added one additional very active cooler star, LQ Hya (HD 82558), to test applicability of our method to faster rotators. Except this early-K dwarf, believed to have a mass of about 0.8 M_⊙, all our targets have masses within ±10% of the solar value (see Table 1 in See et al. 2019).

Table 4 summarises relevant parameters of our targets. The first three columns list the HD number, the commonly used name, and the spectral type adopted from the HIPPARCOS input catalogue (Turon et al. 1993). This is followed by the effective temperature T_eff and surface gravity log g taken mainly from Valenti & Fischer (2005). The stellar ages and rotational periods reported in columns six and seven are adopted primarily from Vidotto et al. (2014) and See et al. (2019), respectively. The eighth column in Table 4 provides Rossby numbers calculated by dividing the rotational period by the convective turnover time. The latter was computed according to the prescription given by Cranmer & Saar (2011). The ninth column reports v_e sin i determined for each star as part of our spectrum synthesis analysis described in Sect. 2.4. The last two columns in Table 4 list two widely used proxies of the stellar magnetic activity. The ratio of the X-ray to bolometric luminosity log L_X/L_bol was adopted from Vidotto et al. (2014) and Wright et al. (2011) or calculated using the X-ray fluxes from Boller et al. (2016). The Ca II H&K chromospheric emission indicator log $R_{\text{HK}}^{\prime }$ is the median value of the measurements found in the catalogue by Boro Saikia et al. (2018).

As discussed above, we used v_mac calculated with the modified calibration of Doyle et al. (2014) and employed the same v_mic  =  0.85 km s⁻¹ for all G dwarfs. For HD 82558 we used v_mac  =  1.5 km s⁻¹ and v_mic  =  0.5 km s⁻¹ (Cole et al. 2015; Flores Soriano & Strassmeier 2017). All model atmospheres required for the line profile synthesis were extracted from the MARCS (Gustafsson et al. 2008)³ model atmosphere grid. The logarithm of surface gravity is close to 4.5 for the majority of stars in our sample. Therefore, we used a pair of MARCS models with log g = 4.5 and T_eff values bracketing the stellar effective temperature. Theoretical spectra were obtained using linear interpolation between SYNMAST calculations with these two atmospheric models. For HD 131156A (ξ Boo A), which has log g = 4.65, we used a bilinear interpolation between theoretical spectra calculated with four model atmospheres. For HD 82558 (LQ Hya) we used a single model atmosphere with T_eff = 5000 K and log g = 4.0 (Cole et al. 2015).

2.6. Observational data

High-resolution archival spectra of target stars were collected from the two sources. We used the PolarBase archive (Petit et al. 2014)⁴ to retrieve 755 observations of 15 stars obtained with the twin spectropolarimeters ESPaDOnS and Narval, installed at the 3.6 m Canada-France-Hawaii Telescope (CFHT) and the 2 m Telescope Bernard Lyot (TBL), respectively. These data were processed with an automatic reduction pipeline (Donati et al. 1997) running at the telescopes and are available from PolarBase in fully reduced format. Each observation covers the 3700−10 500 Å spectral window at the resolving power of R = 65 000. We have also used 107 spectropolarimetric observations of five stars acquired with the HARPSpol instrument mounted at the 3.6 m ESO telescope. These spectra cover the 3780−6910 Å wavelength region at the resolution of R = 110 000. The HARPSpol spectra is a mixture of public data available from the ESO archive⁵ and observations which we have collected during several recent visitor observing runs described elsewhere (Hackman et al. 2016; Rosén et al. 2016; Lehtinen et al. 2020). All HARPSpol observations were reduced with the REDUCE package (Piskunov & Valenti 2002) following the procedure detailed in Makaganiuk et al. (2012) and Rusomarov et al. (2013). Observations from all three instruments were normalised to the continuum with the method described by Rosén et al. (2018).

Most of the data sets considered here were obtained for the purpose of monitoring rotational variation and reconstructing global magnetic field topologies with ZDI. Consequently, each observing epoch is represented by anywhere between two and 44 individual observations obtained over the time span from one night to a few months. On average, there are ten observations taken over 25 days. Our first-look analysis did not reveal any differential variability of the magnetically sensitive Fe I lines relative to Fe I 5434 Å. This indicates that the small-scale magnetic fields investigated in this paper are distributed approximately uniformly over stellar surfaces. Moreover, in all cases the amplitude of rotational modulation due to cool spots was found to be significantly smaller than the differential magnetic intensification signature (e.g. see Fig. 6). This justifies co-adding all spectra obtained within the same observing run to enhance the signal-to-noise ratio and average out profile distortions caused by cool spots. This approach yielded 78 high-quality average spectra for 15 stars, with the largest number of epochs (ten) available for HD 206860 and the smallest number (two) for HD 29615 and HD 175726. Information on individual average spectra employed in our study is given in Table B.1. The first three columns of this table list the observing epoch, the number of individual observations used to calculate the average spectrum, and the facility where these data were acquired.

3. Magnetic fields of active Sun-like stars

3.1. HD 1835 (BE Cet)

This star is part of the group of young solar analogues included in the “Sun in Time” sample (Ribas et al. 2005; Güdel 2007). However, its magnetic activity was studied relatively infrequently in the past compared to other stars in that sample. Here we report three field strength measurements based on the two epochs of HARPSpol data and a pair of spectra obtained with ESPaDOnS. We infer ⟨B⟩ = 0.61−0.75 kG, with the difference between the extreme values being statistically significant. Comparison between one of the observations and the model spectra, shown in Fig. A.1, reveals a clear evidence of both Zeeman broadening and intensification effects. A single ⟨B⟩ = 0.45 kG determination can be found for this star in the literature (Saar 1987). The global magnetic field of HD 1835 was investigated by Rosén et al. (2016) based on the earliest data set analysed here.

3.2. HD 20630 (κ¹ Cet)

This object is considered to be one of the best young Sun proxies (e.g. Ribas et al. 2010; Fichtinger et al. 2017; Lynch et al. 2019) and is intensely studied in this role. We derived four field strength measurements in the range from 0.45 to 0.55 kG using observations spanning from 2012 to 2017. The difference between the extreme ⟨B⟩ determinations is not statistically significant. Figure A.2 shows an example of the fit to the four Fe I lines with and without magnetic field. Multiple previous estimates of ⟨B⟩ are available for HD 20630, with values ranging between 0.32 and 0.52 kG (Saar 1987; Saar & Baliunas 1992). These determinations are generally consistent with our results. The global field topology of HD 20630 was independently mapped by Rosén et al. (2016) and do Nascimento et al. (2016) using the circular polarisation spectra from the two earlier data sets analysed here.

3.3. HD 29615

Magnetic activity of this very young rapidly rotating (v_e sin i = 20 km s⁻¹, P_rot = 2.34 d) star was studied by Waite et al. (2015) and Hackman et al. (2016). These authors focused on mapping distribution of cool spots using high-resolution spectra and reconstructed global magnetic field topology with ZDI. Despite a significant rotational broadening, we detect an unambiguous signature of magnetic intensification (Fig. A.3). The strengthening of magnetically sensitive lines does not depend on rotational phase and is much stronger than the line profile variability caused by cool spots. We derive a mean field strength of 1.30−1.38 kG, with about 50% of the stellar surface covered by ≈ 2.7 kG field. Our analysis was based on the two sets of HARPSpol spectra collected in 2013 and 2017.

3.4. HD 39587 (χ¹ Ori)

This star is another frequently studied young solar analogue. We were able to obtain eight individual ⟨B⟩ measurements using spectra recorded in the time interval from 2007 to 2017. An example of the field detection based on the earliest data set is presented in Fig. A.4. Most of our field strength determinations cluster around ⟨B⟩ = 0.46 kG. The difference between the extreme values (0.38 and 0.50 kG) is marginally significant, hinting at a long-term trend of decreasing average photospheric magnetic field strength. All our field strength measurements appear to be weaker than the single estimate ⟨B⟩ = 0.60 kG available in the literature (Saar 1987). The global field topology of HD 39587 was studied by Rosén et al. (2016) using 4 out of 8 observing epochs considered here.

3.5. HD 56124

This star is frequently included in studies based on ZDI analyses of cool stars (e.g. Vidotto et al. 2014; See et al. 2019), although a detailed account of its magnetic mapping is yet to be published (Petit et al., in prep.). HD 56124 is the least active object in our study in terms of its log L_X/L_bol and log $R_{\text{HK}}^{\prime }$ indices (Table 4). Magnetic field effects on its line profiles are very subtle (see Fig. A.5). Nevertheless, we were able to determine a consistent and formally significant ⟨B⟩ of about 0.22 kG from the Narval spectra corresponding to four different observing epochs in 2008−2012. The mean magnetic field of HD 56124 is the weakest among the stars studied here and is formally compatible with the solar average field strength inferred with our method in Sect. 2.4. HD 56124 is also the only object other than the Sun for which our analysis yields a magnetic field filling factor below 10%.

3.6. HD 72905 (π¹ UMa)

This is another very frequently studied young solar twin with an age, rotation rate, and activity indices very similar to HD 39587. We derived five mean field strength measurements using Narval spectra taken in 2007−2016. A comparison of one of the observed spectra with the best fitting magnetic model calculation is presented in Fig. A.6. All our ⟨B⟩ determinations fall in a narrow range around 0.59 kG, with an insignificant scatter. These results show that HD 72905 possess a stronger average magnetic field than HD 39587. The same difference was also found in the ZDI analysis of the global fields of these two stars (Rosén et al. 2016). No previous Zeeman broadening estimates of magnetic field strength are available for HD 72905.

3.7. HD 73350 (V401 Hya)

This object was included in the ZDI study of four solar twin stars with different rotation rates (Petit et al. 2008). Here we analysed four epochs of Narval observations, including the data used by Petit et al. (2008). We inferred ⟨B⟩ to be in the 0.43−0.52 kG interval and found no conclusive evidence of the field strength variation from one epoch to the next. One of our field strength measurements is illustrated in Fig. A.7.

3.8. HD 76151

This star was also part of the ZDI study by Petit et al. (2008). It is an old star with a relatively low activity level and the longest rotational period (P_rot = 20.5 d) among the stars studied here. We derived five ⟨B⟩ measurements based on the spectropolarimetric observations collected with Narval in 2007−2015. One of these measurements is shown in Fig. A.8. All field strength determinations obtained for HD 76151 are consistent, within error bars, with ⟨B⟩ = 0.41 kG. This indicates that HD 76151 is a more active object compared to HD 56124 despite having a slightly longer rotational period. The tomographic mapping of the global magnetic field topology with ZDI points to a similar disparity between magnetism of these two stars (Petit et al. 2008; See et al. 2019).

3.9. HD 82558 (LQ Hya)

LQ Hya is one of the most frequently studied very active, young, rapidly rotating late-type dwarf stars. It is a popular target for Doppler mapping of the surface distributions of cool spots and global magnetic field (Donati 1999; Donati et al. 2003; Kovári et al. 2004; Cole et al. 2015; Flores Soriano & Strassmeier 2017). This star is the fastest rotator (v_e sin i = 28 km s⁻¹, P_rot  =  1.60 d) and shows the strongest X-ray and chromospheric emission in our sample. It is also somewhat cooler and less massive than the rest of the stars studied here. Due to its rapid rotation, we had to slightly change the field strength measurement methodology by excluding Fe I 5506.78 Å from the group of magnetically sensitive lines. This was motivated by the difficulty of modelling the blending of the blue wing of this line by Mn I 5505.87 Å and Fe I 5505.68 Å. Nevertheless, as illustrated by Fig. A.9, the evidence of a magnetic intensification in the two remaining magnetic diagnostic lines is unambiguous and cannot be confused with the rotational profile variations. Simultaneous fit of the three Fe I lines requires ⟨B⟩ of about 2 kG. This mean field strength was consistently obtained from three epochs of HARPSpol observations and one data set obtained with ESPaDOnS. The same four epochs of spectropolarimetric data were analysed with ZDI by Lehtinen et al. (2020), allowing us to make a direct comparison of the global and total magnetic fields for this star.

Considering previous Zeeman broadening studies of HD 82558, Saar (1996) reported a single B ⋅ f measurement of 2.45 kG using near-infrared spectra, which is in good agreement with the outcome of our study. On the other hand, Saar et al. (1992) investigated the possibility of deriving a spatially resolved distribution of magnetic field strength for HD 82558 by combining temperature DI with a magnetic intensification analysis in the optical. Their preliminary study, based on a small number of low-quality spectra, suggested a high-contrast map of B ⋅ f with the extremes at 0.1 and 2.5 kG and a surface-averaged field strength of only ≈1.0 kG. Such a large field strength variation across the stellar surface appears to contradict our observations of the magnetically sensitive lines in the spectrum of HD 82558 because these features do not exhibit any noticeable additional variability compared to the reference magnetic null line (see Fig. A.9).

3.10. HD 129333 (EK Dra)

EK Dra is the most active object among the well-established young solar twins. It is frequently studied with DI (Strassmeier & Rice 1998; Järvinen et al. 2007, 2018) and ZDI (Rosén et al. 2016; Waite et al. 2017) inversion techniques, providing a key reference point for investigations of different activity-rotation-age relationships among young Suns. Here we analysed six epochs of observing data obtained with the Narval and ESPaDOnS spectropolarimeters. Four of these data sets were previously used to produce global magnetic field maps. A detailed discussion of the derivation of ⟨B⟩, B, and f for one of the average spectra of HD 129333 was given in Sect. 2.4 and illustrated by Figs. 5 and 6. Considering results for all six epochs, we found ⟨B⟩ to be in the 1.36−1.48 kG range for this star. The difference between these extremes is not statistically significant. However, there is a trend of increasing ⟨B⟩ from the earliest observation in 2006 to the latest data set taken in 2016.

3.11. HD 131156A (ξ Boo A)

This star is a moderately active late-G dwarf in a wide binary system with an active K dwarf. Both components of this system are often used for benchmarking relations between different magnetic activity indicators (e.g. Wood & Linsky 2010; Finley et al. 2019). The global magnetic field of HD 131156A was studied in detail by Morgenthaler et al. (2012). Here we used the same six epochs as were analysed in that paper as well as two more recent Narval data sets. These observations cover the time interval from 2005 to 2015. We determined the mean field strength of 0.78 kG, with a formally insignificant scatter around this value for individual epochs. There is, however, a hint of a long-term trend in the ⟨B⟩ data as the field strength is systematically increasing from 2008 to 2010 until the latest observing epochs. An example of the fit to observations corresponding to one of the recent epochs is shown in Fig. A.10.

The separate determination of B and f, albeit rather uncertain, yields anomalous results for HD 131156A. All other stars in our sample exhibit B around 3 kG with a magnetic filling factor gradually increasing with the stellar activity level up to f ≲ 50%. In contrast to this behaviour, HD 131156A shows a weaker field covering 69 ± 28% of the stellar surface. This is the largest magnetic filling factor derived in our study. There are no obvious reasons, besides a lower mass and thus a somewhat thicker convection zone, for HD 131156A to show a different surface field strength distribution compared to earlier G dwarfs.

3.12. HD 166435

This star is a moderately active solar analogue with a rotational period of 3.43 d. We derived ⟨B⟩ = 0.69 kG with essentially no scatter in the mean field strength values corresponding to four observing epochs. One of the field strength determinations is illustrated in Fig. A.11. All observations of HD 166435 used in our study were obtained with Narval in 2010−2016. The corresponding ZDI analysis (Petit et al., in prep.) is not published yet, but the summary of the global field mapping results is available in the literature (Vidotto et al. 2014; See et al. 2019).

3.13. HD 175726

This object is similar to HD 166435 in terms of the rotational period and the fact that only a summary of ZDI results has been published. At the same time, HD 175726 exhibits a systematically weaker magnetic activity according to the coronal and chromospheric emission indices. It is therefore not surprising that we obtained nearly a factor of two weaker mean field, ⟨B⟩ = 0.37 kG, for this star compared to HD 166435. These measurements relied on the average spectra corresponding to two Narval data sets collected in 2008 and 2012. Figure A.12 shows a comparison between theoretical model and observations for the former epoch.

3.14. HD 190771

This star was included in the ZDI study of the four solar twins by Petit et al. (2008). Subsequently, Petit et al. (2009) reported a polarity reversal of the global magnetic field and published ZDI results for three separate epochs. Here we analysed observations obtained at nine epochs, covering ten years starting from 2007 and ending in 2016. This is the second largest data set (after HD 206860) investigated in this study. Our field strength measurements indicate ⟨B⟩ = 0.59 kG without any evidence of an epoch-to-epoch scatter or a long-term variation. An example of the field strength determination for this star is shown in Fig. A.13.

3.15. HD 206860 (HN Peg)

This active, young Sun-like star was studied with ZDI by Boro Saikia et al. (2015) and Rosén et al. (2016). Combining the data analysed in these studies with newer observing material yields ten epochs spread between 2007 and 2016. This is the largest collection of spectra, in terms of the number of epochs, analysed in our study. The modelling of the Fe I magnetically sensitive lines suggests ⟨B⟩ around 0.45 kG, with individual field strength determinations ranging from 0.40 to 0.55 kG. One of our field strength measurements is shown in Fig. A.14. Given the formal error bars, the ⟨B⟩ scatter appears to be significant. A long-term behaviour of the mean field strength suggests a quasi-periodic variation, with apparent magnetic minima in 2008 and 2014.

4. Discussion

4.1. Correlation with stellar parameters

The overall results of the Zeeman intensification magnetic field measurements carried out in this study are summarised in Table 5. Columns 2−4 of this table list the time-averaged values of the magnetic field strength B, the filling factor f, and the mean field strength ⟨B⟩. These parameters, and the corresponding asymmetric error bars, were obtained by calculating median values of individual measurements in Table B.1 for targets with three or more observing epochs and using mean values otherwise. In this section, we use these time-averaged magnetic field measurements to assess correlations between photospheric magnetic field characteristics and fundamental parameters, rotation, and magnetic activity indicators of the target stars.

Figure 7a illustrates dependence of the mean magnetic field strength on the stellar age. The latter was obtained from diverse literature sources and has different reliability, depending on the star. Typically, the age is known relatively well for young stars (age ≲ 1 Gyr) thanks to memberships in open clusters and young moving groups. Conversely, the age of older objects in our sample (HD 56124, HD 76151, HD 166435, HD 190771) is constrained with far lesser precision by comparing stellar spectroscopic parameters with theoretical isochrones (e.g. Valenti & Fischer 2005). Despite this caveat, Fig. 7a shows a clear overall decline of the mean field strength with age. The targets studied here can be broadly separated into three age groups with different magnetic field characteristics. The youngest group (age ≤ 120 Myr), represented by HD 29615, HD 82558, and HD 129333, has ⟨B⟩ = 1.3−2.0 kG. The intermediate group (ages between 200 and 600 Myr), comprising eight stars, shows fields in the 0.4−0.8 kG range. The four oldest stars (age ≥ 2.7 Gyr) seemingly exhibit ⟨B⟩ in a wide range of 0.2−0.7 kG and lack any correlation with age.

Fig. 7. Correlation of the mean magnetic field strength with the stellar age (panel a), rotational period (panel b), and Rossby number (panel c). Individual targets are identified with Latin letters. The vertical dotted line at Ro = 0.1 in panel c indicates the saturation limit. The dashed line shows the best fitting relation ⟨B⟩∝Ro−0.58 ± 0.07.

Next, we examine correlation between the mean field strength and stellar rotation. Figure 7b shows a general decrease of ⟨B⟩ with increasing P_rot. This trend is arguably dominated by the three rapid rotators with the strongest fields (HD 29615, HD 82558, HD 129333) and becomes less pronounced when considering the rest of the sample. The magnetic field-rotation relation can be recast in terms of the dependence on Rossby number (Fig. 7c), which is known to reduce the scatter. As established by many studies (e.g. Noyes et al. 1984; Wright et al. 2011; Douglas et al. 2014; Vidotto et al. 2014; Folsom et al. 2016), both indirect magnetic activity indicators and direct field measurements correlate with Ro until stellar dynamo reaches a saturated state at Ro ≈ 0.1. All our targets except HD 82558 (LQ Hya) have Ro > 0.1 and therefore are expected to be in the unsaturated dynamo regime. A weighted least-squares power law fit of ⟨B⟩ as a function of Ro, excluding HD 82558, yields

$$\begin{array}{c}\hfill log⟨B⟩=(2.65\pm 0.05)-(0.67\pm 0.11)·log\mathrm{Ro}.\end{array}$$(3)

One can note that the same targets tend to deviate from the general trends in all three panels of Fig. 7. For example, HD 56124 and HD 175726 appear to have magnetic fields that are too weak for their age and rotation. It cannot be excluded that parameters other than those considered here (e.g. inclination of stellar rotational axis, phase in a long-term activity cycle, etc.) contribute to the scatter in magnetic field-rotation-age relationships.

Finally, we consider separate determinations of the local magnetic field strength B and filling factor f. As discussed above, the relative precision with which these parameters can be constrained individually is much worse than that of their product ⟨B⟩ = B ⋅ f. Nevertheless, our error determination procedure accounts for the anti-correlation between B and f, allowing us to derive useful constraints on these parameters in many cases. For example, it can be established that a homogeneous one-component (f = 1) model is excluded at more than a 99.7% confidence level for all stars except HD 29615, HD 82558, and HD 131156A.

An interesting result emerges when we look at the plot of B and f as a function of ⟨B⟩ (Fig. 8). The upper panel of Fig. 8 demonstrates that all stars except HD 131156A (ξ Boo A) exhibit essentially the same local field intensity B = 3.2 ± 0.6 kG. This means that any variation of the mean field strength ⟨B⟩ is due to changes of the magnetic filling factor f, as conclusively demonstrated by Fig. 8b. Excluding HD 131156A, the relation between f and ⟨B⟩ in the unsaturated dynamo regime can be represented by

$$\begin{array}{c}\hfill logf=-(3.12\pm 0.20)+(0.86\pm 0.07)·log⟨B⟩.\end{array}$$(4)

Fig. 8. Correlations of the magnetic field strength B (panel a) and magnetic filling factor f (panel b) with the average field strength ⟨B⟩ = Bf. Individual targets are identified with Latin letters. The dashed line in the lower panel shows the best fitting relation log f ∝ ⟨B⟩0.86.

Thus, our results give further strong support, now based on a set of accurate field strength measurements for a sample of structurally similar G dwarf stars, that magnetic regions exhibit approximately the same local field intensity and that stellar magnetic activity is primarily modulated by the fraction of stellar surface occupied by these regions (Saar & Linsky 1986a; Montesinos & Jordan 1993; Cranmer & Saar 2011).

The common local field strength of 3.2 kG derived in our study noticeably exceeds the thermal equipartition field strength B_eq = 1.7−1.9 kG that can be calculated with $B_{\mathrm{eq}}=\sqrt{8\pi P_{\mathrm{gas}}}$ for the gas pressure at τ₅₀₀₀  =  1 in T_eff = 5750−6000 K, log g = 4.5 MARCS model atmospheres. A marginal evidence of superequipartition fields was found for previous heterogenous samples of, mostly K, dwarf stars (e.g. Cranmer & Saar 2011). The presence of such fields is now definitively established for active M dwarfs (Shulyak et al. 2017; Kochukhov & Lavail 2017; Kochukhov & Shulyak 2019) and T Tauri stars (Sokal et al. 2020). This suggests that local field strength is not limited by the confinement of magnetic flux tubes by the photospheric gas pressure as was repeatedly assumed in the past (Saar & Linsky 1986a; Cranmer 2017; See et al. 2019).

4.2. Correlation with activity proxies

The coronal (X-ray) emission is a frequently considered indirect proxy of magnetic activity (e.g. Pevtsov et al. 2003; Vidotto et al. 2014; Wright et al. 2011). Using our new magnetic field measurements and the ratio of X-ray to bolometric luminosity reported for the target stars in Table 4, we are able to assess a correlation between these quantities for the first time for a sample of active solar analogues. The correlation between ⟨B⟩ and log L_X/L_bol is shown in Fig. 9a. It should be noted that log L_X/L_bol may exhibit a large variation over stellar activity cycles. Based on the solar X-ray variation Vidotto et al. (2014) estimated the representative cyclic changes to be as large as 0.65 dex in log L_X/L_bol. Despite this, Fig. 9a suggests a relatively tight correlation between magnetic field strength and X-ray luminosity. We derive

Fig. 9. Panel a: correlation between the total magnetic field strength and X-ray-to-bolometric luminosity ratio log LX/Lbol. Individual targets are identified with Latin letters. The dashed line shows the best fitting power-law relation ⟨B⟩∝(LX/Lbol)0.37. The horizontal bar at the top illustrates typical variation of log LX/Lbol over stellar activity cycles. Panel b: correlation between the total magnetic field strength and the Ca HK emission measure log $R_{\text{HK}}^{\prime }$. The dashed line corresponds to the best fitting relation $⟨B⟩\propto {R{\prime }_{\mathrm{HK}}}^{1.09}$.

$$\begin{array}{c}\hfill log⟨B⟩=(4.47\pm 0.09)+(0.37\pm 0.02)·log{L}_{\mathrm{X}}/L_{\mathrm{bol}}\end{array}$$(5)

by fitting a power law relation with ⟨B⟩ measured in G.

Equation (5) is equivalent to L_X/L_bol ∝ ⟨B⟩^{2.68 ± 0.16}. In comparison, Vidotto et al. (2014) deduced L_X/L_bol ∝ ⟨B⟩^2.25 through indirect means, by combining correlations of ⟨B⟩ and L_X/L_bol with Ro published by Saar (2001) and Wright et al. (2011) respectively. Large uncertainties involved in this estimate of ⟨B⟩ vs. L_X/L_bol dependence precluded Vidotto et al. (2014) from making a meaningful comparison with the corresponding relation for the global field ⟨B_V⟩ measured by ZDI studies, L_X/L_bol ∝ ⟨B_V⟩^{1.61 ± 0.15}. Our analysis shows that the latter relation is significantly flatter than the dependence of X-ray to bolometric luminosity on the total field strength ⟨B⟩. This may be taken as an indication that small-scale fields provide a dominant contribution to the coronal emission of active Sun-like stars.

The Ca H&K chromospheric emission is another popular magnetic activity proxy that loosely correlates with the global field strength (Marsden et al. 2014) and is coupled to X-ray emission (Mittag et al. 2018). The relation between log ⟨B⟩ and the chromospheric emission measure log $R_{\text{HK}}^{\prime }$ for our stellar sample is presented in Fig. 9b. HD 29615 is missing in this plot since no log $R_{\text{HK}}^{\prime }$ measurements can be found for this star in the catalogue by Boro Saikia et al. (2018) or elsewhere in the literature. A weighted least-squares power law fit yields

$$\begin{array}{c}\hfill log⟨B⟩=(7.59\pm 0.54)+(1.09\pm 0.12)·logR{\prime }_{\mathrm{HK}}\end{array}$$(6)

with ⟨B⟩ in G.

The steep ramping up of X-ray activity with increasing ⟨B⟩, as given by Eq. (5), indicates that the coronal X-ray emission cannot itself be taken as a direct measure of the efficiency of magnetic field generation. This means that the canonically steep slope L_X/L_bol∝ Ro⁻² of the X-ray rotation-activity relation (e.g. Pizzolato et al. 2003) does not represent correctly the relation between dynamo efficiency and the stellar rotation rate, nor can observed deviations from this slope be directly interpreted in terms of predictions from dynamo models (Wright et al. 2011). In order to understand the correct scaling of the dynamo efficiency in stars, it is necessary to look at the scaling of the magnetic field itself (Eq. (3)), which indicates a rather shallow rotation dependence. Alternatively, the chromospheric Ca H&K emission also provides a decent proxy for the magnetic field strength, since the two follow a nearly linear relation with each other (Eq. (6)).

4.3. Comparison with ZDI results

The sample of active Sun-like stars investigated here was chosen among the stars monitored with high-resolution spectropolarimetry and analysed with ZDI. This gives us a unique opportunity to directly compare characteristics of the global magnetic field reconstructed using ZDI with the true total magnetic field strength inferred by our Zeeman intensification analysis. Previously similar comparisons were done using large heterogeneous samples of active stars based on historic Zeeman broadening measurements collated from diverse literature sources spanning several decades (Vidotto et al. 2014; See et al. 2019). In contrast, in this study we are able to intercompare parameters of the global and total magnetic fields obtained with modern techniques from the same observational data.

We use the mean field strength ⟨B_V⟩ as a measure of the strength of global magnetic field distributions recovered with ZDI. This is the same global field strength parameter as the one considered by See et al. (2019). On the other hand, it differs from the mean unsigned radial field ⟨|B_r|⟩ analysed by Vidotto et al. (2014). The main difference between these two characteristics, besides the fact that ⟨B_V⟩ is more readily available from original ZDI publications, is that ⟨B_V⟩ includes contributions from both poloidal and toroidal global field structures whereas ⟨|B_r|⟩ is insensitive to toroidal field since the latter does not have a radial component.

Columns 5 and 6 in Table 5 provide ⟨B_V⟩ for all stars in our sample along with the references to ZDI studies where magnetic maps were published. In the case of HD 73350, HD 76151, HD 190771 (Petit et al. 2008, 2009), and HD 131156A (Morgenthaler et al. 2012) the mean field strengths reported in the ZDI papers are incompatible with the actual magnetic maps presented in those studies. For these three stars we adopted ⟨B_V⟩ from the compilation by See et al. (2019). The same paper was used as a source of ⟨B_V⟩ values corresponding to the unpublished ZDI analyses of HD 56124, HD 166435, and HD 175726 (Petit et al., in prep.). Whenever multiple ZDI maps were available for the same star, we calculated the median global field strength value and adopted the error bars corresponding to the full range covered by individual ⟨B_V⟩ measurements.

Figure 10a presents ⟨B_V⟩ as a function of ⟨B⟩. A weighted least-squares fit with a power law function yields the following relation between these two quantities

Fig. 10. Panel a: strength of the global magnetic field estimated with ZDI as a function of the total magnetic field strength derived in this study. Individual targets are identified with Latin letters. The dashed line shows the best fitting power-law relation ⟨B⟩∝⟨BV⟩0.48. Panel b: ratio ⟨BV⟩/⟨B⟩ as a function of the total field strength. Panel c: same ratio as a function of ve sin i.

$$\begin{array}{c}\hfill log⟨B⟩=(2.19\pm 0.08)+(0.48\pm 0.05)·log⟨B_V⟩\end{array}$$(7)

for the field strengths measured in G. This relation is somewhat flatter and considerably more precise than ⟨B⟩∝⟨B_V⟩^{0.78 ± 0.12} derived by See et al. (2019) for the M > 0.5 M_⊙ part of their heterogeneous active star sample.

The ratio ⟨B_V⟩/⟨B⟩ is plotted as a function of ⟨B⟩ in Fig. 10b. As expected, the mean global field strength obtained with ZDI is vastly weaker than the total field strength diagnosed from Stokes I. The relation between ⟨B_V⟩ and ⟨B⟩ cannot be described by a constant scaling factor (Cranmer 2017). Instead, the global-to-total field strength ratio clearly grows with ⟨B⟩. The highest ⟨B_V⟩/⟨B⟩ of ∼10% is found for HD 82558 (LQ Hya). This ratio drops to about 6% for HD 29615 and HD 129333 and lies in the 1−5% range for the remaining targets. In other words, ZDI recovers merely 0.01−1% of the total magnetic field energy (approximated by square of the mean field strength). This means that the total magnetic field measured in our study is dominated by a small-scale magnetic component invisible to ZDI. Nevertheless, the tightness of the relation described by Eq. (7) suggests that the large- and small-scale fields are closely coupled and are likely to be produced by the same underlying dynamo mechanism.

The systematic increase of the fraction of the total magnetic field recovered by ZDI studies with ⟨B⟩ evident in Fig. 10b can be due to an intrinsic shift of the magnetic energy to larger spatial scales in more active stars. Alternatively, this effect can be a consequence of a better surface resolution achieved by ZDI studies of faster rotating stars due to a larger Doppler broadening of their line profiles. A correlation of ⟨B_V⟩/⟨B⟩ with v_e sin i seen in Fig. 10c supports the latter hypothesis.

4.4. Effect of magnetic field on abundance determination

The strengthening of spectral lines due to magnetic intensification is usually ignored in stellar abundance analyses, especially in the context of modern spectroscopic studies dealing with large stellar samples. Based on the outcome of our in-depth investigation of magnetic fields of selected active solar analogues, we are able to perform a quantitative assessment of the likely errors incurred by neglecting magnetic field. To this end, we calculated a set of synthetic profiles of the three magnetically sensitive Fe I line for fixed B = 3.2 kG and changed the filling factor to get B ⋅ f in the interval from 0.05 to 2.5 kG. We then adjusted Fe abundance in the non-magnetic spectrum synthesis model until it matched the equivalent widths measured in magnetic spectra. All calculations adopted the same model atmosphere with T_eff = 5750 K, log g = 4.5, and assumed v_mic = 0.85 km s⁻¹.

The resulting error of the Fe abundance determination is illustrated in Fig. 11. The curve corresponding to the Fe I 5497.5 Å line represents the upper limit of possible errors of an optical abundance analysis based on a single absorption feature. The curves corresponding to the other two lines are closer to average multi-line abundance determination bias. This plot indicates that an abundance error of up to 0.4−0.5 dex can be encountered for the most active stars (⟨B⟩ = 1.5 − 2.0 kG), such as HD 82558. The errors are still significant (0.15−0.30 dex) for stars with a more modest activity levels characterised by ⟨B⟩ = 0.5 − 1.0 kG. This analysis thus suggests that metallicities and individual abundances reported for active stars (e.g. Valenti & Fischer 2005; Brewer et al. 2016) may have been systematically overestimated. An unusually high line-to-line abundance scatter or a trend of abundance with wavelength might be a symptom of unaccounted Zeeman intensification.

Fig. 11. Error of Fe abundance determination incurred by ignoring magnetic intensification of Fe I lines.

5. Conclusions

In this paper we have developed a new methodology for measuring average unsigned magnetic field strength in the atmospheres of Sun-like stars. This method takes advantage of the differential Zeeman intensification of Fe I lines in the optical wavelength region. The new magnetic field diagnostic procedure is less restrictive compared to traditional Zeeman broadening analysis. In particular, it can be applied to moderately fast rotators and does not require spectroscopic data of exceptionally high quality. The technique can be implemented using typical high-resolution optical echelle spectra that are currently provided by many different night-time astronomical facilities around the world.

We have performed magnetic field measurements for a sample of 14 G dwarf stars and one early K dwarf, thereby considerably extending the number of solar analogues at different ages and activity levels with direct observational constraints on the total surface magnetic field. We derived 78 individual mean field strength (⟨B⟩=B ⋅ f) measurements and, in most cases, were able to provide meaningful constraints on the local field strength B and the fraction f of the stellar surface covered by a magnetic field. These data greatly improve our understanding of different evolutionary phases of the solar-like magnetic dynamo and can be used to estimate magnetic characteristics of the Sun during the first ∼1 Gyr of its life. Main conclusions of our study can be summarised as follows.

• We have measured mean magnetic field strengths in the range from 0.2 to 2.0 kG and detected a systematic decline of the field intensity from B ⋅ f >  1 kG in stars younger than ∼100 Myr to weaker fields in older stars.

• A clear anti-correlation of the mean field strength and the stellar rotational period, or Rossby number, is established and a new empirical calibration of B ⋅ f as a function of Ro is obtained.

• Separate analysis of the local field strength B and magnetic filling factor f suggests that the increase of B ⋅ f associated with the transition from less active to more active stars is mostly due to a systematic growth of f from ∼10% to > 50% of the stellar surface. The local field strength B remains at the level of ∼3 kG in all stars.

• The mean magnetic field strength determined in our study exhibits clear correlations with the chromospheric and coronal emission indicators. This allowed us to derive new calibrations of the widely used activity indices log L_X/L_bol and log $R_{\text{HK}}^{\prime }$ in terms of the surface magnetic field strength.

• We have compared our magnetic field strength measurements with results of ZDI analyses of the global magnetic field topologies, which in most cases relied on the same observational data as used in our study. We found a tight correlation of the ratio of the global and total magnetic field strengths, ⟨B_V⟩/⟨B⟩, with the stellar activity level. In terms of this ratio, the fraction of magnetic field detected by spectropolarimetric studies varies from ∼10% in the most active stars to less than 1% in the least active objects.

• Ignoring magnetic intensification in an optical spectroscopic abundance analysis of active Sun-like stars may lead to overestimation of element concentrations by up to 0.4–0.5 dex if this analysis is based on lines with a large magnetic sensitivity.

The most important correlations derived in our study are schematically illustrated in Fig. 12. This plot combines information from Figs. 7–10 to show how the field strength, the filling factor and the ratio of the global to total field strength change with Rossby number and how these trends relate to X-ray emission of the target stars.

Fig. 12. Summary of the magnetic field analysis results. The average field strength is plotted as a function of Rossby number, similar to Fig. 7c. The symbol size increases with the magnetic filling factor f. The inner and outer outlines reflect uncertainty of the filling factor determination. The symbol shape corresponds to the ratio of the global magnetic field strength obtained with ZDI and the total field strength, ⟨BV⟩/⟨B⟩. This ratio ranges from about 4% to less than 0.1%. The symbol colour corresponds to the X-ray-to-bolometric luminosity ratio, log LX/Lbol.

There are numerous further applications of the magnetic diagnostic procedure developed in our paper. For example, this technique can be systematically applied to large volumes of high-resolution spectroscopic data accumulated for Sun-like stars by the radial velocity exoplanet surveys (e.g. Adibekyan et al. 2012; Delgado Mena et al. 2017; Hojjatpanah et al. 2019). This method can be also utilised to investigate the variation of magnetic field strength during stellar activity cycles for active stars with sufficiently long spectroscopic time-series (Alvarado-Gómez et al. 2018; Yana Galarza et al. 2019). One can also envisage incorporating Zeeman intensification of the Fe I lines studied here into Doppler imaging inversions with the goal of simultaneously reconstructing both cool spot maps and distributions of small-scale magnetic field intensity (Saar et al. 1992).

Acknowledgments

O.K. acknowledges support by the Knut and Alice Wallenberg Foundation (project grant “The New Milky Way”), the Swedish Research Council (projects 621-2014-5720, 2019-03548), and the Swedish National Space Board (projects 185/14, 137/17). J.J.L. acknowledges the Academy of Finland “ReSoLVE” Centre of Excellence (project number 307411) and the Max Planck Research Group “SOLSTAR” funding. This work is based on observations made with the HARPSpol instrument on the ESO 3.6 m telescope at La Silla (Chile), under the programmes 091.D-0836 and 0100.D-0176. Also based on archival observations obtained at the Bernard Lyot Telescope (Pic du Midi, France) of the Midi-Pyrénées Observatory, which is operated by the Institut National des Sciences de l’Univers of the Centre National de la Recherche Scientifique of France, and on observations obtained at the Canada-France-Hawaii Telescope, which is operated by the National Research Council of Canada, the Institut National des Sciences de l’Univers of the Centre National de la Recherche Scientifique of France, and the University of Hawaii. This study has made use of the VALD database, operated at Uppsala University, the Institute of Astronomy RAS in Moscow, and the University of Vienna.

References

Appendix A: Additional figures

Fig. A.1. Same as Fig. 6, but for HD 1835 (2013.70 epoch).

Fig. A.2. Same as Fig. 6, but for HD 20630 (2013.70 epoch).

Fig. A.3. Same as Fig. 6, but for HD 29615 (2017.95 epoch).

Fig. A.4. Same as Fig. 6, but for HD 39587 (2007.09 epoch).

Fig. A.5. Same as Fig. 6, but for HD 56124 (2008.09 epoch).

Fig. A.6. Same as Fig. 6, but for HD 72905 (2014.32 epoch).

Fig. A.7. Same as Fig. 6, but for HD 73350 (2011.05 epoch).

Fig. A.8. Same as Fig. 6, but for HD 76151 (2009.95 epoch).

Fig. A.9. Same as Fig. 6, but for HD 82558 (2016.05 epoch).

Fig. A.10. Same as Fig. 6, but for HD 131156A (2013.33 epoch).

Fig. A.11. Same as Fig. 6, but for HD 166435 (2010.53 epoch).

Fig. A.12. Same as Fig. 6, but for HD 175726 (2008.59 epoch).

Fig. A.13. Same as Fig. 6, but for HD 190771 (2010.48 epoch).

Fig. A.14. Same as Fig. 6, but for HD 206860 (2013.70 epoch).

Appendix B: Individual magnetic field measurements

All Tables

All Figures

Fig. 1. Response of the continuum normalised spectrum in the 4000−10 000 Å wavelength range to a 1 kG magnetic field covering the entire stellar surface. Lower panel: synthetic spectra calculated for the solar atmospheric parameters and abundances with (light curve) and without (dark curve) magnetic field. Upper panel: difference between these two spectra. The lines showing the strongest Zeeman response are identified according to the numbering adopted in Table 1.

In the text

Fig. 2. Theoretical profiles of the four Fe I spectral lines studied in this paper. Calculations for the solar parameters and magnetic field strengths ranging from 0 to 4 kG are shown for representative rotational, macroturbulent and instrumental broadening (solid lines) and without broadening (dotted lines). The bar plots above line profiles schematically show the Zeeman splitting patterns for a field strength of 4 kG.

In the text

	Fig. 3. Equivalent width change of the Fe I lines λ 5497.5, 5501.5, and 5506.8 Å relative to the Fe Iλ 5435.5 Å line in response to the variation (from top to bottom) of the effective temperature, surface gravity, microturbulent velocity, Fe abundance, and magnetic field strength.
In the text

	Fig. 4. Solar flux spectrum obtained with the HARPS spectrograph (symbols connected by thin solid lines) compared to best fitting synthetic profiles (solid lines) for the four diagnostic Fe I lines studied in this paper.
In the text

Fig. 5. Illustration of the determination of magnetic field strength B and filling factor f for HD 129333 (2009.02 epoch). The greyscale corresponds to the χ2 of the fit to the three magnetically sensitive Fe I lines, as quantified by the side bar. The best fitting parameters, B = 3.5 kG and f = 0.4, are indicated with the cross. Solid lines show 68.3, 95.5, and 99.7% confidence limits. The dashed line corresponds to the B ⋅ f = const curve.

In the text

Fig. 6. Magnetic field measurement for HD 129333 (2009.02 epoch). The average observed profiles of the four Fe I lines are shown with black symbols connected by thin solid lines. The grey curve in the background corresponds to ± twice the standard deviation for each pixel of the observed spectra. The thick red solid line shows theoretical spectrum for the magnetic model yielding the lowest χ2 according to Fig. 5. The dashed blue line shows the corresponding non-magnetic spectrum.

In the text

	Fig. 7. Correlation of the mean magnetic field strength with the stellar age (panel a), rotational period (panel b), and Rossby number (panel c). Individual targets are identified with Latin letters. The vertical dotted line at Ro = 0.1 in panel c indicates the saturation limit. The dashed line shows the best fitting relation ⟨B⟩∝Ro−0.58 ± 0.07.
In the text

	Fig. 8. Correlations of the magnetic field strength B (panel a) and magnetic filling factor f (panel b) with the average field strength ⟨B⟩ = Bf. Individual targets are identified with Latin letters. The dashed line in the lower panel shows the best fitting relation log f ∝ ⟨B⟩0.86.
In the text

Fig. 9. Panel a: correlation between the total magnetic field strength and X-ray-to-bolometric luminosity ratio log LX/Lbol. Individual targets are identified with Latin letters. The dashed line shows the best fitting power-law relation ⟨B⟩∝(LX/Lbol)0.37. The horizontal bar at the top illustrates typical variation of log LX/Lbol over stellar activity cycles. Panel b: correlation between the total magnetic field strength and the Ca HK emission measure log $R_{\text{HK}}^{\prime }$. The dashed line corresponds to the best fitting relation $⟨B⟩\propto {R{\prime }_{\mathrm{HK}}}^{1.09}$.

In the text

	Fig. 10. Panel a: strength of the global magnetic field estimated with ZDI as a function of the total magnetic field strength derived in this study. Individual targets are identified with Latin letters. The dashed line shows the best fitting power-law relation ⟨B⟩∝⟨BV⟩0.48. Panel b: ratio ⟨BV⟩/⟨B⟩ as a function of the total field strength. Panel c: same ratio as a function of ve sin i.
In the text

	Fig. 11. Error of Fe abundance determination incurred by ignoring magnetic intensification of Fe I lines.
In the text

Fig. 12. Summary of the magnetic field analysis results. The average field strength is plotted as a function of Rossby number, similar to Fig. 7c. The symbol size increases with the magnetic filling factor f. The inner and outer outlines reflect uncertainty of the filling factor determination. The symbol shape corresponds to the ratio of the global magnetic field strength obtained with ZDI and the total field strength, ⟨BV⟩/⟨B⟩. This ratio ranges from about 4% to less than 0.1%. The symbol colour corresponds to the X-ray-to-bolometric luminosity ratio, log LX/Lbol.

In the text

	Fig. A.1. Same as Fig. 6, but for HD 1835 (2013.70 epoch).
In the text

	Fig. A.2. Same as Fig. 6, but for HD 20630 (2013.70 epoch).
In the text

	Fig. A.3. Same as Fig. 6, but for HD 29615 (2017.95 epoch).
In the text

	Fig. A.4. Same as Fig. 6, but for HD 39587 (2007.09 epoch).
In the text

	Fig. A.5. Same as Fig. 6, but for HD 56124 (2008.09 epoch).
In the text

	Fig. A.6. Same as Fig. 6, but for HD 72905 (2014.32 epoch).
In the text

	Fig. A.7. Same as Fig. 6, but for HD 73350 (2011.05 epoch).
In the text

	Fig. A.8. Same as Fig. 6, but for HD 76151 (2009.95 epoch).
In the text

	Fig. A.9. Same as Fig. 6, but for HD 82558 (2016.05 epoch).
In the text

	Fig. A.10. Same as Fig. 6, but for HD 131156A (2013.33 epoch).
In the text

	Fig. A.11. Same as Fig. 6, but for HD 166435 (2010.53 epoch).
In the text

	Fig. A.12. Same as Fig. 6, but for HD 175726 (2008.59 epoch).
In the text

	Fig. A.13. Same as Fig. 6, but for HD 190771 (2010.48 epoch).
In the text

	Fig. A.14. Same as Fig. 6, but for HD 206860 (2013.70 epoch).
In the text