# Validation of the method for electric field inversions using a magnetic flux
emergence simulation

Andrey N. Afanasyev  DKIST Ambassador  Laboratory for Atmospheric and Space
Physics, University of Colorado Boulder
1234 Innovation Drive, Boulder, CO 80303, USA  National Solar Observatory,
University of Colorado Boulder, Boulder, CO, USA  Institute of Solar-
Terrestrial Physics of SB RAS, Irkutsk, Russia  Maria D. Kazachenko
Department of Astrophysical and Planetary Sciences, University of Colorado
Boulder, Boulder, CO, USA  National Solar Observatory, University of Colorado
Boulder, Boulder, CO, USA  Yuhong Fan  High Altitude Observatory, National
Center for Atmospheric Research, Boulder, CO, USA  George H. Fisher  Space
Sciences Laboratory, University of California Berkeley, Berkeley, CA, USA
Benoit Tremblay  Laboratory for Atmospheric and Space Physics, University of
Colorado Boulder
1234 Innovation Drive, Boulder, CO 80303, USA  National Solar Observatory,
University of Colorado Boulder, Boulder, CO, USA  Department of Astrophysical
and Planetary Sciences, University of Colorado Boulder, Boulder, CO, USA

###### Abstract

Knowledge of electric fields in the photosphere is required to calculate the
electromagnetic energy flux through the photosphere and set up boundary
conditions for data-driven magnetohydrodynamic (MHD) simulations of solar
eruptions. Recently, the method for inversions of electric fields from a
sequence of vector magnetograms and Doppler velocity measurements was improved
to incorporate spherical geometry and a staggered-grid description of
variables. The method was previously validated using synthetic data from
anelastic MHD (ANMHD) simulations. In this paper, we further validate the
method, using approximately one–hour long MHD simulation data of magnetic flux
emergence from the upper convection zone into the solar atmosphere. We
reconstruct photospheric electric fields and calculate the Poynting flux, and
compare those to the actual values from the simulations. We find that the
accuracy of the reconstruction is quite good during the emergence phase of the
simulated ephemeral active region evolution and decreases during the shearing
phase. Analysing our results, we conclude that the more complex nature of the
evolution (compared to the previously studied ANMHD case) that includes the
shearing evolution phase is responsible for the obtained accuracy decrease.

Solar magnetic flux emergence — Solar active region magnetic fields — Solar
active region velocity fields — Solar active regions — Solar photosphere —
Solar atmosphere

††journal: ApJ††software: (Fisher et al., 2020), MFE (Fan, 2017), FLCT (Fisher
& Welsch, 2008)

##  1 Introduction

Although it is the magnetic field and magnetic energy stored in various
topological structures that are believed to cause all the phenomena of solar
transient activity, the information on both magnetic and electric fields is
required to estimate electromagnetic energy fluxes through the solar
atmosphere. Those fluxes provide valuable information on the energy
transported from the solar interior to the atmosphere, contributing to the
understanding of magnetic heating of the higher solar atmosphere (e.g. Welsch
2015) as well as to the understanding of physical processes in solar eruptions
and flares. In particular, electric fields at the photosphere are used to set
up boundary conditions for data-driven magnetohydrodynamic (MHD) simulations
of the emergence of a new magnetic flux, formation of flux ropes, and
subsequent eruptions (e.g. Fan & Gibson, 2004; Fan, 2010; Toriumi et al.,
2020; Hoeksema et al., 2020, see Cheung & Isobe, 2014 for a review).

Because direct observations of electric fields based on the Stark effect are
quite difficult to achieve due to the low sensitivity of measurements (e.g.
Foukal & Behr, 1995), much effort was put into development of the techniques
to estimate photospheric electric fields using other measured quantities.
Several methods (Kusano et al., 2002; Fisher & Welsch, 2008; Schuck, 2008, see
Welsch et al., 2007 for comparison) estimated the plasma velocities in the
photosphere to obtain the electric field from the ideal Ohm’s law
($\mathbf{E}\propto\mathbf{V}\times\mathbf{B}$). These methods generally track
changes in the magnetic field from a sequence of magnetograms, although some
of them also used the vertical component of Faraday’s law
($\partial\mathbf{B}/\partial t\propto\nabla\times\mathbf{E}$) as an
additional constraint. The weakness of such approaches lay in extracting the
information about apparent (or optical) flows, which could be affected by the
temporal evolution of magnetic fluxes, e.g. their emergence or submergence,
and therefore did not correspond to real plasma flows in the photosphere. That
affected the quality of the reconstructed electric field, leaving the problem
unresolved. The recent development of techniques using neural networks,
however, promises a significant improvement in plasma velocity estimating
methods, which could provide accurate estimates of the photospheric electric
fields. In particular, neural networks such as DeepVel (Asensio Ramos et al.,
2017; Tremblay et al., 2018) and DeepVelU (Tremblay & Attie, 2020) reconstruct
plasma velocities from surface observations by emulating flows in realistic
MHD simulations of the photosphere. Note that flows estimated using DeepVel
and DeepVelU are neither optical nor physical as the neural networks currently
do not solve an optical flow equation or physical equations.

A different approach to estimate electric fields was proposed by Fisher et al.
2010. The idea was to uncurl Faraday’s law, using the poloidal-toroidal
decomposition (PTD) for the magnetic field time derivative. The solution of
the obtained set of equations for the electric field was not unique, however,
and the gradient of any scalar function could be added into the inductive
electric field, with the resulting field still satisfying the induction
condition. In particular, the pure inductive electric field was not
perpendicular to the magnetic field, as required by the ideal Ohm’s law. That
fact motivated Fisher et al. 2012 to include the additional observational
information about the line-of-sight (LOS) component of the plasma velocity
from photospheric Dopplergrams as well as about horizontal optical flows
obtained using the Fourier Local Correlation Tracking (FLCT) technique (Fisher
& Welsch, 2008) on a sequence of LOS magnetograms. Those data provided
contributions into the potential component of the electric field. Yeates
(2017) and Pomoell et al. (2019) also pointed out the importance of the non-
inductive electric fields. The perpendicularity between the resulting electric
field and magnetic field in the method was finally achieved with an iterative
technique eliminating the electric field component parallel to the magnetic
field vector. On the other hand, the iterative technique appeared to provide
the correct use of optical horizontal flows as an available proxy, instead of
real horizontal plasma velocities. In addition, the iterative technique
allowed other missing observational information to be taken into account
phenomenologically.

