Superstrong photospheric magnetic fields in sunspot penumbrae

Recently, there have been some reports of unusually strong photospheric magnetic fields (which can reach values of over 7 kG) inferred from Hinode SOT/SP sunspot observations within penumbral regions. These superstrong penumbral fields are even larger than the strongest umbral fields on record and appear to be associated with supersonic downflows. The finding of such fields has been controversial since they seem to show up only when spatially coupled inversions are performed. Here, we investigate and discuss the reliability of those findings by studying in detail observed spectra associated with particularly strong magnetic fields at the inner edge of the penumbra of active region 10930. We apply classical diagnostic methods and various inversions with different model atmospheres to the observed Stokes profiles in two selected pixels with superstrong magnetic fields, and compare the results with a magnetohydrodynamic simulation of a sunspot whose penumbra contains localized regions with strong fields (nearly 5 kG at $\tau=1$) associated with supersonic downflows...


Introduction
Sunspots are the largest concentrations of magnetic flux on the solar surface and, due to the significant suppression of the convective motions caused by the strong fields, they are seen as dark regions on the photosphere.The various brightness levels observed in sunspots indicate differences in temperature, and therefore different magnetic regimes.Such a relation between brightness/temperature and magnetic field strength in the photosphere has been extensively studied (e.g., Lites et al. 1990;Solanki 1993;Keppens & Martínez Pillet 1996;Mathew et al. 2003Mathew et al. , 2004;;Tiwari et al. 2015).
Thus, the strongest magnetic fields in sunspots are generally found within their central dark regions or umbrae, where the field is closely vertical with typical strengths between 2.5-4 kG (e.g., Livingston 2002).The largest field strength ever recorded within an umbral region is nearly 6.1 kG (Livingston et al. 2006), but an even larger value (above 6.2 kG) was recently reported by Okamoto & Sakurai (2018), and was observed in a light bridge.Moreover, field strengths in excess of 7 kG have been found within sunspot penumbrae (van Noort et al. 2013;Siu-Tapia et al. 2017).Penumbrae are the less dark regions partially or completely surrounding umbrae and where the magnetic field is filamentary, strongly inclined, and generally weaker than in the umbrae (typically the field varies from about 1.5-2.5 kG at the inner penumbral boundary to 0.5-1 kG towards the outer boundary) (e.g., Lites et al. 1990;Skumanich et al. 1994;Westendorp Plaza et al. 2001).
Sunspot penumbrae are additionally characterized by a gas outflow with speeds of several kilometers per second, the socalled Evershed flow (EF; Evershed 1909), which is directed G −1 ], for the pair of Fe I absorption lines at 6301.5 and 6302.5 Å in the atmosphere inferred by the SPINOR 2D inversions for the penumbral pixel marked with a blue '+' in Fig. 3a.Some physical parameters for such an atmosphere are presented in Table 2.The horizontal dashed lines indicate the selected nodes for the SPINOR 2D inversion: log(τ) = 0, −0.8 and −2.
along the penumbral filaments, i.e. along the bright elongated channels with more inclined fields within the penumbra itself.Tiwari et al. (2013), in a detailed analysis of the penumbral fine structure, found that the EF has its sources in hot upflows that occur at the inner endpoints of the penumbral filaments (heads) and part of it sinks in concentrated cooler dowflows that occur at the outer endpoints (tails) where the magnetic field bends over vertically and displays a strengthening of about 1.5-2.5 kG on average, so that it can reach up to 3.5 kG in individual tails.
The superstrong penumbral fields reported by van Noort et al. (2013) (reaching up to 7.5 kG) were observed in some tails of complex penumbral filaments, i.e. those with a single tail but with more than one head, near the penumbral periphery in supersonic downflow areas (estimated line-of-sight velocities of up to 22 km s −1 ).Such values were inferred by using a sophisticated inversion technique that takes into account the spatial coupling between the pixels of the observed image when considering the instrumental effects that cause the image degradation (SPINOR 2D, van Noort 2012;van Noort et al. 2013).
In contrast, the superstrong fields reported by Siu-Tapia et al. (2017) (> 7 kG based on SPINOR 2D inversions) were observed in the tails of inverted penumbral filaments carrying a counter-EF (CEF), i.e. a gas inflow towards the sunspot umbra at photospheric heights, near the inner penumbral boundary.Unlike in van Noort et al. (2013) and Tiwari et al. (2013), Siu-Tapia et al. (2017) considered the outer endpoint of the CEF-carrying filaments to be heads, since at such places the sources of the CEF were found; and the inner endpoints of the CEF-carrying filaments to be tails, where the sinks of such flow were observed.Therefore, they found on average field strength of ∼4.5 kG (when excluding all those fields in excess of 7 kG) at the tails of the CEF-carrying filaments.Such superstrong penumbral fields were also found to be associated with the supersonic downflows occurring at the sinks of the CEF.
The superstrong penumbral fields reported by both van Noort et al. (2013) and Siu-Tapia et al. (2017) are unusual and are even stronger than the strongest umbral fields found by Livingston et al. (2006) of ∼ 6.1 kG.Moreover, in both studies, spatially coupled inversions were employed, and therefore, the reality of these extraordinary field strengths needs to be confirmed with other techniques.The existence of field strengths in excess of 7 kG possibly has an impact on the approximations made by the inversion code, which can lead to errors in the density stratifications if the non-vertical field components are significant, and therefore in the atmospheric stratifications of the physical parameters inferred by the code.Hence, the reliability of the inversion results for such peculiar pixels needs to be tested.
This work is aimed at discussing the reliability of those penumbral fields in excess of 7 kG reported by Siu-Tapia et al. (2017) concentrated at the inner boundary of the penumbra of the main sunspot in active region (AR) 10930.Here we examine the reality of such findings by applying some classical diagnostic methods and various inversion techniques with different model atmospheres to the Stokes profiles observed with the spectropolarimeter (SP) of the Solar Optical Telescope (SOT) on board Hinode.Finally, since the mere possibility of such superstrong fields in sunspots leads to the question about the necessary physical conditions for their appearance, we analyze the physical structure of a sunspot MHD simulation with CEFs which was reported by Rempel (2015) (see also Siu-Tapia et al. 2018).
This work is organized as follows: In section 2 we describe our data and inversion technique.In section 3, some classical diagnostic methods are applied to peculiar pixels with extreme spectra.In section 4, various inversions with different model atmospheres are considered.The results are presented in section 5.In section 6, we compare the observations with MHD sunspot simulations and study the mechanisms that can amplify the magnetic field in the penumbra.Finally, in section 7 we discuss our results and draw our conclusions.

