A Comparison between Physics-based and Polytropic MHD Models for Stellar Coronae and Stellar Winds of Solar Analogs

The Zeeman–Doppler Imaging (ZDI) technique (Donati & Semel 1990) has provided, for the first time, synoptic observations of surface magnetic fields of low-mass stars (mostly M, K, and G types). Despite some criticism of the uncertainty of these stellar "magnetograms" (e.g., Kochukhov et al. 2010; Reiners 2012), the growth in availability of stellar magnetograms (see Vidotto et al. 2014, for a summary of these observations) has led the stellar astrophysics community to adopt techniques that were used by the solar community to extrapolate the three-dimensional coronal magnetic field. In particular, ZDI data have been used to drive three-dimensional models for stellar coronae and stellar winds (in a similar manner to the way in which solar magnetograms are used to drive models for the solar corona and solar wind).

The first approach adopted from solar physics is the so-called "potential field source surface" (PFSS) method (Altschuler & Newkirk 1969). In the PFSS model, the three-dimensional magnetic field is assumed to be static (i.e., there is no forcing on the field by electric currents), and as such, the field can be described as a gradient of a scalar potential. This scalar potential can be obtained by solving Laplace's equation, assuming that the field is purely radial above a spherical surface—the "source surface"—which is set at a random distance. The potentiality of the solar coronal field has been debated for some time among the solar community (see, e.g., Riley et al. 2006). In particular, the validity of a spherical source surface and the choice of its distance have been challenged frequently (see Cohen 2015, for a recent review on the PFSS issues). Currently, the PFSS method is mainly used to calculate the initial state of the magnetic field in magnetohydrodynamic (MHD) models or to obtain a tentative description of the coronal field. Because of its limitations, the solar community is currently transitioning from using the PFSS method to more sophisticated field extrapolation methods, such as the nonlinear force-free technique (see the review by Wiegelmann & Sakurai 2012). In contrast, the PFSS method is used more frequently to extrapolate stellar coronal fields based on the ZDI maps (see, e.g., Jardine et al. 1999).

Following the pioneering coronal modeling work by Pneuman & Kopp (1971), modern, multi-dimensional MHD models for the solar corona have first applied a thermally driven, polytropic Parker wind (Parker 1958) on top of a potential field (Linker et al. 1990; Usmanov 1993; Linker & Mikic 1995). The Parker wind solution is obtained analytically assuming a (nearly) isothermal, radially expanding flow. A steady state is obtained when a pressure balance between the Parker wind and the potential field is achieved. Despite confirmation from spacecraft of the existence of a super-Alfvénic solar wind, the Parker wind and its implementation in these older MHD models have failed to reproduce the observation. In particular, they could not produce the observed fast solar wind. Thus, we should account for additional acceleration missing from the Parker model.

Holzer (1977) and Usmanov et al. (2000) have pointed out that in order to properly reproduce the solar corona and the solar wind, MHD models need to include additional momentum and energy terms beyond the set of ideal MHD equations. While the physical interpretation of these terms is still under debate, they account for the observed coronal heating and wind acceleration. In an intermediate stage, the polytropic models have been extended to use a varying polytropic index, γ, as a function of some local properties of the gas (e.g., Roussev et al. 2003; Cohen et al. 2007; Feng et al. 2010; Jacobs & Poedts 2011). In particular, models were developed to relate the local value of γ to the observed empirical relation between the solar wind speed and the expansion geometry of the magnetic flux tube (Wang & Sheeley 1990; McComas 2007). These models successfully reproduced the observed solar wind. However, due to the lack of detailed thermodynamics, they did not provide good agreement with observations of the solar corona itself (i.e., comparison with EUV/X-ray line-of-sight (LOS) observations). Recent MHD models for the solar corona and solar wind adopt a physics-based approach and include additional energy and momentum sources, as well as thermodynamic terms (heating and cooling), and a multi-fluid approach to include the separation between ions and electrons (e.g., most recent models by Lionello et al. 2014a, 2014b; van der Holst et al. 2014; Usmanov et al. 2016).