The new method was later referred to as the PDFI method for the techniques and
data used for electric field estimates: PTD – Doppler – FLCT – Ideal, where
“Ideal” stood for the iterative technique to get perpendicularity of the
electric and magnetic fields, coming from Ohm’s law for ideal plasmas.
Kazachenko et al. (2014) significantly expanded the scope of the PDFI method
by including non-normal viewing angles for observations, improved performance
with the use of the FISHPACK library, and performed a test to validate the
method. More specifically, the results of the electric field inversions based
on synthetic magnetograms from anelastic magnetohydrodynamic (ANMHD)
simulations (Abbett et al., 2004) were compared to the actual ANMHD electric
fields. Kazachenko et al. (2014) found the best values for free parameters in
the PDFI method and applied it to real solar data to study the evolution of
the energy and magnetic helicity fluxes in a solar active region (Kazachenko
et al., 2015). Fisher et al. (2020) further improved the PDFI method by
including spherical geometry as well as switching to a more conservative
staggered-grid description of variables, hence the new method name became .

So far, the ANMHD simulation data were the only synthetic test data used for
validation of the software. However, as mentioned in Kazachenko et al. (2014),
the vertical emergence of the bipole magnetic flux was probably exaggerated in
the synthetic data from the ANMHD simulations, because the simulations are
limited to modelling the rise of the flux tube in the convection zone without
modelling the emergence process into the solar atmosphere. Both observations
(e.g. Liu & Schuck, 2012) and MHD simulations of flux emergence into the solar
atmosphere (e.g. Manchester et al., 2004; Fan, 2009) have shown that
photospheric shearing motions in the emerging region contribute significantly
to the free-energy and helicity accumulation in the corona, which are largely
absent in the synthetic test data from the ANMHD simulations. In addition, the
ANMHD data were smooth in time and space, thus not revealing small-scale
structures or small-scale flow patterns that could potentially affect the
reconstructed electric field in a significant way. Recently, Lumme et al.
(2019) performed an analysis of the dependence of the reconstructed energy and
magnetic helicity fluxes on the temporal cadence of the input data, using
observations of a solar active region, and found that the results were stable
for a range of the sampled dataset cadences from several minutes up to two
hours.

In this paper, we test the electric field inversion method, using synthetic
photospheric data obtained from an MHD simulation of a twisted magnetic flux
tube emerging from the top layer of the convection zone into the atmosphere
and the corona, similar to the simulation in Fan (2009), where photospheric
shearing and rotational motions of the polarities develop to drive the build-
up of a corona flux rope. We apply the method to reconstruct the photospheric
electric field from a sequence of the synthetic photospheric magnetograms and
vertical velocity data extracted from the simulation, and compare the results
to the actual electric field in the MHD simulation.

The paper is organised as follows. In Section 2, we first give a model
description of the flux emergence simulation, and then analyse the synthetic
photospheric data extracted from the simulation in Section 3. In Section 4, we
perform FLCT reconstructions of the horizontal optical flows from simulated
magnetograms of the vertical component of the photospheric magnetic field. In
Section 5, we show reconstructions of the photospheric electric field and
discuss our results. Section 6 presents our conclusions.

##  2 Magnetic flux emergence simulation

Here we give a description of the interior-to-corona simulation of the
emergence of a twisted magnetic flux tube, of which the photospheric magnetic
and velocity field data are used in this work for the validation of the
inversion method. The detailed results of the flux emergence simulation will
be presented in a separate paper.

The simulation uses the “Magnetic Flux Eruption” (MFE) code, which solves the
following semi-relativistic magnetohydrodynamic equations (Fan, 2017):

$\frac{\partial\rho}{\partial t}=-\nabla\cdot(\rho{\bf v}),$ (1)
$\displaystyle\frac{\partial({\rho\bf v})}{\partial t}$ $\displaystyle=$
$\displaystyle-\nabla\cdot\left(\rho{\bf v}{\bf v}\right)-\nabla p+\rho{\bf
g}+\frac{1}{4\pi}(\nabla\times{\bf B})\times{\bf B}$
$\displaystyle+\frac{{v_{A}}^{2}/c^{2}}{1+{v_{A}}^{2}/c^{2}}\left[{\cal
I}-{\hat{\bf b}}{\hat{\bf b}}\right]$ $\displaystyle\cdot\left[(\rho{\bf
v}\cdot\nabla){\bf v}+\nabla p-\rho{\bf g}-\frac{1}{4\pi}(\nabla\times{\bf
B})\times{\bf B}\right],$ $\frac{\partial{\bf B}}{\partial
t}=\nabla\times({\bf v}\times{\bf B}),$ (3) $\nabla\cdot{\bf B}=0,$ (4)
$\frac{\partial e}{\partial t}=-\nabla\cdot\left({\bf
v}e\right)-p\nabla\cdot{\bf v}-\nabla\cdot{\bf q}_{s}+Q_{\rm rad}+H_{\rm
em}+H_{\rm num},$ (5) $p=\frac{\rho RT}{\mu},$ (6)

where

$e={\frac{p}{\gamma-1}},$ (7) $v_{A}=\frac{B}{\sqrt{4\pi\rho}}.$ (8)

