Formulating Mass-loss Rates for Sun-like Stars: A Hybrid Model Approach

We observe an enhanced stellar wind mass-loss rate from low-mass stars exhibiting higher X-ray flux. This trend, however, does not align with the Sun, where no evident correlation between X-ray flux and mass-loss rate is present. To reconcile these observations, we propose a hybrid model for the stellar wind from solar-type stars, incorporating both Alfvén wave dynamics and flux emergence-driven interchange reconnection, an increasingly studied concept guided by the latest heliospheric observations. For establishing a mass-loss rate scaling law, we perform a series of magnetohydrodynamic simulations across varied magnetic activities. Through a parameter survey concerning the surface (unsigned) magnetic flux (Φsurf) and the open-to-surface magnetic flux ratio (ξ open = Φopen/Φsurf), we derive a scaling law of the mass-loss rate given by Ṁw/Ṁw,⊙=Φsurf/Φ⊙surf0.52ξopen/ξ⊙open0.86 , where Ṁw,⊙=2.0×10−14M⊙yr−1 , Φ⊙surf=3.0×1023Mx , and ξ⊙open=0.2 . By comparing cases with and without flux emergence, we find that the increase in the mass-loss rate with the surface magnetic flux can be attributed to the influence of flux emergence. Our scaling law demonstrates an agreement with solar wind observations spanning 40 yr, exhibiting superior performance when compared to X-ray-based estimations. Our findings suggest that flux emergence may play a significant role in the stellar winds of low-mass stars, particularly those originating from magnetically active stars.