Despite this evolution in modeling the solar corona and wind (and perhaps due to some possible disconnection of the astrophysics community from solar physics literature) recent models for stellar coronae of solar analogs have returned to using the spherical polytropic Parker wind imposed on a potential magnetic field. These models have used available ZDI maps or low-order idealized magnetic fields to study and scale stellar mass-loss rates and stellar spindown with stellar age and rotation periods, as well as to characterize the coronae and winds of planet-hosting stellar systems (e.g., Matt et al. 2012; Vidotto et al. 2014; Matsakos et al. 2015; Réville et al. 2015; Strugarek et al. 2015), while they were not calibrated by solar input parameters and solar magnetogram data. Instead, some of the models have adopted an idealized "fiducial" solar case, which is represented by solar mass, radius, rotation period, and a dipole magnetic field.

The goal of this paper is to demonstrate that the polytropic Parker wind imposed on a potential field struggles to reproduce solar observations. We use high-resolution solar magnetograms to drive both a polytropic Parker wind MHD model and a physics-based MHD model for the solar corona and solar wind, and compare the output from each model with solar observations. In the next section we describe the models and the chosen input parameters, we present and discuss the results in Section 3, and conclude our findings in Section 5.

In this study, we compare the solutions of two types of MHD models, both performed using the BATS-R-US MHD code (Powell et al. 1999; Tóth et al. 2012).

2.1. Thermally Driven Polytropic MHD Model

2.2. The AWSOM Model

2.3. Model Parameters

In order to keep the two models as consistent with each other as possible we set them both with the same magnetogram input data, the same initial conditions for the three-dimensional magnetic field (a PFSS extrapolation), the same initial density structure, and the same initial condition for the wind speed (a Parker solution with T = 3 MK). We also specify the boundary conditions to be the same in both models, with the coronal number density n = 10⁸ [cm⁻³] and the coronal temperature T = 3 MK. We use this high coronal temperature to obtain an upper limit to the wind speed in the polytropic model. Figure 1 shows the initial conditions for the number density and wind speed.

Zoom In Zoom Out Reset image size

Despite the different grid geometries used in the two settings, we trust that the comparison presented here is valid because we are interested in comparing the overall wind acceleration via global properties, such as the maximum wind velocity and total mass-loss rate, and because the grid size near both the inner and outer boundaries is comparable in both settings. For the polytropic setting, we test three cases with γ = 1.01, 1.05, and 1.1. We also test how the coronal base density affects the solution by performing one case with γ = 1.05 and n = 10⁷ [cm⁻³].

2.4. Input Data

Figure 2 shows the results for the wind speed and temperature from the polytropic model with γ = 1.01, 1.05, 1.1, and from AWSOM. The wind speed in the polytropic solutions is quite uniform, with speeds ranging from 100 km s⁻¹ to about 500 km s⁻¹ for γ = 1.01, 100 km s⁻¹ to about 400 km s⁻¹ for γ = 1.05, and 100 km s⁻¹ to about 300 km s⁻¹ for γ = 1.1. The fast/slow wind contrast is much more visible in the AWSOM solutions, with clear regions of slow wind with a speed of about 300–400 km s⁻¹, and regions of fast wind with a speed above 600 km s⁻¹.

Zoom In Zoom Out Reset image size

The MHD solutions presented here use single-fluid plasma, which assumes that the ion temperature, T_i, and the electron temperature, T_e, are equal, and that the plasma temperature is ${T}_{p}=({T}_{i}+{T}_{e})/2=p/2\rho$ (where p and ρ are the simulated pressure and density, respectively). The plasma temperature is half the sum of the electron and ion temperatures in the initial Parker wind and the coronal base temperature. Here we show this reduced temperature as it represents a more realistic coronal temperature (especially in the context of EUV/X-ray observations described below).

The temperature plots show, as expected from an isothermal solution, an almost uniform coronal temperature of 1.5 MK for the γ = 1.01 case. For the γ = 1.05 and γ = 1.1 cases, the temperature contrast is more notable and it follows the density and magnetic field structure. In particular, the temperature contrast clearly follows the helmet streamers (the closed field regions) in the solar minimum case. In the AWSOM solutions, the temperature contrast is very clear, and the temperature is much higher (over 2 MK) in the helmet streamers. The temperature does not exceed the prescribed 1.5 MK in the polytropic solutions.

