© ESO, 2015

1. Introduction

The Extreme ultraviolet Imaging Spectrometer (EIS) on board Hinode provides detailed observations of active regions, showing the existence of localised regions with systematic, relative low velocity upflows, that last from hours to days. However, neither the source region nor the mechanism which drives these upflows have been identified yet, although several theories have been put forward. There has been a general consensus to use the expression “outflows” for the phenomenon, with the indirect assumption, that these flows are related to open field lines, that allow the upflows to contribute directly to the slow solar wind. As most models have been based on local field extrapolations, it is not certain these field lines are really open and allow the plasma to escape the sun; therefore, we refer to these flows as upflows in this paper. For a detailed review of the observational contributions to this field see Vanninathan et al. (2015, hereafter Paper I).

Over time different theories have been suggested to explain the mechanisms, that drive the upflows. It has been suggested, that three- and five-minutes oscillations can drive upflows (Del Zanna 2008; Marsch et al. 2004). These oscillations have been identified in the upflow regions by Guo et al. (2010) and Ugarte-Urra & Warren (2011), who propose the underlying processes to be sporadic in both space and time. From a potential magnetic field extrapolation, Harra et al. (2008) found, that reconnection on one side of the investigated active region generated long loops which connect across the region. They concluded, that the observed upflows were the result of an expansion of the newly formed loops as they move up toward a new equilibrium.

In a different investigation Marsch et al. (2008) identified circulation of plasma through magnetic funnels and closed loops as a natural state of the plasma, that is continuously in motion. They suggested photospheric convection as the ultimate driving mechanism of the plasma motions in the corona. Using potential field extrapolations of a region containing two active regions Boutry et al. (2012) found, that about 18% of the mass upflow from one region was accounted for as downflow in the neighbouring region. Del Zanna et al. (2011) found a relation between radio noise storms and coronal upflows, indicating, that magnetic reconnection may be initiating the process. From magnetic field extrapolations they find magnetic null points and associated separatrices in the corona which could be associated with magnetic reconnection. A study of an active region-coronal hole (AR-CH) boundary by Van Driel-Gesztelyi et al. (2012) supports the hypotheses, that upflows are caused by magnetic reconnection along null separatrices. They identify a magnetic null point in the corona and propose this to be the site of continuous reconnection. Del Zanna (2008) speculates, that the observed blue-shifts could be the result of chromospheric evaporation driven by magnetic reconnection. The author hypothesises, that this mechanism is similar to the one observed during the gradual phase of two ribbon flares. Here, the upflows are seen from the footpoints of the newly formed loops as a result of the indirect heating from the reconnection point high above the photosphere. Srivastava et al. (2014) made a simple magnetohydrodynamic (MHD) simulation of a stratified loop and imposed an excess pressure at the photospheric region. Their experiment shows, that upflow velocities within the observed range can be generated by this mechanism. Baker et al. (2009) combined EIS observations and constant alpha force free magnetic field extrapolation based on Helioseismic and Magnetic Imager (HMI) data to study the magnetic field structures relation to the upflow regions. They identify quasi separatrix layers (QSLs; Priest & Démoulin 1995) between short closed and longer “open” field lines to be located in the vicinity of the upflow region and propose reconnection to take place along or near the QSLs. This magnetic field evolution can result in smooth extended upflow region as long as the QSLs are continuously stressed. Depending on the driver and field structure this model could generate a host of different effects such as chromospheric evaporation, reconnection jets, siphon flows along closed loops and upflows due to pressure gradients generated after the reconnection between the two different loops.

By studying an AR-CH complex through numerical simulations Murray et al. (2010) was able to model the upflows as a result of a new magnetic flux emergence. The new flux will push its way into the coronal hole region where it initiates compression waves on the open field lines. This short duration event is a precursor for the build up of a local current sheet where reconnection will be initiated. Following this Harra et al. (2010) found new flux emergence to be associated with a decaying active region that displayed upflows. Harra et al. (2012) investigated the role of a new emerging flux in producing upflows through MHD simulations where a new bipole emerges into an open magnetic field configuration. They conclude that the observed blue-shifts are either reconnection jets or pressure-driven upflows. The two MHD models (Murray et al. 2010; Harra et al. 2012) used in these investigations are identical and only the extent of the interpretation of the results differ slightly.

Although several indirect observations (McIntosh & De Pontieu 2009; He et al. 2010; Nishizuka & Hara 2011) have indicated a correlation between jets in the chromosphere and coronal upflows, in Paper I it has been shown through chromospheric spectroscopic imaging that the upflows persists despite the absence of chromospheric jets, indicating that the driving should be located higher in the atmosphere.

There has been a general assumption, that the upflow takes place along open magnetic field lines, allowing the plasma to directly contribute to the slow solar wind. An analysis of an active region by Culhane et al. (2014) shows this is not a necessary case. Their global potential field source surface (PFSS) model (which is the only investigation taking into account the large scale magnetic field) suggests, that a QSL related reconnection takes place close to the active region populating the field lines with active region FIP enhanced plasma. These field lines are further required to pass through a second reconnection with a distant 3D null, experiencing a spine-fan reconnection (Priest & Pontin 2009; Pontin et al. 2013; Galsgaard & Pontin 2011) process before the plasma finally enters an open field region. This requires a very precise line up of the different field line systems, that seems unlikely for the typical upflow regions and even a change from a potential to a non-linear force free field (NLFFF) model would most likely remove this alignment.

In this paper we investigate how restricted photospheric convective motions (only horizontal flow components) may assist in creating a systematic upflow in an active region by driving a localised reconnection process. A data driven MHD model is used in the investigation, where a time sequence of magnetograms are used to provide an assumed boundary driving flow, and the initial magnetic field structure is provided by a potential field extrapolation from a HMI magnetogram. The layout of the paper is as follows. In Sect. 2 the initial conditions and the numerical approach are discussed. The results are discussed in Sect. 3. Section 4 presents the results relative to the observations. The conclusions are given in Sect. 5.

2. Initial conditions

