What is the gas temerature expected within a spot if only a reduction in gas temperature is accounted for, owing to the presence of a magnetic pressure term?
We'll begin under the assumption that the gas is subject to a polytropic equation of state such that
\begin{equation} P_{\rm gas} = K \rho^{\gamma}, \end{equation}where $P_{\rm gas}$ is the thermal gas pressure, $K$ is a constant of proportionality, $\rho$ is the gas density, and $\gamma$ is the ratio of specific heats. We will also assume, for simplicity, that the magnetic field is in equipartition with the gas, meaning that the internal energy of the gas, $U_{\rm gas}$, is equal to the energy of the magnetic field, $U_{\rm mag}$. Finally, if we assume the gas can be described as ideal, then we know that $\gamma = 5/3$, $P_{\rm gas} = nkT$, and $U_{\rm gas} = 3nkT/2 = 3P_{\rm gas}/2$.
It's possible to re-write the ideal gas equation of state by substituting the gas (mass) density for the gas number density using the mean molecular weight, $\mu$. We then have that $n = \rho/(\mu m_H)$, where $m_p$ is the mass of a proton. Thus,
\begin{equation} P_{\rm gas} = \frac{\rho}{\mu m_p}kT. \end{equation}Utilizing Equation (1), however, we find that
\begin{equation} P_{\rm gas}^{(\gamma - 1)/\gamma} = \frac{k}{K \mu m_p}T = \mathbb{K} T. \end{equation}where we have collected all constants to define a new constant $\mathbb{K}$. Note that the mean molecular weight is here taken to be constant with $\mu = 1$ (pure hydrogen gas), for simplicity.
Neglecting magnetic tension forces, the magnetic pressure is equal to the magnetic energy and is spatially isotropic, $P_{\rm mag} = U_{\rm mag} = B^2 / 8\pi$. Since the magnetic field is in energy equipartition with the gas, we have
\begin{equation} U_{\rm mag} = \frac{B^2}{8\pi} = \frac{3}{2}P_{\rm gas}, \end{equation}which provides a convenient estimate of the magnetic field strength at any given point in the gas. If we assume that the total pressure at a given point must be the same, regardless of the presence of a magnetic field, then we can write that
\begin{equation} P_{\rm tot} = P_{\rm gas} + \frac{B^2}{8\pi} = \frac{5}{2} P_{\rm gas}, \end{equation}which must be equivalent to the total pressure when no magnetic field is present (i.e., the gas pressure in the absence of a magnetic field),
\begin{equation} P_{\rm tot} \equiv P_{\rm gas,\, 0} = \frac{5}{2}P_{\rm gas} \end{equation}where $P_{\rm gas,\, 0}$ is the gas pressure in the absence of a magnetic field.
Under these approximations, we can estimate the change in gas temperature caused by the presence of a magnetic field, neglecting effects related to the transport of energy. Using Equation (3), we can write
\begin{equation} \frac{T_{\rm gas}}{T_{\rm gas,\, 0}} = \left(\frac{2}{5}\right)^{(\gamma - 1)/\gamma} = \left(\frac{2}{5}\right)^{2/5} \approx 0.693. \end{equation}In the case of the Sun, if we take the background photospheric temperature to be approximately 5779 K, then a rough estimate for the temperature within a sunspot would be 4000 K, or 1779 K cooler than the background photosphere. Quite surprisingly, this agrees with rough estimates for the temperatures within starspot umbra (Solanki 2003).
If, instead, the magnetic pressure were equal to the gas pressure within a spot, the umbral temperature ratio would be equal to
\begin{equation} \frac{T_{\rm gas}}{T_{\rm gas,\, 0}} = \left(\frac{1}{2}\right)^{(\gamma - 1)/\gamma} = \left(\frac{1}{2}\right)^{2/5} \approx 0.758, \end{equation}which implies spot temperatures of order 4300 K. Again, consistent with observations of umbral temperatures.
One can also estimate the change in density resulting from the decrease in gas pressure and temperature. From Equation (1),
\begin{equation} \frac{\rho_{\rm gas}}{\rho_{\rm gas\, 0}} = \left(\frac{P_{\rm gas}}{P_{\rm gas\, 0}}\right)^{1/\gamma}, \end{equation}which yields a density ratio of approximately 50% for an ideal gas with $\gamma = 5/3$. When typical values for the density in the solar photosphere ($\log\rho \sim 6.4$ near $\tau = 1$) are used, this leads to a density change of approximately $2\times10^{-7}$ g cm$^{-2}$. This is consistent with density changes observed in 3D radiation MHD simulations by Kitiashvili et al. (2010).
Given agreement between estimated this simple estimate of umbral temperatures and estimates from starspot observations, what further consequences might lead to the exclusion of this hypothesis for the cooler nature of sunspots?
Development from the presence of a strong concentration of magnetic flux to the appearance of sunspot pores, if goverend by the thermal evolution of the gas, must occur on a timescale related to the sound crossing time. It is only over this timescale that the gas can effectively communicate information about the presence of the magnetic field. Typical sunspot pores are on the order of 1 Mm in size (Bray & Loughhead 1964) and the adiabatic speed of sound at the solar surface is on the order of 8 km/s. The latter is not exact, but provides a starting point.
\begin{equation} \tau = \frac{R_{\rm pore}}{c_s} \sim \frac{10^3}{10^1} \textrm{ s} = 10^2 \textrm{ s}, \end{equation}or approximately 2 minutes. The sound crossing time is, therefore, quite small compared to the pore formation timescale, which is estimated to be between several hours and several days. Limitations imposed by the sound crossing time are therefore not significant.
It is possible, however, that the process of magnetic suppression of convection and the cooling of the photospheric layers are intertwined. Suppression of convection is the result of the interaction of convective flows with the magnetic field via Lorentz forces. These forces are the same forces that lead to the isotropic pressure that could potentially be responsible for the cooling of the layers by offsetting the gas pressure required to maintain hydrostatic equilibrium.
The presence of a magnetic field in a plasma would cause a general cooling of the surface layers, in the event that the system was left to equilibrate. Since the sound crossing time is on the order of 2 minutes, this is expected to occur quite rapidly. However, the surface is not static, as convective updrafts are constantly supplying fresh, warm material to the surface. Convection in the near surface layers has an overturn time of order several minutes (solar granulation timescale), of the same order as the time required for the material to arrive in equilibrium with the magnetic field. Therefore, it may be unlikely for the magnetic field to efficiently cool the gas until the system is unable to supply warm material from deeper layers. This is particularly important since, as a glob of plasma may cool, it will have a tendency to sink while convection occurs. Suppression of convection starves the upper layers of warm material and prevents cooler material from traveling inward, permitting the gas in the near surface layers to equilibrate to a cooler temperature, as dictated by the additional pressure contribution to the equation of state.
Needs validation.
Using the above theoretical development, we find the temperature contrast on other stars is equal to that for sunspots. However, measurements of starspot temperature contrasts find a correlation between starspot temperature contrast, and the effective temperature of the star (Berdyugina 2005).
In [1]:
import numpy as np
In [2]:
# confirm sunspot estimate above
Bz = np.sqrt(12.*np.pi*10**4.92)
print "Sunspot Bz (G): {:8.3e}".format(Bz)
print "Sunspot umbral temperature (K): {:6.1f}".format(5779.*0.4**0.4)
While the temperature ratio appears fixed in Equations (7) and (8), variation in the derived values can be estimated from the variation in $\gamma = c_p / c_v$. There will also be variation in the specific relation between the gas internal energy and the gas pressure, but this is more difficult to quantify from stellar evolution model output.
To address this problem, a small grid of models was run with masses in the range $0.1$ — $0.9 M_{\odot}$ with a mass resolution of $0.1 M_{\odot}$. Output from the code was modified to yield $\gamma$ directly as a function of radius, temperature, pressure, and density.
In [3]:
%matplotlib inline
import matplotlib.pyplot as plt
In [4]:
Teffs = np.arange(3000., 6100., 250.)
Tspot_fixed_gamma = Teffs*(0.4)**0.4
# gamma where tau ~ 1 (note, no 0.2 Msun point)
Gammas = np.array([1.22, 1.29, 1.30, 1.29, 1.36, 1.55, 1.63, 1.65])
ModelT = np.array([3.51, 3.57, 3.59, 3.61, 3.65, 3.71, 3.74, 3.77]) # tau = 1
ModelT = 10**ModelT
ModelTeff = 10**np.array([3.47, 3.53, 3.55, 3.57, 3.60, 3.65, 3.69, 3.73]) # T = Teff
Tratio = 0.4**((Gammas - 1.0)/Gammas)
Tspot_physical_gamma = np.array([ModelT[i]*0.4**((Gammas[i] - 1.0)/Gammas[i]) for i in range(len(ModelT))])
# smoothed curve
from scipy.interpolate import interp1d
icurve = interp1d(ModelT, Tspot_physical_gamma, kind='cubic')
Tphot_smoothed = np.arange(3240., 5880., 20.)
Tspot_smoothed = icurve(Tphot_smoothed)
# approximate Berdyugina data
DeltaT = np.array([ 350., 450., 700., 1000., 1300., 1650., 1850.])
BerdyT = np.array([3300., 3500., 4000., 4500., 5000., 5500., 5800.])
BSpotT = BerdyT - DeltaT
print Tratio
In [5]:
fig, ax = plt.subplots(1, 1, figsize=(8,4))
ax.set_xlabel('Photospheric Temperature (K)', fontsize=20.)
ax.set_ylabel('Spot Temperature (K)', fontsize=20.)
ax.tick_params(which='major', axis='both', length=10., labelsize=16.)
ax.grid(True)
ax.plot(Teffs, Tspot_fixed_gamma, '--', lw=2, dashes=(10., 10.), markersize=9.0, c='#1e90ff')
ax.plot(BerdyT, BSpotT, '-', lw=2 , dashes=(25., 15.), c='#000080')
ax.fill_between(BerdyT, BSpotT - 200., BSpotT + 200., facecolor='#000080', alpha=0.1, edgecolor='#eeeeee')
ax.plot(ModelT, Tspot_physical_gamma, 'o', lw=2, markersize=9.0, c='#800000')
ax.plot(Tphot_smoothed, Tspot_smoothed, '-', lw=3, c='#800000')
Out[5]:
In [6]:
fig, ax = plt.subplots(1, 1, figsize=(8,4))
ax.set_xlabel('Photospheric Temperature (K)', fontsize=20.)
ax.set_ylabel('T(phot) - T(spot) (K)', fontsize=20.)
ax.tick_params(which='major', axis='both', length=10., labelsize=16.)
ax.grid(True)
ax.plot(Teffs, Teffs - Tspot_fixed_gamma, '--', lw=2, dashes=(10., 10.), markersize=9.0, c='#1e90ff')
ax.plot(BerdyT, DeltaT, '-', lw=2 , dashes=(25., 15.), c='#000080')
ax.fill_between(BerdyT, DeltaT - 200., DeltaT + 200., facecolor='#000080', alpha=0.1, edgecolor='#eeeeee')
ax.plot(ModelT, ModelT - Tspot_physical_gamma, 'o', lw=2, markersize=9.0, c='#800000')
ax.plot(Tphot_smoothed, Tphot_smoothed - Tspot_smoothed, '-', lw=3, c='#800000')
Out[6]:
Results of this simple model indicate that equilibration of a polytropic gas within a magnetic structure located near the photosphere ($\tau_{\rm ross} = 1$) provides a reasonable approximation to observed spot temperatures from low-mass M dwarfs up to solar-type stars. Above 5400 K, the gas is sufficiently ideal that the model predicted relationship (red line) is asymptotic to the case of a purely ideal gas (small-dashed light-blue line). Below that temperature, the simple model traces the relationship provided by Berdyugina (2005). Difficulties below 4000 K may be the result of either model inaccuracies, either stemming from atmospheric structure or the simple approximation of energy equipartition, or observational complications that arise from measuring M dwarf photospheric temperatures.
We can also estimate umbral magnetic field strengths by extracting the model gas pressure at the same optical depth ($\tau_{\rm ross} = 1$),
In [7]:
# log(Pressure)
p_gas = np.array([6.37, 5.90, 5.80, 5.70, 5.60, 5.45, 5.30, 5.15])
B_field_Eeq = np.sqrt(12.*np.pi*0.4*10**p_gas)/1.0e3 # in kG
B_field_Peq = np.sqrt( 8.*np.pi*10**p_gas)/1.0e3
# smooth curves
icurve = interp1d(ModelT, B_field_Eeq, kind='cubic')
B_field_Eeq_smooth = icurve(np.arange(3240., 5880., 20.))
icurve = interp1d(ModelT, B_field_Peq, kind='cubic')
B_field_Peq_smooth = icurve(np.arange(3240., 5880., 20.))
In [8]:
fig, ax = plt.subplots(1, 1, figsize=(8,4))
ax.set_xlabel('Photospheric Temperature (K)', fontsize=20.)
ax.set_ylabel('Equipartition Magnetic Field (kG)', fontsize=20.)
ax.tick_params(which='major', axis='both', length=10., labelsize=16.)
ax.grid(True)
ax.plot(ModelT, B_field_Eeq, 'o', lw=2, markersize=9.0, c='#800000')
ax.plot(np.arange(3240., 5880., 20.), B_field_Eeq_smooth, '-', lw=2, c='#800000')
ax.plot(ModelT, B_field_Peq, 'o', lw=2, markersize=9.0, c='#1e90ff')
ax.plot(np.arange(3240., 5880., 20.), B_field_Peq_smooth, '-', lw=2, dashes=(20., 5.), c='#1e90ff')
Out[8]:
The two curves represent two different approximations, one is energy equipartition (red curve) and the other that the magnetic pressure is precisely equal to the gas pressure (blue curve). These values do not represent surface averaged magnetic field strengths, but the strengths of local concentrations of magnetic flux. Based on energy equipartition, we do not expect spot magnetic field strengths to be considerably larger than those estimated from the red curve.
Finally, we can estimate a "curve of cooling" relating the strength of a magnetic field to the temperature within the magnetic sturcture. Since stars the properties of the photospheric layers are so different for stars as a function of effeictve temperature, it'll be helpful to compute curves at several effective temperatures characteristic of low-mass M dwarf, a late K dwarf, and an early K/G dwarf.
In [14]:
B_field_strengths = np.arange(0.1, 8.1, 0.1)*1.0e3
log_P_phot = np.array([5.15, 5.60, 5.90])
Gamma_phot = np.array([1.65, 1.36, 1.29])
Exponents = (Gamma_phot - 1.0)/Gamma_phot
print Exponents
fig, ax = plt.subplots(3, 1, figsize=(8, 12), sharex=True)
i = 0
ax[2].set_xlabel('Spot Magnetic Field (kG)', fontsize=20.)
for axis in ax:
B_field_eq = np.sqrt(8.*np.pi*(0.4*10**log_P_phot[i]))
axis.grid(True)
axis.set_ylabel('$T_{\\rm spot} / T_{\\rm phot}$', fontsize=20.)
axis.tick_params(which='major', axis='both', length=10., labelsize=16.)
axis.plot(B_field_strengths/1.0e3, (1.0 - B_field_strengths**2/(8.*np.pi*10**log_P_phot[i]))**Exponents[i],
'-', lw=3, c='#800000')
axis.plot(B_field_eq/1.0e3, (1.0 - B_field_eq**2/(8.0*np.pi*10**log_P_phot[i]))**Exponents[i],
'o', markersize=12.0, c='#555555')
i += 1
Shulyak et al. (2010, 2014) have measured distribution of magnetic fields on the surfaces of M dwarfs by modeling FeH bands in high resolution, high S/N Stokes I sectra. They find M dwarf spectra are typically best fit by a uniform 1 kG magnetic field everywhere on the star, but with the addition of local concentrations across the surface. These local concentrations can reach upward of 7 – 8 kG. We now ask, do the stars that possess these field lie close to the locus defined by models?
In [10]:
Shulyak_max_B = np.array([[3400., 100., 6.5, 0.5],
[3400., 100., 6.0, 0.5],
[3300., 100., 7.5, 0.5],
[3100., 100., 6.5, 0.5],
[3100., 50., 1.0, 0.5]])
In [11]:
fig, ax = plt.subplots(1, 1, figsize=(8,4))
ax.set_xlabel('Photospheric Temperature (K)', fontsize=20.)
ax.set_ylabel('Equipartition Magnetic Field (kG)', fontsize=20.)
ax.tick_params(which='major', axis='both', length=10., labelsize=16.)
ax.grid(True)
ax.fill_between(np.arange(3240., 5880., 20.), B_field_Peq_smooth, B_field_Eeq_smooth, facecolor='#555555', alpha=0.3, edgecolor='#eeeeee')
ax.errorbar(Shulyak_max_B[:,0], Shulyak_max_B[:,2], xerr=Shulyak_max_B[:,1], yerr=Shulyak_max_B[:,3], lw=2, fmt='o', c='#555555')
Out[11]:
We find that the data from Shulyak et al. (2014) for the maximum mangetic field strengths (in local concentrations) of active stars are on the order of those expected from either energy of pressure equipartition. The fact that two have weaker magnetic fields only indicates that the there were no strong local concentrations, similar to the surface of the quiet Sun.
Note, however, that there is some uncertainty in this comparison. Shulyak et al. (2014) quote the effective temperature, which for M dwarfs is not necessarily equal to the photospheric temperature. Model atmospheres predict that the "effective temperature" for an M dwarf occurs in the optically thin layers above the opaque photosphere. Thus, it is more likely that the photospheric temperature for the Shulyak sample is greater than those quoted above.
If we instead convert the quoted effective temperatures to photospheric temperatures using stellar model atmospheres, we find
In [12]:
Shulyak_max_B = np.array([[3685., 100., 6.5, 0.5],
[3685., 100., 6.0, 0.5],
[3578., 100., 7.5, 0.5],
[3379., 100., 6.5, 0.5],
[3379., 50., 1.0, 0.5]])
fig, ax = plt.subplots(1, 1, figsize=(8,4))
ax.set_xlabel('Photospheric Temperature (K)', fontsize=20.)
ax.set_ylabel('Equipartition Magnetic Field (kG)', fontsize=20.)
ax.tick_params(which='major', axis='both', length=10., labelsize=16.)
ax.grid(True)
ax.fill_between(np.arange(3240., 5880., 20.), B_field_Peq_smooth, B_field_Eeq_smooth, facecolor='#555555', alpha=0.3, edgecolor='#eeeeee')
ax.errorbar(Shulyak_max_B[:,0] + 100., Shulyak_max_B[:,2], xerr=Shulyak_max_B[:,1], yerr=Shulyak_max_B[:,3], lw=2, fmt='o', c='#555555')
Out[12]:
which shows that localized concentrations of magnetic fields exceed those estimated by energy or pressure equipartition.
In general, these comparisons are further complicated by uncertainties regarding the formation height of molecules and atomic features with respect to the background plasma. Molecules may reveal the magnetic field strength at an optical depth (in the ambient plasma) below $\tau = 1$, in which case the gas pressure used to compute the equipartition magnetic field would be larger.