© ESO 2020

1. Introduction

The study of the sunspot number as well as that of solar variability have been good indicators of magnetic activity (Eddy 1976; Stuiver & Quay 1980; Kivelson & Russell 1995). Continuous records of the sunspot number have existed since the first observation by Galileo at the beginning of the seventeenth century¹. The solar activity from observational data exhibits a period of about 11 years with amplitude variations that remain an open issue in modern astrophysics (Thomas & Weiss 2012).

One of the most intriguing phenomena observed in the sunspot counting survey is the long period phase, posteriorly called the Maunder Minimum (MM hereafter) when almost no sunspots were observed on the solar surface (Eddy 1976). It lasted for about 70 yr between the sixteenth and seventeenth centuries. There is still no complete explanation for this period of extremely low activity (Wright 2004; Lubin et al. 2010), however, it is historically known that it coincided with a short glacial era on Earth (Eddy 1976).

The fundamental question about how singular the magnetic activity evolution of the Sun has motivated long-term surveys of stellar activity of solar-type stars. Chromospheric stellar activity basically consists of the temporal variations of the stellar magnetic field caused by mass loss and radiation outside the star (Hall 2008). These variations have a connection with the flux of Ca II H&K lines, which have been widely used as a spectroscopical indicator of stellar activity (Auer & Mihalas 1969; Skumanich 1970; Wilson 1968). Regarding the Sun, for instance, Wilson (1963, 1968); Wilson & Skumanich (1964) established the connection between solar activity and the strength of the Ca II H&K lines.

Strong efforts have been made over the last decades seeking to increase the spectroscopic data (Schwarzschild 1958; Auer & Mihalas 1969). In his seminal work, Skumanich (1972) showed that activity, as given by the Ca II H&K lines, rotation, and lithium abundance, decay as the inverse square root of time. This paper was a milestone in stellar chromospheric activity studies. The first observational program to continuously measure the activity of main sequence stars, the HK project at the Mount Wilson Observatory² (MWO hereafter), was described in detail by Vaughan et al. (1978). Improved bases of the project and its results of spectroscopical measurements were related by Wilson (1978). Later on, Duncan et al. (1991) presented a large set of time series spectroscopical measurements from the Mount Wilson survey containing data for 1296 stars observed between 1966 and 1983. Baliunas et al. (1995) presented enhanced and extended results of Mount Wilson stars and flagged the activity of the stars as cyclic and flat.

From these observational results, it is accepted that young stars have strong activity levels (do Nascimento et al. 2016). However, some young stars present extremely low activity over many years. This unexpected phase is a puzzle for theoreticians and observing astronomers. Many authors have proposed using the solar chromosphere as a pattern to study the activity evolution of other solar-type stars. It could be an approach to investigate the occurrence of MM in stars similar to the Sun (Schrijver 1987; Schrijver et al. 1989, 1992; Schröder et al. 2012).

In this work, we investigated the low spectroscopic activity of the solar analog HD 43587 over 50 years of measurements. This object has been widely referred to as a flat activity star and an MM candidate (Baliunas et al. 1995; Wright 2004; Lubin et al. 2010; Schröder et al. 2012). We used the fluxes in the Ca II H&K lines, through the Mount Wilson S-index (S_MW), whose photometric measurements were taken with the 1.5 m telescope at MWO. Since the MWO survey, this object has been observed by several others surveys (detailed below in the Sect. 2.1), as well as being described by Wright (2004), Morel et al. (2013), and Radick et al. (2018). We focused on the continuous low activity level presented from spectroscopic measurements, showed by the time-averaged value ⟨S_MW⟩ = 0.154. In addition, the observations from the CoRoT satellite provide a photometric activity proxy S_ph = 102 ppm as described by Mathur et al. (2014), sustaining the flat activity status.

This paper is organized as follows; in Sect. 2 we briefly present the observational data and reduction used along this work, and in Sect. 3 we present our results. Section 4 shows an attempt to simulate the star HD 43587 by using the EULAG-MHD code. Finally, we present our conclusions in Sect. 5.

2. Observational data

In this section, we briefly present the characterization of the solar analog HD 43587, and explain all of the spectroscopic data used, as well as the S_MW calibration for each instrument. In the last subsection, we explain the photometric data from the CoRoT satellite.

2.1. The solar analog HD 43587

Our object HD 43587 (more precisely HD 43587Aa) belongs to a quadruple system composed of two distant main-sequence visual binaries and localized in the Orion constellation. It is a bright sun-like star, and its fundamental parameters are presented in Table 1. This system has been largely observed by many authors (Duncan et al. 1991; Baliunas et al. 1995; Catala et al. 2006; Morel et al. 2013). This star was also monitored by using the CFHT (Canadian-France Hawai Telescope), for four seasons: January 2002 (one night), June 2002 (two nights), September 2004 (one night) and September 2005. In all cases, H₂ and [Fe II] filters (Catala et al. 2006) were used.

HD 43587 was a CoRoT primary target in the seismology program with ID CoRoT 3474, observed during the LRa03 for 145 days. Complementary observations were done with the spectrograph High Accuracy Radial velocity Planet Search (HARPS), mounted on the 3.6 m telescope in La Silla Observatory/ESO (Boumier et al. 2014). In addition, it has been monitored by several spectroscopical long-term observational surveys dedicated to activity, such as the HK project of MWO (Duncan et al. 1991; Baliunas et al. 1995) between 1966 and 1990, and the California and Carnegie Planet Search, at the Keck Observatory between 1995 and 2000 (Wright 2004). This object was also continuously monitored by Automatic Photometric Telescopes of the Tennessee State University at the Fairborn Observatory, and the Solar-Stellar Spectrograph at the Lowell Observatory between 1995 and 2016 (Radick et al. 2018). It was also observed by the NARVAL spectrograph, mounted on the Bernard Lyot Telescope at the Observatoire du Pic du Midi between 2009 and 2010³ (Petit et al. 2014). Spectroscopic analysis of this star performed by the McDonald Observatory 2.1 m Telescope and Sandiford Cassegrain Echelle Spectrograph derived a projected rotation velocity v sin i = 4.9 km s⁻¹ (Luck 2017).