Here the equations are solved in Cartesian geometry. In the above equations,
symbols have their usual meanings, where ${\bf v}$ is the velocity field,
${\bf B}$ is the magnetic field and $B$ is its magnitude, $\rho$, $p$, and $T$
are respectively the plasma density, pressure and temperature, $e$ is the
internal energy density, $c$ is the (reduced) speed of light, ${\cal I}$ is
the unit tensor, ${\hat{\bf b}}={\bf B}/B$ is the unit vector in the magnetic
field direction, ${\bf g}=-g{\hat{\bf z}}$ is the constant gravitational
acceleration, $R$, $\mu$, and $\gamma$, are respectively the gas constant, the
mean molecular weight, and the adiabatic index of the perfect gas. We assume
the adiabatic index ${\gamma}=5/3$. The internal energy equation (Eq. 5) takes
into account the following non-adiabatic effects: i) the field-aligned heat
conduction flux, ${\bf q}_{s}$, which we compute using the hyperbolic heat
conduction approach (Rempel, 2017); ii) the optically thin radiative cooling,
$Q_{\rm rad}$:

$Q_{\rm rad}=N^{2}\Lambda(T),$ (9)

where $N=\rho/m_{p}$ is the proton number density assuming a fully ionised
hydrogen gas, with $m_{p}$ being the proton mass, and the radiative loss
function $\Lambda(T)$ is as given in Athay (1986), with the cooling suppressed
for $T$ below $2\times 10^{4}$ K; iii) an empirical coronal heating, $H_{\rm
em}$, that is a function of height $z$:

$H_{\rm em}=\frac{F}{L_{H}}\exp\left[-(z-z_{0})/L_{H}\right],$ (10)

which is non-zero for $z>z_{0}$, with $F=5\times 10^{5}\;{\rm
erg\;s^{-1}cm^{-2}}$, $L_{H}=11\;{\rm Mm}$, and $z_{0}=3\;{\rm Mm}$ above the
photosphere. Finally, we also evaluate the effective heating resulting from
numerical diffusion of velocity and magnetic fields in the momentum and
induction equations, and added it to the internal energy in Equation (5) as
$H_{\rm num}$.

The setup of the simulation is as follows. We use a Cartesian domain as shown
in Figure 1 with the horizontal sizes: $x\in[-26.5,26.5]\;{\rm Mm}$,
$y\in[-26.5,26.5]\;{\rm Mm}$, and the vertical size: $z\in[-6.1,49.8]\;{\rm
Mm}$, where $z=0$ corresponds to the photosphere (the grey horizontal plane in
Figure 1).

Figure 1: The initial setup of the simulation domain. The grey plane is the
photosphere. The colour image in the vertical cross-section shows the
temperature of the initial background equilibrium plane-parallel atmosphere
permeated with a uniform open magnetic field represented by the field lines in
the cross-section. A twisted magnetic flux tube (as represented by the twisted
field lines) oriented along the $x$-axis is inserted into the interior layer
below the photosphere. See text for detailed description.

The domain is resolved with a non-uniform grid of $480(x)\times 480(y)\times
512(z)$ cells, with a peak vertical resolution of $\Delta z=31$ km at the
photosphere and a peak horizontal resolution of $\Delta x=\Delta y=46$ km
centred on the emerging region. We first initialise a plane-parallel,
vertically stratified atmosphere, as shown in Figure 2, which consists of an
adiabatically stratified polytropic layer representing the top of the
convection zone, an approximately isothermal layer that represents the
photosphere and chromosphere, and a high temperature isothermal layer
representing the corona, permeated with a uniform slanted open magnetic field
of 6 G (illustrated by the field lines shown in the vertical cross-section in
Figure 1).

Figure 2: The temperature (black curve) and density (blue curve) profiles of
the initial equilibrium plane-parallel atmosphere.

This initial one-dimensional atmosphere with a uniform open magnetic field has
been numerically relaxed to an approximate hydrostatic and thermal equilibrium
along the magnetic field, with the transition region and hot corona maintained
by the empirical coronal heating (Eq. 10), balanced by the optically thin
radiative cooling and the field-aligned heat conduction. A horizontal twisted
magnetic flux tube is inserted into the convection zone layer below the
photosphere (see Figure 1), whose magnetic field is given as follows:

${\bf B}=B_{x}(r)\,\hat{\bf x}+B_{\theta}(r)\,\hat{\bf\theta},$ (11)

where

$B_{x}(r)=B_{0}\exp(-r^{2}/a^{2}),$ (12) $B_{\theta}(r)=qrB_{x}(r),$ (13)

$\hat{\bf x}$ is the tube axial direction, $\hat{\bf{\theta}}$ denotes the
azimuthal direction in the tube cross-section, $r$ denotes the radial distance
to the tube axis which is located at $(y,z)=(0,-1.69)\;{\rm Mm}$, $a=0.38$ Mm
is the radius of the tube, $q$ is the rate of the field-line twist and is set
to be $-a^{-1}$, and $B_{0}=5677$ G is the field strength at the initial axis
of the flux tube. The total flux of the subsurface flux tube is $2.63\times
10^{19}$ Mx, much smaller than the flux of typical active regions, but
corresponds to the flux scale of ephemeral regions (e.g. Hagenaar et al.,
2008). Inside the flux tube, the gas pressure differs from that of the initial
background plane parallel atmosphere $p_{0}(z)$ by:

$p_{1}(r)=-\frac{B_{0}^{2}}{2}\exp(-2r^{2}/a^{2})\left[1-\frac{1}{2}q^{2}a^{2}\left(1-\frac{2r^{2}}{a^{2}}\right)\right],$
(14)

such that the pressure gradient from $p_{1}$ balances the Lorentz force of the
twisted flux tube. Furthermore, the middle segment of the flux tube is made
buoyant by introducing a density change $\rho_{1}$ within the flux tube
relative to the density $\rho_{0}(z)$ of the background plane parallel
atmosphere:

$\rho_{1}=-\rho_{0}(z)\frac{B_{x}^{2}(r)}{2p_{0}(z)}\left[(1+\epsilon)\exp(-x^{2}/\lambda^{2})-\epsilon\right]$
(15)

where the length scale of the buoyant segment $\lambda=1.54$ Mm, and
$\epsilon=0.1$.