Figure 3 shows a comparison between full-disk images of the Sun from the Extreme Ultraviolet Imaging Telescope (EIT) on board the Solar and Heliospheric Observatory (SOHO)³ and from the Yohkoh Soft X-ray Telescope (SXT),⁴ and synthetic images produced by AWSOM. The EIT bands are 171, 195, and 284 Å, and the X-ray images are integrated over the range 2.4–32 Å. The synthetic images are obtained by integrating the square of the electron density, n_e, multiplied by a response function for a particular temperature bin, Λ(T), along the LOS through the three-dimensional solution. The response functions are calculated from the CHIANTI atomic database (see, e.g., Dere et al. 1997; Landi et al. 2013) in a similar way to that in which the actual observed images are produced, taking into account the particular instrument calibration (see Downs et al. 2010; Oran et al. 2015, for more details about the production of synthetic LOS images in BATS-R-US).

Zoom In Zoom Out Reset image size

The comparison in Figure 3 shows a very good agreement between the observations and the images produced from AWSOM. The location and overall structure of the active regions is well reproduced in the modeled images, as well as the location and size of the coronal holes (the dark regions associated with the lower density in the open field regions). In the solar minimum case, the active region at the center of the disk is reproduced, along with the active region close to the right limb. The active region in the center is blank in the 171 Å band, probably because the model overheated the area with a temperature above the one that that particular line responds to. For solar maximum, the overall structure of the belts in the active regions is reasonably reproduced, along with the coronal hole boundaries. We performed a similar procedure to produce synthetic EIT/SXT images from the polytropic models. However, due to the lower coronal plasma temperature and the low temperature contrast all the LOS images are blank and do not show any features.

The usage of ZDI maps to drive MHD models for the coronae and winds of cool stars is important in the context of stellar spindown and stellar evolution. Table 1 summarizes the maximum wind speed, the total mass-loss rate, and the total angular momentum loss rate of each of the solutions.

Table 1.  Simulations: Global Parameters of the Solar Corona

Case	CR 1916 Parker	CR 1916 Parker	CR 1916 Parker	CR 1962 Parker	CR 1962 Parker	CR 1962 Parker	CR 1916 Parker	CR 1916 Parker	CR 1916 AWSOM	CR 1962 AWSOM
	γ = 1.01	γ = 1.05	γ = 1.1	γ = 1.01	γ = 1.05	γ = 1.1	γ = 1.05 n = 107	γ = 1.05 n = 5 × 107
Maximum Speed (km s−1)	460	350	275	518	440	520	1400	367	680	770
$\dot{M}\,({10}^{-14}\,{M}_{\odot }\,{\mathrm{yr}}^{-1})$	9.1	5.3	1.5	6.8	5.2	0.9	0.4	2.0	1.3	1.6
$\dot{J}\,({10}^{28}\,{\rm{g}}\,{\mathrm{cm}}^{2}\,{{\rm{s}}}^{-2})$	9.5	5.5	1.8	6.8	5.2	2.4	0.5	2.1	1.1	2.7

Download table as:  ASCIITypeset image