The primary focus of this study is the physical properties of the stellar wind, which play significant roles in the evolution of the stellar system as follows.The magnetized stellar wind transports angular momentum outward through magnetic braking (Kraft 1967;Weber & Davis 1967;Mestel 1968;Sakurai 1985;Kawaler 1988), resulting in the slowing down of the stellar rotation (Skumanich 1972;Barnes 2003;Irwin & Bouvier 2009;van Saders et al. 2016).The stellar wind also influences planetary evolution due to its non-negligible impact on planetary outflows (Bourrier & Lecavelier des Etangs 2013;Kislyakova et al. 2013;Cohen et al. 2015;Bourrier et al. 2016;Kislyakova et al. 2019;McCann et al. 2019;Vidotto & Cleary 2020;Mitani et al. 2022).It has been suggested that a portion of Earth's water may have originated from the solar wind (Daly et al. 2021).To quantitatively understand these processes, the mass-loss rate of the stellar wind (denoted as Ṁw ) is an essential parameter.Consequently, the scaling law of Ṁw serves as a crucial building block in the study of the stellar wind and its effects on the evolution of the whole stellar system.
The difficulty in measuring mass-loss rates has resulted in a limited number of observed values, which are indirectly obtained through methods such as astrospheric Lyman-alpha absorption (Wood et al. 2002(Wood et al. , 2005(Wood et al. , 2014(Wood et al. , 2021)), slingshot prominences (Jardine & Collier Cameron 2019;Jardine et al. 2020), and planetary responses to the stellar wind (Vidotto & Bourrier 2017;Kavanagh et al. 2019).It is important to note that direct observations using radio waves (Lim & White 1996;Lim et al. 1996;Gaidos et al. 2000;Fichtinger et al. 2017;Vidotto & Donati 2017) and X-rays (Wargelin & Drake 2002) provide only the upper limit of Ṁw .The observed mass-loss rate roughly exhibits a power-law dependence on the X-ray flux (measured in the ROSAT band, Wood et al. 2021;Vidotto 2021).However, the mass-loss rate of the solar wind remains nearly constant over the activity cycle, despite a factor of 20 variation in the X-ray flux (Johnstone & Güdel 2015).This discrepancy between solar and stellar observations of Ṁw presents a major concern when employing the observational scaling laws.
The purpose of this work is to establish a theoretical scaling law for Ṁw (Schröder & Cuntz 2005;Cranmer & Saar 2011;Suzuki et al. 2013;Suzuki 2018;Chebly et al. 2023) that aligns with both solar and stellar observations.Solar wind models generally rely on either Alfvén-wave acceleration (Ofman & Davila 1995, 1998;Matthaeus et al. 1999;Dmitruk et al. 2002;Ofman 2004;Suzuki & Inutsuka 2005;Verdini & Velli 2007;Verdini et al. 2010;Perez & Chandran 2013;van der Holst et al. 2014;van Ballegooijen & Asgari-Targhi 2016;Usmanov et al. 2018;Shoda et al. 2019;Réville et al. 2020;Matsumoto 2021;Magyar & Nakariakov 2021;Schleich et al. 2023) or loop-opening acceleration (Fisk et al. 1999;Fisk 2003;Rappazzo et al. 2012;Lionello et al. 2016;Bale et al. 2023;Drake et al. 2023;Kumar et al. 2023;Raouafi et al. 2023;Chitta et al. 2023a,b).Recent observations of magnetic switchbacks (Bale et al. 2019;Kasper et al. 2019;Bale et al. 2021;Telloni et al. 2022) appear to support the loop-opening scenario (Fisk & Kasper 2020;Tenerani et al. 2020;Zank et al. 2020;Drake et al. 2021;Magyar et al. 2021), though alternative or hybrid models for switchback generation have been put forward (Squire et al. 2020;Schwadron & McComas 2021;Shoda et al. 2021;Wyper et al. 2022).It has been proposed that a hybrid model (incorporating both Alfvén wave and loopopening effects) is required to accurately model the solar wind (Wang 2020;Iijima et al. 2023).With this in mind, we aim to derive the scaling law of Ṁw based on the hybrid model of solar (and stellar) wind and assess its compatibility with solar and stellar observations.The rest of this paper is organized as follows.Section 2 provides a comprehensive description of the proposed hybrid model, detailing the key aspects and components that contribute to its development.In Section 3, we present the simulation results, illustrating the performance of the hybrid model and its agreement with solar and stellar observations.Finally, Section 4 offers a summary of the paper, along with a discussion of potential applications and limitations of the study, highlighting areas for future research and improvement.

Model overview and assumption
In this study, we focus on solar-type stars, namely stars possessing the same luminosity, mass, radius, and metallicity as the Sun.The only stellar parameter that we consider to exhibit variability is the rotation period.Altering stellar rotation has two effects.First, wind characteristics change due to centrifugal force when equatorial rotation exceeds a small fraction of the breakup speed (Matt et al. 2012).Second, the rotation period directly impacts the stellar dynamo, affecting the properties of the stellar magnetic field (Saar 1996(Saar , 2001;;Reiners et al. 2009;See et al. 2019;Reiners et al. 2022).In this study, we solely consider changes in the stellar magnetic field for simplicity.Our assumption is validated for rotation rates less than ten times solar, or equivalently, in the unsaturated regime in terms of X-ray emission (Pizzolato et al. 2003;Wright & Drake 2016).
We adopt the assumption that non-magnetic photospheric parameters -temperature, density, scale height, and correlation length -align with those of the Sun, as detailed in Table 1.The rationale is that non-magnetic processes (nuclear reaction, radiative transfer, and convection) presumably shape stellar internal structures (Paxton et al. 2011).Nonetheless, the potential influence of magnetic fields on the photosphere deserves further investigation in the future.
In this work, B r, * symbolizes the field intensity of localized magnetic patches located within intergranular lanes.Following Cranmer & Saar (2011), we assume that these fields exhibit a magnetic pressure equivalent to the surrounding gas pressure, referred to as the equipartition field.Specifically, an assertion consistent with solar observations (Tsuneta et al. 2008).Generally, the equipartition field substantially surpasses the surface-averaged magnetic field in magnitude, and the proportion between these two field intensities is often interpreted as the filling factor of the magnetic field on the surface (Saar 1996(Saar , 2001;;Cranmer 2017).Although the assumption of equipartition is commonly used in the literature, it is important to note that the magnetic field could be larger than the equipartition on the photosphere (Okamoto & Sakurai 2018;Kochukhov et al. 2020).We also note that the employed value of B r, * is slightly larger than the typical field strength of emerging flux (Centeno et al. 2007).

Basic equations
In this study, we solve the MHD equations in a onedimensional system along a super-radially expanding flux tube.For simplicity, we assume the flux tube to be vertically aligned.Let r be the coordinate along the axis of the flux tube and x and y be the transverse components.The basic equations are then written as follows (Shoda & Takasao 2021). where and the total pressure p T and total energy density e are expressed as where e int is the internal energy density per unit volume.We employ the equation of state for partially ionized hydrogen gas to determine e int as a function of p and ρ (Vögler et al. 2005;Shoda & Takasao 2021).f open denotes the filling factor of the open magnetic flux.The terms D v x,y and D b x,y represent the phenomenological terms of turbulent dissipation.The heating rates by thermal conduction and radiation are denoted by Q cnd and Q rad , respectively.The heating rate of interchange reconnection driven by flux emergence is represented by We assert that although turbulent dissipation components do not feature in the energy equation, our model accounts for turbulent heating.This is clarified by elucidating the conservation law of mechanical energy (e mech = e − e int ), as shown below: where represents the turbulent dissipation rate of Alfvén wave.We note that the divergence of a vector field X is given by where X r is the r component of X.
The internal energy equation is obtained by subtracting Eq.( 7) from the equation for e.The result is This equation includes the turbulent dissipation rate as a heating term.Thus, by solving Eq.s (2)-( 4), the heating effect of turbulent dissipation is fully incorporated.The same arguments holds for numerical dissipation and heating.

Flux-tube parameters
f open regulates the shape of the open field lines and is in reality determined by the MHD force balance.However, since our model does not account for force balance in the cross-field direction, we need to prescribe f open as a function of r.In this study, we formulate f open (r) using its surface value f open * and the flux-tube expansion parameter (η exp ) as follows where The functions f 1 (r) and f 2 (r) represent the flux-tube expansion in the chromosphere (Ishikawa et al. 2021) and corona (Kopp & Holzer 1976), respectively.The parameter η exp represents the flux-tube expansion factor in the chromosphere.The expansion length scale is set to where Φ surf and Φ open represent the unsigned surface and open magnetic fluxes, respectively.We infer η exp , presuming equivalence between the unsigned magnetic flux at the photospheric and coronal base levels.This implies the absence of magnetic loops closing within the photosphere or chromosphere-an assumption that may overstate reality given that a subset of these loops peak at heights less than 1000 km (Wiegelmann & Solanki 2004;Cranmer & van Ballegooijen 2010).The product With the formulation of f open (r), we can deduce the radial component of the magnetic field B r from the magnetic-flux conservation of the flux tube:

Alfvén-wave turbulence
Alfvén-wave turbulence is likely responsible for converting wave energy into thermal energy.Since Alfvénwave turbulence is fundamentally a multi-dimensional process, with energy cascading primarily in the crossfield direction, the effect of turbulent dissipation must be phenomenologically incorporated in one-dimensional models.In this work, we employ the conventional onepoint closure model for the phenomenology of Alfvénwave turbulence (Hossain et al. 1995;Matthaeus et al. 1999;Dmitruk et al. 2002).It is important to note, however, that there are improved models available, such as the nearly-incompressible MHD turbulence model and the MHD shell model.
Following Shoda et al. (2018), we express D v x,y and D b x,y in terms of Elsässer variables (Elsasser 1950) as follows: where In this study, instead of solving the radial transport of the correlation length (Matthaeus et al. 1994;Breech et al. 2008;Cranmer & van Ballegooijen 2012;Zank et al. 2018), we prescribe λ ⊥ as a function of r.It is often assumed that λ ⊥ scales with the radius of the flux tube (Hollweg 1986;Verdini & Velli 2007;van der Holst et al. 2014), which is likely satisfied in the magneticallydominated region.However, in the lower atmosphere, the spatial scale of fluctuation is possibly regulated by granulation rather than the magnetic field expansion.Indeed, the granular pattern of velocity appears to persist up to the low chromosphere (Salucci et al. 1994;Espagnet et al. 1995;Kostik et al. 2009).Based on these arguments, we assume that λ ⊥ undergoes radial expansion until the flux tube completes its expansion in the chromosphere, and thereafter follows field-aligned expansion.This can be expressed as: where the correlation length in the photosphere is set as λ ⊥, * = 2000 km (Table 1).

Thermal conduction and radiation
The heat generated due to thermal conduction can be expressed in terms of the heat flux, denoted by F cnd , as follows: For the calculation of F cnd , we utilize the Spitzer-Härm type flux (Spitzer & Härm 1953) incorporating quenching in low-density regions (Shoda et al. 2020): Introducing quenching enables us to mitigate the stringent time-step constraint arising from thermal conduction.However, it is crucial to ensure that conduction is not quenched within the middle corona, where it significantly influences the energy balance.In this study, we use ρ cnd = 10 −20 g cm −3 , leading to conduction quenching above r/R ⊙ ≈ 10 − 20, where the wind is supersonic.This results in a considerably reduced conductive flux compared to observations at 1 au (Salem et al. 2003;Bale et al. 2013); nevertheless, it is not expected to substantially impact the mass-loss rate, as quenching starts to take effect in the supersonic region.
The radiative heating term, Q rad , consists of two components: the optically thick (Q thck ) and optically thin (Q thin ) contributions.These are combined as follows (Shoda & Takasao 2021): where the parameter p rad /p * is set to 0.1.Q thck is approximated by an exponential cooling process (Gudiksen & Nordlund 2005), which drives the internal energy toward a reference value on a specified timescale: where the timescale τ thck is defined as: and e int ref denotes the reference internal energy, which corresponds to a temperature equivalent to the surface temperature.Consequently, the optically thick cooling mechanism serves to bring the atmospheric temperature closer to the surface value.
The optically-thin cooling rate, denoted as Q thin , is expressed in terms of the radiative loss function (Λ) as follows: where n H and n e represent the number densities of hydrogen atoms and electrons, respectively.The function Λ(T ) is computed using the CHIANTI atomic database version 7 (Dere et al. 1997;Landi et al. 2012) with the photospheric abundance and extended to accommodate lower-temperature regimes (T < 10 4 K) as provided by Goodman & Judge (2012), following the methodology outlined in the literature (Iijima 2016;Shoda & Takasao 2021).

Phenomenology of flux emergence and interchange reconnection
In Eq. ( 4), the term Q FE denotes the plasma heating rate attributed to interchange reconnection.In this study, we consider interchange reconnection as a source of additional heating and disregard changes in flux tube shape and coronal-base parameters due to altered field connectivity.The energy released by this process is presumed to be conveyed through flux emergence, thereby linking Q FE with the energy flux of flux emergence, represented as F FE .The relationship between these variables can be expressed as follows: Considering that flux emergence and interchange reconnection are inherently multi-dimensional processes, we offer a phenomenological depiction of F FE , subject to the ensuing assumptions.
Primarily, we hypothesize that F FE remains temporally invariant.Although the time-dependent nature of flux emergence frequently engenders various transient phenomena, our primary interest resides in the timeaveraged properties of the stellar wind.As a result, we exclusively investigate the time-averaged energy transport mediated by flux emergence.
Secondly, we contend that the energy flux associated with flux emergence bears a direct proportionality to the coronal background magnetic field, as substantiated by the analysis of the coronal magnetic field's recycling timescale (Wang 2020(Wang , 2022)).Specifically, we propose that the energy flux discerned at the coronal base is represented by where X cb represents the value of X measured at the base of the corona, specifically defined as the location where the time-averaged temperature reaches 4.0×10 5 K (see also Section 3.1).We remark that the F FE cb value for B r,cb = 10 G is one-third of the initial value presented by Wang (2020) (2.9 × 10 5 erg cm −2 s −1 ).This deviation arises due to different estimates of the loop height h ER , typically associated with the total magnetic flux Φ ER as follows (Cranmer & van Ballegooijen 2010): In Wang (2020), the mean loop height ⟨h ER ⟩ = 6.4 × 10 3 km corresponds to the average magnetic flux ⟨Φ ER ⟩ = 1.0 × 10 19 Mx.Given the higher prevalence of flux emergence with smaller magnetic flux (Thornton & Parnell 2011), the mean magnetic flux may decrease, implying an average loop height less than 6.4 × 10 3 km (Wiegelmann & Solanki 2004;Cranmer & van Ballegooijen 2010).In addition, considering magnetic loops lower than the transition region do not directly contribute to coronal heating, emerging flux with total magnetic flux less than 1.0×10 18 Mx should be disregarded.Based on these premises, we estimate ⟨h ER ⟩ = 2.0 × 10 3 km.As the energy flux of interchange reconnection is assumed to vary with ⟨h ER ⟩ (Wang 2020), we lower the F FE cb value by a factor of 6.4/2.0 ≈ 3.
Considering the aforementioned premises, the functional form of F FE is delineated as follows: where we adopt the values r FE /r * = 1.2 and σ FE /r * = 0.02 to ensure that the heating due to interchange reconnection is confined to the low (subsonic) corona.We have investigated the dependence of Ṁw on r FE and σ FE , finding it largely unchanged as long as heating by flux emergence is confined to the subsonic region (see Appendix A for detail).
It is noteworthy that a previous study (Wang 1994) incorporated an energy flux similar to Eq. (31).A comparison between our flux-emergence energy flux (solid line, given by Eq.( 31)) and the basal-heating energy flux from Wang (1994) for both interplume (IPL, dashed line) and plume (PL1, dotted line) cases is shown in Figure 2 The energy flux in our model finds its position between the IPL and PL1 models, with respect to the basal value and the damping length.Given the successful modeling of plume and interplume in Wang (1994), it is plausible that our basal heating likely mirrors the actual one.

Boundary condition
The simulation domain encompasses the region from the stellar surface, constituting the inner boundary, to beyond the super-Alfvénic region, which serves as the F FE (Eq.( 25)) model (PL1) in Wang (1994) model (IPL) in Wang (1994) Figure 2. Comparison between F FE (solid line, as given in Eq. ( 31)) and the basal-heating energy fluxes utilized in Wang (1994).Dotted and dashed lines signify the plume (PL1) and interplume (IPL) models in Wang (1994), respectively.
outer boundary.We impose the free boundary condition at the outer boundary, as all MHD waves exhibit outward propagation.For the inner boundary, we utilize the fixed boundary condition for ρ, p, and B r , as follows: where X in represents the value of X at the inner boundary.Additionally, we enforce a zero vertical velocity at the inner boundary: v r,in = 0.
To impose the horizontal velocity and magnetic field such that 1. the reflected (inward) Alfvén waves penetrate the lower boundary, and 2. a (nearly) constant net Poynting flux is injected from the surface, we apply a free boundary condition for the inward Elsässer variables: where Subsequently, we adjust the amplitude of z + ⊥,in to maintain a fixed time-averaged wave energy flux (Shoda et al. 2020), the value of which is taken from the literature (Cranmer & van Ballegooijen 2005): where and X represents the time-average of X.
The spectra for z + ⊥,in are provided as follows: where ϕ i,x,y is a random phase function, and N + 1 designates the total mode count.Furthermore, define the least and greatest wave frequencies.In our study, N is set to 20.This implementation results in the outward Elsässer energy exhibiting a Kolmogorov-type spectrum (z + ⊥,in 2 ∝ f −5/3 ), as observed in the spectrum of the coronal transverse waves (Morton et al. 2019).
The outer boundary of the simulation domain, denoted r out , is adjusted to situate in the super Alfvénic region (in a quasi-steady state).Specifically, r out /r * = 44.9 for Φ surf /Φ surf The iterative computation of grid spacing is carried out as delineated below: where r i and ∆r i represent the radial distance of the i-th grid center and the i-th grid size, respectively.
The parameter r exp signifies the radial distance beyond which grid expansion occurs, and ε denotes the ratio of the neighboring grid size in the grid-expansion region.∆r min is held constant at 20 km, whereas ∆r max varies within the range of 2000 − 6000 km contingent upon the value of Φ surf .In this work, we establish ∆r 0 = ∆r min , r 0 = r * , r exp /r * = 1.0144, and ε = 0.0158.In numerically solving Eqs. ( 2)-( 4), we employ the HLLD solver (Miyoshi & Kusano 2005) in conjunction with the third-order SSP Runge-Kutta integration (Shu & Osher 1988;Gottlieb et al. 2001).We adopt the crosssection-weighted variables (Shoda & Takasao 2021) to render the form of the fundamental equations analogous to the conventional one-dimensional MHD equation.To guarantee elevated spatial accuracy within r < r exp , the fifth-order MP5 reconstruction (Suresh & Huynh 1997) is implemented, whereas the second-order MUSCL reconstruction ( van Leer 1979) is utilized in the grid-expansion region.In the area where ∆r = ∆r max , the MP5 reconstruction is applied to transverse velocity and magnetic field, while the MUSCL reconstruction is employed for the remaining variables.Thermal conduction is solved separately using the super-timestepping method (Meyer et al. 2012(Meyer et al. , 2014)).The influence of flux emergence and interchange reconnection on temperature profiles is distinctly evident.For Φ surf /Φ surf ⊙ ≲ 4, the temperature attains its maximum value at r/r * ≈ 3.5, where the solar wind undergoes efficient acceleration.Conversely, when Φ surf /Φ surf ⊙ ≳ 4, the temperature reaches its peak at r/r * ≈ 1.2, the location at which the majority of the emerging-flux energy is transformed into heat.This increased heat deposition due to flux emergence could potentially result in a two-step acceleration process for the stellar wind, as observed in the case of Φ surf /Φ surf ⊙ = 21.72.This acceleration occurs first in close proximity to the Sun, at r/r * ≈ 1.1, and is driven by the gas-pressure force (Parker 1958;Kopp & Holzer 1976), followed by a second acceleration in the supersonic region, at r/r * ≈ 10, propelled by the wave-pressure force (Alazraki & Couturier 1971;Belcher 1971).It is im-portant to note that this initial phase of wind acceleration may manifest as the active-region outflow observed within the solar corona (Sakao et al. 2007;Harra et al. 2008;Brooks et al. 2015).We find a comparable nonmonotonic trend in v r as a function of r in the previous models (Cranmer et al. 2007;Cranmer 2010;Panasenco et al. 2019), as identified in the case where the flux-tube expansion varied non-monotonically with r.
In the profiles of v r , it is evident that the stellar-wind velocity at the outer boundary remains relatively constant, despite the considerable diversity in open magnetic flux.Solar wind observations suggest that the governing factor for wind velocity is the flux-tube expansion factor (Wang & Sheeley 1990;Arge & Pizzo 2000;Riley et al. 2015), rather than the open magnetic flux (Wang 2020).In the proposed model, this expansion factor is equivalent to 1/ξ open and remains fixed, as illustrated in Figure 3. Consequently, our model offers a coherent explanation for the observed behavior of the solar wind.
The mass density of the solar wind exhibits a monotonic increase as a function of Φ surf .This can be attributed to the fact that the majority of emerging-flux heating transpires below the subsonic point, thereby leading to an enhancement in mass density at the critical point and further along (Leer & Holzer 1980;Pneuman 1980;Hansteen & Velli 2012).Consequently, the massloss rate, denoted by Ṁw , displays a monotonic increase with Φ surf , provided that ξ open ≡ Φ open /Φ surf remains constant.A comprehensive discussion on the behavior of Ṁw in relation to magnetic parameters can be found in Section 3.2.
The substantial diversity seen in Figure 3 indicates a primary reliance of the wind-driving mechanisms on Φ surf .We verify this by computing the energy fluxes due to flux emergence (F FE , see Eq. ( 31)) and Alfvén waves (F AW , see Eq. ( 34)) at the coronal base, r = r cb , where the average temperature is 0.4 million Kelvin (T = 4.0 × 10 5 K).The top panel in Figure 4 plots F FE cb (blue diamonds) and F AW cb (red circles) against Φ surf /Φ surf ⊙ .Flux-emergence energy flux proportionally rises with Φ surf , as shown by Eq. ( 29), while Alfvén-wave energy flux gently ascends and plateaus at high Φ surf , consistent with Shoda et al. (2020).Consequently, ≈ 20.This does not mean, however, that the influence of flux emergence is negligible in Φ surf /Φ surf ⊙ < 20.Given that mass loss is affected by the heat deposited beneath the sonic critical point (Leer & Holzer 1980;Pneuman 1980;Hansteen & Velli 2012), the heating rate in the subsonic corona is crucial.We thus calculate the averaged heating rates from Alfvén waves (H AW ) and flux emergence (H FE ) in the subsonic corona as follows: where r sc is the sonic critical point's radial distance.
In the above calculation, we assume that the loss in the energy flux is attributed to heating.The detailed radial profiles of heating and cooling rates are shown in Appendix B.
The lower panel of Figure 4 plots H FE (blue diamonds) and H AW (red circles) against Φ surf /Φ surf ⊙ .Unlike energy flux, the flux-emergence effect overtakes the Alfvén-wave effect at Φ surf /Φ surf ⊙ ≈ 3.This finding aligns with the observed regime transition in Figure 3. Thus, flux emergence could significantly impact the mass-loss rate, even when it plays a minor role in overall stellar-wind energetics.

Theoretical scaling law of Ṁw
The bottom panel of Figure 5 illustrates the comparison between the scaling law in Eq. ( 41) and the simulation data.The results demonstrate the efficacy of Eq.( 41) in accurately representing the simulation data.
The correlation between stellar surface magnetic flux, denoted as Φ surf , and stellar X-ray luminosity, represented by L X , is typified by power-law relations (Pevtsov et al. 2003;Kochukhov et al. 2020;Zhuleku et al. 2020).As a result, the scaling law expressed in Eq. ( 41) can be recast in terms of L X .For example, Toriumi & Airapetian (2022) and Toriumi et al. (2022) illustrate that By integrating Eq.s ( 41) and ( 42), the resulting equation is which displays congruity with the observation-based scaling law given by (Vidotto 2021): provided that ξ open does not diminish substantially with L X .In order to elucidate the impact of flux emergence on the scaling law of the mass loss rate ( Ṁw ), we present in Figure 6 a comparison of the dependence of Ṁw on the surface magnetic flux (Φ surf ) for scenarios incorporating flux emergence and those excluding it (F FE = 0).For the purpose of this analysis, we maintain a constant open-to-surface flux ratio of ξ open = 0.1.In the absence of flux emergence, it is observed that Ṁw exhibits a nearly constant relationship with Φ surf , in stark contrast to the monotonically increasing trend observed when flux emergence is present.It is noteworthy that the near-constant behavior of Ṁw aligns with the results of previous stellar-wind simulations (Shoda et al. 2020).Thus, the observed enhancement in Ṁw corresponding to an increase in activity level (Wood et al. 2021; Vidotto 2021) be ascribed to the basal heating engendered by flux emergence.

Comparison with the solar wind observations
To substantiate the theoretical scaling law presented in Eq. ( 41), we examine its compatibility with solar observations.Our objective is to juxtapose the long-term trends of mass-loss rates ascertained through Eq. ( 41) and those directly measured in the solar wind.
Solar-disk observations spanning several decades provide data on the total unsigned magnetic flux of the Sun.While the open magnetic flux is not directly observable, it can be inferred either through magneticfield extrapolation techniques (Wang & Sheeley 2002) or via in-situ measurements in the solar wind, presuming that the r 2 B r is nearly uniform in latitude (Smith & Balogh 1995).The discrepancy between these two estimation methods is currently acknowledged as the openflux problem (Linker et al. 2017(Linker et al. , 2021;;Arge et al. 2023).In this study, we procure the observed data for Φ surf and Φ open from the literature (Cranmer 2017), which have been averaged over 0.25-year intervals.
In addition to the theoretical scaling law in Eq. ( 41), we also evaluate the compatibility of the stellar-based empirical scaling law in Eq. ( 44) with solar observations.Since L X in Eq. ( 44) represents the X-ray luminosity measured in the ROSAT band (5.2 − 124 Å, Truemper 1982;Güdel 2007), we employ the daily spectral irradiance data from the XUV Photometer System (XPS, Woods et al. 2005;Woods & Rottman 2005) onboard the Solar Radiation and Climate Experiment (SORCE) satellite to compute the solar L X .In addition, we also estimate the solar L X using an empirical conversion from the GOES band (1−8 Å) the ROSAT band (Judge et al. 2003;Cranmer et al. 2021).The solar X-ray luminosity is set to L X,⊙ = 1.0 × 10 27 erg s −1 (Peres et al. 2000;Johnstone & Güdel 2015).
Figure 7 illustrates the time evolution of three-monthaveraged parameters pertinent to the solar magnetic field and solar wind.The upper panel displays the filling factors of both the unsigned surface magnetic flux and open magnetic flux, which is congruent with Figure 2a of Cranmer (2017).In the central panel, observed massloss rates are represented by lines, while estimated rates are denoted by symbols.The solid and dashed lines represent the observed mass-loss rate as determined from 1-au (OMNI) and Ulysses data, respectively, exhibiting analogous behavior.The red circles indicate theoretical rates derived from Eq. ( 41), and blue crosses represent the empirical mass-loss rate derived from data over five solar cycles (Finley et al. 2019): The lower panel is similar to the central panel, juxtaposing observations with the empirical X-ray-based estimations (indicated by pluses) from Eq. ( 44).Here, red and blue pluses correspond to L X data ascertained from GOES and SORCE/XPS, respectively.
Our theoretical prediction better represents the longterm trend of Ṁw than the X-ray empirical law.Furthermore, it closely matches the empirical equation presented in Eq. ( 45).Both our model and fitted relation demonstrate an enhanced agreement with observational data in the later phase (year ≥ 2000), which may be attributed to the improved precision of magneticfield measurements (refer to Cranmer (2017) for an indepth explanation).Notably, neither approach captures the observed Ṁw peak around 1991-1992, suggesting an anomalous phase of the solar wind during this period.
Based on this comparative analysis, we deduce that the hybrid model utilized in this study offers a suitable representation of solar (and potentially stellar) wind properties.Furthermore, we assert that Eq. ( 41) provides a valuable tool for estimating the mass-loss rate of solar-type stars.Investigating the applicability of our scaling law in the stellar context is beyond the scope of this study.However, it is pertinent to assess if our law offers any advantages over the X-ray-based scaling law or if both are analogous.

DISCUSSION
In the present investigation, we introduce an innovative stellar-wind model founded upon the hybrid scenario, which incorporates both Alfvén waves and interchange reconnection processes.A comprehensive parameter survey leads to the derivation of a scaling relationship for Ṁw , as presented in Eq. ( 41).This relationship aligns with the observed increase in Ṁw as a function of X-ray flux in low-mass stars and the nearconstant nature of solar-wind mass flux, thereby offering potential utility for estimating mass-loss rates in solartype stars.
The application of Eq. ( 41) is complicated by the necessity of calculating the open magnetic flux, a task currently unattainable even for our Sun (Linker et al. 2017).However, the latest improvements in modeling stellar large-scale magnetic field structures, via global MHD simulations and Zeeman-Doppler Imaging (Folsom et al. 2016(Folsom et al. , 2018;;See et al. 2019), indicate potential feasibility for deducing open magnetic flux ( Ó Fionnagáin et al. 2019;Airapetian et al. 2021;Evensberget et al. 2021Evensberget et al. , 2022Evensberget et al. , 2023)).A comprehensive study of the relation between the Sun's surface and open magnetic flux (Yoshida et al. 2023) could also be crucial in eliminating Φ open from Eq. ( 41).
While our model demonstrates efficacy, it is essential to acknowledge the presence of numerous free parameters that could potentially influence our conclusions.A subset of these parameters is linked to the formulation of F FE cb , which represents the energy flux discharged by interchange reconnection in the corona, as delineated in Section 2.6.Additionally, there exists an uncertainty in the net energy flux of Alfvén waves, as denoted by Eq. ( 33).Given the centrality of these parameters to our model, it is imperative that their values are validated or adjusted through three-dimensional simulations in future research.
Despite the aforementioned inherent limitations, it can be reasonably postulated that the influence of flux emergence may play a significant role in modulating mass-loss rates.Consequently, it is meritorious to reexamine extant models concerning wind mass-loss rates predicated on the Alfvén-wave framework (Cranmer & Saar 2011;Shoda et al. 2020).Future stellar wind models should concurrently integrate effects of Alfvén wave and flux emergence, alongside stellar rotation (Weber & Davis 1967;Belcher & MacGregor 1976;Holzwarth & Jardine 2007;Shoda et al. 2020).
Finally, we address the implications of our work for the magnetic braking problem in low-mass stars.Low-mass stars with Rossby numbers near the solar value are found to halt their spin-down (van Saders et al. 2016(van Saders et al. , 2019;;Metcalfe et al. 2023), a phenomenon under active investigation in stellar physics.A plausible explanation is a shift in magnetic morphology, where higher-order multipolar components become dominant (Garraffo et al. 2016(Garraffo et al. , 2018)).This leads to reduced open magnetic flux (Yoshida et al. 2023) and mass-loss rate as per Eq. ( 41), subsequently lowering magnetic braking efficiency (Matt et al. 2012;Réville et al. 2015;Finley & Matt 2018).Our findings therefore support this hypothesis, although we need further investigations such as deriving a formula for open magnetic flux as a function of surface magnetic field and mass-loss rate.
Numerical computations were carried out on the Cray XC50 at the Center for Computational Astrophysics (CfCA), National Astronomical Observatory of Japan.This work made use of matplotlib, a Python library for publication quality graphics (Hunter 2007), andNumPy (van der Walt et al. 2011).This work was supported by JSPS KAKENHI Grant Nos.JP20KK0072 (PI: S. Toriumi), JP21H01124 (PI: T. Yokoyama), JP21H04492 (PI: K. Kusano), and JP22K14077 (PI: M. Shoda).findings.This appendix explores the impact of these parameters on Ṁw .
In the original F FE formulation given by Eq. ( 31), we assume exp[(r FE − r * )/σ FE ] ≫ 1.This holds in our reference case where r FE /r * = 1.2 and σ FE /r * = 0.02.However, for smaller r FE or larger σ FE values, this may not hold, leading to a reduced F FE at the photosphere.To address this, during a parameter survey over r FE and σ FE , we modify F FE as follows.The top panel of Figure 8 illustrates how the massloss rate varies with r FE .The value of σ FE is fixed to the reference value (σ FE /r * = 0.02).The horizontal dashed line marks the mass-loss rate without flux emergence.The mass-loss rate peaks at r FE = 0.1, and the reason for this trend is explained as follows.As r FE increases, more heat is deposited into the supersonic region, which boosts wind velocity rather than mass flux (Leer & Holzer 1980;Pneuman 1980).This leads to a decreasing mass-loss rate with larger r FE , approaching the no-flux-emergence value.The small drop at r FE < 0.1 likely stems from heating in the chromosphere, which reduces the energy flux transmitted into the corona.Since the mass-loss rate is nearly proportional to the coronal energy flux (Cranmer & Saar 2011;Shoda et al. 2020), heat deposition in the chromosphere results in a smaller mass-loss rate.
The bottom panel of Figure 8 shows the effect of varying σ FE , keeping r FE /r * = 1.2 constant.The mass-loss rate remains constant for σ FE /r * ≤ 0.1 and declines for σ FE /r * ≥ 0.1.As discussed above, this decrease occurs because heat deposition in the supersonic region becomes larger with increasing σ FE .These results confirm the robustness of our results as long as flux emergence heating is limited to the subsonic region.

B. RADIAL PROFILES OF HEATING RATES
In this appendix, we detail the radial profiles of heating and cooling rates in stellar winds.Figure 9 displays the radial heating and cooling rates in stellar winds.Panels differ by surface magnetic flux: top (Φ surf /Φ surf ⊙ = 1.36), middle (Φ surf /Φ surf ⊙ = 5.43), and bottom (Φ surf /Φ surf ⊙ = 21.72).Thick and thin lines represent positive (heating) and negative (cooling) rates, respectively.Colors signify mechanisms: red for interchange reconnection (Q FE /ρ), orange for Alfvénwave turbulence (Q turb /ρ), blue for thermal conduction (Q cnd /ρ), and green for radiation (Q rad /ρ).We note that Alfvén-wave turbulence is the likely dominant dissipation mechanism for Alfvén waves in solar wind (Shoda et al. 2018).A grey vertical line in each panel marks the sonic point.
Regardless of the magnitude of surface magnetic flux, heating by interchange reconnection is localized to the subsonic region, as per Eq.(31).In contrast, heating by Alfvén-wave turbulence is widely distributed and peaks in the supersonic region.It is important to clarify that the heating discussed here is in terms of rate per unit mass.For larger surface magnetic flux, heating by Alfvén-wave turbulence diminishes, whereas heating by interchange reconnection intensifies.These trends align with the results in the bottom panel of Figure 4.It is also noteworthy that heating by interchange reconnection is balanced by cooling via thermal conduction.Consequently, the lower atmosphere (comprising the coronal base, transition region, and chromosphere) experiences increased conductive heating, particularly when the surface magnetic flux is large, as illustrated in Figure 9.This conductive heating in the chromosphere contributes to coronal mass enhancement through chromospheric evaporation (Yokoyama & Shibata 2001), making flux emergence and interchange reconnection key factors in elevating solar wind density and mass flux.
where H * is the pressure scale height on the photosphere.f open * and η exp are connected to the unsigned surface magnetic flux (Φ surf ) and open magnetic flux (Φ open ) through magnetic-flux conservation.
signifies the ratio of open to surface magnetic flux and is denoted as ξ open : and λ ⊥ is the perpendicular correlation length of turbulence.Previous threedimensional reduced MHD simulations(Perez & Chandran 2013;van Ballegooijen & Asgari-Targhi 2016, 2017;Chandran & Perez 2019) and shell-model calculations(Verdini et al. 2019) have suggested that c d is significantly smaller than unity.The inferred best-fit parameter is c d = 0.1 (van Ballegooijen & Asgari-Targhi 2016), which is used in this work.
Figure 3 presents a comparison of simulation outcomes for varying surface magnetic flux magnitudes.A constant open-to-surface magnetic flux ratio is maintained (ξ open = 0.2).Each subplot illustrates the timeaveraged radial profiles for temperature (upper), velocity (central), and density (lower).Lines exhibiting a more intense red hue signify a greater surface magnetic flux.Circles situated in the upper subplot represent the point of maximum temperature.The position of the sonic point (defined as v r = γp/ρ) is indicated by diamonds in the central subplot.The influence of flux emergence and interchange reconnection on temperature profiles is distinctly evident.For Φ surf /Φ surf

ΦΦΦ
Figure 3. Time-averaged radial profiles of stellar-wind parameters: mass density (top panel), temperature (middle panel), and radial velocity (bottom panel).Line color signifies the corresponding surface magnetic flux.Temperaturemaximum points are denoted by circles in the top panel, while sonic points are represented by diamonds in the middle panel.The open-to-surface magnetic flux ratio is maintained at a constant value of ξ open = 0.2.

Figure 5
Figure 5 displays the simulated data for Ṁw across a range of Φ surf and ξ open .The three parameters are normalized using solar reference values: Ṁw,⊙ = 2.0 × 10 −14 M ⊙ yr −1 , Φ surf ⊙ = 3.0 × 10 23 Mx, and ξ open ⊙ = 0.2.The top and middle panels indicate that Ṁw increases with Φ surf and ξ open in a manner approximating a power-law relationship.Inspired by these trends, we fit the simulation data employing double power-law relations with respect to Φ surf and ξ open .Specifically, the fitting results yield the following scaling law: Ṁw Ṁw,⊙ = Φ surf Φ surf ⊙

Figure 4 .
Figure 4. Upper panel: Coronal-base energy flux plotted against Φ surf /Φ surf ⊙ .The blue diamonds denote the energy flux of the emerging flux (F FE cb ), while the red circles indicate the Alfvén wave energy flux (F AW cb ).Lower panel: heating rate in the subsonic corona plotted against Φ surf /Φ surf ⊙ .Blue diamonds and red circles represent heating rates resulting from the emerging flux (H FE ) and Alfvén wave (H AW ), respectively.

Figure 6 .
Figure 6.Comparative analysis of the Φ surf -Ṁw relations incorporating and excluding the influence of flux emergence, maintaining a constant open-to-surface flux ratio at ξ open = 0.1.

Figure 7 .
Figure 7.Long-term trends (quarterly averages) of solar magnetic field and solar wind parameters.Upper panel: Surface (Φ surf , displayed in black) and open (Φ open , in red) solar magnetic flux.Central panel: The solar wind mass-loss rate inferred from OMNI (solid line) and Ulysses (dotted line) datasets, juxtaposed with the open-flux-based fitting (blue crosses, refer to the text for detail) and theoretical massloss rate (red circles).Lower panel: Similar to the central panel, however, observations are counterpointed with an Xray-based empirical scaling law (denoted by pluses), where different colors specify different X-ray data sources (refer to the text for detail).

Figure 8 .
Figure 8. Dependence of mass-loss rate on flux-emergence heating parameters.Top: Variation with radial distance of heat deposition (r FE ).Bottom: Variation with spatial extent of heat deposition (σ FE ).Reference cases are indicated by red symbols.The surface magnetic flux and open-to-surface flux ratio are set at Φ surf /Φ surf ⊙ = 10.86 and ξ open = 0.1, respectively.
the survey, we fix the surface magnetic flux to Φ surf /Φ surf ⊙ = 10.86 and the open-to-surface magnetic flux ratio to ξ open = 0.1.

Table 1 .
Surface parameters employed in our model.