Observations and spatially coupled inversions
For this study we use the spectropolarimetric observations of AR NOAA 10930 from the Hinode SOT/SP instrument (Lites et al. 2001(Lites et al. , 2013)), which simultaneously measures the full Stokes profiles of a pair of absorption lines of Fe I that are formed in the lower solar photosphere with central wavelengths at 6301.5 and 6302.5 Å (having effective Landé factors g L =1.67 and 2.49, respectively).
The SP instrument scanned the main sunspot in AR 10930 (whose umbra displayed a negative magnetic polarity) on Dec 8 2006 at an heliocentric angle of ∼47 • while operating in normal mode, i.e., using an exposure time of 4.8 seconds per slit position and a spatial sampling of 0.16 .These observations were corrected for dark current, flat field, orbital drift and instrumental cross-talk by reducing the row data with IDL routines of the Solar-Soft package (Lites & Ichimoto 2013).
As described by Siu-Tapia et al. (2017), the atmospheric properties of the main sunspot in AR 10930 were derived by inverting the observed Stokes profiles with the SPINOR 2D inversion code (van Noort 2012;van Noort et al. 2013), which is the spatially coupled version of the SPINOR inversion code (Frutiger 2000), based on the STOPRO routines (Solanki 1987) that solve the radiative transfer equations for polarized light.
The method is able to invert 2-D maps of spectropolarimetric data that have been degraded spatially in a known way (information that is contained in the point-spread function or PSF of the telescopes).Then, the code is able to use the information contained in the spectral dimension and the known spatial degradation properties to constrain a parameterization of the atmosphere over the whole FOV.The image degradation is applied to the solution rather than by deconvolving the original data themselves, which is the clasical approach, while the code performs a coupled inversion of all the pixels simultaneously.
According to van Noort (2012), the spatially coupled inversion method is stable to oversampled data and produces an inversion result with a resolution up to the resolution limit of the telescope.To be able to reach the diffraction limit of Hinode SOT, we oversampled all Stokes maps by a factor two, to 0.08 , following van Noort et al. (2013).Thus, by considering a spatial grid that is denser than the original data and by employing a single-component atmospheric model per pixel, the inverted atmospheric parameters returned by the code correspond to the best fits to the Stokes profiles once the blurring effect of the telescope's PSF has been taken into account.This significantly improves the spatial resolution, allowing structures at the diffraction limit of the telescope to be properly resolved.For the SPINOR 2D inversions used in this work (cf.Siu-Tapia et al. 2017) all the atmospheric free parameters were initially defined at three optical depth nodes, placed at log(τ) = −2, −0.8, and 0. This set of optical depth nodes was shown to work well for the inversion of Hinode/SP observations of sunspots by Siu-Tapia et al. (2017).At each of the three chosen nodes the temperature T , magnetic field strength B, field inclination relative to the line-of-sight γ LOS , field azimuth φ, line-of-sight velocity v LOS , and a microturbulent velocity v MIC , were fitted, leading to 18 free parameters in total.
In Figure 1, the response functions (RFs) of Stokes I and V to the magnetic field B are plotted as functions of wavelength and optical depth for the pair of Fe I absorption lines at 6301.5 and 6302.5 Å.To compute these RFs we have used the atmospheric parameters retrieved by SPINOR 2D inversions (see Table 2) in a penumbral pixel with superstrong field (labeled as LFP 1 and marked with a blue '+' in Fig. 3a).These plots show that the selected nodes at log(τ) = 0, −0.8 and −2 for performing our inversions all formally lie within the formation region of the lines.In particular, the lowest node at log(τ) = 0 is also sensitive to the magnetic field changes, which means that information about the magnetic field at such node is available for the analyzed wavelength range.

Sunspot's magnetic and velocity field as inferred by SPINOR 2D
As reported by Siu-Tapia et al. ( 2017), the SPINOR 2D inversions return very large magnetic field strengths, B ≥ 5 kG at all three height nodes (log(τ) = −2.0,−0.8 and 0) in the inner center-side penumbra of the main sunspot in AR 10930 (see field strength map in Fig. 2a for the height node at log(τ) = 0 only).Such strong magnetic fields are located at places that coincide with supersonic sinks of the counter-EF (CEF, see Fig. 2b and Figs. 3 and 5 in Siu-Tapia et al. (2017)).These are among the largest field strengths ever observed in penumbral environments (see also van Noort et al. 2013).
In particular, field strengths above 7 kG (whiter regions near the inner penumbral boundary in Fig. 2a) are even larger than the largest umbral fields found by Livingston et al. (2006) of ∼ 6.1 kG, who also found that only a very small fraction of sunspots (around the 0.2% in a 9-decade record of ∼ 32000 sunspots) have umbral fields stronger than 4 kG.
The scatter-plot in Figure 3b shows the SPINOR 2D inversion results for the magnetic field inclination angle in the localreference-frame (LRF) after the azimuthal disambiguation was resolved with the Non-Potential Magnetic Field Computation method (NPFC, Georgoulis 2005) versus the line-of-sight (LOS) velocity at log(τ) = 0, for all the 226 pixels where the SPINOR 2D inversions return B ≥ 7 kG (hereafter large field pixels or LFPs, yellow markers in Fig. 3a).
Most of the LFPs (see main population in Fig. 3b) appear associated with supersonic LOS velocities, i.e. around v LOS = 10 km s −1 (the sound speed is around C s ∼ 5 − 8 km s −1 in penumbrae, e.g., van Noort et al. 2013) and even larger than 20 km s −1 in some pixels.The peak of the LFPs distribution occurs near γ LRF = 140 • .
Therefore, according to the SPINOR 2D inversions, the LFPs are of umbral polarity (which is negative) with a large longitudinal component and, under the assumption of field-aligned flows, they contain supersonic downflows at log(τ) = 0. See Siu-Tapia et al. (2017) for more details on their brightness and thermal structure.
Fig. 4: Observed Stokes profiles (dashed curves) in two pixels (left: pixel 1, right: pixel 2) located at the inner penumbra of the CEF region (see blue crosses on maps displayed in Fig. 3a).From top to bottom: Stokes I, Q, U and V.The vertical lines in Stokes I and V panels indicate the position of the λ − (blue) and λ + (red) used to calculate the magnetic field strength with three different methods (see the main text for a description of the methods): the solid lines in the Stokes V panels show the splitting derived from method 1; the solid lines in the Stokes I panels correspond to the line splitting obtained with method 2; and, the dashed lines in both I and V panels indicate the splitting derived from method 3. The vertical dashed green line shows the wavelength position used to separate the two Fe I lines.7).