This star has been considered by various authors as a strong MM candidate (Baliunas et al. 1995; Lubin et al. 2010; Schröder et al. 2012; Morel et al. 2013). Its long-time activity series shows a long flat activity lasting about 50 years. Its time-averaged activity index is ⟨S_MW⟩ = 0.154 ± 0.004, while for the Sun, ⟨S_MW⟩_⊙ = 0.167 ± 0.006. The observational spectroscopic data considered here were obtained in part from already published long-term surveys, as well as from HARPS and NARVAL spectra. These data are detailed as follows:

• NSO/MWO: measurements from MWO obtained between 1966 and 1983 (Duncan et al. 1991). Additionally, we mainly used a dataset updated from National Solar Observatory (NSO) with S-index measurements until 1995⁴;

• NARVAL: 50 high-resolution spectra obtained in the archive available from Narval/Polarbase (Petit et al. 2014);

• HARPS: three high-resolution spectra from HARPS/ESO (Morel et al. 2013);

• CCPS: S-index activity measurements from the California and Carnegie Planet Search (CCPS) program at Keck Observatory, available in Wright (2004).

• LO: a dataset of chromospheric Ca II H+K emission observations between 1994 and 2016 from the Lowell Observatory (LO), available in Radick et al. (2018).

Recently, Castro et al. (2020) presented an evolution status analysis of this star using an optimization approach, by modelling through two different evolution codes, TGEC (Hui-Bon-Hoa 2008) and CESTAM (Marques et al. 2013), and by using spectroscopic, photometric, and seismic constraints. Results of the optimization are presented in Table 1. These results confirm the solar analog status of HD 43587, and indicate a star that is slightly more massive and older than the Sun, which is in agreement with its observed low activity level.

2.2. The Mount Wilson S-index

Long-term observational programs at the MWO have been relevant to the detection of anomalous behaviour in stellar activity like the extreme low-activity cases, for instance, young solar-like stars probably in MM phase. In fact, the MM has been investigated by many authors since the 70s (Eddy 1976; Duncan et al. 1991; Baliunas et al. 1995; Wright 2004; Lubin et al. 2010; Schröder et al. 2012; Radick et al. 2018), initially in the MW project and afterward in other surveys. Over the years several indexes to measure stellar activity have been proposed. Here, we adopt the well-known MW S-index (see Eq. (1)), used to measure the stellar chromospheric activity and to investigate the stellar activity evolution.

Stellar spectra provide valuable information about stellar chromospheres. To measure stellar chromospheric activity, we used the flux in the Ca II H&K lines, to derive the S_MW index, which was initially proposed in the MW project (Duncan et al. 1991; Baliunas & Jastrow 1990; Baliunas et al. 1995).

The S_MW is computed as follows:

$$\begin{array}{c}\hfill S_{\mathrm{MW}}=\alpha \frac{H+K}{R+V},\end{array}$$(1)

where α is an observational constant defined by Duncan et al. (1991), H and K mean the flow strengths in triangular bandpasses with FWHM = 2.18 and centralized at 3933.633 Å (K) and 3968.469 Å (H), and R and V are continuous rectangular bands centralized at 3901.07 Å (V) and 4001.07 Å (R).

Since 1960, many authors have investigated the spectroscopic activity using S-index time series (Duncan et al. 1991; Baliunas & Jastrow 1990; Wright 2004). Few stars with low activity or flat activity, compared to the Sun at the minimum, have been reported. The value of the time average ⟨S_MW⟩ for these objects is slightly smaller than ⟨S_MW⟩_⊙. Subsequently, some studies were performed to understand the nature of flat-activity stars, among them, Schröder et al. (2012) compared the solar chromospheric activity as a pattern in the context of solar-like stars. In spite of these efforts, it has been difficult to precisely determine any star in the MM phase.

Since the first spectroscopic activity measurements, the solar S-index has varied between 0.16 (minimum) and 0.22 (maximum). The S-index time series of the Sun from Egeland et al. (2017) is shown in Fig. 1. According to Schröder et al. (2012), the low activity levels of the Sun are a good pattern for investigate the activity evolution of other solar-like stars under a long flat-activity profile. These time-lapses of lower activity are useful to study the chromospheric flat-activity phase in solar-type objects like HD 43587. This may be possible evidence of flat-activity stars at the MM phase, as proposed by Baliunas et al. (1995).

Fig. 1. Sun’s chromospheric activity over 50 years according to Egeland et al. (2017). Upper: complete time series measurements made over 50 years depicting the SMW periodic profile. Bottom: both the Lomb-Scargle (left), and the autocorrelation function (ACF, see Brockwell 2002) (right) periodograms showing the Pcyc ≈ 11 yr solar activity cycle.

2.3. S-index from NARVAL spectra

We also used spectra obtained with the NARVAL spectrograph, mounted on the Telescope Bernard Lyot (TBL) located at Pic du Midi, France. The spectra are available in the PolarBase⁵ archive (Petit et al. 2014). We used 50 high-resolution spectra collected in September 2009, January 2010, and February 2010. It consists of five spectra per night, totalizing 10 nights. The spectra corresponding to the maximum and minimum of the activity observed during this campaign are shown in Fig. 2.

Fig. 2. HD 43587 spectra from NARVAL spectrograph. The magenta spectrum corresponds to the maximum SMW value of the dataset, observed on February 13, 2010. The blue spectrum corresponds to the minimum SMW value, observed on January 18, 2010. Both spectra were normalized from 3889 Å to 4013 Å. An upward shift in flux of 0.05 was applied in the blue spectra to better visualization. The inset panel shows the region between 3930 Å and 3938 Å in detail.

Since different spectrometers have different bandpasses, it is necessary to use a different calibration to derive the S_MW-index for each instrument. For NARVAL spectra, S_MW is derived following the method proposed by Marsden et al. (2014), where a non-renormalized spectrum was used, similar to that used by Wright (2004). The activity index is computed as follows:

$$\begin{array}{c}\hfill S_{\mathrm{MW}}=\frac{12.873H+2.502K}{8.877R+4.271V}+1.183\times {10}^{-3}.\end{array}$$(2)

The coefficients in Eq. (2) were obtained using a least-square fit (see Fig. 4 in Marsden et al. 2014) to the 113 stars observed by NARVAL (Marsden et al. 2014) and by the California Planet Searcher (Wright 2004).