The data driven experiment presented in this paper is based on the observational study of NOAA 11123 analysed in Paper I. This investigation shows the southern edge of the active region to host a systematic upflow with a characteristic speed of 5–20 km s^-1 but also a second and much faster component of 105 km s^-1. The upflow lasts for at least the 3.5 h period the observations cover and it does not seem to be initiated by magnetic flux emergence, but more likely by some effect, that acts well above the chromosphere, but at the same time allows for ample diffusion of photospheric magnetic field (see Paper I). Having detailed observations by the Solar Dynamical Observatory (SDO) Atmospheric Imaging Assembly (AIA) and HMI makes the event an ideal target for a data driven investigation, aiming at understanding the physical implications of the photospheric driving on the systematic upflow.

To preform a data driven 3D MHD model, several independent steps are required. The line of sight HMI observations provide a time series of magnetograms with a time resolution of 45 s. These data allow for either a potential or a constant alpha force free 3D model extrapolations of the magnetic field structure. From the available full disk HMI data a suitable region around the active region was selected to limit edge effects, that may influence the extrapolated magnetic field structure. A NLFFF extrapolation was also attempted, but due to the location of the active region and its strong flux imbalance, the result was unsatisfactory for further analysis (private discussion with Thomas Wiegelmann). The selected region is shown in the top panel of Fig. 1, where the line of sight magnetic field component at the edges of the region is smoothed to zero over a short distance. This allows for a horizontally 2D periodic 3D magnetic field extrapolation. The HMI magnetograms are de-rotated to the same physical time, making a direct comparison between the time dependent structure of the region easier (see Paper I for more details).

Potential magnetic field extrapolations are obtained using a fast Fourier transform (FFT) approach. The first of these field extrapolations was used as the initial magnetic field for the MHD experiment. Before it can be used in the MHD code, it has to be further manipulated, making sure it is potential and divergence free when applying the high order stencil adopted in the code. A representation of the 3D magnetic field structure is shown in the bottom panel of Fig. 1. This shows selected field lines, traced both from the vicinity of the most important 3D null points in the domain, and from randomly scattered points on the photospheric surface.

A full discussion of the field topology is provided below. Combining model field-line traces with the information about the location of the upflow region identified in Paper I, shows, that the upflow region originates from the interface region between the open and closed flux region bounded by the green, orange and blue field lines (see Fig. 5 below). Here, open field lines are only to be taken as indicative as the extrapolation is conducted over a smaller region of the Sun not taking into account active regions further away.

Fig. 1 Top: initial HMI magnetogram used to derive the potential magnetic field. The dynamic range is ±1400 G. Bottom: 3D field line structure of the closed region areas with field lines traces from the vicinity of the dominating null points (green, orange, red, blue and yellow). These are combined with field lines of the “open” field (red) region located around the active region. The orientation is rotated slightly relative to the magnetogram, with the left, right and up down edges being the same in the two representations. The size of the magnetogram is 256′′×256′′.

A potential field model is not the most appropriate starting condition for the simulation as there is no free magnetic energy and this has to be build up before the magnetic field can start driving the upflow. When comparing the extrapolated field line structure with EUV observations (AIA 171 Å and 193 Å), only minor variations between the two extrapolations are seen, and redoing the model calculations with a constant alpha force free magnetic field does not increase significantly the amount of free magnetic energy. This also has the problem, that this type of extrapolation requires flux balance on the photospheric boundary to be performed which in the present example would significantly change the magnetic field configuration.

To stress the magnetic field configuration, a boundary flow is derived using the local correlation tracking (LCT) algorithm developed by Fisher & Welsch (2008). This algorithm utilises the information provided by the changes in the positions of the magnetic fragments between two consecutive HMI magnetograms to obtain representative horizontal velocity vectors for, that time interval. Repeating this on the whole dataset provides a time dependent 2D horizontal velocity flow, that to a given accuracy represents the change in the magnetic field positions as a function of time. In applying the LCT approach one needs to choose a stencil width over which the correlation is done. This value is critical for two reasons. First, its size determines the physical size of the flow structures, that may be determined by the algorithm. Second, choosing a very small stencil size allows small scale flow structures to be detected, while at the same time the algorithm becomes critically dependent on noise in the data. This can provide large fluctuations in the flow field within very short distances. This high amplitude noise subsequently has to be smoothed/filtered out before using the data for driving the model. The increasing noise level with a decreasing sigma value puts a limit on how small length scales the velocity structures determined by this method may have. Compared with the underlying HMI data, it is clear, that a lot of small scale motions exist in these data. These motions are responsible for the rapid advection of small fragments, representing typically a diffusion-like processes of the magnetic field structure. In the present approach we, naively, assume small scale motions are unimportant for the general evolution of the magnetic field, and therefore we apply the velocity field determined by the standard sigma = 19 value. These small scale motions naturally have implication on introducing small scale “noise” on top of the large scale advection of the magnetic footpoints. If they are important for stressing the global system, we are missing this energy input in this simulation. We use the full cadence of the HMI data, i.e. 45 s. The small scale convective motion does not change dramatically on this timescale (45 s) as the typical lifetime of the granulation is 5 min, i.e. around 6 times longer. One would, therefore, expect the LCT applied on a full cadence HMI dataset to provide a velocity profile, that slowly changes with time. The length scale, that is resolved with the applied sigma value, is much larger than the granular length scale and does represent an even longer timescale in the data. This makes it surprising to find a velocity profile, that changes dramatically between each dataset and, that this result seems to be independent on the choice of sigma smoothing of the LCT algorithm. An example of the time dependent changes in a single velocity component for a single pixel is shown in Fig. 2. This may be the result of the noise in the independent datasets, that gives rise to problematic velocity information. The critical is, that we still find similar fluctuations in the flow determination when a longer time base between observations is used in deriving a velocity profile.

Fig. 2 Change in one velocity component (in arbitrary units) as a function of HMI snapshot number (time) for a randomly chosen pixel in the plane.

Before the driving profile is used in the experiment, a simple smoothing of the data has been done over time to limit the largest and most rapid changes in the flow profile. Further to this, the flow velocity, like the magnetic field, is damped toward zero at the physical boundaries of the domain to allow it to be consistently used in a 2D periodic domain. These velocity data are used as input for the imposed boundary driving for the MHD experiment described below. The driver information is defined at discreet time intervals of 45 s, a time interval, that is much larger than the time stepping used in the MHD code. To provide a smooth change in the driving velocity as a function of time a third order interpolation algorithm is applied to the time varying LCT velocity to get the actual flow vector at a given time.