The values from Table 1 show that the polytropic models cannot produce consistent agreement with the observed properties of the Sun. The solutions with γ = 1.05 and 1.1 cannot produce the fast (above 600 km s⁻¹) wind, while the γ = 1.01 solution does produce a faster wind but it consistently overestimates the observed solar mass-loss rate of about ${\dot{M}}_{\odot }\approx (2\mbox{--}3)\times {10}^{-14}\,{M}_{\odot }\,{\mathrm{yr}}^{-1}$ (Cohen 2011). Reducing the boundary density of the γ = 1.05 solution from 10⁸ to 5 × 10⁷ [cm⁻³] does not produce fast wind either, while the overestimation of the mass-loss rate is reduced. Reducing the boundary density even further to 10⁷ [cm⁻³] leads to a much faster wind, but too low a mass-loss rate of 0.4 × 10⁻¹⁴ M_⊙ yr⁻¹. The AWSOM solution produces fast winds above 600 km s⁻¹, and mass-loss rates of (1.3–1.6) × 10⁻¹⁴ M_⊙ yr⁻¹, within a factor of 2 or so of the observed one. The two AWSOM solutions for the two solar epochs produce consistent values for the mass-loss and angular momentum loss rates, while the polytropic solutions differ slightly between solar minimum and solar maximum. Surprisingly, the polytropic solutions overestimate the mass-loss rate but underestimate the angular momentum loss rate for γ = 1.01 and 1.05 compared to the AWSOM solutions. The mass-loss rate is comparable between the γ = 1.1 and the AWSOM solutions, but the angular momentum loss rate is overestimated by the former.

The trends mentioned above can be explained by the difference in size of the Alfvén surfaces in the different solutions (shown in Figure 2). For the γ = 1.01 and 1.05 cases, the Alfvén surface is smaller than that in the AWSOM solution. Thus, the overall mass loss is higher due to the higher density on the surface (closer to the solar surface), but the angular momentum loss rate is smaller due to the shorter lever arm that applies a torque on the star. The Alfvén surface is slightly bigger in the γ = 1.1 solutions than in the AWSOM solutions. This explains the similar mass-loss rate but slightly larger angular momentum loss rate. Nevertheless, the wind speeds and temperatures of the γ = 1.1 solutions do not agree with solar observations.

An explanation for these trends can also be found in Figure 4, where we show the wind's radial speed, number density, and temperature extracted along an open field line with the same footpoint location for all the different solutions. Figure 4 shows that the wind is accelerated much faster in the AWSOM solution to higher values, whereas the polytropic speeds do not exceed 400 km s⁻¹. It is only in the AWSOM solution that the temperature first rises and then falls adiabatically, while in the polytropic solutions it only falls. Finally, the density drops much faster in the AWSOM solutions than in the polytropic solutions. Figure 5 shows a similar extraction along the same field line for the Alfvén Mach number, M_A, and the local mass-loss rate, ${{dM}}_{i}={\rho }_{i}{u}_{i}4\pi {r}_{i}^{2}$, where the index i stands for a particular point along the field line. The plots show that the wind exceeds the Alfvén speed in the polytropic solutions with γ = 1.01 and 1.05 at a height 4–5 solar radii lower than in the AWSOM solutions, where M_A = 1 around r = 14–15 R_⊙. The wind exceeds the Alfvén speed at higher radii for the γ = 1.1 polytropic solution than for the AWSOM solution. The trends of the local mass-loss rate are similar to those of the density, where dM is 5–10 times higher at the location of M_A = 1 for the polytropic solutions than for the AWSOM solution. This explains the overestimation of the mass-loss rate for γ = 1.01 and 1.05 compared to AWSOM. This trend of a slower decline of the local mass-loss rate is compensated by the further Alfvén point in the γ = 1.1 solution. Therefore, this solution does not overestimate the mass-loss rate like the the solutions with lower values of γ do, because it accounts for a lower mass-loss rate at a greater distance.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Leer & Holzer (1980) found that if the energy deposition occurs below the Alfvén point it affects mostly the mass flux, without affecting the final wind speed much. Conversely, they found that energy deposition above the Alfvén point affects mostly the final wind speed, without affecting the mass flux much. In our simulations, the Alfvén point seems to be closer to the surface in the polytropic model than in the AWSOM model. This result might be due to the fact that the polytropic model is set uniformly, and it is driven and constrained exclusively by the boundary conditions. Thus, the polytropic model is set by a small number of global constraints. The AWSOM model, on the other hand, accounts for the heating and acceleration at each point of the domain individually, and the sources and sinks of energy and momentum are defined in a local manner, while also taking into account more detailed processes. This difference seems to make a significant difference to the location of the Alfvén point. Following Leer & Holzer (1980), it is possible that, since the Alfvén point in the polytropic model is rather low, the mass flux cannot be regulated much above it, resulting in an overestimation of the mass-loss rate. The resulting relatively high density above the Alfvén point (compared to the AWSOM density at similar heights) prevents the acceleration of the plasma to the high speeds obtained by the AWSOM model.