Analysis
As a very simple alternative estimation of the field strength in the selected LFPs, we compute the magnetic field strength directly from the splitting of the two observed Fe I line profiles at 6301.5 Å (line 1) and 6302.5 Å (line 2), using the Zeeman splitting formula: where λ 0 is the central wavelength of the line, λ ± are the centroids of the right and left circularly polarized line components (σ-components), m and e are the electron mass and charge, g L is the effective Landé factor of the transition (g L = 1.67 for line 1 and g L = 2.5 for line 2), and c is the speed of light.Given the huge field strengths inferred by the SPINOR 2D inversions, the Zeeman splitting should be complete, so that this approach is expected to prove the inversion results in case the super-strong fields are real.However, the Zeeman splitting approach can be misleading if there are strong gradients with height or multiple components in the resolution element.A strong magnetic field near optical depth unity which rapidly decreases with height might not display a complete Zeeman splitting in the profiles but very broad wings of the Stokes I and V profiles.
The critical point about using Equation 1 is determining λ − and λ + due to the large asymmetries observed in all four Stokes profiles from the LPF.We use three ways: Method 1) Selecting λ − where Stokes V takes on the largest negative value in the blue wing of the corresponding line, and λ + where Stokes V is largest in its red wing.Method 2) Placing λ − and λ + where the two deepest minima of Stokes I away from the central wavelength are found (using line 2 only, since line 1 is insufficiently split).Method 3) Applying the center of gravity method or COG method (Semel 1967(Semel , 1970;;Rees & Semel 1979;Cauzzi et al. 1993), in which λ ± = λ ± COG , where λ ± COG are the center of gravity wavelengths of the centroids of the right and left circularly polarized components, respectively, of the corresponding line, i.e.: The position of λ ± COG and the wavelengths delimiting the integration intervals in Equation 2for each of the Fe I lines, are indicated in Figure 4 for the observed spectra in the two selected LFPs (dashed vertical lines).The resultant field strengths obtained with the methods described above are listed in Table 1, for both LFPs.
The B values computed with methods 1 and 2 are indeed very large, lying between ∼ 7 kG (obtained from line 2) and ∼ 9.4 kG (from line 1) in both LFPs.On the contrary, the COG method provides considerably smaller field strength values in both LFPs: B ∼ 2.6 kG from line 1 and 3.6 − 3.9 kG from line 2. The difference between the COG and direct splitting methods likely stems from the fact that the profiles are not simple Gaussians, but show complex shapes indicating a range of field strengths (of which methods 1 and 2 sample only the largest), but is also partly due to the non-longitudinal direction of the field (methods 1 and 2 determine field strength, while method 3 only gives the LOS component).
Scatter plots in Figure 5 show the results of all three methods applied to all 226 LFPs found in the penumbra.Method 1 (top panels) returns mainly two different solutions in the two Fe I lines, one of which is centered at ∼ 3.5 kG in both lines and a second solution that is centered at ∼ 9 kG for line 1 and at ∼ 7 kG for line 2. The fact that method 1 senses two very different field strengths in both cases is a consequence of the complex shape of the multi-lobed Stokes V profiles in the observed spectra, see e.g. Figure 6 which displays four different shapes of Stokes V profiles that prevail among the 226 LFPs.For profiles like the one shown in Figure 6a, method 1 selects the innermost lobes in both lines (given that in this kind of profiles the inner lobes on the red wing of the lines satisfy the criterion used by method 1 to select λ + ), therefore returning field strengths that are almost half of the SPINOR 2D inversion result (B S P2 ∼ 7.3 kG).Similar situations occur in profiles with shapes as shown in Figures 6b  and 6d, where the inner lobes on the red wings for line 1 and line 2, respectively, are the largest.In contrast, in profiles like those shown in Figures 6b (line 2), 6c (both lines), and 6d (line 1), the method selects the most external lobes, therefore computing very large field strengths as displayed in Figure 5 (top panels).
In contrast, method 2 (middle panel in Fig. 5) finds a field strength distribution centered at ∼ 7 kG; whilst method 3 (bottom panels) mainly senses field strengths of the order of 2.5 kG for line 1 and around 3 kG for line 2.
The correlation between the field strengths from the methods and the results from SPINOR 2D are very low, particularly for method 3.However, there is some level of correlation between SPINOR 2D and the cloud of pixels displaying the strongest field solution in method 1.Likewise, there is a good correlation between the clouds of pixels displaying the weaker field solution in method 1 with the solutions in method 3.
The large discrepancy between the methods has to do with the fact that only the LFPs have been plotted.Ideally, to check if the methods are consistent at lower field values but depart from each other only at strong fields, scatter-plots with a larger sample of pixels covering a broader range of magnetic field values in the x axis are needed.However, it is not possible to apply the direct Zeeman splitting method to profiles that do not show distinguishable sigma components, i. e. to most of the penumbral pixels, except for the LFPs.
For the COG method, differences are expected due to velocity gradients, vertical field gradients, temperature, etc., that are not taken into account by the method (the gradients in B and v LOS are not taken into account in Eqs. 1 and 2).The presence of such gradients is indicated by the asymmetries of the Stokes profiles.Moreover, if there are spatially unresolved areas containing field inhomogeneities (e.g.different field strengths and/or different field polarities located next to each other within the same resolution element), then the resultant line profiles from such pixels will be the summation of all the different magnetic components, so that the number of lobes in the observed Stokes V profiles will depend on the number of unresolved components (e.g., Stenflo 1993).Most Stokes V profiles from the 226 LFPs display more than two lobes, suggesting the presence of different magnetic field components in each of those pixels.Nonetheless, according to methods 1 and 2, the very large wavelength separation between the most external lobes observed in the Stokes V profiles could still reflect that one of the field components is particularly strong.In such a scenario, an important aspect to consider is if the Paschen-Back effect would play an important role under the presence of such strong fields, i.e. if the splitting of atomic levels caused by such strong external magnetic fields would dominate over the LS -coupling.
The Paschen-Back regime occurs when ∆E ik µ B g L MB, where ∆E ik is the energy difference in the atom between the terms of the multiplet structure, µ B is the Bohr magneton, g L is the Landé factor, M is the magnetic quantum number, and B is the magnetic field strength.The Fe I 5 P 2 − 5 D 2 λ = 6301.5Åand the Fe I 5 P 1 − 5 D 0 λ = 6302.5Ålines belong to the same multiplet No. 816 (Moore 1945).The minimum energy difference of the lower levels of these lines is ∆E ik ≈ 0.032 eV and leads to a magnetic field 'threshold value' for the Paschen-Back effect of the order of 10 3 kG, a value that is well above the field strength values inferred by methods 1, 2, and the SPINOR 2D Occurrence (%)  (red) used to calculate the magnetic field strength with method 1.The headers of the plots indicate the magnetic field strengths as derived from: SPINOR 2D inversions (B S P2 ), method 1 applied to the Fe I line at λ = 6301.5Å(B Z1 ), and method 1 applied to the Fe I line at λ = 6302.5Å(B Z2 ).
inversions.Nonetheless, according with laboratory experiments the Paschen-Back effect actually takes place under the presence of magnetic fields with strengths 10 − 100 kG when the multiplet splitting energy difference is so small that the two adjacent lines of the same multiplet are separated by a distance of less than 1 Å (Frisch 1963).In the present case, since the analyzed pair of Fe I lines are separated by around 1Å and the Zeeman splittings being discussed are seemingly very large (corresponding to field strengths approaching 10 kG in some LFPs according to method 1), it is possible that the Paschen-Back effect starts playing some role in the magnetic splitting for those cases.However, even within the Paschen-Back regime, the magnetic field could still be measured with sufficient accuracy according to the laboratory results of Moore (1945).
Another plausible scenario is that the observed multi-lobed Stokes V profiles are produced by two unresolved atmospheric components that display large differences in their Doppler velocity (e.g., Solanki & Montavon 1993;Martínez Pillet 2000;Borrero & Bellot Rubio 2002;Schlichenmaier & Collados 2002;Bellot Rubio et al. 2004).One of the components could be associated with the umbral magnetic field in the sunspot (i.e.nearly at rest) and the second one with the filamentary CEF penumbra (strongly redshifted).This possibility is considered in the following section based on the fact that all the LFPs appear to be mostly located near, or at the umbral/penumbral edge (see Figs. 2a and 3a).

Inversions
In order to gain more insight into the reliability of the large magnetic field strengths returned by the SPINOR 2D inversion code, we now apply 5 additional inversions, considering two atmospheric components in some of them.It is noteworthy that the inclusion of a second atmospheric component causes the number of free parameters, n, to increase to almost twice that in a single component model used in the SPINOR 2D inversions.While this can and should lead to a much better fit for a complex Stokes profile (lower χ 2 ), it also increases the risk of obtaining artificial (unphysical) results.
The merit functions, χ 2 , are not computed in the same way by the different inversion codes; and hence, they are in principle not comparable among each other.To perform a valid comparison between the different fits, in the following, the minimum χ 2 is chosen to be the sum over the squared differences between the observed and the synthetic profiles resulting from the best fits in each case.
In Tables 2 and 3, we show some of the parameters corresponding to the best fits of the Stokes profiles in the two selected LFPs, obtained from all 6 different inversions.These inversions can be classified into two categories:

Height-dependent inversions
The first category includes height dependent inversions (Table 2), i.e. inversions in which all the parameters are allowed to vary with optical depth (with nodes being set at log(τ) = −2.0,−0.8 and 0), such as (a) SPINOR 2D inversions, with 18 free parameters; (b) SPINOR 1D (Frutiger et al. 2000) 2-component inversions, applied to the observed Stokes profiles (best fits are shown by green curves in Fig. 7), with 37 free parameters; (c) SPINOR 1D single-component, applied to the deconvolved Stokes profiles, with 18 free parameters and (d) SPINOR 1D 2-component inversions, applied to the deconvolved Stokes profiles, with 37 free parameters.
Inversions (c) and (d) are intended to resemble the SPINOR 2D technique as far as accounting of the instrumental effects over the observed Stokes profiles is concerned, although the treatment is not as consistent as that by SPINOR 2D.We retrieve the deconvolved Stokes profiles from the spatially degraded observed Stokes profiles by using an effective point-spread function (PSF) (Danilovic et al. 2008) constructed from the pupil function of the 50-cm Hinode SOT (Suematsu et al. 2008) and applying the Richardson-Lucy deconvolution method (see Richardson (1972) and Lucy (1974) for details).The resultant deconvolved Stokes profiles in the selected pixels are displayed in Figure 8 (black dashed curves) together with their best fits obtained by SPINOR 1D (c) and (d) (orange and green curves in Fig. 8, respectively).

Height-independent inversions
The second category corresponds to height independent inversions (Table 3), i.e., we assume no variation with optical depth of atmospheric parameters by using (e) a 2-components Milne-Eddington inversion (Skumanich & Lites 1987;Lagg et al. 2004Lagg et al. , 2009;;Borrero et al. 2011) applied to the observed Stokes profiles (best fits are shown by purple curves in Fig. 7), with 15 free parameters; and (f) 2-component 1-node SPINOR 1D inversions applied to the deconvolved Stokes profiles (best fits are shown by purple curves in Fig. 8), with 17 free parameters.

Results
The SPINOR 2D best-fits give B ∼ 8.3 kG at log(τ) = 0 for the profiles in pixel 1, while B = 8 kG is obtained for both pixels at log(τ) = −0.8(see Table 2), with χ 2 = 14 and 35 for each fit, respectively.Note that these fits do not succeed in perfectly reproducing the reversed central lobes of the Stokes V profiles observed in both pixels.The deconvolved spectra displayed in Figure 8 also show the reversed central lobes in Stokes V, which suggests that they are not a result of the mixing of signals due to instrumental effects.Nonetheless, inversions (c) also fail in reproducing the central reversed lobes, providing a much poorer fit than SPINOR 2D (χ 2 = 44 and 47 for each pixel, respectively, when taking into account an increased noise of ∼ 5σ in the deconvolved Stokes profiles caused by the deconvolution itself) and featuring B ∼ 7 kG at log(τ) = 0 in both pixels.Thus, in both cases, the 1-component inversions cannot fit these four-lobe While the height-dependent 2-component inversions produce better fits to the four-lobe profiles than the 1-component inversions, χ 2 = 8 and 12 for each pixel respectively in inversions (b), and χ 2 = 14 and 18 respectively in inversions (d), they also use a much larger number of free parameters than the 1-component inversions so that a comparison is not straightforward.The components 1 in inversions (b) and (d) give B ∼ 4−4.2 kG and v LOS ∼ 1 km s −1 at log(τ) = 0 in both pixels, roughly consistent with the umbral environment in the vicinity of those pixels (B ∼ 4.2 kG at log(τ) = 0, according with SPINOR 2D).Furthermore, the v LOS in the components 1 are relatively small at all three atmospheric layers, as expected for umbral environments.In contrast, the components 2 in inversions (b) and (d) give B ∼ 4.4 and ∼ 4.9 kG, respectively, for pixel 1 at log(τ) = 0; and B ∼ 5 and 5.8 kG, respectively, for pixel 2; with v LOS 16 km s −1 and with the filling factors α, i.e. the mixing ratio between the 2 components, being slightly larger for the second components.The second component in both pixels could then correspond to the tails of the penumbral filaments harboring the CEF, with a large redshift and stronger fields than in the surrounding umbra; which is qualitatively compatible with the results from SPINOR 2D, but quantitatively suggests lower penumbral field strengths of the order of 4 − 5 kG (although approaching 6 kG in one pixel).
The SPINOR 1D inversion code can find a solution involving 2 components, one of which is strongly wavelength shifted to mimic the seemingly very strongly split spectral line.This nearly halves the field strength, although even in this case, we get B values reaching up to nearly 6 kG.These are still very large field strengths and are atypical for penumbral environments.They are close to the record measurement of 6.1 kG in sunspot umbrae (Livingston et al. 2006).However, due to the large number of free parameters involved, it is difficult to judge if the results they provide are more reliable than those from SPINOR 2D inversions.
A formal approach to compare the different results obtained from inversions that consider different model assumptions would be through the comparison of their error bars.Unfortunately, specifying the uncertainties in the fitted atmospheric parameters that take into account the possible degeneracies between parameters is an intrinsic difficulty facing inversions.Especially in the case of the spatially coupled inversions, the changes in the parameters of a single pixel severely affect the result, and therefore the uncertainties of the parameters, of the neighboring pixels.This fact makes the computation of formal errors for a single pixel principally impossible, and therefore the SPINOR 2D inversion code does not provide errors.
As a simple proxy for the quality of the fits, we compare the models by using the Bayesian Information Criterion (BIC; Schwarz 1978), which is based on the crude approximation of Gaussianity of the posterior with respect to the model parame- ters: where χ 2 min is the merit function of the best-fits to the Stokes profiles in each model, n is the number of free parameters and N is the number of observed points.The computed values of the BIC for each fit are shown in Tables 2 and 3.The model with the smallest value of the BIC is the preferred one.The heightdependent model preferred by the BIC is SPINOR 2D in both pixels.However, one of the fundamental problems of this criterion is that it penalizes all parameters equally, not taking into account situations in which data do not constrain some parameters (see e.g.Asensio Ramos et al. 2012).
The height-independent 2-component inversions (e) and (f) give nearly identical results to each other in both pixels, with lower field strengths of around 3.5 kG in component 1, which is almost at rest (v LOS < 1 km s −1 ) compared to component 2 in which B ∼ 4 kG and v LOS ∼ 16 km s −1 .However, the resultant values of B and v LOS given by the two height-independent inversions (e) and (f) generally resemble the results from the heightdependent 2-component inversions (b) and (d) at log(τ) = −0.8.This is not surprising since the sensitivity to v LOS and B pertur- bations is higher at log(τ) = −0.8than at log(τ) = 0 for the Fe I 6302.5 Å line (see, e.g., response functions computed by Cabrera Solana et al. ( 2005) and Fig. 1).Even if the absorption line is not formed at a single depth, the height-independent inversions mainly provide information on the physical conditions prevailing at depths at which the line is more sensitive.Thus, if stronger magnetic fields are present in deeper layers of the solar atmosphere (e.g. at log(τ) = 0), as suggested by the results of inversions (a), (b), (c) and (d), they cannot be retrieved by inverting the Stokes profiles of the current wavelength range with a height-independent inversion technique only.
The BIC values obtained for inversions (e) and (f) are only slightly better than the ones obtained for SPINOR 2D in both pixels.Even if they both succeed in capturing many relevant aspects of all four Stokes profiles (despite the increased noise in the deconvolved spectra), it is very unlikely that the physical parameters in the pixels of interest display no gradients with height.Generally in sunspots, one would rather expect large gradients with height of the physical parameters (particularly of the magnetic field) as supported by MHD simulations.In such cases, height-dependent inversions provide a more appropriate model atmosphere.The SPINOR 2D inversions are the most reliable model in this sense, since they take into account the height stratification of the physical parameters while keeping a good balance between the quality of the fit and the number of free parameters in the model, according to the obtained BIC values in the two studied pixels with peculiar spectra.

Strong photospheric penumbral fields in a MURaM MHD simulation
We now use the 3D high-resolution sunspot simulation by Rempel (2015) (see also Siu-Tapia et al. 2018) with a pixel resolution of 48 km in the horizontal direction and 24 km in the vertical direction to investigate the possible origin of super-strong penumbral magnetic fields associated with supersonic downflows in CEF-carrying filaments, and to compare their synthetic spectra with the observed Stokes profiles reported in the previous sections.
As reported by Rempel (2015) and Siu-Tapia et al. ( 2018), the sunspot simulation covers a time-span of 100 solar hours and after t ∼ 50 hours, it displays radially aligned penumbral filaments with fast Evershed ouflows along them; but in some regions of the penumbra, the filaments carry instead a CEF, i.e. radial inflows directed toward the sunspot umbra and strong downflows (v z <−8 km s −1 ) at the end of such filaments, i.e. in their inner endpoints where the magnetic field is noticeably enhanced, up to values of around 5 kG and γ ∼ 40 • near the local τ = 1 level (see for example Fig. 9, cf. Fig. 2 in Siu-Tapia et al. (2018)).
In a rather simple and quick attempt to quantify the effect of the 3D atmospheric structure and magnetic field on the profiles of the analyzed spectral lines, we employed the forward part of the SPINOR code (STOPRO) to solve the radiative transfer equation in the MURaM cube analyzed by Siu-Tapia et al. (2018), which was obtained with non-gray radiative transfer.Figure 10a shows the emergent synthetic spectra from a vertical column located close to the inner edge of the simulated penumbra and which, at the log(τ) = 0 level, intersects with the tail of a CEF-carrying filament and contains supersonic downflows at that height.Similarly to the observed Stokes profiles for the pair of Fe I lines from the LFP in AR 10930, the synthetic Stokes profiles from the MHD simulations are highly asymmetric, display large redshifts, and show multi-lobed Stokes V profiles.
In Figure 10b, we display degraded Stokes profiles which were obtained by convolving the synthetic spectra in Fig. 10a with an effective PSF=0.16(Danilovic et al. 2008) constructed from the pupil function of the 50-cm Hinode SOT (e. g., CEF sinks (inner penumbra) Fig. 12: Average contributions to the induction equation from advection (solid black), field-line stretching (solid red), flow divergence (solid blue) and a residual term (solid green) to the vertical field component at the sinks of the NEF (left) and at the sinks of the CEF (right).The averages have been focused on localized regions that have fairly strong fields (B z < −2 kG for the sinks of the NEF and B z > 2 kG for the sinks of the CEF) and supersonic downflows at the local log(τ) = 0 level.The residual term represents the numerical magnetic diffusivity as an approximated magnitude that is calculated by the negative sum of the advective, stretching and divergent terms; but it also contains potential contributions from time variation.The dashed and dotted lines represent approximated terms that have a major contribution to the advective term (black) and to the stretching term (red), whose general expressions are indicated in the lower labels.The dash-dotted blue line represents the contribution to the divergence term due to flows perpendicular to the vertical field component, i.e. radial and azimuthal flow components.
Such complex shapes in the MHD case are mainly the result of the large vertical gradients in all the atmospheric quantities and the magnetic field structure.In particular, the stratification of the field strength (Figure 11a), the flow velocity (Figure 11b), and field inclination (Figure 11c), qualitatively resemble the SPINOR 2D results in the LFPs (cf.Table 2), i.e., while the field strength and downflow speeds increase with depth, the field inclination increases with height.Another important aspect that might contribute to the complexity of the emergent spectra is the presence of a shock which is seen as the sudden transition of the downflow speed from subsonic to supersonic near log(τ) = −1, and as the sudden variation of all the physical parameters shown in Figure 11.
The presence of strong magnetic fields at the local log(τ) = 0 level (in excess of 5 kG) in the simulations are mainly due to the influence of the neighboring umbral field and the highly depressed surfaces of constant optical depth as reported by Siu-Tapia et al. (2018).Such surface depression is of the order of 400−600 km in the analyzed case (see Figure 11f), and is related to a strong decrease of the gas density beneath log(τ) = −1 (see Figure 11e).Thus, given the similarity of the physical scenarios that the observations and the simulations present, the determination of the physical processes involved in the maintenance and amplification of the field in the simulations can give us insight into the possible origin of super-strong photospheric magnetic fields observed in sunspot penumbrae.

Induction equation
In order to study how the vertical field is maintained in the penumbra, we evaluate the different terms of the induction equation in the simulation box: We analyze the different terms in eq. 4 during the time-period from 60 to 70 solar hours (the range of time during which the counter-EF are found in the simulations) of the total of 100 hours simulations run (see Siu-Tapia et al. (2018) for details), using the transformation to cylindrical coordinates to separate the direction along and perpendicular to the penumbral filaments, i.e. r, θ, and z coordinates; and by separating outflows (v r >0) from inflows (v r <0), and upflows (v z > 0) from downflows (v z < 0) in the simulated penumbra.
As reported by Siu-Tapia et al. ( 2018), the vertical field component is noticeably enhanced at the tails of the penumbral filaments, i.e. at the filament end-points hosting sinks in both, those filaments carrying a normal-EF (NEF, i.e. outflows) and those carrying a counter-EF (CEF, i.e. inflows), which are mainly located close to the outer and inner penumbral boundary, respectively (see e.g.Fig. 3 in Siu-Tapia et al. (2018)).Therefore, we explore the mechanisms that can lead to the amplification of the magnetic field strength at those places.We use different masks in order to separate the sinks (downflows) that occur at the tails of the NEF-carrying filaments (regions where v z < 0 and v r > 0 in the outer penumbra) from those sinks that happen at the tails Fig. 13: Force balance in the radial (left) and vertical (right) directions.The force terms have been horizontally and temporally averaged over the places with strongest fields and supersonic downflows of the CEF in the inner penumbra.Black: pressure forces, red: Lorentz forces, blue: acceleration forces, green: residual forces.The vertical gray dashed lines are placed at the average height of the log(τ) = 0 level in the selected regions and are located nearly 400 km beneath the average height of the log(τ) = 0 surface in the quiet sun (i.e.z = 0 km, black dashed line).
of the CEF-carrying filaments (regions where v z < 0 and v r < 0 in the inner penumbra).
Figure 12 displays the contributions of each term in equation 4 to the vertical field component, horizontally averaged over the sinks of the NEF (left plot) and over the sinks of the CEF (right plot), and focusing the averages on localized regions that have fairly strong fields (B z < −2 kG for the sinks of the NEF and B z > 2 kG for the sinks of the CEF) as well as supersonic downflows at log(τ) = 0.In both cases, NEF and CEF sinks, the contributions from stretching, advection and divergence are mostly in balance, implying that the residual terms, which have potential contributions from the numerical magnetic diffusivity and from time variations (green lines), do not play a significant role in shaping the magnetic structure of the penumbra in these simulations during the analyzed time period.
The roles of advection and divergence in the vertical induction (black and blue solid lines respectively) are opposed to each other in both penumbral regions.Besides, they appear with a sign swap in the outer penumbra (left panel) compared to the inner penumbra (right panel).However, the roles of advection and divergence for maintaining the vertical field component are the same in both regions given that B z < 0 at the sinks of the NEF and B z > 0 at the sinks of the CEF, which causes the sign swap in the vertical induction.
Thus, at the sinks of the NEF (outer penumbra, left plot), there is an opposed contribution from the advective term to the maintenance of the (negative) vertical field component at all heights of the analyzed z−range (where z = 0 corresponds to the average height of the log(τ) = 0 surface in quiet sun).In contrast, the stretching term behaves as a source for the (negative) vertical field component in the outer penumbra, above z ∼ −200 km.The major contribution to this term comes from vertical stretching (red dotted line), i.e.B z ∂v z ∂z < 0 due to a strong downward transition of the downflow speeds (from subsonic to supersonic) that leads to the steepening of ∂v z ∂z near z ∼ −200 km (generally, the supersonic NEF sinks are shallower than in the CEF case, for an example see the white contours in the v z panel of Fig. 3 in Siu-Tapia et al. (2018) which enclose regions where v z < −8 km s −1 ).Deeper down, below z ∼ −200 km, there is a radial shear profile (red dashed line) due to a radial outward increase of the downflow speeds, B r ∂v z ∂r < 0, that also contributes to the maintenance of the (negative) vertical field component in the outer penumbra.However, the major source comes from the divergence term (blue line) due to the converging aspect of the downflows, i.e. ∇ • v < 0. The dominant contribution to this term in the near-surface layers is given by flows that are perpendicular to the vertical field component, v r and v θ , i.e. by a horizontal convergence of the downflowing material (dash-dotted blue line).The remainder is due to the term −B z ∂v z ∂z , which becomes negative below z ∼ −200 km.
Similarly, advection plays an opposite role for the maintenance of the (positive) vertical field component at the sinks of the CEF in the inner penumbra (right plot in Fig. 12).However, the role of stretching is not significant above z ∼ −400 km in this case (height at which most of the CEF sinks become supersonic, see an example in Figs.11b and 11f).Notwithstanding, due to the strong downward acceleration of the gas at the sinks of the CEF, vertical stretching acts as a source for the (positive) vertical field above z ∼ −400 km, i.e.B z ∂v z ∂z > 0. Likewise, a radial inward increase of the downflow speed at the sinks of the CEF leads to a radial shear term that contributes positively below z ∼ −400 km, i.e.B r ∂v z ∂r > 0. The dominant positive contribution to the (positive) vertical field component is given by the divergence term (i.e.∇ • v < 0, since B z is positive in the inner penumbra), which means that, similarly to the NEF sinks, the CEF sinks are also convergent downflows, which implies compression and amplification of the magnetic field.

Force balance
In order to investigate how the strongest fields in the penumbra are balanced, we follow the analysis performed by Rempel (2011)  Fig. 14: Radial (solid) and vertical (dashed) Lorentz force terms separated for the magnetic pressure (red) and the magnetic tension (magenta) forces at the places with strongest fields and supersonic downflows of the CEF in the inner penumbra.The radial (black solid) and vertical (black dashed) gas pressure forces are overplotted for comparison.Same format as in Fig. 13.
Each term in the above equation is then separated into a radial and a vertical component; and a residual force term is introduced as the negative sum of the pressure, the acceleration, and the Lorentz force contributions in the corresponding direction given that the viscosity force terms are not explicitly calculated.
Figure 13 shows the average radial and vertical force balance at the strong field regions of the inner penumbra with supersonic downflows of the CEFs.In both directions the forces are mostly in balance and the residual forces are nearly zero.These plots show that acceleration forces (blue lines) are very significant in the radial direction but almost negligible in the vertical one, which means that the system is close to magnetohydrostatic in the vertical direction given that the Lorentz force (red line) nearly completely balances with the pressure forces (black line).
The average height of the log(τ) = 0 level over the selected regions lies nearly 400 km below its average height in the quiet sun.Such average height is considerably deeper than in the whole penumbra (approximately at z = −200 km) and even deeper than in the inner penumbra (z ∼ −250 km).
Figure 14 displays the Lorentz force separately for the magnetic pressure term (−∇[B 2 /8π], red lines) and the magnetic tension term ([B • ∇]B/4π, magenta lines) in the radial (solid) and vertical (dashed) directions, which have been averaged over the strong field regions of the inner penumbra with supersonic downflows (same mask as in Fig. 13).Looking at the individual components of the Lorentz force, we see in the vertical direction that the strongest fields in the simulation are less force-free in the deeper domain, where the gas pressure (black dashed line) is strong enough to balance, but they become mostly force-free near the observable photosphere and in the layers above.In contrast, in the radial direction, the magnetic pressure force dominates over the magnetic tension force at most heights.However, the gas pressure force term (black solid line) is also large enough in the radial direction to considerably contribute to keep the balance.In addition, we have seen that the contribution of the radial acceleration forces (blue line in Fig. 13) plays an important role for maintaining the force balance in the radial direction.

Discussion and conclusion
Inversion techniques are currently the most powerful tools to infer the physical properties of the solar atmosphere from polarization line profiles, being able to provide reliable and robust results from many types of Stokes profiles according to numerical tests (e.g., Ruiz Cobo 2007).There are several different inversion techniques in the literature, each of them with its own advantages and shortcomings, which largely depend on the addressed problem.Certainly, after any Stokes inversion the results need always to be validated and one needs to be aware that the resulting model atmosphere is not necessarily the real one since the solution might not be unique or the model underlying the inversion not appropriate to the actual solar situation.Nonetheless, as stated by Sabatier (2000), by means of Stokes inversions it is generally possible to retrieve as much information as possible for a model which is proposed to represent the system in the real world.
The existence of B > 7 kG in the inner penumbra of a sunspot would require an unusually deep Wilson depression to be consistent with idealized magnetohydrostatic models of sunspots (Livingston et al. 2006).However, besides the non negligible dynamical effects of the studied penumbra, the magnetic field does not have to be in maximally non-force-free state as usually assumed for the photosphere.We have found in the MHD simulations that the strongest fields in the penumbra (∼ 5 kG) are vertically close to force-free in the observable photosphere and the gas pressure is sufficient to reach a force balance in the deeper layers.Although these fields are less force-free in the radial direction, the radial gas pressure force provides a sufficient balance to keep the system in near equilibrium, but is only with the contribution of radial acceleration forces that an almost complete balance can be reached.In this sense, the existence of 7 kG magnetic fields in the photosphere would be possible in a highly dynamical environment, such as that inferred by the SPINOR 2D inversions, i.e. with supersonic counter-Evershed flows sinking supersonically in the inner penumbra.Such superstrong field concentrations would likely fan out significantly with height and remain close to potential in the observable photosphere.
Field strengths larger than 7 kG in penumbral environments have previously been reported by van Noort et al. (2013).They obtain such large field strengths in supersonic downflow regions in the peripheral penumbra of a sunspot, also by using SPINOR 2D inversions.Nonetheless, they observe B > 7 kG only in the deepest layer (log(τ) = 0) and obtain a good agreement between the inversion results and a MURaM simulation of a sunspot by Rempel (2012).They propose a scenario in which the high magnetic field values are the result of a field intensification in the deep photosphere due to the interaction of the supersonic downflows with an external magnetic barrier (e.g., with a plage region).Such a scenario could also be valid for the present observations, with the umbral field playing the role of the magnetic barrier.
According to the SPINOR 2D results (e.g.Fig. 3b and Table 2), most of the magnetic fields whose strength is above 7 kG (LFPs) are nearly vertical and have the same polarity as the umbra, which is negative.In addition, they are associated to supersonic downflows, even exceeding v LOS = 20 km s −1 at log(τ) = 0 in some LFPs.Moreover, the regions in the sunspot where B > 5 kG are also associated with supersonic LOS flow velocities (see, e.g., Fig. 5a in Siu-Tapia et al. (2017)).As displayed in Figure 3a, the LFPs with B > 7 kG (yellow pixels) are surrounded by those weaker fields, which are still in excess of 5 kG (red pixels), in a supersonic flow environment (see also Fig. 3).A very low density of the supersonic downflowing material could also explain the observation of unusually strong magnetic fields in the penumbra, since it would cause the optical depth layers to be strongly depressed.In addition, their close vicinity to the umbral field also plays a role.This is in agreement with MHD simulations of counter-Evershed flows (Siu-Tapia et al. 2018, and Section 6), which show a ∼ 400 km depressed τ = 1 surface in average (and up to 600 km, with respect to its average height in the quiet sun) at the tails of the CEF-carrying filaments, where the supersonic downflowing material becomes a very low density gas (see also Fig. 11e).As a consequence, deep-lying field strengths of the order of 5 kG in the inner penumbra (near the umbra) are visible at the τ = 1 level in those simulations.
The resultant synthetic Stokes profiles associated to photospheric fields of the order of 5 kG and downflow speeds of 10−12 km s −1 in the simulated penumbra display large asymmetries, redshifts and multiple lobes in their Stokes V profiles, in agreement with the observed spectra in the LFPs.This result supports the possibility that the observed Stokes profiles associated with the LFPs in AR10930, which display even more extreme characteristics, are produced by larger fields and stronger downflows as inferred by the SPINOR 2D inversions (Fig. 3b).
The strong penumbral magnetic fields in the simulations (nearly 5 kG at the local τ = 1 level) are mainly due to the influence of the neighboring umbral field, the highly depressed surfaces of constant optical depth, and the formation of shocks by the transition of the downflow speeds from subsonic to supersonic; but there is also a local intensification of the field that can be associated to the converging aspect of the supersonic downflows, which lead to a compression and intensification of the vertical magnetic field component in the inner penumbra.
A similar mechanism amplifies the negative vertical field component at the sinks of the NEF in the outer penumbra in the simulations, where the field bends over and the downflows are convergent.There, the negative vertical field component can additionally be intensified by vertical stretching due to the strong downward acceleration of the gas up to supersonic speeds.These results are in agreement with the proposed scenario by van Noort et al. (2013) to explain the possible observation of superstrong fields associated with peripheral downflows in the penumbra of the leading spot of NOAA AR 10933.
According to the Bayesian Information Criterion (BIC, Schwarz 1978), which assesses a best fit based on the balance between the quality of the fit and the number of free parameters considered by the model, the preferred height-dependent model atmospheres for the LFPs are those provided by the SPINOR 2D inversions.Furthermore, given the extreme characteristics of the observed Stokes profiles associated with the sinks of the CEF in the LFPs and their resemblance to the synthetic spectra derived from the MHD simulations, we finally cannot easily discard the possibility that we are dealing with actual observations of B ∼ 7 kG and even larger in regions of supersonic downflowing material (up to v LOS ∼ 35 km s −1 ) with very low densities, where similar mechanisms to those occurring in the analyzed simulations might explain their origin.
The major strength of SPINOR 2D lies on its simultaneous coupled inversion of all the pixels to self-consistently take into account the influence of straylight from neighboring pixels.This approach is able to reproduce complex multi-lobed profiles with a simple, one-component atmosphere per pixel, keeping an acceptable number of free parameters, which significantly enhance the reliability and the robustness of the inversion results.However, we have seen that the highly complex observed Stokes pro-files cannot be perfectly reproduced with any of the presented inversion techniques without almost doubling the number of free parameters.The inherent complexity of these profiles may involve physical aspects that are not considered within the assumptions and approximations made by the inversion codes, which could lead to significant errors in the returned values.These assumptions include hydrostatic equilibrium, which is unlikely to be satisfied in such a dynamic environment.
Finally, the field estimations performed by means of the Zeeman splitting and the COG method (methods 1, 2, and 3) provide results that are not entirely consistent among each other.Unfortunately, all these methods are unable to take into account the errors from the instrumental effects.Moreover, methods 1 and 2 are reliable for ideal cases only, as when they are applied to single-component profiles produced by homogeneous magnetic fields, which is clearly not the case for the LFPs with highly asymmetric and multi-lobed Stokes profiles.As a consequence, the possibility of two (horizontally) unresolved structures with a large velocity difference cannot be ruled out either, in spite of the shortcomings of the 2-component height-dependent and heightindependent inversions.Unusually strong penumbral magnetic fields (5-6 kG) also show up as the more plausible physical solution in some of the 2 components models.However, the existence of two Doppler shifted components (one carrying nearly 15 km s −1 and another one with gas almost at rest) coexisting over several neighboring pixels would require an unresolved fine structure with sub-resolution canals of two types over an extended area.These extreme gradients in velocity required to be present in the one resolution element to produce the observed spectra are considerably less plausible.A solution where the pixels are smoothly connected (as in the 2D inversions), with a height gradient within the line forming region could represent a more plausible scenario and is in agreement with MHD simulations.
For future work, it would be interesting to investigate how the presence of superstrong magnetic fields in the photosphere affects the shape of the observed Stokes profiles in the pair of Fe I lines when considering the Paschen-Back effect.This will be useful to determine whether or not such effects need to be taken into account by the inversion codes when dealing with photospheric field strengths of the order of 7 kG or larger.
Table 2: Parameters resulting from 4 different height dependent inversions which were applied to the two sets of observed Stokes profiles ((a) and (b)) and to their corresponding deconvolved Stokes profiles ((c) and (d)) from the two selected LFPs.From left to right: number of free parameters n, pixel identification number P, atmospheric component C, optical depth node log(τ), field strength B, field inclination γ, LOS velocity v LOS , merit function χ 2 , filling factor α for each component, and Bayesian Information Criterion (BIC) value.Table 3: Parameters resulting from two height independent inversions: (e) two-components Milne-Eddington inversions applied to the two sets of observed Stokes profiles (black dashed lines in Fig. 7) and, (f) two-components SPINOR 1D inversions applied to the corresponding deconvolved Stokes profiles (black dashed lines in Fig. 8) of the two selected LFPs (blue markers on Fig. 3a).

Height-dependent inversions
Columns are in the same format as in Table 2.

Fig. 1 :
Fig. 1: Normalized response function of Stokes I (left) and Stokes V (right) to the magnetic field B, multiplied by 1/∆τ C [10 −5G −1 ], for the pair of Fe I absorption lines at 6301.5 and 6302.5 Å in the atmosphere inferred by the SPINOR 2D inversions for the penumbral pixel marked with a blue '+' in Fig.3a.Some physical parameters for such an atmosphere are presented in Table2.The horizontal dashed lines indicate the selected nodes for the SPINOR 2D inversion: log(τ) = 0, −0.8 and −2.
Fig. 2: A portion of the (center-side) penumbra of the main sunspot in AR 10930 (cf.Fig. 5 in Siu-Tapia et al. (2017)).(a) Map of the magnetic field strength at log(τ) = 0, as inferred with SPINOR 2D inversions.The arrow points towards the solar disk-center.(b) Line-of-sight velocity as inferred with SPINOR 2D inversions at log(τ) = 0. Positive values (red-to-yellow colors) indicate plasma motions away from the observer, while negative values (blue colors) indicate motions towards the observer.Solid and dashed black contour lines were placed at I c /I QS = 0.26 and I c /I QS = 0.94, respectively.The large penumbral sector where red-to-yellow colors dominate harbors a fast counter-EF which supersonically sinks (with speeds larger than 12 km s −1 ) near the inner penumbral boundary (Siu-Tapia et al. 2017).Notice that the color-bar scales have been saturated.
Fig. 3: (a) Normalized continuum intensity map as observed by the Hinode SOT/SP in AR 10930.The red markers indicate the location of 845 pixels where the SPINOR 2D inversions return B ≥ 5 kG at log(τ) = 0, of which 226 harbor fields larger than 7 kG at log(τ) = 0 (yellow markers, also named in the text such as large-field-pixels or LFPs).The blue markers highlight the location of two LFPs ('+': LFP 1, 'x': LFP 2) selected to test the robustness of our inversion solutions (see text).The arrow points towards the solar disk-center.(b) Scatter plot of the magnetic field inclination in the local-reference-frame γ LRF (after removal of the field azimuth ambiguity through the NPFC method) versus the line-of-sight velocity v LOS , from SPINOR 2D at log(τ) = 0, for the 226 LFPs.

Fig. 5 :
Fig. 5: Scatter plots of the 226 LFPs where SPINOR 2D returns B > 7 kG at log(τ) = 0 (B S P2 ).The y axis indicates magnetic field values obtained with each of the 3 alternative methods (direct Zeeman splitting and COG methods, B Z and B COG respectively) as described in the text and displayed in Table 1.From top to bottom: Methods 1, 2, and 3 for line 1 at λ = 6301.5Å(subscript 1, left plots) and for line 2 at λ = 6302.5Å(subscript 2, right plots).Dashed lines represent expectation values if both methods give identical results (white).
Fig. 6: Four examples of observed Stokes V profiles in the LFPs.The vertical lines indicate the position of the λ − (blue) and λ +(red) used to calculate the magnetic field strength with method 1.The headers of the plots indicate the magnetic field strengths as derived from: SPINOR 2D inversions (B S P2 ), method 1 applied to the Fe I line at λ = 6301.5Å(B Z1 ), and method 1 applied to the Fe I line at λ = 6302.5Å(B Z2 ).

Fig. 7 :
Fig. 7: Observed Stokes profiles (dashed curves), and their best fits returned by SPINOR 2D (orange curves), SPINOR 1D 2component height-dependent inversions (green curves) and Milne-Eddington 2-component inversion (purple curves) in the two selected LFPs.Plots are in the same format as plots in Figure 4.

Fig. 8 :
Fig.8: Deconvolved Stokes profiles (dashed curve), and their best fits returned by SPINOR 1D inversions, using a single-component height dependent model atmosphere (orange curves), a 2-component height-dependent model (green) and a 2-component height independent model atmosphere (purple curves) for the two selected LFPs.Plots are in the same format as plots in Figure4.
Fig. 9: A portion of the inner penumbra in the MURaM sunspot simulation by Rempel (2015) with some filaments hosting a counter-EF (see also Siu-Tapia et al. 2018).The maps show: (a) the magnetic field strength B [G]; (b) the field inclination with respect to the vertical γ [ • ], i.e. γ=0 • represents a vertical field of umbral polarity, γ=90 • a horizontal field, and γ=180 • a vertical field of opposite polarity to the umbra; (c) radial flow velocity v r [km s −1 ]; and (d) the vertical flow velocity v z [km s −1 ].Negative v r and v z values (red-to-yellow colors) indicate inflows and downflows, respectively.This sign convention differs from the one used in observational studies, where negative values denote flows moving towards the observer along the line-of-sight.The black contour lines where placed at I c /I QS < 0.45 near the umbra(left)-penumbra(right) boundary.All maps show the corresponding physical parameters at log(τ) = 0.
Fig. 10: (a) A set of emergent synthetic Stokes profiles in the MURaM sunspot simulation at the location of a supersonic downflow in the tail of a CEF-carrying filament.The spectra were synthesized by using the SPINOR code (STOPRO routines).(b) Degraded Stokes profiles obtained from the convolution of the synthetic spectra with a point-spread-function (PSF) of 0.16 arcseconds, similar to the case of the Hinode telescope.

Fig. 11 :
Fig. 11: Vertical profiles, in a log(τ)-scale, of the MHD physical parameters leading to the emergent synthetic Stokes profiles shown in Fig. 10, at the location of a supersonic downflow in the tail of a CEF-carrying filament: (a) magnetic field strength, (b) vertical flow velocity (negative values indicate downflows), (c) temperature, (d) field inclination, (e) gas density, and (f) geometrical height at the location of the grid-cell containing the supersonic downflows from the MURaM simulation.Vertical dashed lines delimit the approximate τ−range where the lines show a significant response, i.e., between log(τ) = −2 and 0.
and bySiu-Tapia et al. (2018) and use the following force balance equation which is derived from the momentum equation by assuming stationarity:ρg − ∇p Pressure + j × B Lorentz −ρ(v • ∇)v Acceleration

Table 1 :
Results of direct measurement of Zeeman splitting (Methods 1 and 2) and center-of-gravity (COG) method in two LFPs.