As a supplement to the potential magnetic field a simple hydrostatic 1D model of the solar atmosphere is applied. The domain covers the region from the photosphere and into the corona. The model is defined by specifying the temperature profile with height and then integrating the hydrostatic equation using a simple integration algorithm. Subsequently, this profile is iterated to fit the hydrostatic criteria using the chosen grid spacing and the high order finite difference algorithm applied in the MHD code. For simplicity the plasma is assumed to follow the ideal gas law, and optical thin radiation and heat conduction are ignored in the energy equation. These terms are ignored in the present investigation, as we are testing the idea, that a magnetic reconnection process directly provide the reason for the continuous upflow. The main aim is, therefore, to look for these processes, that to a very large extent develop independently of the accuracy of the energy equation. When starting from a potential magnetic field and a simple 1D hydrostatic equilibrium, one needs the magnetic field to be stressed for a period of time before it reaches a level where it is possible to ignore an artificial background heating source to maintain a hot corona. From the analysis in Paper I it is clear, that the upflow does not originate from below the lower transition region, where most evaporation flows are initiated. Including these terms in the investigation is therefore a direct contradiction to what we are trying to do, namely a direct data simulation where the magnetic field is to be stressed in a way, that to a feasible degree represents the driving in the observations and hopefully ends driving a reconnection process. Figure 3 shows a 1D representation of the field quantities imposed on the initial conditions. Bourdin et al. (2014) did a similar experiment with a different active region. There they show, that to include heat conduction and radiation from this type of initial conditions takes about an 1 hour solar time to reach an energetic steady state. After this the calculations become expensive due to the restriction on the time step from the heat conduction. How this influences the magnetic structure and its ability to dissipate magnetic energy through reconnection has not been investigated.

Fig. 3 1D hydrostatic model atmosphere parameters in code units as functions of height (Mm), showing the density (solid line), the temperature (dashed line) and the gas pressure (dot-dashed line). For comparison, the average magnetic field strength is shown as the triple-dotted-dashed line with the long-dashed line representing the average plasma beta of the model.

The experiment is evolved in time by solving the non-ideal MHD equations (as given in the paper by Baumann et al. 2013) using the Copenhagen Stagger Code which contains hyper-diffusive viscous and ohmic terms self-consistent in the momentum, energy and induction equations (Nordlund & Galsgaard 1997). The numerical box has 512 × 512 × 511 points in the (x,y,z) directions defined on a uniform grid. The model represents a domain of size 184.3 × 182.3 × 91.4 Mm. This gives a grid resolution of 0.36 Mm in the horizontal direction and 0.18 Mm in the vertical direction. The scaling between the code units and real units are given by setting t₀ = 100 s, L₀ = 10⁶ m and ρ₀ = 10^-4 kg/m³ from which it follows, that T₀ = 10⁴ K and B₀ = 1121 Gauss. The boundary conditions of the domain are periodic in the two horizontal directions (x,y) for simplicity. The domain is closed in the vertical direction (z). This prevents mass flow in/out of the domain and limits the driving velocity on the lower boundary to be only the horizontal flows defined by the LCT approach. As the LCT approach provide a proxy for the flow speed and structure, we take the freedom in the experiment, namely to allow changes in the imposed driving amplitude, and the time interval between the different datasets. Such changes allow for the investigation of different characteristics in the flow patterns impact on the dynamical evolution of the magnetic field configuration and is a typical way to speedup the slow driving time relative the coronal transit time. This is often used to save computing time and to prevent diffusion releasing the imposed free energy before it has reached a level where more dynamical phenomena can take place. This is a common issue when conducting numerical experiments where the diffusion is many times larger than in the real solar environment. Table 1 contains the driver parameters used for the various models discussed in the paper. Additional to the zero normal velocity component imposed on the top and bottom boundaries, a dedicated damping layer is implemented in a smaller part close to the top of the domain. This prevents upward propagating waves to be reflected downward by the imposed zero normal boundary velocity. The damping is implemented using a Newton-cooling type algorithm, that brings the plasma density, the internal energy and the momentum back toward their initial values on a timescale, that increases exponentially with distance from the top boundary. The damping timescale is chosen such, that a short wavelength amplitude pulse is fully damped on its way through this region (the characteristic Alfvén crossing time). This approach was first used in Moreno-Insertis & Galsgaard (2013) where it proved to work excellent even in the situation where a high velocity jet impacts the top boundary for a long period of time. A side effect of the layer is, that it acts as a sink/source for the plasma density and thermal energy when the region below the damping zone changes its characteristic parameters.

3. Results

This section discusses the general evolution of the experiments, using Exp-1 as the reference case. Significant deviations found in the other experiments are discussed when appropriate. The time evolution of the general magnetic field topology, the identified velocity patterns, the typical Poynting flux, and the mass evolution of the active region are discussed below.

3.1. Magnetic topology evolution

To understand the dynamical evolution of the experiments an insight of the magnetic field topology is required. This provides information about the spacial locations where interaction between different magnetically connected domains may be able to drive systematic flows for longer periods of time. That is if the involved boundary regions are exposed to the correct systematic stress. Understanding the magnetic skeleton defining the topology of the magnetic field requires knowledge of the magnetic null point locations (Bungey et al. 1996) and their local structure, spine axis and fan plane (Parnell et al. 1996), the presence of separatrix lines connecting individual null points and the possible existence of QSLs (Priest & Démoulin 1995; Titov 2007), that are independent of the null points.

To investigate the skeleton, magnetic null points are identified using the method developed by Haynes & Parnell (2007). This approach assumes the magnetic field to be well represented by a trilinear interpolation between the grid-points. For the initial magnetic field configuration, defined below, this assumption is only challenged close to the photospheric boundary where the fastest changes in the magnetic field exist. Null points found just a few grid points above this surface easily fulfils this requirement due to the exponential decay in complexity of the field with height.

The initial potential magnetic field is extrapolated from the HMI data taken on November 13, at 13:02  UT, as described above. The size of the domain is chosen such, that the field line structure around the active region complex is not critically influenced by the imposed 2D periodic boundary conditions. In contrast to this the field lines connectivity close to the periodic boundaries does not represent the real solar magnetic field and a large fraction of these would most likely be closed over a larger distance. Figure 4 shows the horizontal positions of the null points in the domain, that are located at least one grid point above the bottom boundary. Their colour indicates the height of the nulls. Two things are to be noticed. First, a large fraction of the nulls are scattered around the edges of the 2D horizontal domain, where they are introduced by the artificial damping of the normal component of the magnetic field, and also by the rapid changing of the weak field structure. These nulls are, therefore, not considered in the following investigation. Second, a number of nulls are located in the close vicinity of the active region. These will be discussed below, as they are important for understanding the topological structure of the magnetic field around the observed upflow region.