2.4. S-index from HARPS spectra

We used three high-resolution spectra observed between December 2010 and January 2011, obtained with the spectrograph HARPS, in the high-efficiency mode (for detailed information, see Sect. 2 of Morel et al. 2013). These spectra are also available in the SISMA database (Rainer et al. 2016).

We computed the instrumental S-index from HARPS spectra following the same calibration as prescribed by Lovis et al. (2011):

$$\begin{array}{c}\hfill S_{\mathrm{HARPS}}=\alpha ·8·\frac{1.09\mathrm{\AA }}{20\mathrm{\AA }}·\frac{H+K}{R+V}·\end{array}$$(3)

The calibration to the MW scale was done using a linear calibration between both indexes obtained with standard stars. According to Lovis et al. (2011), the activity index transformation is

$$\begin{array}{c}\hfill S_{\mathrm{MW}}=1.111·S_{\mathrm{HARPS}}+0.0153.\end{array}$$(4)

2.5. Photometric observation from CoRoT mission

The CoRoT satellite observed the object HD 43587 over 145 days (Boumier et al. 2014), during the third long run (LRa03). The LRa03 was in the direction of the Galactic anticenter and was performed from October 2009 to March 2010. We made all reductions from the light curve extracted from the level 2 official products (Samadi et al. 2007). The possibility that the clear absence of a rotational modulation is due to the fact that the star is seen pole on, is ruled out by the determination of the projected rotation velocity by Luck (2017) and mentioned in Sect. 2.1.

3. Analysis and discussion

This section is divided in three parts. The first one presents the analysis of the photometric data collected by the CoRoT satellite. The second part presents the long-term time series analysis of the S_MW and the inferred activity cycle. In the third part, we discuss about the indirect determination of the HD 43587 rotation period.

3.1. The HD 43587‘s CoRoT light curve

Based on the photometric data from the CoRoT satellite, we investigated the variability of the light curve of HD 43587. Aiming to determine the rotational period, we computed the Generalized Lomb-Scargle (GLS) periodogram (Zechmeister & Kürster 2009) as well as the Autocorrelation Function (Brockwell 2002) (see Fig. 3). We used the “DATEJD” as the time array and the “RAWFLUX” in electrons per second as the raw flux to create the light curve. The error relative to the raw flux, “RAWFLUXDEV” in electrons per second, was used to compute S_ph, as explained below. We also took into account the “RAWSTATUS” flag to make sure we used valid data. The resulting valid data exhibits a significant jump in the flux in the first 4.65 days. This discontinuity creates a false tendency on the overall aspect of the light curve, and, even when corrected, this discontinuity does not have the same dispersion in flux as the rest of the light curve. Thus, we decided to completely remove the first 4.65 days of the data. Then, we followed two approaches. First, we tried to recover the GLS and ACF periodograms (red lines in the lower panels of Fig. 3) from the raw light curve (after removing the jump). We initially obtained the periodicities $P=93.{98}_{-43.49}^{101.31}$ days with the GLS and P = 86.81 days with the ACF, as indicated by the vertical red dashed lines in the lower panels of Fig. 3. The modulation of ∼90 days is due to the accumulation of discontinuities in the flux data, as can be seen in the upper-right panel (red dots), which creates a false tendency. The second approach was to detrend the entire light curve in each discontinuity, thus eliminating the long-term modulation. We then obtained periodicity peaks at ∼7, ∼10, ∼14, and ∼21 days with the GLS, and a peak at 16.96 days with the ACF. However, these GLS peaks do not have any statistical significance as they lay far below the 5% false alarm probability, as presented in the middle-left panel of Fig. 3.

Fig. 3. CoRoT light curve of HD 43587. Top-left panel: raw light curve including discontinuities (red dots) and the detrended light curve (black dots). Top-right panel: detail of the discontinuities at the end of the light curve (red dots) and the flat aspect of the curve after the reduction (black dots). Middle-left panel: GLS periodogram for the detrended light curve with false alarm probability of 5% and 1% (green dots and blue dashed line, respectively). Middle-right panel: GLS periodogram for raw light curve with false alarm probability the same as the detrended analysis. Bottom panel: ACF periodogram for the same two reductions.

In accordance with do Nascimento et al. (2014), such stars should present a rotation period around the solar value, ranging between 14 and 30 days. Even though the rotation period inferred by the ACF method belongs to this interval, the peak corresponding in the ACF curve is not strong enough, compared to the rest, to be considered as a reliable result. Thus, even after eliminating the false long-term tendency, it was impossible to determine a valid modulation period for HD 43587 from the CoRoT light curve.

We also derived the photometric activity proxy, by computing the standard deviation of the light curve corrected after subtracting the photon noise S_ph (Mathur et al. 2014). We obtained S_ph = 102 ppm, which returns a rotational period of 43.38 d. This value of P_rot is characteristic of subgiant and RGB stars for this stellar mass, according to Mathur et al. (2014). On the other hand, the weak variability of this light curve, represented by S_ph = 102 ppm over a 50-year time span, reinforces the flat-activity profile of this star. However, the lack of modulation could also be an indication of a dynamo shutdown. Therefore, these preliminary photometric results strengthen the argument in favor of HD 43587 being under a long flat-activity phase.

Finally, we investigated the possibility of revisiting the photometric rotational modulation of this object through a combination of CoRoT and TESS data, in order to enhance the accuracy of the measurements. Unfortunately, after locating this star in the TESS field of view, we realized that HD 43587 was situated between two CCDs in one of the TESS field sectors.

3.2. S-index time-series analysis

We computed the S-index for the different spectra obtained with the spectrographs NARVAL and HARPS. We used our own tools developed according to the prescriptions of Wright (2004) and Marsden et al. (2014). These tools extract the activity index from the Ca II H&K lines emission fluxes, performing a discrete counting from the spectral lines H, K, R, and V. We used Eq. (2) to derive the S_MW directly from NARVAL spectra. From HARPS spectra, we computed the instrumental S_HARPS following Eq. (3) and calibrated it to the MW scale S_MW using Eq. (4).