The results show that the AWSOM model agrees with all elements of the observed solar properties—the wind speed, coronal temperature and density structure, and total mass-loss rate (the total observed angular momentum loss rate is harder to determine). None of the polytropic solutions can provide such a consistency with solar observations, and some of the solutions might even be considered "unphysical," due to a maximum wind speed that is lower than the minimum observed solar wind speed. The choice of the number density, n, at the coronal base in the models is crucial to get agreement with both the observed total mass-loss rate and the observed three-dimensional coronal density. Our work suggests that this base density should be about n = 10⁸ cm⁻³, which is lower than the choice in other polytropic models for stellar coronae (e.g., Réville et al. 2015; Vidotto et al. 2015).

In order to extend this work to the stellar context, we perform similar simulations to the planet-hosting star HD189733 using a ZDI map which is taken from the data published by Fares et al. (2010). In the simulation of HD189733, we use the same parameter setting as in the solar runs. Figure 6 shows the solution for HD189733, which shows overall a similar range of wind speeds and coronal temperatures to the solar results. The figure also shows the three-dimensional coronal magnetic field, colored with temperature contours with hotter closed field lines (shown in red) and colder open field lines (shown in green/yellow). We also show synthetic EIT/SXT images of HD189733, where a good agreement is obtained between the location of the colder open field lines and the coronal holes in the LOS images. Table 2 shows the global parameters of HD189733 using the γ = 1.01 polytropic model and AWSOM. The mass-loss rate is overestimated by the polytropic model compared to AWSOM, but the angular momentum loss rate is comparable. However, just as in the solar case, the polytropic model cannot produce a fast wind. The mass-loss rate for HD189733 is similar to the solar one in the AWSOM model, but the angular momentum loss rate is much higher, probably due to the difference in rotation period (about 12 days compared to the solar rotation period of 25 days).

Zoom In Zoom Out Reset image size

Table 2.  Simulations: Global Parameters of HD189733

Case	Parker γ = 1.01	AWSOM
Maximum Speed (km s−1)	500	946
$\dot{M}\,({10}^{-14}\,{M}_{\odot }\,{\mathrm{yr}}^{-1})$	4.1	1.6
$\dot{J}\,({10}^{28}\,{\rm{g}}\,{\mathrm{cm}}^{2}\,{{\rm{s}}}^{-2})$	4.1	1.7

Download table as:  ASCIITypeset image

We perform a comparative study between a polytropic wind MHD model and the physics-based AWSOM model in order to test which model agrees better with solar observations. The AWSOM model produces a clear bimodal solar wind, with fast wind above 600 km s⁻¹, a good agreement (within a factor of 2) with the observed total solar mass-loss rate, and a good agreement with the coronal temperature and density structure (validated by solar full-disk LOS observations). The polytropic wind model, with three choices of the polytropic index γ, fails to produce solutions that are consistent with solar observations. In particular, the polytropic model cannot produce the fast solar wind, and the coronal density drops too slow in this type of model. This leads to inconsistency with the observed total mass-loss rate, the solar wind speed, and the coronal base density.

Our recommendation is that polytropic models, which are used to estimate mass-loss rates and other parameters of solar analogs, must be first calibrated with solar observations. Alternatively, these models can be calibrated with models that capture more detailed physics of the solar corona (such as the AWSOM model) and that can reproduce solar observations in a consistent manner. Without such a calibration, the results of the polytropic models cannot be validated, but they can be wrongly used by others.

The work presented here was funded by a NASA Living with a Star grant number NNX16AC11G. Simulation results were obtained using the Space Weather Modeling Framework, developed by the Center for Space Environment Modeling, at the University of Michigan with funding support from NASA ESS, NASA ESTO-CT, NSF KDI, and DoD MURI. The simulations were performed on the Smithsonian Institute HYDRA cluster.