Fig. 4 Horizontal positions of the null points in the initial potential magnetic field relative to the frame centre. The colour coding of the points indicate the heights of the nulls points. The square-root of the height is shown to compress the scale. The contour lines (full and dashed) show the locations of the strong magnetic flux concentrations (±500 G). Top: 74 null points above 0.54 Mm. Bottom: the 13 central null points listed in Table 2.

Fig. 5 Representation of the initial magnetic field line topology. The different coloured field lines are traced from the vicinity of different null points, above the active region. These indicate the characteristic connectivity of the magnetic flux. The spatial orientation and size are comparable to Fig. 1. The base is the real magnetogram scaled to ±500 G. The white contour line in front of the footpoints of the magnetic field lines represents the footprints of the outflow region marked for a more detailed investigations in Fig. 6 in Paper I.

Fig. 6 Top: quasi separator layer map indicating the locations of rapid changes in the field line mapping across the core region of the domain using the invariant measure defined by Titov (2007). The data represent ±45 Mm in the x and y around the centre of the domain. The data are logarithmically scaled over values from 2 to 1000. White areas have low values. Bottom: QSL map and its relation to the null points close to the active region.

Figure 5 shows the skeleton of the magnetic field in the vicinity of the active region. The data slice in the bottom of the frame represents the magnetogram on to which the region of the footpoints of the upflow defined in Paper I is marked with a white line. The different coloured field lines represent the structure of the magnetic field close to the most important null points, clearly showing the spine axis and fan planes, that separate the magnetic field line connectivity into independent flux regions. The nulls represented by blue and green field lines both have their fan planes extending down toward the photosphere, rooted in the negative (red) polarity flux. The spine axis of these nulls are, therefore, connecting either to the positive (purple) flux in the photosphere or to the top boundary of the domain. The three other coloured field line regions are more complicated and can, from Table 2, be seen to contain more than one null point. The spine axis of the dominating null points (yellow and red) are typically in the horizontal direction, which is clearly seen for the yellow null region. The fan planes of these nulls, on one side, reach toward the spine axis of the blue and green nulls, and therefore intersect with their fan plane. This implies, that several separator lines connecting between the different nulls, providing the division of the magnetic field line connectivity into several independent flux domains. Depending on where the field lines are traced from the photosphere they either belong to a closed or an open loop system with regard to the extrapolation of the limited field region. Especially, it is noticed, that all field lines starting from the negative flux in the active region are closed field lines connecting to the positive flux regions due to the general excess of positive flux in the active region. The negative flux regions therefore contain both closed and open field lines, where the open field lines connect mainly to the top boundary of the domain. This is seen in the bottom frame of Fig. 1, where the red field lines outside the null dome regions connect from the photosphere to the top boundary. This open field line structure is an artefact of the limited region used for the field extrapolation. The position of the thirteen null points are given in Table 2 for the initial condition in Exp-1. Investigations of magnetic field extrapolations have often used the concept of QSL to identity locations in space where the field line mapping changes particularly fast and reconnection therefore may be initiated. A QSL does not need to contain a mathematical discontinuity in the field line mapping, as implied by the presence of 3D null points, while a null point can be seen to be embedded in a QSL region (Pontin et al. 2013). Quasi separator layer regions, like null points, fan planes and separator lines, are found to be locations in space where current easily accumulates when the magnetic field configuration is appropriately stressed. Figure 6 shows a QSL map for the central region of the domain. It is seen, that the main QSL regions are located around the locations where the fan planes of the overlying null points intersect the photosphere. Other QSL regions are also seen further away from these null defined regions. These may be associated with the many null points located closer to the domain boundary (see top frame of Fig. 4). This indicates, that the magnetic field contain lots of small-scale variations due to the actual mapping outside the active region. It is clear from the analysis of the magnetic skeleton, that the underlying structure is very complex, which allows many different reactions to the imposed boundary driving.

Fig. 7 Comparison of the magnetograms from the HMI observations and the numerical experiment. The top frame represents HMI data 193.5 min after the initial state seen in the top frame of Fig. 1. The bottom panel represents the flux distribution after the same time derived by the MHD code. The dynamical scale in the frames represent ± 1500 G.

The imposed driving does not introduce large changes to the flux distribution on the model photosphere over the duration of the experiment and it does not allow for flux to emerge or submerge. This is both an advantage and a limitation with the imposed flow pattern and the general handling of the boundary conditions. Comparisons between the HMI observations and the evolved model show clear differences and, that these increase with time. An example is shown in Fig. 7, that represents one time in the evolution. A more detailed comparison shows, that most of the differences are due to small scale dynamical evolution in the observed magnetograms, that is not represented by the driver obtained from the LCT.

3.2. Velocity evolution

To identify the mechanism behind the upflows, regions with systematic positive vertical flow velocities need to be located in the experiment. For this, three different visualisations of the data are used.

Fig. 8 Top: representation of a vertical slice (y, z-plane) in the observed upflow region of the domain, showing the vertical flow speeds of ±8 km s-1. Middle: the same but for a horizontal plane (x, y plane inside the corona). Bottom: time slice of the density variation in height, showing how structures are moving upward in the atmosphere with a speed of roughly 150 km s-1. The plot also show the period of the oscillations is of the order of 3.7 min. The snapshots used for the top and middle frames represent t = 132.5 min.

First, we investigate the vertical velocity in a 2D vertical plane, that in its horizontal direction passes through the upflow region identified by the observations (north-south orientation). This shows complicated flow patterns covering the whole plane, with an alternating velocity structure, that is advected upward in time. An example of this is shown in the top panel of Fig. 8. This flow pattern represents velocities with a peak of the order of 10 km s^-1. The alternating flow patterns are wave motions, that represent internal compression/decompression speeds in the vertical direction, while the wave patterns move upward along the field lines with a much higher speed.

Second, by only investigating a single vertical plane one can easily miss different flow behaviours in other parts of the domain. This can be helped by investigating the vertical flow speed in a horizontal plane, at the base of the coronal region, where one can easily look for systematic flow patterns. An example is shown in the middle panel of Fig. 8. This shows the whole plane is covered with alternating velocity patterns. Following these structures in time, the patterns are found to be constantly changing without showing any localised regions containing a systematic flow pattern.