Afterward, we combined the whole computed S_MW time series with the previous published data as described in Sect. 2. The complete dataset contains 1524 measurements obtained over 50 years between 1966 and 2016. This approach allows us to build a long activity time series calibrated to the MW scale, which is presented in the upper panel of Fig. 4. The standing flat-activity profile, ⟨S_MW⟩ = 0.154 and $log\langle R_{\text{HK}}^{\prime }\rangle =-4.97$ (see Sect. 3.3), is reinforced by the very low variability of the S-index, quantified by the standard deviation σ_S = 0.0043.

Fig. 4. Entire spectroscopic chromospheric activity measurements with the Sindex calibrated to the Mt. Wilson Scale for HD 43587 between 1966 and 2016. Upper panel: all measurements from Duncan et al. (1991), Wright (2004), Hall (2008), and Radick et al. (2018), combined with the computed SMW from the NARVAl and HARPS spectra archives. The red solid line is the sinusoidal curve fitting for 10.436 yr, the blue solid line is the long-term trend found to be larger than 50 yr. The yellow dashed line shows the mean SMW for the Sun. Bottom left panel: GLS periodogram of the whole SMW time series (solid gray line), and removing the long trend of over 50 years (solid black line). The Gaussian fit (solid red line) of this second periodogram indicates an activity cycle of 10.436 years (vertical black dashed line). Bottom-right panel: phase of the SMW (black circles) and the folded fit with the found period (red circles).

From this complete activity time series, we investigate the periodicity of a reliable activity cycle of HD 43587 consistent with a solar analog star by using the GLS and the Bayesian generalized Lomb-Scargle periodogram with a trend (BGLST) as proposed by Olspert et al. (2018). The results obtained with both methods are similar and are detailed below.

3.2.1. Generalized Lomb-Scargle analysis

We used the GLS analysis to compute the activity cycle period, as well as to construct a phase curve for the best determination. The GLS periodogram of the whole S_MW time series is represented by the solid gray line in the bottom left panel of Fig. 4. We found two peaks, one around 10 years, and another one, representing a long-term tendency of over 50 years. This latter trend is represented in the upper panel by the solid blue line. We removed this tendency from the time series and generated a new GLS periodogram represented by the solid black line in the bottom-left panel. We used a Gaussian fit (solid red line) to find the cycle period without the long term trend. We thus inferred a cycle period, P_cyc = 10.44 ± 3.03 yr. This result is represented by a sinusoidal fitting in the upper panel (red solid line) along with the blue long term trend. The bottom-right panel shows the phase of the folded time series with the sinusoidal model plotted over it. This period of the activity cycle is in good agreement with other solar analog stars (see e.g., Metcalfe et al. 2016). This result revises the activity cycle length proposed by Egeland (2017) who, using a shorter dataset, found P_cyc ≈ 36 yr.

We also found a long-term trend that might be a longer cycle with a length exceeding 50 years. This possible longer period could be similar in nature to the period of the solar Gleissberg cycle, which is about 88 yr. We note that HD 43587 could be around the minimum of this cycle as indicated by the blue line in Fig. 4. A considerable observational difficulty lies in the fact that, for HD 43587, but also for other analog stars, this secondary cycle may be too long to be determined by using the currently available data. Nevertheless, in the case of HD 43587, this result is consistent with the possibility of a primary and a secondary activity cycle (Brandenburg et al. 2017).

3.2.2. Bayesian generalized Lomb-Scargle periodogram with trend (BGLST) analysis

We also performed a periodogram analysis based on the Bayesian generalized Lomb-Scargle method with linear trend (BGLST), as presented in Olspert et al. (2018). When applied to our data, we found a similar period of ∼10 years. The BGLST method found a linear trend from the beginning of the data until the middle and tried to adjust the same tendency until the final fit. This “linear fit” moves the final fit a bit upwards. Thus, the BGLST method returns a result similar to the GLS analysis presented in the previous section. However, instead of removing a long cyclical trend, it results in a periodic signal on top of a linear trend as shown by the magenta line in Fig. 5. In addition, we performed a Gaussian processes (GP) regression with quasi-periodic cycles by using the GPy code (GPy 2014). We did not find any meaningful discrepancies on the overall period found by GLS alone. Figure 5 shows the analysis made with GLS (red line), GP (blue line and shaded blue region), and BGLST (magenta line). The green solid line corresponds to the long trend fit.

Fig. 5. Comparison between the analysis made with the GLS, GP, and BGLST methods. The green line shows the long trend fit, the red line is a sinusoidal representative fit created with GLS, the magenta line is the BGLST fit. The blue line represents the GP mean, and the blue shaded region is the confidence region from the GP fit.

3.3. Stellar rotation

Since the photometric analysis using the CoRoT light curve does not provide a reliable and consistent value for the rotational period of HD 43587, in this section we explore three different alternatives to compute this parameter. First, we directly used the long-term spectroscopic measurements of S_MW over 50 yr to look for a rotation period in the GLS periodogram. We analyzed the time series from 46340.99 to 46525.65 days (JD-2400000), containing 195 data points spread over 184 days, which represents the denser region of data in the time series. The period values found in the GLS were below the 10% false alarm probability, and thus discarded as true rotation period.

$$\begin{array}{c}\hfill P_{\mathrm{rot}}=f(B-V)·g(t).\end{array}$$(5)

Alternatively, we estimated the rotation period using gyrochronology as proposed by Barnes (2009) and Meibom et al. (2009), presented in Eq. (5). It consists of computing both the age-dependence function, g(t)  =  tⁿ (with the age t in Myr) and the mass-dependence function, f(B − V)  =  a(B − V − c)^b. The constant terms are a  =  0.770 ± 0.014, b  =  0.553 ± 0.052, c  =  0.472 ± 0.027 and n  =  0.519 ± 0.007. We used the two ages determined by Castro et al. (2020), inferred from two different evolution codes, as presented in Table 1 and described in detail in Sect. 2.1. The color index used is B − V = 0.61 (Van Leeuwen 2007). Therefore, we found two rotation periods: $P_{\mathrm{rot}}^{\mathrm{T}}=25.0\pm 4.4$ d with TGEC age, and $P_{\mathrm{rot}}^{\mathrm{C}}=22.9\pm 3.7$ d with CESTAM age. From these results, we can deduce a theoretical Rossby number, Ro, which is an important parameter of hydromagnetic dynamo theory. It is defined as the ratio of the rotation period, P_rot, to the convective turnover time, τ_c, which, in turn, is computed as the ratio of the mixing length to the turbulent velocity one pressure height scale above the base of the convection zone (Gilman 1980). From the inferred rotation period with the TGEC model’s age and convective turnover time τ_c = 12.8 ± 1.0 d, we deduce a Rossby number Ro^T = 1.96 ± 0.38. Unfortunately, we had no access to the convective turnover time in CESTAM models and could not infer a Rossby number from them.