Given the initial state described above, the subsurface flux tube develops an
$\Omega$-shape, with its central buoyant segment rising towards the
photosphere. When the apex of the flux tube reaches the stably stratified
photosphere, the magnetic buoyancy instability sets in, which causes the upper
boundary of the flux tube to expand into the atmosphere and begins the flux
emergence process, similar to what has been found in many previous simulations
(see e.g. Manchester et al., 2004; Fan, 2009; Archontis & Hood, 2013). In this
case, as a result of the flux emergence, a twisted coronal flux rope is built
up in the corona as shown in Figure 3.

Figure 3: 3D magnetic field at $t=45$ min from the flux emergence simulation,
more specifically during the shearing phase. The field lines are coloured by
temperature. The horizontal cross-section corresponds to the photosphere with
grey-scale image showing the vertical magnetic field.

The flux rope reconnects with the ambient open field, leading to the
development of a rotating “blow-out” jet with untwisting field lines in the
jet column (see Figure 3). The detailed result of the coronal magnetic field
evolution resulting from the flux emergence will be described in a separate
paper (Manek et al. in preparation). It is found that the formation of the
coronal flux rope is not due to the bodily emergence of the subsurface flux
rope, but is mainly due to the shearing and rotational motions of the polarity
concentrations at the photosphere during the flux emergence, which transport
twist from the subsurface flux tube into the corona (e.g. Fan, 2009). In the
following section, we extract the time sequence of the photospheric vector
magnetic field and velocity field evolution from the flux emergence simulation
as synthetic data, and apply the method to the synthetic magnetograms to
examine how well it infers the electric field at the photosphere.

Figure 4: Slices of the solar surface in the MFE simulation representing the
vertical component of the magnetic field and the three components of the
plasma velocity at three different times: 14.99 min (left column), 28.01 min
(middle), and 40.24 min (right). The black contours are drawn for values of
$\pm 200$ G of the vertical component of the magnetic field.

Figure 5: Analysis of the MFE simulated emergence of a new magnetic flux into
the solar atmosphere. The left panel shows the evolution of the magnetic flux
through the computational domain and its change rate. The solid lines
represent the unsigned magnetic flux (black), signed (green), positive (blue),
and absolute value of the negative magnetic flux (red). The dashed lines
represent the change rates for the unsigned (black), and positive and negative
(blue) magnetic fluxes. The right panel shows the evolution of the total
energy flux $F$ through the computational domain and the injected energy $W$.
The solid lines represent the shear component (blue), emergence component
(red), and full energy flux (black). The purple and green solid lines
demonstrate the upward and downward contributions to the full energy flux,
respectively. The dashed lines represent the amount of the injected energy as
calculated from the full energy flux (black), and the emergence (red) and
shear (blue) components of the energy flux. The axes for dashed lines are
located on the right side of each panel. The shaded areas in both panels
represent the emergence phase (10–22 min) of the simulated active region’s
evolution. The second dashed vertical line of each plot marks the beginning of
the shearing phase ($t>22$ min).  Figure 6: Comparison between the horizontal
$x$\- and $y$-components (top and bottom rows, respectively) of the plasma
velocity reconstructed using FLCT (left column) and computed by the MFE
simulation (middle column) at $t=14.99$ min, i.e. during the emergence phase
of the simulation. The right panels show the scatter plots, linear fit (red
dashed line), Spearman correlation coefficient and slope of the fit, and
reference one-to-one line (black dashed).  Figure 7: The same as Figure 6 at
$t=28.01$ min, i.e. during the shearing phase of the simulation.

##  3 Analysis of the photospheric synthetic data

We use simulation data for the magnetic and velocity fields at the
photospheric level and the ideal Ohm’s law to calculate our reference electric
fields for the validation of . The dataset spans approximately one hour in
time. After a 10–min initial stage of the flux emergence, the two polarities
of the bipolar region start moving in the opposite directions. Figure 4 shows
the typical patterns for the vertical component of the photospheric magnetic
field as well as for the components of the plasma velocity at three different
times. The times have been selected to demonstrate different stages of the
emerging magnetic flux evolution. An important feature of the simulated data
is the shear structures seen in the 28 and 40 min snapshots for both
horizontal components of the plasma velocity in the regions corresponding to
the magnetic polarities. The vertical magnetic field does not reveal such
shearing structures, although showing some fine swirling dynamics. These shear
structures appear to be associated with rotational motions of the twisted
emerging flux, which are detected in solar white-light images as rotating
sunspots (e.g. Kazachenko et al., 2009). We also note that vertical (Doppler)
velocities decrease with time. The late stage of the simulated evolution
contains downflow regions in the middle of the active region, which are caused
by the draining and tilting motions of U-shaped, anchored field lines of the
partially emerged flux rope.

Figure 5 analyses the evolution of the flux and energy characteristics of the
MFE simulated active region. We calculate the magnetic flux and its change
rate through the photospheric cross-section of the MFE computational domain.
The photospheric magnetic flux is a combination of the emerging sub-surface
flux tube and a uniform slanted open magnetic field permeating the simulation
domain. The open field flux is of positive polarity on the photosphere. We
note that the (absolute) value of the negative polarity flux on the
photosphere (due to flux emergence) increases from zero to a peak value of
about $3.9\times 10^{19}$Mx (red solid curve in left panel of Figure 5),
significantly exceeding the total axial flux of the subsurface flux tube
($2.63\times 10^{19}$ Mx). This is because some of the twisted field lines of
the partially emerged flux tube thread through the photosphere with more than
one stitch, and thus the unsigned emerged flux of each polarity on the
photosphere can exceed the total axial flux of the subsurface twisted flux
tube.

In the left panel of Figure 5, we can distinguish two phases of the magnetic
flux evolution. The first one (around 10–22 min) is characterised by a rapid
change in the magnetic flux due to the flux emergence, which is also observed
in the magnetic flux change rate. After 22 min, the magnetic flux does not
increase so significantly. Throughout the simulation, the signed flux stays
constant because of the symmetry in the flux variations at both polarities,
with its non-zero value being due to the constant background magnetic field of
6 G.

Analysing the energy fluxes, we break the vertical component of the Poynting
flux, $S_{z}$, into two terms,