Third, to get an estimate of the upward motion of the flow patterns, one can plot the relative changes in the plasma density from a single column in height as a function of time. This shows the upward propagation of density structures. The inclination of the patterns represents the propagation speed, bottom panel of Fig. 8. Using this to estimate the propagation velocity in the coronal part of the domain, between 20 and 80 Mm in height, one gets a velocity of the order of 150 km s^-1. This simple visualisation does not take into account the magnetic field line structure, but at this hight in the atmosphere its direction is dominantly in the vertical direction.

Fig. 9 Characteristic speeds of the Alfvén velocity and the sound speed. The full line shows the Alfvén velocity for the initial setup. The dotted line shows the comparable profile at t = 132.5 min. The same is shown for the sound speed with the dashed line showing the initial profile while the dot-dashed represents the later time.

To check the speed of the upward propagating oscillations with the characteristics flow speeds in the atmosphere, Fig. 9 shows the Alfvén and sound speeds for t = 0 and t = 132.5 min. The overall speed scales between 0 and 485 km s^-1. The typical sound speed is 150–170 km s^-1 in the coronal region, while the Alfvén velocity is higher. This shows, that a propagating wave in this coronal model will obtain rather different speeds depending on where in the domain it is located and which wave type it represents. The typical wave speed observed in the density perturbation shown in Fig. 8 indicates a flow speed of around 150 km s^-1, which strongly indicates the waves to be a slow magnetosonic wave with a speed close to the sound speed.

Figure 10 shows a 3D representation of the field line structure combined with blobs of the vertical positive/negative flow velocity. Inspecting such visualisations of the data shows the oscillating flow patterns to be space filling and upward propagating in time.

Fig. 10 Magnetic field line structure is shown in the same way as in Fig. 5. Additionally, the red/purple blobs represent the vertical down/up flow velocities indicating the compression/decompression patterns of the propagating waves. (See the attached movie for a time-dependent evolution of this visualisation, it represents the time period from 83.3–325.0 min.)

Repeating the experiment with different driving parameters (Table 1), one finds the same characteristic behaviour, showing, that the typical dynamical evolution is independent of the imposed driving velocity and characteristic timescale of the driver, and even of the driving pattern. The change in the driving time influences the characteristic timescale of the oscillations and wavelength of the travelling perturbations.

3.3. Poynting flux evolution

To measure the magnetic energy flux through different heights of the corona the Poynting flux is calculated: (1)The first contribution represents the advection of magnetic energy through the surface, while the second term is the work done by advecting the field lines in the horizontal direction in the plane. Time series of Pf for different heights in the corona show, that the time resolution of the saved data from the experiments is too low to directly follow the evolution of the Poynting flux patterns in the surface, this was also seen for the velocity components, and especially the horizontal components. It is clear, that the patterns fluctuate significantly in both space and time. Initially this fluctuation is mainly concentrated above the active region, while the influence spreads to cover the whole plane as time progresses. To smooth out the influence of the low time resolution in the time series (50 s between datasets), one can instead look at the surface integrated Poynting flux as a function of time. This quantity indicates whether magnetic flux, on average, is transported through the plane as a function of time. The result is shown in Fig. 11, where the top frame shows a snapshot of the Poynting flux through a vertical plane and the bottom frame shows the integrated Poynting flux through the plane as a function of time. This shows a small positive integrated Poynting flux through this plane, that is about 2 orders of magnitude less than the peak values, indicating, that magnetic energy, on average, is advected up through the plane to the volume above. Compared to the existing magnetic energy above this hight, the magnetic energy input is totally insignificant for the evolution of the magnetic energy over the time period covered by the experiment.

Fig. 11 Top: snapshots of the Poynting flux at a constant height in the corona. The plane is from t = 107.5 min at a height of 35.9 Mm. Bottom: time variation of the surface integrated Poynting flux as a function of time for the same height in the model.

Driving with a slower speed shows comparable spatial and time structures of the Poynting flux, with the space and time variations being even closer to zero.

3.4. Mass evolution

The initial density profile is shown as part of Fig. 3. Due to the large variations with hight, the simplest way to quantify density changes is to measure it relative to the initial profile. Doing this, two characteristic changes are seen in the data. First, it is found, that the density of the medium transition region (height 3−6 Mm) is nearly doubled as a time average of the full experiment. This significant increase in mass arises from two regions, one is mainly from an expansion of the layer below 3 Mm. This is caused by the imposed boundary stress which locally increases the magnetic and gas pressure of this region, making it expand in height. The second is a minor contribution from a general decrease in the density from the region above 6 Mm. The coronal mass is slowly being depleted. Large local variations from this behaviour is seen, showing fluctuation in both space and time. These variations depend on both the magnetic field line structure and the imposed flow profile. An example is shown in Fig. 12 where the percentage change of the density for a vertical plane is shown. In Fig. 13 the average percentage change for a horizontal plane is presented for both Exp-1 and Exp-2. This shows a constant decrease of the mass in the coronal region, with the actual rate being dependent on the speed of the imposed driving velocity.

Fig. 12 Percentage change of the density profile in a vertical plane between 80% and 120% of the initial value. The plane is the same as shown in Fig. 8.

Fig. 13 Time variation of the density of a plane at a constant height in the corona – data from different heights in the corona show the same time evolution. The full line represents experiment Exp-1, and the dashed line Exp-2.

The massflux (m·u_z) through a horizontal surface has also been investigated. It shows the decrease to be fixed in time, and, that the loss rate depends on the height in the atmosphere. The panels in Fig. 14 show snapshots of both a vertical and a horizontal plane at t = 132.5 min. Averaging the values in the horizontal plane gives a time dependent evolution which shows a near constant depletion of mass with time and a decreasing magnitude as a function of height in the atmosphere. The panels in Fig. 14 show large variations in the mass flux in space (and time), following a pattern which resamples the previous discussion of the mass density and vertical velocity. Supporting the presence of only wave phenomena in the experiment. In the table of Fig. 5 of Paper I the density along a presumed loop is measured. Using the potential magnetic field model it is possible to estimate the magnetic field strength in these region and from this estimate the area change or density change assuming, that the same amount of plasma is moving long the field lines all the time. The result of this exercise is shown in Table 3. This indicates, that the observed change in density to first order can be explained by the general expansion of the magnetic field under the assumption, that the boxes are represented by nearly the same magnetic field line.