In addition, we inferred a rotation period by using the semiempirical calibrations from Noyes et al. (1984), using the time-average observational S_MW detailed in the Sect. 2. We derived the chromospheric contribution, $R_{\text{HK}}^{\prime }$, from the S_MW measurements. It is computed as $R_{\text{HK}}^{\prime }=R_{\text{HK}}-R_{\text{phot}}$, where log R_phot = −4.02 − 1.40(B − V) and R_HK = 1.340C_cf ⋅ S_MW, where C_cf is a color-dependent factor proposed by Middelkoop (1982). Noyes et al. (1984) found a cubic fit relation between the Rossby number and the activity index, $R_{\text{HK}}^{\prime }$ (see Eq. (3) in Noyes et al. 1984). We applied this calculation for each value of S_MW in our dataset, and determined Ro for each of them. The mean Rossby number with its standard deviation is Ro^N = 2.05 ± 0.05. We then determined the convective turnover time by using the empirical relation presented in Eq. (4) of Noyes et al. (1984), which depends on (B − V). With the resulting value, τ_c = 9.67 d, and the definition of the Rossby number Ro = P_rot/τ_c, we derived a mean $P_{\mathrm{rot}}^{\mathrm{N}}=19.8\pm 0.7$ d.

With these results, we infer an average rotation period, ${\overline{P}}_{\mathrm{rot}}=22.6\pm 1.9$ d for HD 43587, which is consistent with the rotation of a solar analog. The good estimate of ${\overline{P}}_{\mathrm{rot}}$ found is reinforced when we examine the ratio of the cycle frequency to the rotation frequency given by ω_cyc/Ω = P_rot/P_cyc as suggested by Brandenburg et al. (2017). In the case of the primary cycle, P_cyc = 10.44 yr, we have $log({\overline{P}}_{\mathrm{rot}}/P_{\mathrm{cyc}})=-2.23$, while the solar value is log(P_rot/P_cyc)_⊙ = −2.20. In the case of the possible secondary cycle, P_cyc ≥ 50 yr, we found $log({\overline{P}}_{\mathrm{rot}}/P_{\mathrm{cyc}})\le -2.87$. For the Sun, the long Gleissberg cycle (GC) gives a ratio log(P_rot/P_cyc)_⊙, GC = −3.10. If we examine the position of HD 43587 in Fig. 5 of Brandenburg et al. (2017), which plots this ratio as a function of $log\langle R_{\text{HK}}^{\prime }\rangle$, we can see that HD 43587 lies on the inactive branch region, highlighting that the $log\langle R_{\text{HK}}^{\prime }\rangle =-4.97$ value is below the Vaughan-Preston gap ($log\langle R_{\text{HK}}^{\prime }\rangle =-4.75$).

Finally, another interesting fact concerning HD 43587 is that its Rossby number, as determined above, is close to the critical Rossby number, Ro ∼ 2, beyond which the efficiency of magnetic braking is dramatically reduced (Van Saders et al. 2016; Metcalfe et al. 2016). This region corresponds to an activity level $log{R}_{\text{HK}}^{\prime }~-4.95$. This suggests that HD 43587 lies in a zone below the Vaughan-Preston gap where the differential rotation, key mechanism of the dynamo process, becomes inefficient (Metcalfe et al. 2016).

4. Simulating HD 43587 with the EULAG-MHD code

A relevant tool currently available for studying stellar magnetism is global numerical modeling. Several recent works have considered this approach to study the magnetic activity of solar-like stars rotating at different rates (e.g., Brown et al. 2011; Augustson et al. 2013; Strugarek et al. 2017; Warnecke 2018; Viviani et al. 2018; Guerrero et al. 2019b). The results have provided a new understanding of the relation between stellar rotation and magnetic activity. To substantiate our observational analysis we also perform numerical simulations of an atmosphere that mimics the stratification of HD 43587. Our goal is not to make a detailed analysis of the simulations, but to get some insights about the internal activity of this star providing support to the observational conclusions.

4.1. Description of the models

We simulated the internal dynamics of the star HD 43587 with the EULAG-MHD code (Smolarkiewicz & Charbonneau 2013). The fundamental difference between the models presented here and previous solar global dynamo simulations is the adoption of a fiducial fit for its stratification. Thus, structure parameters as the location of the tachocline or the velocity field in the convection zone are properly captured. The code solves the anelastic equations of Lipps & Hemler (1982) extended to the MHD case and in the inviscid regime. One relevant point distinguishing EULAG-MHD simulations from global dynamo modeling with other codes is the use of the implicit subgrid scale (SGS) method MPDATA. It allows us to keep the numerical viscosity sufficiently low in such a way that it is possible to resolve the transition between convective and radiative zones without much transfer of angular momentum. Within the convection zone, the numerical dissipation introduced by the MPDATA scheme plays the role of the SGS contribution to the large scale motions and magnetic field.

We consider a spherical shell, 0 ≤ φ ≤ 2π, 0 ≤ θ ≤ π, radially spanning from r = 0.6R to r = 0.96R, and a mesh resolution of 128 × 64 × 48 grid points in φ, θ and r, respectively. The full set of equations are described in Guerrero et al. (2016). The boundary conditions for the velocity field are impermeable and stress free at both boundaries. For the potential temperature perturbations we consider zero radial derivative of the radial convective flux at the bottom, and zero convective radial flux at the top. The magnetic boundary conditions correspond to a perfect conductor at the bottom and pseudo-vacuum at the top. The simulations start with white noise perturbations in the potential temperature, the radial velocity, and the azimuthal magnetic field. Since the main input parameter of the model is the stellar stratification, we present details of its setup below.