$S_{z}=\frac{c}{4\pi}\left[\mathbf{E}\times\mathbf{B}\right]_{z}=\frac{B_{h}^{2}}{4\pi}v_{z}-\frac{B_{z}}{4\pi}\left(\mathbf{B}_{h}\cdot\mathbf{v}_{h}\right),$
(16)

where $\mathbf{E}$ is the electric field vector, and the subscripts stand for
the horizontal (h) and vertical (z) vector components (e.g. Abbett & Fisher,
2012). The emergence term determined by the vertical component of the plasma
velocity corresponds to the first term in Equation (16), while the shearing
term is determined by the horizontal flows and corresponds to the second term
in the equation. The right panel in Figure 5 shows that the first, rapid phase
in the magnetic flux evolution corresponds to the deposit of energy into the
solar atmosphere predominantly due to the emergence component of the energy
flux (the red solid line), while the quiet phase corresponds to the shearing
evolution of the active region from the viewpoint of the deposited energy flux
(the blue solid line). Therefore, throughout the paper, we will refer to the
first phase of the simulated active region evolution (10–22 min) as the
emergence phase, while the subsequent evolution of the active region (after 22
min) will be referred to as the shearing phase (see Figure 5). The shearing
phase starts once the contribution to the vertical electromagnetic energy flux
from the horizontal (shearing) motions dominates that from the vertical plasma
motions. Physically, this is due to the decrease of vertical plasma flows
during the course of simulated active region evolution because of the
completion of the flux emergence process. The analysis of the upward and
downward contributions into the full vertical energy flux demonstrates the
dominant energy inflow to the atmosphere from below. We also note that
cumulative contributions to the injected energy from the emergence and
shearing energy flux components in the simulation are equal to each other.

Figure 8: Evolution of the Spearman correlation coefficient (blue) and slope
parameter (yellow) of the horizontal velocity FLCT reconstructions for the
entire duration of the simulation. Upper (lower) panel shows results for the
$x$($y$)-component of the flows. The shaded areas represent the emergence
phase (10–22 min) of the simulated active region evolution. The second dashed
vertical line of each plot marks the beginning of the shearing phase ($t>22$
min) The horizontal dashed lines mark the limits for the strong, moderate, and
weak correlation.

##  4 FLCT inversions of horizontal flows

To test the software, we first apply the FLCT technique (Fisher & Welsch,
2008) to the MFE $B_{z}$ magnetograms to obtain estimates of horizontal plasma
flows in the photosphere. For each timestep, we use a pair of magnetograms
that are one timestep behind and ahead of the given time.

For the FLCT technique, we have to specify free parameters. To find the best
value for the windowing parameter $\sigma$, we performed several FLCT
reconstructions of horizontal velocities with $\sigma=$ 5, 10, 15, 20, 25, 50
pixels, and selected the optimal value $\sigma=20$ pixels by computing the
correlation between the reconstructed and reference MFE velocities (cf.
$\sigma=5$ in Fisher et al., 2012, $\sigma=15$ in Kazachenko et al., 2014). We
chose the threshold value of 14 G for the magnitude of the minimum magnetic
field, which is about 1% of the maximum value of $B_{z}$ during the entire
simulation. A threshold value is often set up when working with observational
data to eliminate noisy data in the weak field regions. Here we introduce the
threshold to reduce the computational time and also for consistency with the
previous studies. Moreover, magnetic field masks are used by the software when
computing the Doppler contribution to the electric field (Fisher et al.,
2020). We performed several FLCT reconstructions with other threshold values
(even one order of magnitude larger) and have not found a significant effect
on the results.

Figure 6 and Figure 7 compare the FLCT reconstructions of horizontal optical
flows to the MFE simulation flows at two different times. We use the Spearman
correlation coefficient, $c_{S}$, and the slope angle of the linear fit as
statistical metrics, both determined from the scatter plots, and find that the
quality of inversions is worse for the $y$-component of the plasma velocity
than for the $x$-component for both times and that quality of both components
is worse at $t=28.01$ min than at $t=14.99$ min. These discrepancies are most
likely due to the more complex evolution of the magnetic flux emergence, in
particular, due to the fact that the shear (rotational) structures are not
detectable in the simulated magnetograms (see Section 3 and Figure 4) and thus
cannot be tracked by the FLCT algorithm in both directions. Indeed, the FLCT
technique tracks changes in the structure of the $B_{z}$-component and
provides information on the optical horizontal flows. In the case of study, we
have two kinds of motions in the photospheric level: i) translational motions
of two regions of the opposite magnetic polarities, and ii) their rotational
motions. However, unlike the $y$-component of the horizontal flow, its
$x$-component has a significant contribution from the translational motions of
the polarities, which is captured by the FLCT technique. That is why the
reconstruction of the $x$-component of horizontal photospheric flows is better
during the entire course of the simulation.

In Figure 8, we show the evolution of the correlation coefficient between the
reconstructed and simulation velocities and slope parameter over the total
duration of the simulation. Throughout the paper, we use the following
correlation scores: strong correlation ($c_{S}>0.8$), moderate
($0.5<c_{S}<0.8$), and weak ($0.2<c_{S}<0.5$). The FLCT reconstruction of the
$y$-component of horizontal flows shows strong correlation with MFE velocities
during the emergence phase of the simulated magnetic flux evolution (see
Section 3), and a moderate-to-weak correlation during the shearing phase which
is consistent with the results from Figure 7 taken during this phase. The
metrics for the reconstruction of the $x$-component of horizontal flows
decrease when transitioning from the emergence phase to the shearing phase,
but per our criteria the correlation remains strong throughout the entire
simulation duration. We should also keep in mind that the FLCT technique
tracks apparent flows rather than physical plasma velocities, which could
justify the lower values obtained for the slope parameter.