Fig. 14 Top: mass flux in a horizontal plane. Bottom: mass flux in a vertical plane. The data is scaled to ± 10-7 kg m/s.

4. Discussion

It is important to find the physical mechanism, that drives the systematic upflows from the edges of active regions. When this mechanism has been identified we may also have found the mechanism, that drives the slow solar wind from the boundary region of coronal holes. For this study we used a numerical 3D MHD model based on an observational magnetic field where HMI data are used as the basis for a potential field extrapolation. To stress the field configuration, a boundary driving derived from a LCT approach is imposed. The results discussed above show, that the experiment does not provide a region containing a systematic upflow. This indicates, that there is no continuous physical process, in this model scenario, that can drive such a flow, i.e. magnetic reconnection or expanding loop systems. Instead, magneto-acoustic waves are continuously being launched from the lower atmosphere and propagate into the corona.

The experiment is expected to include most of the important physical conditions for reproducing upflow if these are directly driven by a reconnection process, so why is it failing in doing so? There may be several reasons for this. Here, we discuss the key building blocks and their possible impact on the experiment.

4.1. Potential field model

The numerical experiment starts from a potential magnetic field model extrapolated from HMI data, while the solar coronal magnetic field is more likely to be closer to a NLFFF with small local alpha values. A NLFFF model requires high quality vector magnetogram observations. Vector magnetogram observations exist for the analysed region, but deriving a NLFFF model was not possible due to poor data quality and a significant imbalance in the magnetic flux at the source surface for this particular region.

Would a NLFFF model be better than a potential model? First, a NLFFF model contains free magnetic energy, which implies, that much less imposed boundary stress may drive the magnetic configuration into a situation where some fraction of the free magnetic energy can be released through a local long duration reconnection process. Second, a NLFFF would typically reduce the characteristic height of the closed magnetic field system. A more compact field would provide for a higher Poynting flux input for the same boundary motion, leading to a larger stress of the field and more free energy to be used for driving the upflow.

Would it make a significant difference to the result? Creating a long duration systematic upflow in a NLFFF system, compared to a potential system, is unlikely to make a significant difference as the magnetic connectivity of the two magnetic fields are related, implying, that the boundaries between closed and open field regions cannot be significantly different as both extrapolations will contain magnetic null points approximately at the same locations and of the same type. With other words, the topology of the two fields will host comparable connectivity regions for energy accumulation and release.

Recent numerical simulations of magneto convection (Nordlund et al. 2009) show the photospheric magnetic field to contain much smaller length-scales than at present can be resolved in satellite observations. These smaller scales introduce both shorter characteristic dynamical time- and length-scales of the 3D magnetic field configurations. This enhanced salt and pepper distribution of the magnetic field on small scales leads to a more complicated local field-line structure. The small scale convection allows for an additional energy input to the lower regions of the atmosphere through flux braiding (Parker 1972; Priest et al. 2002) and, therefore, a higher transferee of energy into plasma heating through dissipation of current concentrations than can be inferred from our modelling. How this energy transfers into higher heights on the other hand is not clear, but it may be transferred further up in the atmosphere by MHD-waves, that may be dissipated by, for instance, phase mixing (Heyvaerts & Priest 1983) providing additional heating of the coronal domain, that may be important when including conduction and radiation in the model.

4.2. Driving mode

The LCT approach has a number of limitations, it only measures horizontal velocities based on the tracking of magnetic features. This works fine as long as the magnetic elements are small independent units, that move “randomly” around in the 2D plane and the time resolution in the data is sufficiently high to resolve the motions. When the flux concentrations become too large and their structures are too smoothly distributed, the algorithm have problems deriving the correct flow velocities. This is particularly a problem in the active region area, where the derived flow speeds are clearly suppressed compared to the rest of the area.

Due to the kernel size (19 pixels) used for the detection of flows, small scale motions, which seems to be very important for the local evolution of the magnetic fragments when inspecting the HMI data, cannot be detected. Ignoring the small scale evolution of the flow field clearly limits the evolution of the magnetic field, and possible also the general energy flux into the corona (van Ballegooijen et al. 2014). Smaller kernel sizes do provide smaller flow structures, but also introduce a significant increase in the number of locations where excessive flow speeds are derived due to noise/errors in the two compared magnetograms. The LCT is, therefore, not representing the real small scale motions in the photosphere, that may be important for driving high frequency perturbations into the magnetic configuration.

A physical process missing in the LCT approach, is the emergence and submergence of magnetic flux, which requires a vertical flow velocity. This process continuously takes place on different length scales and is important for the evolution of the coronal magnetic field. Observations of this active region (see Paper I) indicate, that no significant flux changes due to large events of emergence or submergence take place around the expected footpoint region of the field lines showing upflows, but this region shows a general decrease in the unsigned flux over time. It can, therefore, not be ruled out, that the change in flux is important for driving the upflows and signals, that some important effect is taking place somewhere above the photosphere. To get an idea of the change in the boundary structure between the HMI data and the model one can make direct time comparisons. This shows differences in the flux distribution and, that this continues to increases with time. The HMI data show a decrease in the unsigned flux with time, while it remains close to constant over time for the numerical experiment. This implies, that the two field representations diverge with time and, that these differences are most pronounced on small and medium length scales (Fig. 7). It, therefore, cannot be expected, that the simulation result will represent the evolution of the smaller scale magnetic field. As the upflow last for hours to days, to us, this indicates, that the small scale motion cannot be too important for the general driving of these phenomena if we assume it to be driven by a systematic reconnection process. The proposed driving for the wave motion as mechanism of the upflows also does not seem to directly initiate upflows either.

The LCT approach has a natural limitation on the velocity structures on scales smaller than the kernel size. This will limit the energy flux into the magnetic field and therefore only provide a lower limit to the local energy flux.

4.3. Atmospheric model

The initial atmosphere model is a simple 1D representation, which contains the basic properties of the solar atmosphere. Imposing a different hydrostatic model would not have a significant impact on the dynamical result of the experiment. The major difference would be to change the slope of the transition region and by this the wave speed through this region. In terms of creating a systematic upflow, it is unlikely to have a significant impact on the results, as the evolution of the atmosphere is mostly controlled by the structure and stressing of the magnetic field.