As described in Guerrero et al. (2016) and Cossette et al. (2017), two states are necessary in the formulation of convective simulations with EULAG-MHD: the isentropic reference state, and the ambient state. The later allows adiabatic displacements of fluid parcels in the super-adiabatic layers, and inhibits convective motions where the ambient is sub-adiabatic. Both states are described by hydrostatic equilibrium equations:

$$\begin{array}{cc}& \frac{\partial T_{\mathrm{i}}}{\partial r}=-\frac{g}{R_{\mathrm{g}}(m_{\mathrm{i}}+1)},\hfill \end{array}$$(6)

$$\begin{array}{cc}& \frac{\partial {\rho }_{\mathrm{i}}}{\partial r}=-\frac{{\rho }_{\mathrm{i}}}{T_{\mathrm{i}}}(\frac{g}{R_{\mathrm{g}}}-\frac{\partial T_{\mathrm{i}}}{\partial r}),\hfill \end{array}$$(7)

where the index i stands either by s, for the reference state, or by e, for the ambient state. In the reference state, the polytropic index m_s is equal to 1.5. The ambient stratification is a bi-polytropic fit to the mixing length (ML) stratification model obtained with the evolutionary code TGEC as described by Castro et al. (2020).

Our fit takes into account the recent results of Guerrero et al. (2019b), indicating that the interface between the radiative and convective zones might play a relevant role in the dynamo generating mechanism. Their results also show that the magnetic field is stored in the stable layer where it is prone to magnetoshear instabilities (see e.g., Miesch et al. 2007), in turn developing helical non-axisymmetric motions and currents (source of magnetic field via the so-called α-effect). As demonstrated in Guerrero et al. (2019a), these instabilities are highly dependent on the stratification of the atmosphere. For this reason we calibrate the ambient state to reproduce the Brunt–Väisälä frequency profile of TGEC in the stable stratified layer. Since the motions in the convection zone define the development of the differential rotation, in this region we adjust the polytropic index allowing for a radial profile of the turbulent velocity compatible with that of the TGEC model (see Fig. 7). Thus, the quantity m_e in Eq. (6) is given by

$$\begin{array}{c}\hfill m_{\mathrm{e}}(r)=m_{\mathrm{rz}}+\frac{\mathrm{\Delta }m}{2}[1+\mathrm{erf}\left(\frac{r-r_{\mathrm{tac}}}{w_{\mathrm{t}}}\right)]\end{array}$$(8)

where m_rz = 2.0 is the polytropic index of the radiative zone, Δm = m_cz − m_rz, with m_cz = 1.49989 the polytropic index in the convective zone, and w_t = 0.1R the width of transition between both regions. The radiative zone extends until r_tac = 0.76R, the convective zone corresponds to the upper envelope. The profiles of temperature and density are shown in Fig. 6. The yellow solid and dashed lines depict the ML stratification from TGEC model, and the blue and red solid (dashed) lines correspond to the fitted density (temperature) profile for the ambient and reference states, respectively.

Fig. 6. Radial profiles of the temperature, T (dashed lines), and density, ρ (solid lines). The blue and red lines correspond to the ambient and isentropic states, respectively. The yellow lines show the TGEC stratification.

We performed three simulations with the same stratification and different rotational period. They are labeled as P29, P25, and P21 for 29, 25, and 21 days, respectively. The rotation periods cover the uncertainties in the determination of P_rot from the analysis presented in Sect. 3.3. In Fig. 7, we compare the convective velocity from TGEC model with the turbulent velocity, u_rms obtained in the three global simulations performed with EULAG-MHD. Since rotation constrains the convective motions, the amplitude of u_rms decreases with the decrease of P_rot. All three velocity profiles exhibit, however, the same trend. The values of u_rms volume averaged in the convection zone are presented in Table 2. The Rossby number presented in the same table was computed as it was in Gilman (1980), with τ_c computed from the pressure scale height and the resulting u_rms one pressure scale height above the star’s tachocline. At this radial point, the turbulent velocity in simulations P25 and P21 is smaller than the one predicted by the mixing length model. Thus, Ro is smaller than the one obtained from the TGEC model for the same rotational period, but still within the determination error of Ro^T.

Fig. 7. Radial profile of the turbulent velocity for simulations P29 (red), P25 (yellow) and P21 (green). The blue line shows the TGEC stratification model.

Also in contrast with TGEC, the turbulent velocity in the simulations is different from zero in the radiative zone. These motions have two origins, overshooting of the downward flows coming from the convection zone, and non-axisymmetric motions developed in the stable layer due to MHD instabilities. As we present below, the realistic stratification of HD 43587 captured in the ambient state allows us to reproduce some aspects of the magnetism of this star which are compatible with the observational results obtained in the previous sections.

4.2. Simulation results

Figure 8 shows the differential rotation profiles resulting from the simulations (a) P29 to (c) P21. The mean angular velocity, $\overline{\mathrm{\Omega }}$, is an average on longitude and time (along the last two magnetic cycles of each simulation). All the cases show a clear formation of a radial shear layer at the base of the convection zone. For the three rotational periods studied, the resulting differential rotation is solar-like, with slow poles and a fast equator. This kind of differential rotation appears when the timescale of the convective flow is similar to the rotational one, meaning rotationally constrained motions. In this case, the transport of the angular momentum is directed towards the equator (see e.g., Featherstone & Miesch 2015; Guerrero 2020).

Fig. 8. Mean differential rotation in the meridional plane from simulations (a) P29, (b) P25 and (c) P21. The color scale shows the angular velocity with respect to the rotation frame. The dashed lines depict the tachocline location.

In the simulations P29−P21, we notice that the latitudinal gradient of angular velocity diminishes as the rotational period decreases (Figs. 8a–c) This is quantitatively verified in the third column of Table 2 showing the relative differential rotation, ΔΩ/Ω_eq, where ΔΩ = Ω_eq − Ω₄₅. In all the cases, this value is smaller than the solar one, ∼0.1, and about one order of magnitude smaller than in other stars where it was inferred by astroseismology (Benomar et al. 2018).