In order to understand better the influence of rotational plasma motions, we
plot the temporal evolution of the absolute vertical vorticity, namely the
$z$-component of the curl of the photospheric flows, averaged throughout the
photospheric regions with the applied vertical magnetic field mask. Figure 9
shows the results for different mask values. We use the mask to exclude the
contributions to the average value of the vorticity from the periphery,
however, one should be careful with that because the number of contributing
pixels is different at each time, and therefore the result can be highly
affected by the emerging nature of the magnetic flux evolution. In particular,
we believe that this is the cause for the significant initial drop of the
average vertical vorticity seen in Figure 9 before $t=15$ min. The snapshots
of the vertical vorticity show that at the early stages, there are pixels with
high vorticity values, but their total masked area (used for averaging) is
still small, leading to the high values seen in the vorticity plot. As time is
progressing and flux emerges, the area of masked pixels quickly increases,
giving the obtained vorticity drop by $t=15$ min. Also, the analysis of the
vertical magnetic field evolution shows that during the time interval before
$t=15$ min, the magnetic polarities experience significant expansion in all
directions without rotations. This can explain very good correlations for both
horizontal flow components during that interval. After $t=15$ min, the
vorticity grows, or in other words, rotational motions become more pronounced,
which affects the FLCT reconstructions. As seen from the plots for the 100 G
and 200 G masks, the vorticity saturates by $t=30$ min, that is why the
reconstruction metrics stop decreasing, staying almost constant. After around
45 min (for the blue curve, and earlier for the yellow one), the vorticity
decreases, so we can see an increase in the correlation coefficient both for
$V_{x}$ and $V_{y}$ components. Concluding this Section, we however stress
that this analysis of the dependence of the FLCT results on the vertical
vorticity may not be considered as a strict explanation, because the vorticity
plots presented in Figure 9 depend on the used strong-field mask.

Figure 9: Evolution of the average absolute vertical vorticity calculated
using the MFE simulation horizontal velocities for different magnetic field
mask values: 100 G (blue solid line), 200 G (yellow dashed), 50 G (green dash-
dotted). The shaded area represents the emergence phase (10–22 min) of the
simulated active region evolution. The second dashed vertical line marks the
beginning of the shearing phase ($t>22$ min).

##  5 electric field inversions

In this section, we use the method to reconstruct electric fields in the
photosphere from the MFE simulated magnetograms and Dopplergrams, and FLCT
estimates of the horizontal flows. Figures 10, 11, and 12 show the obtained
electric fields and vertical component of the Poynting flux at three different
times around 15, 28, and 40 min as well as the scatter plots comparing the
reconstructions to actual MFE simulated values. We calculate the Spearman
correlation coefficient and slope angle of the linear fit from the scatter
plots. We find good agreement of the spatial structure of the inverted
electric fields with the simulation data, although at the snapshots
corresponding to the shearing phase ($t=28$ and 40 min snapshots), the
$x$-component of the reconstructed electric field is highly underestimated.

Figure 13 shows the evolution of the correlation coefficient and slope
parameter for the three components of the electric field and the vertical
component of the Poynting flux over the entire duration of the simulation.
During the emergence phase (about 10–22 min, see Section 3), the results have
strong-to-moderate correlation (0.5–1.0) with the MFE electric fields. The
slope parameter is within 40% difference with the ideal reconstruction unity
value. During the shearing phase, the results have moderate-to-weak
correlation (0.5–0.8), the slope is within 50% difference for the $y$\- and
$z$-components, whereas the $x$-component of the electric field is
significantly underestimated. We find low values for the correlation
coefficient for the vertical component of the electric field at about 40–50
min. While the slope parameter is quite high (around 0.6–0.7), and the plotted
maps show all the spatial structures seen in the actual MHD data (see e.g.
Figure 12), that correlation is unexpectedly low.

The reconstruction of the vertical Poynting flux is good during the emergence
phase, with strong-to-moderate correlation and a slope within 20% of the ideal
reconstruction. During the shearing phase, the method still adequately
reconstructs the vertical Poynting flux with moderate correlation except for a
short time interval close to 30 min, and a slope within 60% of the ideal unity
value. We note that despite the significant underestimation and very low value
of the slope for $E_{x}$, the vertical Poynting flux reconstruction quality is
good. In addition, we calculate the total vertical energy flux, integrating
spatially the PDFI_SS-inverted vertical Poynting flux, as shown in Figure 14.
The reconstructed total energy flux shows good agreement (within 20%) with the
MFE values.

The poor value of the slope for the $x$-component of the electric field
requires some additional attention. The discrepancy could result either from
the structures in $V_{y}$ which are not detected in $B_{z}$ maps (see Section
4) or from the more sophisticated nature of the active region evolution,
including the shearing phase, which was not pronounced in the ANMHD test. To
shed some light, we have run the software using actual MHD horizontal plasma
velocities instead of the FLCT estimates. The obtained reconstruction (Figure
15) shows significantly better slope values for $E_{x},E_{y}$, and $S_{z}$,
but is still far from the close-to-ideal reconstruction (cf. ANMHD case,
Kazachenko et al., 2014). We therefore conclude that the poor FLCT
reconstruction of the horizontal plasma flows due to the $V_{y}$ shear
structures is not the only source of the discrepancy, although it does affect
the electric field inversions. The shearing nature of the later stage of the
simulated active region evolution comes into play. We speculate that during
the shearing phase some discrepancy is due to the fact that the method misses
the non-inductive (i.e. potential) contribution into the vertical component of
the electric field. Another improvement would be to include additional
observational information capturing the rotation of magnetic polarities in
order to reconstruct the shear velocity structures properly. Indeed, the
simplified simulation of the flux emergence without taking into account the
effect of convection results in a very smooth shape of polarity flux
concentrations. By using other proxies in the FLCT technique or/and other
velocity inference methods, we could reconstruct better the shear velocity
structures in the magnetic polarity regions, and therefore improve the quality
of the reconstructions of the electric field (cf. Figure 15).