Here, differences in the gas pressure between open and closed magnetic field regions may be important for driving plasma flows when magnetic reconnection between them takes place. This is not, initially, included in the present experiment, but may be build up over time as the model is stressed by the imposed boundary motions.

4.4. MHD approach

The MHD approach used in the calculations is simplified, ignoring optically thin radiation and anisotropic heat conduction in the energy equation, two physical contributions, that are important in the solar atmosphere. Including these terms from the start of the experiment would not be simple. Without a direct heating contribution, the upper atmosphere would cool rapidly resulting in an unrealistic coronal model. One could, as others have done (Gudiksen & Nordlund 2002), drive the initial field configuration to a steady state situation, where the energy input is converted into localised heating throughout the model. At this time the optically thin radiation and anisotropic conduction are activated. From Chen et al. (2014) it is seen to require approximately one hour solar time for a closed loop system before a steady state is reached, and possibly much longer for an open field region. This is, therefore, not a first option in the present case, as we are mainly interested in seeing whether magnetic reconnection can drive the outflow as a response to the driving. If one would include this, the experiment should start at least an hour earlier to allow the field to build up enough complexity to form current structures on various scales. This is needed to provide enough heating through dissipation to balance the effects of heat conduction and optically thin radiation. Another problem is the LCT driver, which does not provide a sufficiently correct flow field to represents the actual stressing of the magnetic field. Using only the scalar magnetograms in the process at the same time prevents the method to take into account the continuously ongoing emergence and submergence of flux, that with an increasing likelihood takes place in active regions as the investigated time interval increases.

Would this lead to the formation of a systematic upflow of the order of 10–50 km s^-1 if the conditions were correct? It is found in some reconnection experiments (Yokoyama & Shibata 2001), that energy release in the corona can be channelled down toward the photosphere by heat conduction. This may heat the lower atmosphere resulting in an expansion where the increased plasma pressure creates an upflow creating an atmospheric evaporation scenario. This effect was seen in Srivastava et al. (2014) where they conducted a 2D simulations where they imposed an energy deposition at the lower boundary, resulting in a slow pressure driven upflow. If this is the case, then the experiment needs to be updated to include these energy terms and a systematic coronal energy release would be required. The problem with this type of experiments, is the lack of physical explanations, that therefore only pushes the understanding to some unexplained mechanism(s).

The above limitations of including short length scale magnetic field and boundary driving, point in a direction where one assumes, that more energy may be transferred into the lower atmosphere. On average this would heat all over the solar surface and is, therefore, unlikely to be responsible for local heating processes, that may drive systematic upflows lasting for hours due to an evaporation process. In this context, the energy conduction from the corona to the lower transitions region may locally be important. But this requires a sufficient enhancement in the energy deposition in the corona to allow an additional energy conduction allowing for an evaporation process. In the experiment waves are travelling up through the “open” field region. If inhomogeneities are present in the corona some of this energy may be tapped locally by phase-mixing (Heyvaerts & Priest 1983). This is a more realistic candidate to heat the corona in the open field regions than direct flux braiding (Parker 1972).

4.5. Location of reconnection

Where is the reconnection most likely to take place? Having magnetic null points in the model suggests, that these are the locations where reconnection can/will take place. But most of these are high in the atmosphere and an ongoing reconnection here would very likely result in a local hot region and would also initiate high velocity outflows. Something like this is seen in Moreno-Insertis & Galsgaard (2013) where an emerging flux system interacts with an open field region. The observational characteristics from this type of event clearly does not fit the observations of the slow upflow. The direct reconnection process high in the atmosphere, therefore, seems rather unlikely. This also rules out the pressure driven upflow (atmospheric evaporation) caused as a secondary effect as it is seen in postflare loops. To slow down the process, the reconnection has to take place lower in the atmosphere, where the impact of the plasma heating is much less. Here, the dominating driver of the flows may be the general pressure difference between the two flux regions, where the short closed flux systems typically will have an excess pressure compared to the “open” field regions.

4.6. Effect of heat conduction

Including heat conduction in the calculations leads to an energy transfer from the corona to the lower transition region along the magnetic field line structure. In the general situation most of this energy would be radiated away by either optically thin radiation in the hotter regions or by optically thick radiation in the lower atmosphere. On average, a balance between the different energy terms exists, but when additional energy is released locally in the corona, it will effectively be conducted down along the field lines. This may lead to an evaporation flow, that pushes heated material upward into the corona. Chen et al. (2014) showed this interplay with localised coronal heating on closed field lines leads to an evaporation flow, that heats up an emerging loop system. This scenario may be important for the systematic upflows observed at the edges of active regions. In many observations the upflows typically take place on open field lines, that reaches large heights in local models. An additional localised heating event on these lines are, therefore, required to initiate an evaporation flow. This heating is not easily explained as these field lines contain very little current as they are typically close to being potential in structure – a DC mechanism and localised reconnection, therefore, does not seem likely. The alternative is the conversion of wave energy through phase-mixing. This either requires a sufficient density gradient between different field lines to take place within the length scale of the model, or waves with much shorter wavelength than the ones seen in this experiment.

5. Conclusions

For the first time a direct data driven model scenario has been used to simulate the dynamic evolution of an active region, that was observed to produce a long duration, slow upflow event. The model is based on a potential magnetic field extrapolation using HMI data, and it has a complex internal structure defined by the presence of a high number of magnetic null points above the active region. These magnetic nulls would allow for systematic localised reconnection processes, that, directly or indirectly, may be able to drive an upflow. The magnetic field is stressed by imposing a photospheric flow velocity based on a LCT of magnetic elements in the HMI data series. For the conditions in the model, the experiment was not able to reproduce magnetic reconnection driven systematic upflow found in the observations of NOAA 11123 (see Paper I). Our experiment has a preference for producing wave motions, that propagate upward in the atmosphere with a frequency, that needs to be further investigated. However, similar oscillatory motions are seen in a number of observations of active region upflows with typical frequencies of the order of 3–5 min (Guo et al. 2010; Ugarte-Urra & Warren 2011).