In Fig. 9, the colored contours depict the radial magnetic field averaged over longitude, ${\overline{B}}_{\mathrm{r}}$, and the solid (dashed) contour lines represent the positive (negative) mean toroidal magnetic field, ${\overline{B}}_{\phi }$. In general, the simulations result in oscillatory magnetic fields with polarity reversals. Like in Guerrero et al. (2019b), the results suggest the existence of a α² − Ω dynamo operating inside HD 43587. The magnetic fields are anti-symmetric across the equator. For slow rotation the field seems further and further concentrated near the poles. The averages in volume and time of the toroidal and radial fields are also presented in Table 2. We chose to present $⟨{\overline{B}}_{\mathrm{r}}⟩$ instead of the poloidal field because the radial field is directly observed in the Sun.

Fig. 9. Time latitude at r = 0.90R (upper panels), and time radius at θ = 15° latitude (lower panels), butterfly diagrams comparing the magnetic field evolution between simulations (a) P29, (b) P25, and (c) P21. The colored contours show the radial magnetic field intensity. The solid (dashed) lines show the positive (negative) contours of toroidal magnetic field.

The magnetic cycle periods, P_cyc, were computed by Fourier analysis (see Table 2). The values, between 6.4 and 8.6 years are slightly smaller than that of 10.4 years obtained in the chromospheric activity analysis. Nevertheless, the cases for 21 and 25 days are still within the 1σ error of 3.3 years. In Fig. 10 we show the results of the simulations in the P_rot versus P_cyc plane. The figure contains the observational points of F, G, and K stars as presented in Brandenburg et al. (2017). The rotational and magnetic cycles of HD 43587 are presented with a black star. The periods of the three simulations fall in the inactive branch of activity closer to the position of the Sun. The agreement between observation and simulation is encouraging, and indicates that the simulations are capturing the physics of the interior of the star.

Fig. 10. Correlation between Pcyc and Prot for the simulations presented in Table 2 (green points). The blue and red points correspond to F, G, and K star data as presented in Brandenburg et al. (2017). The black star shows the position of HD 43587 with the values inferred in this work, and the yellow point marks the position of the Sun.

In the recent work of Olspert et al. (2018), the results of novel analysis of the MW data suggest that the active branch presented in Fig. 10 may not exist. Instead, they found a branch merging the active and transitional branches. Although this might be the case, we consider that there are not yet sufficient stars observed with fiduciary determination of their periods. Thus, while the activity-rotation problem remains a subject of discussion (Boro Saikia et al. 2018), we decided to present our results together with the canonical A and I branches.

The amplitude of the toroidal magnetic field in the generating radial shear layer is about 0.4 Tesla and is larger in the simulations with a larger period. The radial field component, on the other hand, is about 0.01 Tesla and is larger in the simulations with a shorter period (fast rotating). Figure 9 indicates that the magnetic field migrates upwards, however, at the upper boundary of the model, r = 0.95 R_⋆, both, the toroidal and the radial field components decay about one order of magnitude. Furthermore, the simulation results do not exhibit relevant dynamo action near the upper boundary. Thus, it is expected that it will decrease further when reaching the star surface, which is not considered in the model.

In summary, the simulations are able to reproduce cyclic dynamo action, including field reversals, with a differential rotation smaller than in the Sun, and a magnetic field, which at the top boundary is similar to the radial field observed at the solar surface (∼10⁻³ Tesla). The 10% difference observed in ⟨S_MW⟩ between the Sun and HD 43587 might be correlated to a diminished differential rotation, which in the Sun is at least twice as large as in simulations P25 and P21. These cases are the ones that better reproduce the cycle period of HD 43587 inferred from the spectroscopic chromospheric activity.

5. Conclusions and summary

The spectral analysis is fundamental to understanding astrophysical processes and stellar evolution. Currently, the stellar magnetic activity is one of the most important observables in the context of stellar interiors, and it might impact exoplanet detection as shown by Aigrain et al. (2016). In this work, from the high-resolution spectra archive and already published data, we investigated the chromospheric activity evolution of the solar analog star HD 43587.

We built the larger S_MW time series already published for this star, using data spanning 50 years from 1966 until 2016. The derived ⟨S_MW⟩ = 0.154 and $log\langle R_{\text{HK}}^{\prime }\rangle =-4.97$ indicate a low chromospheric activity, beyond the Vaughan-Preston gap. From the GLS analysis of the S_MW time series, we found a reliable 10.4 yr cycle period (see Fig. 4), which is similar to the solar value. The indirect determinations of a rotation period give ${\overline{P}}_{\mathrm{rot}}=22.6\pm 1.9$ d. These values of P_cyc and ${\overline{P}}_{\mathrm{rot}}$, are consistent for a solar analog star. The inferred Rossby number, Ro^N = 2.05 ± 0.05 from the calibrations prescribed by Noyes et al. (1984), and Ro^T = 1.96 ± 0.38 from TGEC modeling, point out a rotation period about two times larger than the convective turnover time which might represent a diminished differential rotation (Metcalfe et al. 2016).

Based on the inferred rotational period and the internal structure of HD 43587, we performed global simulations with the EULAG-MHD code. We considered three different rotational periods around the value obtained in Sect. 3.3. For all the cases, the angular velocity is solar-like, i.e., fast equator and slow poles, with relative differential rotation smaller than the solar value. In principle, the smaller shear might be responsible for the low levels of activity of the simulated stars. The radial field at the top of the simulations (at r = 0.96 R_⋆) is about 10⁻² T. This value is similar to the radial field at the solar photosphere.

The magnetic field is oscillatory with magnetic field reversals. The magnetic cycle periods obtained in the simulations range from 6.4 yr for the simulations with a rotational period of 29 days, to 8.6 yr for the simulation with P_rot = 21 days. This cycle period is compatible with the one inferred in Sect. 3.2.

It is worth mentioning that given the difficulties of finding the rotation period from the light curves, the indirect determination of the P_rot is the best approach to to characterize the star HD 43587. The fact that the inferred rotation period fits with a solar analog encourages us to presume that it might be correct. Furthermore, ab initio simulations of the star with initial parameters given by the stellar structure (the same as used to determine the rotational period) as well as the inferred P_rot, result in activity cycles close enough to that independently determined from the S-index. These numerical results substantiate the observational ones.

Even though the possibility that HD 43587 is under an MM cannot be ruled out, all our results seem to indicate that it is a solar analog star, older than the Sun, and undergoing a “natural” decrease in its magnetic activity. Nevertheless, some authors argue that a transition in the stellar differential rotation might be the cause of these low levels of activity (Metcalfe & van Saders 2017). The characterization of an MM would only be valid if we could observe the star coming in and out of this activity minimum (Hall 2008). There is no current observation indicating this pattern.