Figure 10: Components of the electric field and vertical component of the
Poynting flux obtained with method (left column) and the MFE simulation
(middle column) at $t=15.2$ min, i.e. during the emergence phase of the
simulation. The right column shows the scatter plots with the linear fit (red
dashed line), Spearman correlation coefficient and slope of the fit, and
reference one-to-one line (black dashed).  Figure 11: The same as Figure 10 at
$t=28.2$ min, i.e. during the early stages of the shearing phase of the
simulation.  Figure 12: The same as Figure 10 at $t=40.4$ min, i.e. during the
later stages of the shearing phase of the simulation.  Figure 13: Quality of
the reconstruction for the entire duration of the simulation: Evolution of the
Spearman correlation coefficient (blue) and slope parameter (yellow). The
shaded area represents the emergence phase (10–22 min) of the simulated active
region evolution. The second dashed vertical line of each plot marks the
beginning of the shearing phase ($t>22$ min). The horizontal dashed lines mark
the limits for the strong, moderate, and weak correlation.  Figure 14: Total
vertical energy flux as computed by the MFE simulation (blue line) and
reconstructed using (red line). The shaded area represents the emergence phase
(10–22 min) of the simulated active region evolution. The second dashed
vertical line marks the beginning of the shearing phase ($t>22$ min).  Figure
15: Evolution of the Spearman correlation coefficient (blue) and slope
parameter (yellow) of the inversions with actual MHD horizontal velocities
instead of FLCT ones. The shaded area represents the emergence phase (10–22
min) of the simulated active region evolution. The second dashed vertical line
marks the beginning of the shearing phase ($t>22$ min). The horizontal dashed
lines mark the limits for the strong, moderate, and weak correlation.

##  6 Conclusion

In this study, we performed a test to validate the electric field inversion
method (Fisher et al., 2020), using the photospheric data extracted from an
interior-to-corona MHD simulation of the emergence of a twisted flux tube.
Unlike the previous validation tests for the method that used just one
snapshot of a rising flux tube in the solar interior (see, e.g. Kazachenko et
al., 2014), in this study we considered the more complex time sequence of
evolution of the emerging flux region at the photosphere. The evolution of the
simulated ephemeral active region consisted of two phases during which the
vertical energy flux was dominated either by vertical plasma motions (the
emergence phase) or by horizontal motions (the shearing phase). The initial
rapid growth of the magnetic flux through the simulated photosphere
corresponded to the emergence phase.

We used synthetic magnetograms and Dopplergrams at the photospheric level to
calculate the inductive part of the photospheric electric field as well as its
non-inductive (potential) parts due to the contributions from vertical and
horizontal plasma flows, and ideal Ohm’s law. The horizontal optical flows
were estimated with the FLCT method (Fisher & Welsch, 2008) tracking the
evolution of the vertical component of the photospheric magnetic field.

Our main findings are the following:

  * •

During the emergence phase of the simulation, the reconstruction of the
photospheric electric field is very good – the correlation between the
reconstructed and MFE simulation fields varies from strong to moderate. During
the shearing phase, the accuracy of the method decreases – the correlation
between the reconstructed and simulation fields is moderate-to-weak, with the
fields, especially the $E_{x}$ component, being underestimated.

  * •

The method reconstructs the vertical component of the Poynting flux well. The
correlation between the reconstructed and MFE simulation energy fluxes varies
mostly from moderate to strong. The spatially integrated total energy flux is
within 20% of the simulation flux throughout both the emergence and shearing
phases of the active region evolution. We note that reliable inversions of the
energy flux is of special importance for data-driven MHD simulations of solar
eruptions.

  * •

Difficulties in reconstructing the $y$-component of the horizontal velocity
using the FLCT method, lead to a decrease in accuracy for the electric field
reconstruction during the shearing phase of the AR evolution. We conclude that
these difficulties are caused by the lack of spatial structure in symmetric
rotating sunspots as seen in $B_{z}$ magnetograms used as an input for the
FLCT, while the $x$-component of the horizontal velocity is reconstructed
better due to strong contributions from the translational motions of the
magnetic polarities at the opposite directions. We speculate that the electric
field reconstruction accuracy is also affected by the lack of the non-
inductive contribution to the vertical component of the electric field.

  * •

To test whether the method can adequately describe the simulated MHD electric
fields, we have performed the reconstructions using the MFE simulation
horizontal plasma velocities (instead of the FLCT velocities). In this case,
the method yields more accurate electric fields, especially for $E_{x}$.
However, the reconstructed electric fields still differ from MFE data. We
believe that the complex nature of the evolution of the simulated active
region that includes the shearing evolution phase, as compared to the ANMHD
case, does not allow the method to provide the highest quality reconstruction
of the photospheric electric field.

The flux emergence simulation used for this validation is still very idealised
in that it does not include the effect of convection so that the emerging flux
pattern consists of very smooth polarity flux concentrations. This makes the
detection of the shear and rotational motions more difficult than reality and
thus may be an underestimate of the sensitivity of the method. We also note
that the reliable detection of the velocities of rotating structures is a
fundamental problem affecting all inversion and tracking methods that depend
on observed time variability. For instance, a rotating sunspot that has
perfect symmetry about the axis will show no evolution and hence zero electric
and velocity field even if the actual field is substantial. In this case,
however, real observations offer an advantage. Many examples of rotating
sunspots show structures in white light and sometimes in magnetic fields,
which clearly rotate, and correlation tracking methods could capture these
motions (Kazachenko et al., 2009).

To summarise, the accuracy of the reconstructions of the electric fields
decreases for the shearing evolution of the active region, however, the
spatial structure of electric fields is captured correctly, correlations are
still decent, and the energy flux is reconstructed very well. As mentioned
previously, the information on the photospheric electric fields is of special
importance for cutting-edge data-driven MHD simulations. This field of physics
of the solar atmosphere still contains a lot of uncertainties (see, for
instance, the recent comparison study by Toriumi et al. (2020) ) and needs yet
to formulate its general methodology. In this light, the information delivered
by the method is of crucial importance and undoubtedly will be used in data-
driven MHD simulations of solar eruptive processes.

We acknowledge support of NASA ECIP NNH18ZDA001N and NASA LWS NNH17ZDA001N
(A.N.A. and M.D.K.), NASA LWS 80NSSC19K0070 (Y.F.), DKIST Ambassador Program
(A.N.A.), and Program of Basic Research No. II.16.1.6 (A.N.A.).