A systematic stressing of the magnetic field by the convective motions cannot be ruled out to drive a plasma upflow. Future experiments using a different approach to determine the imposed stressing of the system will be tried. Here especially the DAVE code (Fisher et al. 2010, 2012; Kazachenko et al. 2014) approach sounds promising, as its approach will provide a more realistic velocity field, that better reproduces the evolution of the magnetic field. Additional trying to include more terms in the energy equation on longer experiments will show whether this has any influence on the general evolution of the region.

The experiment shows significant wave propagation in the atmosphere directly driven by the imposed boundary driver. As the magnetic setup is complicated, it provides a good opportunity to use a wave analysis approach to see how this can be used as diagnostics for the underlying physical variables. This will be investigated in a future paper.

Appendix A

Fig. A.1 Vertical flow-pattern motion upward along magnetic field lines as a function of time. The base represents the magnetogram; the field lines are traced from a small pillow box around the active region in the photosphere. The field lines are coloured with the vertical flow speed, where red is downflow and purple is upflow. The amplitude of the velocity field is scaled to ±10 km s-1. The movie covers the time period from 83.3 to 325 min after the start.

Acknowledgments

This work was supported by a research grant (VKR023406) from VILLUM FONDEN. Brian Walsh for providing the LCT software. VAPOR softwareImagery produced by VAPOR (www.vapor.ucar.edu), a product of the Computational Information Systems Laboratory at the National Center for Atmospheric Research. Z.H. is supported by the China 973 program (2012CB825601) and the National Natural Science Foundation of China under contract 41404135. M.M. is funded by the Leverhulme trust. Research at Armagh Observatory is grant-aided by the N. Ireland Department of Culture, Arts and Leisure and via grants ST/J001082/1 and ST/M000834/1 from the UK STFC.

References

Online material

Movie of Fig. A.1 (Access here)

All Tables

All Figures

Fig. 1 Top: initial HMI magnetogram used to derive the potential magnetic field. The dynamic range is ±1400 G. Bottom: 3D field line structure of the closed region areas with field lines traces from the vicinity of the dominating null points (green, orange, red, blue and yellow). These are combined with field lines of the “open” field (red) region located around the active region. The orientation is rotated slightly relative to the magnetogram, with the left, right and up down edges being the same in the two representations. The size of the magnetogram is 256′′×256′′.

In the text

	Fig. 2 Change in one velocity component (in arbitrary units) as a function of HMI snapshot number (time) for a randomly chosen pixel in the plane.
In the text

	Fig. 3 1D hydrostatic model atmosphere parameters in code units as functions of height (Mm), showing the density (solid line), the temperature (dashed line) and the gas pressure (dot-dashed line). For comparison, the average magnetic field strength is shown as the triple-dotted-dashed line with the long-dashed line representing the average plasma beta of the model.
In the text

Fig. 4 Horizontal positions of the null points in the initial potential magnetic field relative to the frame centre. The colour coding of the points indicate the heights of the nulls points. The square-root of the height is shown to compress the scale. The contour lines (full and dashed) show the locations of the strong magnetic flux concentrations (±500 G). Top: 74 null points above 0.54 Mm. Bottom: the 13 central null points listed in Table 2.

In the text

Fig. 5 Representation of the initial magnetic field line topology. The different coloured field lines are traced from the vicinity of different null points, above the active region. These indicate the characteristic connectivity of the magnetic flux. The spatial orientation and size are comparable to Fig. 1. The base is the real magnetogram scaled to ±500 G. The white contour line in front of the footpoints of the magnetic field lines represents the footprints of the outflow region marked for a more detailed investigations in Fig. 6 in Paper I.

In the text

Fig. 6 Top: quasi separator layer map indicating the locations of rapid changes in the field line mapping across the core region of the domain using the invariant measure defined by Titov (2007). The data represent ±45 Mm in the x and y around the centre of the domain. The data are logarithmically scaled over values from 2 to 1000. White areas have low values. Bottom: QSL map and its relation to the null points close to the active region.

In the text

	Fig. 7 Comparison of the magnetograms from the HMI observations and the numerical experiment. The top frame represents HMI data 193.5 min after the initial state seen in the top frame of Fig. 1. The bottom panel represents the flux distribution after the same time derived by the MHD code. The dynamical scale in the frames represent ± 1500 G.
In the text

Fig. 8 Top: representation of a vertical slice (y, z-plane) in the observed upflow region of the domain, showing the vertical flow speeds of ±8 km s-1. Middle: the same but for a horizontal plane (x, y plane inside the corona). Bottom: time slice of the density variation in height, showing how structures are moving upward in the atmosphere with a speed of roughly 150 km s-1. The plot also show the period of the oscillations is of the order of 3.7 min. The snapshots used for the top and middle frames represent t = 132.5 min.

In the text

	Fig. 9 Characteristic speeds of the Alfvén velocity and the sound speed. The full line shows the Alfvén velocity for the initial setup. The dotted line shows the comparable profile at t = 132.5 min. The same is shown for the sound speed with the dashed line showing the initial profile while the dot-dashed represents the later time.
In the text

	Fig. 10 Magnetic field line structure is shown in the same way as in Fig. 5. Additionally, the red/purple blobs represent the vertical down/up flow velocities indicating the compression/decompression patterns of the propagating waves. (See the attached movie for a time-dependent evolution of this visualisation, it represents the time period from 83.3–325.0 min.)
In the text

	Fig. 11 Top: snapshots of the Poynting flux at a constant height in the corona. The plane is from t = 107.5 min at a height of 35.9 Mm. Bottom: time variation of the surface integrated Poynting flux as a function of time for the same height in the model.
In the text

	Fig. 12 Percentage change of the density profile in a vertical plane between 80% and 120% of the initial value. The plane is the same as shown in Fig. 8.
In the text

	Fig. 13 Time variation of the density of a plane at a constant height in the corona – data from different heights in the corona show the same time evolution. The full line represents experiment Exp-1, and the dashed line Exp-2.
In the text

	Fig. 14 Top: mass flux in a horizontal plane. Bottom: mass flux in a vertical plane. The data is scaled to ± 10-7 kg m/s.
In the text

Fig. A.1 Vertical flow-pattern motion upward along magnetic field lines as a function of time. The base represents the magnetogram; the field lines are traced from a small pillow box around the active region in the photosphere. The field lines are coloured with the vertical flow speed, where red is downflow and purple is upflow. The amplitude of the velocity field is scaled to ±10 km s-1. The movie covers the time period from 83.3 to 325 min after the start.

In the text