The low activity level observed in HD 43587 might be due to a small differential rotation in the convective zone. We have also found evidences that HD 43587 could be at the minimum of a longer activity cycle of over 50 years. However, in the context of observations of magnetic cycles in solar-like stars, the existence of long activity cycles for old objects is rather controversial. The results presented by Brandenburg et al. (2017) show that, besides the Sun, all confirmed objects that have co-existing long and short cycles are younger than 2.3 Gyr. As a matter of fact, the reality of double periodicities in the observed stellar activity cycles is still questioned within the community (Boro Saikia et al. 2018; Olspert et al. 2018).

The results presented above demonstrate the relevance of studying HD 43587 in the context of stellar evolution. Several authors have investigated its features and evolutionary status, as well as its magnetic activity evolution. Its similarity with the Sun, whether evolutive, spectroscopic, or magnetic, gives insight to the future evolution of the Sun and makes it a good candidate for being an Earth-like exoplanet host. Its low activity level provides an excellent opportunity to attempt a signal of an earth-like planet without the interference of the magnetic activity (Korhonen et al. 2015). In this way, we estimate that it is fundamental to keep monitoring the activity of HD 43587 to confirm the co-existence of a short and long period cycle, as well as to increase the accuracy of the cycle and rotation period determination. We also consider that it should be a key target for the ESA PLATO mission.

Acknowledgments

We thank the anonymous referee for the useful comments and suggestions that improved this paper. Research activities of the Ge³ (http://ge3.dfte.ufrn.br) stellar team at the Federal University of Rio Grande do Norte are supported by continuous grants from Brazilian scientific promotion agencies CAPES and CNPq. R.R.F. and R.B. acknowledge CAPES and CNPq graduate fellowships. J.D.N. and M.C. acknowledge support from CNPq (Bolsa de Produtividade). We also thank T. Morel for providing the HARPS spectra used in this paper. This work made use of the computing facilities of the Laboratory of Astroinformatics (IAG/USP, NAT/Unicsul), whose purchase was made possible by FAPESP (grant 2009/54006-4) and the INCT-A. This work was carried out with the support of the Coordination for the Improvement of Higher Education Personnel, Brazil (CAPES).

References

All Tables

All Figures

	Fig. 1. Sun’s chromospheric activity over 50 years according to Egeland et al. (2017). Upper: complete time series measurements made over 50 years depicting the SMW periodic profile. Bottom: both the Lomb-Scargle (left), and the autocorrelation function (ACF, see Brockwell 2002) (right) periodograms showing the Pcyc ≈ 11 yr solar activity cycle.
In the text

Fig. 2. HD 43587 spectra from NARVAL spectrograph. The magenta spectrum corresponds to the maximum SMW value of the dataset, observed on February 13, 2010. The blue spectrum corresponds to the minimum SMW value, observed on January 18, 2010. Both spectra were normalized from 3889 Å to 4013 Å. An upward shift in flux of 0.05 was applied in the blue spectra to better visualization. The inset panel shows the region between 3930 Å and 3938 Å in detail.

In the text

Fig. 3. CoRoT light curve of HD 43587. Top-left panel: raw light curve including discontinuities (red dots) and the detrended light curve (black dots). Top-right panel: detail of the discontinuities at the end of the light curve (red dots) and the flat aspect of the curve after the reduction (black dots). Middle-left panel: GLS periodogram for the detrended light curve with false alarm probability of 5% and 1% (green dots and blue dashed line, respectively). Middle-right panel: GLS periodogram for raw light curve with false alarm probability the same as the detrended analysis. Bottom panel: ACF periodogram for the same two reductions.

In the text

Fig. 4. Entire spectroscopic chromospheric activity measurements with the Sindex calibrated to the Mt. Wilson Scale for HD 43587 between 1966 and 2016. Upper panel: all measurements from Duncan et al. (1991), Wright (2004), Hall (2008), and Radick et al. (2018), combined with the computed SMW from the NARVAl and HARPS spectra archives. The red solid line is the sinusoidal curve fitting for 10.436 yr, the blue solid line is the long-term trend found to be larger than 50 yr. The yellow dashed line shows the mean SMW for the Sun. Bottom left panel: GLS periodogram of the whole SMW time series (solid gray line), and removing the long trend of over 50 years (solid black line). The Gaussian fit (solid red line) of this second periodogram indicates an activity cycle of 10.436 years (vertical black dashed line). Bottom-right panel: phase of the SMW (black circles) and the folded fit with the found period (red circles).

In the text

	Fig. 5. Comparison between the analysis made with the GLS, GP, and BGLST methods. The green line shows the long trend fit, the red line is a sinusoidal representative fit created with GLS, the magenta line is the BGLST fit. The blue line represents the GP mean, and the blue shaded region is the confidence region from the GP fit.
In the text

	Fig. 6. Radial profiles of the temperature, T (dashed lines), and density, ρ (solid lines). The blue and red lines correspond to the ambient and isentropic states, respectively. The yellow lines show the TGEC stratification.
In the text

	Fig. 7. Radial profile of the turbulent velocity for simulations P29 (red), P25 (yellow) and P21 (green). The blue line shows the TGEC stratification model.
In the text

	Fig. 8. Mean differential rotation in the meridional plane from simulations (a) P29, (b) P25 and (c) P21. The color scale shows the angular velocity with respect to the rotation frame. The dashed lines depict the tachocline location.
In the text

	Fig. 9. Time latitude at r = 0.90R (upper panels), and time radius at θ = 15° latitude (lower panels), butterfly diagrams comparing the magnetic field evolution between simulations (a) P29, (b) P25, and (c) P21. The colored contours show the radial magnetic field intensity. The solid (dashed) lines show the positive (negative) contours of toroidal magnetic field.
In the text

	Fig. 10. Correlation between Pcyc and Prot for the simulations presented in Table 2 (green points). The blue and red points correspond to F, G, and K star data as presented in Brandenburg et al. (2017). The black star shows the position of HD 43587 with the values inferred in this work, and the yellow point marks the position of the Sun.
In the text