## References

  * Abbett & Fisher (2012) Abbett, W. P., & Fisher, G. H. 2012, Sol. Phys., 277, 3, doi: 10.1007/s11207-011-9817-3
  * Abbett et al. (2004) Abbett, W. P., Mikić, Z., Linker, J. A., et al. 2004, Journal of Atmospheric and Solar-Terrestrial Physics, 66, 1257, doi: 10.1016/j.jastp.2004.03.016
  * Archontis & Hood (2013) Archontis, V., & Hood, A. W. 2013, ApJ, 769, L21, doi: 10.1088/2041-8205/769/2/L21
  * Asensio Ramos et al. (2017) Asensio Ramos, A., Requerey, I. S., & Vitas, N. 2017, A&A, 604, A11, doi: 10.1051/0004-6361/201730783
  * Athay (1986) Athay, R. G. 1986, ApJ, 308, 975, doi: 10.1086/164565
  * Cheung & Isobe (2014) Cheung, M. C. M., & Isobe, H. 2014, Living Reviews in Solar Physics, 11, 3, doi: 10.12942/lrsp-2014-3
  * Fan (2009) Fan, Y. 2009, ApJ, 697, 1529, doi: 10.1088/0004-637X/697/2/1529
  * Fan (2010) Fan, Y. 2010, The Astrophysical Journal, 719, 728, doi: 10.1088/0004-637x/719/1/728
  * Fan (2017) Fan, Y. 2017, ApJ, 844, 26, doi: 10.3847/1538-4357/aa7a56
  * Fan & Gibson (2004) Fan, Y., & Gibson, S. E. 2004, ApJ, 609, 1123, doi: 10.1086/421238
  * Fisher et al. (2020) Fisher, G. H., Kazachenko, M. D., Welsch, B. T., et al. 2020, ApJS, 248, 2, doi: 10.3847/1538-4365/ab8303
  * Fisher & Welsch (2008) Fisher, G. H., & Welsch, B. T. 2008, in Astronomical Society of the Pacific Conference Series, Vol. 383, Subsurface and Atmospheric Influences on Solar Activity, ed. R. Howe, R. W. Komm, K. S. Balasubramaniam, & G. J. D. Petrie, 373.  https://arxiv.org/abs/0712.4289
  * Fisher et al. (2012) Fisher, G. H., Welsch, B. T., & Abbett, W. P. 2012, Sol. Phys., 277, 153, doi: 10.1007/s11207-011-9816-4
  * Fisher et al. (2010) Fisher, G. H., Welsch, B. T., Abbett, W. P., & Bercik, D. J. 2010, ApJ, 715, 242, doi: 10.1088/0004-637X/715/1/242
  * Foukal & Behr (1995) Foukal, P. V., & Behr, B. B. 1995, Sol. Phys., 156, 293, doi: 10.1007/BF00670229
  * Hagenaar et al. (2008) Hagenaar, H. J., DeRosa, M. L., & Schrijver, C. J. 2008, ApJ, 678, 541, doi: 10.1086/533497
  * Hoeksema et al. (2020) Hoeksema, J. T., Abbett, W. P., Bercik, D. J., et al. 2020, ApJS, 250, 28, doi: 10.3847/1538-4365/abb3fb
  * Kazachenko et al. (2009) Kazachenko, M. D., Canfield, R. C., Longcope, D. W., et al. 2009, ApJ, 704, 1146, doi: 10.1088/0004-637X/704/2/1146
  * Kazachenko et al. (2014) Kazachenko, M. D., Fisher, G. H., & Welsch, B. T. 2014, ApJ, 795, 17, doi: 10.1088/0004-637X/795/1/17
  * Kazachenko et al. (2015) Kazachenko, M. D., Fisher, G. H., Welsch, B. T., Liu, Y., & Sun, X. 2015, ApJ, 811, 16, doi: 10.1088/0004-637X/811/1/16
  * Kusano et al. (2002) Kusano, K., Maeshiro, T., Yokoyama, T., & Sakurai, T. 2002, ApJ, 577, 501, doi: 10.1086/342171
  * Liu & Schuck (2012) Liu, Y., & Schuck, P. W. 2012, ApJ, 761, 105, doi: 10.1088/0004-637X/761/2/105
  * Lumme et al. (2019) Lumme, E., Kazachenko, M. D., Fisher, G. H., et al. 2019, Sol. Phys., 294, 84, doi: 10.1007/s11207-019-1475-x
  * Manchester et al. (2004) Manchester, IV, W., Gombosi, T., DeZeeuw, D., & Fan, Y. 2004, ApJ, 610, 588, doi: 10.1086/421516
  * Pomoell et al. (2019) Pomoell, J., Lumme, E., & Kilpua, E. 2019, Sol. Phys., 294, 41, doi: 10.1007/s11207-019-1430-x
  * Rempel (2017) Rempel, M. 2017, ApJ, 834, 10, doi: 10.3847/1538-4357/834/1/10
  * Schuck (2008) Schuck, P. W. 2008, ApJ, 683, 1134, doi: 10.1086/589434
  * Toriumi et al. (2020) Toriumi, S., Takasao, S., Cheung, M. C. M., et al. 2020, The Astrophysical Journal, 890, 103, doi: 10.3847/1538-4357/ab6b1f
  * Tremblay & Attie (2020) Tremblay, B., & Attie, R. 2020, Frontiers in Astronomy and Space Sciences, 7, 25, doi: 10.3389/fspas.2020.00025
  * Tremblay et al. (2018) Tremblay, B., Roudier, T., Rieutord, M., & Vincent, A. 2018, Sol. Phys., 293, 57, doi: 10.1007/s11207-018-1276-7
  * Welsch (2015) Welsch, B. T. 2015, Publications of the Astronomical Society of Japan, 67, doi: 10.1093/pasj/psu151
  * Welsch et al. (2007) Welsch, B. T., Abbett, W. P., DeRosa, M. L., et al. 2007, The Astrophysical Journal, 670, 1434, doi: 10.1086/522422
  * Yeates (2017) Yeates, A. R. 2017, The Astrophysical Journal, 836, 131, doi: 10.3847/1538-4357/aa5c84