Exploring the magnetic and thermal evolution of a coronal jet
Sushree S Nayak Samrat Sen Arpit Kumar Shrivastav R. Bhattacharyya P. S. Athiray




Sushree S. Nayak, Samrat Sen
0000-0002-4241-627X]Sushree S. Nayak
Center for Space Plasma & Aeronomic Research,
The University of Alabama in Huntsville,
Huntsville, Alabama 35899, USA
corresponding email ids: s.sangeetanayak93@gmail.com, samratseniitmadras@gmail.com

0000-0003-1546-381X]Samrat Sen
Instituto de Astrofísica de Canarias, 38205 La Laguna, Tenerife, Spain


0000-0001-9035-3245]Arpit Kumar Shrivastav
Aryabhatta Research Institute of Observational Sciences, Nainital, India-263002
Joint Astronomy Programme and Department of Physics, Indian Institute of Science, Bangalore 560012, India


0000-0003-4522-5070]R. Bhattacharyya
Udaipur Solar Observatory, Physical Research Laboratory,
Dewali, Bari Road, Udaipur 313001, India

0000-0002-4454-147X]P.S. Athiray
Center for Space Plasma & Aeronomic Research,
The University of Alabama in Huntsville,
Huntsville, Alabama 35899, USA
NASA Marshall Space Flight Center, ST13, Huntsville, AL 35812, USA



§ ABSTRACT 
Coronal jets are the captivating eruptions which are often found in the solar atmosphere, and primarily formed due to magnetic reconnection. Despite their short-lived nature and lower energy compared to many other eruptive events, e.g. flares and coronal mass ejections, they play an important role in heating the corona and accelerating charged particles. However, their generation in the ambience of non-standard flare regime is not fully understood, and warrant a deeper investigation, in terms of their onset, growth, eruption processes, and thermodynamic evolution. Toward this goal, this paper reports the results of a data-constrained three-dimensional (3D) magnetohydrodynamics (MHD) simulation of an eruptive jet; initialized with a Non-Force-Free-Field (NFFF) extrapolation and carried out in the spirit of Implicit Large Eddy Simulation (ILES). The simulation focuses on the magnetic and dynamical properties of the jet during its onset, and eruption phases, that occurred on February 5, 2015 in an active region NOAA AR12280, associated with a seemingly three-ribbon structure. In order to correlate its thermal evolution with computed energetics, the simulation results are compared with differential emission measurement (DEM) analysis in the vicinity of the jet. Importantly, this combined approach provides an insight to the  onset of reconnection in transients in terms of emission and the corresponding electric current profiles from MHD evolutions. The presented study captures the intricate topological dynamics, finds a close correspondence between the magnetic and thermal evolution in and around the jet location. Overall, it enriches the understanding of the thermal evolution due to MHD processes, which is one of the broader aspects to reveal the coronal heating problem. 









§ INTRODUCTION


Solar coronal jets are one of the fascinating eruptions from the solar atmosphere alongside flares and coronal mass ejections (CMEs). They are collimated plasma flows which are morphologically inverted-Y in shape, short-lived, and energetically 4-5 orders less than flares <cit.>. The study of jets is essential as they play an important role in heating the corona, and accelerating the solar wind <cit.>.

The generation mechanism behind these jets is primarily ascribed to magnetic reconnection which involves the reconfiguration of magnetic field topology while releasing stored magnetic energy into Ohmic heating, kinetic energy and accelerating the charged particles from the site of reconnection. Fully or partially open magnetic field lines near the coronal holes or at the boundary of active regions are observed as the favorable sites for ejecting jets while the triggering may be either through flux cancellation <cit.> or flux emergence <cit.>. In one of the important studies by <cit.>, a plausible scenario is reported for the eruption of coronal jets which are now also extensively applicable to many similar bursts of different energetic scales. According to these studies, a pre-existing minifilament escapes the canopy of surrounding open and closed field lines while reconnecting both at that site as well as at its own footpoints. The plasma material ejects out in that collimated channel and is dubbed as a spire carrying a swirling motion. These type of jets are popularly known as “blowout jets” <cit.>. In addition to that, jets are often associated with an initial bright point, also known as a jet bright point (JBP) <cit.> before achieving the fully grown structure.

Flare ribbons are inherent to solar eruptions. The well-known standard flare model a.k.a. the CSHKP model <cit.> pictures a solar eruption as a consequence of release of a magnetic flux rope above the polarity inversion line (PIL) from the solar atmosphere in two dimensions (2D). When the flux rope overcomes the tension force of overarching magnetic loops, that creates a current sheet below ensuing reconnection and producing two ribbons at the footpoints of the reconnected field lines. While this is the case for standard flares irrespective of their eruptive nature, there exists non-standard flares which produces different kind of ribbons namely circular/quasi-circular flare ribbons <cit.>, X-shape ribbons <cit.>, and three-ribbon types <cit.>.


These non-standard type ribbons or flares do not fit the argument of the standard flare model in explaining the flaring or eruptive process. In fact, in the state-of-art scenario, it is not fully understood about its onset, growth, and eruption processes, and merits attention to study. Besides, if they are associated with jet or CME, it is again challenging to comprehend the eruption process, formation of current sheets, or the reconnection affair in broad to capture in three-dimensional (3D) models. But, several attempts with an ensemble of multiple observations, and theoretical/numerical modeling are being made to understand aforementioned criticalities involved in the whole eruption process in these types of atypical events. Especially, when the responsible magnetic configuration in these events is analyzed to detect a potential reconnection triggering location, several potential configurations are reported as an initiation mechanism. For example, <cit.> has observed the formation of current sheets “due to stressing of spine” of the magnetic null point (ideally, |𝐁| =0) skeleton in their simulation which was associated with a quasi-circular flare ribbon. Likewise, <cit.> has explained the role of magnetic null point and quasi-separatrix layers (QSLs; locations bearing a sharp change in magnetic connectivities; <cit.>) in triggering a circular flare ribbon in the active region AR 12192 via magnetohydrodynamics (MHD) simulation initialized with a non-force-free-field (NFFF) extrapolation. Similarly, <cit.> discussed a three-ribbon flare produced due to reconnection along the coronal null line. In a different scenario of an X-shaped ribbon in <cit.>, a hyperbolic flux tube (HFT; intersection of two QSLs; <cit.>) was found in the non-linear-force-free-field (NLFFF) constructed magnetic topology where the dissipation of current sheets leads to ribbons.

With the above backdrop, in this paper, we focus on understanding the onset of an eruptive jet followed by an atypical flare ribbon type or majorly the reconnection process while analyzing the magnetic and thermal evolutions with a side-by-side comparison with observational dynamics. 
The novelty of this approach lies in its envisaged effectiveness in emphasizing the triggering mechanism, the dissipation of magnetic energy, and decay of current and its implications; on the emission near the jet region.


We study a particular event of jet eruption which occurred on February 5, 2015, in AR12280 as reported in <cit.>. The above study was mainly concentrated to investigate the oscillation of a filament that was hit by an erupting jet. However, the exploration of magnetic field topology and the other magnetic properties (e.g. field line twist, field strength, ohmic heating, etc.) are not explored due to the constraint of spectroscopic imaging. The role of magnetic properties becomes significant at low plasma-β in solar corona, and hence it warrants an importance to study the magnetic properties in and around the jet during the eruption using MHD simulations. It also motivates us to enrich the general understanding of the onset process of transients.



The above discussions emphasize the role of magnetic topology and the sites of reconnection in 3D for any transients. However, direct and routine measurements of the coronal magnetic field are still unavailable. Therefore, magnetic field extrapolation has become an alternative to obtain the coronal magnetic field topology using available photospheric magnetic fields. Myriad numerical attempts have been made to perform realistic modeling of coronal magnetic field such as potential source surface models <cit.>, several force-free models <cit.>, non-force-free models <cit.>, magneto-frictional approaches <cit.>, and magnetohydrostatic models <cit.> which includes both static and time-dependent solutions. For further details on features of different coronal models, we direct the readers to <cit.> and references therein. In this work, we have adapted the non-force-free-field extrapolation model devised by <cit.> to construct the coronal magnetic field topology. Then, we have initiated a data-constrained MHD simulation with the NFFF extrapolated magnetic field using the EULAG (Eulerian/Semi-Langrangian) MHD model <cit.> to track the eruption process.  In our study, we have given effort to deliver a collective analysis of the events using magnetic field extrapolation, MHD modeling, and comparison of the modeled parameters with the derived observables  near the jet, and ribbons.


The EULAG-MHD code is based on an incompressible regime, and not capable of estimating thermodynamic evolution directly. This motivates us to use Differential Emission Measurement (DEM) technique to estimate the thermal evolution of the jet region and its surroundings from observational data. This method is extensively used to evaluate the thermodynamic information for various solar observations in general. Also, particularly for the coronal jet studies, it has been used to explore the thermodynamic evolution of the spire, footprints, and source regions <cit.>.


The remainder of the paper is as follows. Section <ref> outlines the event of an active region jet eruption that occurred. In Section <ref>, we have reported the results obtained from DEM analysis of the event. Section <ref> explains the rationale behind the NFFF extrapolation model and discusses the modeled extrapolated field of AR12280. Section <ref> details the EULAG-MHD model, required set-ups, and simulated results. Lastly, Section <ref> summarizes the key findings from the work, highlights the novelty of the study,  and concludes how this work can be useful in a broader aspect to understand the solar atmosphere in future studies.




§ THE JET OBSERVED ON 2015 FEBRUARY 5 IN THE ACTIVE REGION AR12280



We have studied the jet in the active region AR12280 on 2015 February 5 peaked at ≈ 20:42 UT. In Figure <ref>, we have plotted the jet covering the whole active region in the field of view (FOV) in 304 Å pass-band of Atmospheric Imaging Assembly (AIA; <cit.>) onboard Solar Dynamics Observatory (SDO; <cit.>). Panel (a) corresponds to a time instance where no sign of jet activity is detected. Around ≈ 20:30 UT, we spotted an appearance of a bright point in the east side of the AR as marked by the white box in panel (b), which is referred to as JBP. This is a signature of the onset of the jet, which is developed afterward as shown in panel (c) around ≈ 20:46 UT. The jet base region is highlighted by the white box there. The legible three-ribbon structure can also be noticed. Panel (d) depicts the zoomed view of the jet highlighting its base region by the dotted white circle and the direction of the spire by the white arrow. An animation is provided showing the evolution of the jet in 304 channel. The animation starts at 20:24:07 UT and ends at 21:03:19 UT. We have plotted the B_z component of AR12280 at 20:24 UT obtained from Helioseismic Magnetic Imager (HMI; <cit.>) in Figure <ref>. Similarly, we have marked the jet base region in the black box and the bright point location in the white box respectively on the B_z map. The blue and red arrows indicate the direction of the transverse components on the B_z map plotted in gray scale. The green color contours are the location of polarity inversion lines (PILs) in the active region. The PILs are plotted using the routine developed by <cit.>. According to the routine, firstly, all of the pixels in the magnetogram are convolved with a gaussian kernel. Then, centering each pixel of B_z component, five consecutive pixels both in horizontal and vertical directions are scanned and maximum and minimum values are compared in both the directions of those five pixels. If they are found to be of opposite sign and magnitude of these two values exceeds the noise level of 60 Gauss (by trial and error method) on either side of the two arrays, the pixel in between them is located as PIL. In the right panel, we have plotted the magnetogram closer in time to the onset of the jet i.e. 20:36 UT as marked in Figure <ref>.











§ DEM ANALYSIS NEAR THE JET REGION


To understand the thermodynamic changes in the region of interests (ROIs) (see Figure <ref>: bright point region (white box in panel (b), jet base region (bigger white box in panel (c)), we determine differential emission measure (DEM) distributions using six optically thin emissions in AIA pass-bands. These channels have contribution from ionized states of iron with wavelengths 94  Å (Fe 10, Fe 18), 131 Å (Fe 8, Fe 21), 171 Å (Fe 9), 193 Å (Fe 12, Fe 24), 211 Å (Fe 14), and 335 Å (Fe 16). Their temperature response function peaks at log T[K] = (6.05, 6.85), (5.6, 7.05), 5.85, (6.2, 7.25), 6.3, and 6.45 respectively <cit.>. AIA channels are sensitive to a wide range of plasma temperatures. It is also established that the bi-modal temperature response of the hot AIA channels (94 Å , 131 Å , 193 Å ) can introduce systematic overestimation of emission measure at high temperatures <cit.>. We perform DEM analysis using 16 temperature bins from 5.4 ≤ log T ≤ 7.0, with δ log T = 0.1.

We have used the latest version of sparse inversion code <cit.> to obtain DEM, which accounts for the amount of plasma emission in particular temperature intervals integrated along our line of sight.



    DEM (T) =  n_e^2dl/dT ,

where l is the path length along the line of sight and n_e denotes the electron density.

The total emission measure can be calculated from DEM(T),

    EM = ∫ DEM(T) × dT   .

 The DEM(T) in every temperature bin can be utilized to obtain emission measure weighted temperature as,

    T_EM = ∫ T DEM(T) dT/∫ DEM(T) dT  .

 We have plotted the evolution of the EM weighted temperature extracted through DEM analysis in Figure <ref> for the whole active region. Panel (a) is the temperature distribution on the whole region. The bright point location is marked in panel (b). Later in panel (c), it can be seen that the hotter regions are exhibiting near the location of filament eruption. For further inspection, in Figure <ref>, we plotted the emission map in logarithmic scale in panel (a) highlighting the ROIs in black and red boxes.
 In panel (b), we have plotted the average EM weighted temperature inside those two boxes till the decay phase of transients.  
Important is the higher and a sharp ascent in the peak of the temperature near the bright point region than the jet base.










§ INITIAL MAGNETIC CONFIGURATION OF AR12280





 §.§ The Non-Force-Free-Field Extrapolation Model


To investigate the onset of the jet eruption and the formation of multiple ribbon structures, we first obtained the initial magnetic field topology of the AR12280 using the Non-Force-Free-Field (NFFF) extrapolation model <cit.>. The NFFF model is theorized based on the minimum dissipation rate principle (MDR) where the idea uses total dissipation rate as a minimizer keeping the generalized helicity constant for a two-fluid description of plasma <cit.>. The relaxed state carries a non-zero Lorentz force which is used to drive the plasma in the MHD simulation later in our study. The extrapolation solves an inhomogeneous double-curl
Beltrami equation for the magnetic field 𝐁 (<cit.>, and references therein),


    ∇×∇×𝐁 + a_1∇×𝐁 + b_1𝐁 = ∇ψ,

where a_1 and b_1 are constants. The solenoidality of 𝐁 enforces
the scalar function ψ to obey Laplace's equation. The modified
vector 𝐁' = 𝐁 - ∇ψ satisfies the corresponding homogeneous
equation, which represents a two-fluid MHD steady state <cit.> having a solution <cit.>

    𝐁' = ∑_i=1,2𝐁_i,

where each 𝐁 is a linear force-free field satisfying

    ∇×𝐁_i = α_i𝐁_i,

in usual notations. The two sets of constants are related by
a_1 = - (α_1+α_2) and b_1 = α_1α_2. The magnetic field is then

    𝐁 = 𝐁_3 + ∑_i=1,2𝐁_i,

where 𝐁_3 = ∇ψ is a potential field.
We apply the technique described in <cit.>. Briefly, an optimal pair of α is
computed for 𝐁_3 = 0 by minimizing the average normalized
deviation of the magnetogram transverse field B_t from its
extrapolated value b_t, given by

    E_n =  ( ∑_i=1^M |𝐁_t,i - 𝐛_t,i| × |𝐁_t,i|) /(∑_i=1^M |𝐁_t,i|^2 ),

where M = N × N is the total number of grid points on the
transverse plane. Additional minimization of E_n is done by
using 𝐁_3 = ∇ψ as a corrector field for the obtained pair of α's.
Noting a superposition of potential fields results in potential
field only, the above procedure is iteratively repeated until 
E_n plotted against the number of iterations asymptotically
saturates to a minimum. Also
important is to recognize that α_1 and α_2 , alone or combined,
are not the only parameters to determine field line twist <cit.>. With the
twist τ being related to the field-aligned current density

    τ = 𝐉·𝐁/|𝐁|^2 = (α_1𝐁_1 + α_2𝐁_2) · (𝐁_1 + 𝐁_2 + 𝐁_3)/|𝐁|^2,

additional to any modifications in α_1 and α_2, τ can also vary
because of changes in the component fields—  including the
potential field 𝐁_3, an advantage of the NFFF extrapolations.
We note that the values of α_1 and α_2 must be bounded above to
ensure a monotonous decay of 𝐁 with height <cit.>. If required, extra twist can be accommodated in the
extrapolated 𝐁 by varying 𝐁_1, 𝐁_2 and 𝐁_3 to better match the
magnetogram while not exceeding the maximal α_i's. 


The justification of using the NFFF model lies in the following analysis. As per <cit.> on his study of plasma beta, β = p/p_mag (where, p and p_mag represent the plasma and magnetic pressure respectively) variation over the active regions, the force-free approximations on the photosphere does not hold true. There in Figure-3, only the mid corona, a sandwiched region between chromosphere and the upper corona/solar wind accelerated region is shown to have a β value of less than unity. One of the important arguments for the non-force-free-field approximation lies in the use of magnetic field. <cit.> have shown that the magnetic field is not force-free at the photosphere while becomes the same at a height around 400 km above the photosphere. <cit.> have studied the nature of the photopshere in 12 flare producing active regions and found it not to be much deviant from force-free nature. However this may not hold true always if there are twisted structures. In <cit.>, they also found similar arguments but in umbral or inner penumbral regions using high spatial resolution magnetograms from the Solar Optical Telescope/Spectro-Polarimeter on
board Hinode. Besides, for the available force-free extrapolation models based on optimization and MHD relaxation approaches, there can remain a finite residual Lorentz force during the minimization to an equilibrium state. As a remedy to this issue, a preprocessing technique is adapted in the algorithm where the photospheric magnetic field is adjusted. Whereas, the NFFF algorithm operates on observed magnetograms without any apriori preprocessing keeping the original nature of force at the photosphere. Again, the NFFF model does not claim to provide the whole active region to be non-force free. Furthermore, one rationale can be provided in terms of the strength of the Lorentz force, which depends on the magnetic field strength as well as the associated current density (which depends on the spatial gradient of the magnetic field components). Also, the orientation between B and J contribute to strength of the Lorentz force. Therefore, we may not always conclude the force-free or non-force-freeness of a magnetic field configuration based only on the magnetic field strength.
  



 §.§ The NFFF Extrapolated Magnetic Field Topology

To extrapolate the magnetic field over AR12280, we have used the s magnetogram series of HMI. The data in this series consists of magnetic field strength (B), azimuth, and inclination in every 720 seconds with a full disk FOV.  We have reproduced an HMI SHARP-like active region patch from B, azimuth, and inclination using the  available in the  package <cit.> with a slight modification in the coordinate system to it. The objective is to cover a field of view that encloses the base part of the jet sufficiently for the extrapolation while maintaining a flux balance over the cutout. The computational box has 768 uniform grid points or ≈ 276 Mm of physical extent in the x-direction and 512 grid points or ≈ 184 Mm in both y- and z- directions. The E_n has a value of ≈ .36  which amounts to ≈ 36% error in the reconstruction of the transverse field whereas the B_Z is the same observed magnetic field as stated in the Section <ref>. The Figure <ref> depicts the variation of magnetic field (|B|), current density (|J|), and Lorentz force (|J×B|) over the geometrical height normalized to their maximum value in logarithmic scale over the whole computational domain. We have also calculated β inside the jet base region for a pixel having the maximum total B of value ≈2790 G and found to be ≈ .44 at the photospheric level and for a field strength of ≈ 500 G, we found it to be 14. In conclusion, the active region does have a range of β and may not hold the force-freeness condition everywhere. Importantly, the structures in our extrapolated field lies over the PILs where a sharp gradient of the magnetic field connectivity is observed.




In Figure <ref>, we have plotted the magnetic field lines using the Visualisation and Analysis Platform for Ocean, Atmosphere, and Solar Researchers (VAPOR) software <cit.>. The field lines are plotted employing a highly accurate Field-line integration function of VAPOR,  which relies on the adaptive line integration with a fourth-order Runge–Kutta scheme and tri-linear interpolation of field values over cells in the grid; described in Section 2.3.3 of <cit.>. In the extrapolated field, we find a flux rope near the jet base region, plotted in cyan color in the panel (a)  of Figure <ref>, also highlighted in the chartreuse color box. To confirm the structure, we have plotted a direct volume rendering of the twist parameter (T_w) near the flux rope. For this purpose, we have used the code by <cit.>, also available at (http://staff.ustc.edu.cn/ rliu/qfactor.html). The  twist value is measurement of two infintesimally close field lines winding about each other <cit.> and can be casted in the form according to <cit.>;

    T_w = ∫_L𝐉·𝐁/|𝐁|^2 dl.

This flux rope here is a realization of the mini-filament found in case of jets as reported in <cit.>. We will discuss further dynamics in the MHD evolution in Section <ref>. Near the flux rope, we find the maximum twist to be of ≈ 1.33, shown in panel (c). The flux rope is found to be lying close to the bottom boundary. We have provided an inset there highlighting the location of the PIL, indicated by a magenta arrow, where the flux rope is detected. We find a bald patch type topology to the east of the flux rope, highlighted in the black box in panel (a). In panel (b), we have presented the 171 Å channel for the same field of view and at same time instance as panel (a) to mark the resemblance of the overall topology of the active region. We then compute the squashing factor or the Q-value near these regions utilizing the routine of <cit.>. We find the presence of QSLs in the vicinity of flux rope and bald patch. The footpoints of the bald patch are found to trace the high Q-value regions which is evident from the panel (d). The inset in panel (d) shows a flipped side view of the orientation of field lines in the bald patch.









§ FRAMEWORK OF THE EULAG-MHD MODEL AND THE NUMERICAL SET-UP FOR MHD EVOLUTION



To track the evolution of the field lines near the jet base region and its effect on the eruption process, we have performed an MHD simulation using the EULAG-MHD model. The model solves the incompressible Navier-Stokes MHD equations under the assumption of thermal homogeneity (or thermally inactive) and perfect electrical conductivity <cit.>. The governing equations are as follows:


    ∂𝐯/∂ t 
    + (𝐯·∇)  𝐯 =-∇ p
    +(∇×𝐁) ×𝐁+τ_a/τ_ν∇^2𝐯,
          ∇·𝐯=0, 
          ∂𝐁/∂ t=∇×(𝐯×𝐁), 
          ∇·𝐁=0,

written in usual notations. The variables in the MHD equations are normalized as follows 

    𝐁⟶𝐁/B_0,  𝐯⟶𝐯/v_a,  
     L ⟶L/L_0,   t ⟶t/τ_a,  
     p  ⟶p/ρv_a^2.

The constants B_0 and L_0 are generally arbitrary, but they can be fixed using the average magnetic field strength and size of the system. Here, v_a≡ B_0/√(4πρ_0) is the Alfvén speed and ρ_0 is the constant mass density. The constants τ_a and τ_ν represent the Alfvénic transit time (τ_a=L_0/v_a) and viscous dissipation time scale (τ_ν= L_0^2/ν), respectively, with ν being the kinematic viscosity. Utilizing the discretized incompressibility constraint, the pressure perturbation, denoted by p,  satisfies an elliptic boundary-value problem on the discrete integral form of the momentum equation (Equation <ref>) <cit.> and the references within. 


Here we discuss only essential features of the EULAG-MHD and refer readers to <cit.> and references therein. The model is based on the spatio-temporally second-order accurate non-oscillatory forward-in-time multidimensional positive definite advection transport algorithm MPDATA <cit.>.  
Importantly, MPDATA has proven the dissipative property which, intermittently and adaptively, regularizes the under-resolved scales by simulating magnetic reconnections
and mimicking the action of explicit subgrid-scale turbulence models <cit.> in the spirit of
Implicit Large Eddy Simulations (ILES)
<cit.> scheme. Arguably, the residual numerical dissipation is then negligible
everywhere but at the sites of MRs. Moreover, this dissipation being intermittent in time and space, a quantification of it is meaningful only in
the spectral space where, analogous to the eddy viscosity of
explicit subgrid-scale models for turbulent flows, it only acts on the shortest modes admissible on the grid; in particular, in the vicinity of steep gradients in
simulated fields. Such ILESs conducted with the model have already been successfully utilized to simulate reconnections to understand their role in the coronal dynamics <cit.>. In this work, the presented computations continue to rely on the effectiveness of ILES in regularizing the under-resolved scales by the commencement of magnetic reconnections.


The initial magnetic field is supplemented from the NFFF extrapolation and the initial velocity field is set to 𝐯=0. The lateral boundaries (x, and y) are kept open aiming that the net magnetic flux should be conserved. While at the bottom boundary, the z-components of 𝐁 and  𝐯 are chosen to be fixed to their initial values (also termed as line-tied boundary condition) as the flux change during the transient activities is found to be minimum. The top boundary follows the same condition as bottom only exception of it not being fixed throughout the evolution. Also, except the bottom boundary, all variables are calculated by linearly extrapolating from the immediate neighborhood cell values. Notably, the field and the corresponding Lorentz force values at such height become extremely small compared to the counterparts at the lower boundaries <cit.>. As stated earlier in Section <ref>, the simulation is initially driven by the non-zero Lorentz force associated with the extrapolated magnetic field, and the flow is  primarily generated by it. The resulting flow is however made incompressible following the Equation <ref>, an assumption also adapted by <cit.>. Since our focus is to understand the onset of the jet through the topological changes, the assumption seems to be justifiable in the tenuous coronal medium. One reminder is that the density is set to unity and being constant over both space and time. The computational box extension is of the same size as that of the extrapolation, but with ≈ 2 times reduced scaling of δ x. The spatial unit length δ x is .0052 and the time step δ t is set to 1× 10^-3 while satisfying the CFL condition <cit.>. The dimensionless coefficient or the kinematic viscosity τ_a/τ_ν≈ 2 × 10^-4 in the simulation is roughly ≈15 times larger than its coronal value <cit.>. To note, the parameter τ_a/τ_ν is controlled by the spatial resolution and the time step  whilst satisfying the von Neumann stability criteria <cit.>. The larger τ_a/τ_ν, however, only expedites the evolution without an effect on the corresponding changes in the magnetic topology <cit.>. The total simulation time is 5400δ t. To compare with the observational time, we multiplied the total simulation time by 15. Then, it corresponds to a ≈ 40.5 minutes of the observational period.



§ RESULTS AND DISCUSSIONS




 §.§ Dynamics in the field line topologies near all ROIs

We plotted the overall evolution of magnetic topologies near the periphery of the jet and ribbons in Figure <ref> vis-à-vis with the observed dynamics extracted from 131 Å channel of AIA for the whole simulation period. The animation associated with Figure <ref> is available online, that belongs from t = 20:24 to 21:03 UT (rendered with 13 s in the video). The timestamps mentioned on the 131 Å channel are near-cotemporal to the simulation snapshots. The first column with panels (a), (c), (e), (g), and (i) represents the evolution of field lines neighboring the jet and ribbons whereas the second column with panels (b), (d), (f), (h), and (j) are snapshots of progress of transients in 131 Å channel. Initially, the flux rope (yellow color) starts to untwist due to the action of Lorentz force. The ambient loops in red color close to the flux rope are seen to rise toward the eastern part of the active region along with the flux rope. A close correspondence between the loop dynamics and 131 channel can be noticed in panel (e) and (f) where the flux rope opens up ejecting the materials outward, similar to 131 channel. Afterward, the nearly potential yellow color loops in panel (g) and (h) show the post-reconnected loops. According to the standard picture of a jet in <cit.>, there are two types of reconnections takes place in a jetting process, the external reconnection which creates a passage for the flux rope to erupt and the internal reconnection involves the reconnection at its footpoints. Here, the escape of the materials is captured while the internal reconnection not so clearly may be owing to lack of spatial resolution and the location of the flux rope being close to the boundary. Although, we see some kinks in the flux rope, as marked by a circle in panel (f) of Figure <ref>. However, it is important to note that the simulation is magnetically driven whereas in reality, the jets could be due to the combined effect of both magnetically and thermal changes in the surrounding.



Elaborating the discussion on the untwist of the flux rope and the surrounding arcade further, an exclusive evolution of the flux rope is shown in Figure <ref> with a comparison of the filament evolution, overplotted with 304 Å channel on the background till the formation of ribbons. The untwisting of the flux rope (pink color) and footpoints of subsequently formed less potential loops which is also evident from the variation of mean of total current density (|𝐉|)  over the simulation period plotted in Figure <ref>. They are also seen to follow the ribbons in the later phase of the simulation as shown in the panel (c) of Figure <ref>. The directions of Lorentz force and the flow vector are plotted in panels (b), (d), (f) of Figure <ref> and denoted in blue and yellow color arrows respectively. Absence of the plasma flow in the first panel is due to the setting to zero in the first time step in the MHD model as it was purely driven by the Lorentz force. Initially downward and inclined Lorentz force pushes the toroidal fluxes leading to untwisting of the flux rope. To focus it more, we investigate the orientation of the Lorentz force and flow vector at a fixed location within the flux rope as shown in Figure <ref>  for different stages of the flux rope evolution. Similar to Figure <ref>, the blue arrow depicts the Lorentz force and the pink arrow depicts the flow vector. We appreciate the evolution of the directionality of the Lorentz force (at a fixed spatial location) more clearly as shown in Figure <ref>(a) to (d). This shows that the force remains nearly parallel to the flux rope axis during the untwisting phase of field lines. We also notice that the flow vector remains aligned nearly parallel to the flux rope axis and normal to the magnetic field lines at that fixed location during the untwisting phase of the field lines of the flux rope.  These consistent directionalities during the untwisting phase of the evolution represent that the force and velocity flows in the vicinity of the flux rope play the role in the untwisting of the flux rope field lines. We also estimate the spatial average of twist density (|T_w|) over a selected region that contains the flux rope as shown in Figure <ref>.  This is also supporting the trend of twist parameter evolution in Figure <ref> in Section <ref>.
Further supplementing the argument quantitatively, we have plotted the angles between the flow vector (v),and magnetic field (B) in solid line, and the Lorentz force (J × B) in dashed line respectively in Figure <ref>. The angles are computed as cos^-1(v·B)/(|v)||B|) and cos^-1(v·J×B)/(|v)||J×B|). In both cases, the angles are approaching 90^∘ toward the end of the evolution which is in congruent with the topological evolution seen in Figure <ref>. This simultaneously explains the loss of twist and the apparent rise of the field lines in the flux rope as well as the overarching loops there.










In Figure <ref>, we have plotted the evolution of the bald patch with 304 Å channel in the background in panels (a)-(d). There, yellow and blue contours in every panel represent positive and negative polarities respectively, and highlight the footpoint connectivities of the loops of the bald patch. In the background, we have plotted the squashing factor or Q-value to trace the QSL dynamics in panels (a1)-(d1) for the same time stamps as the first row. In panel (b), we have marked the initial bright point location with the green  color box. The continuous reconnections near this bald patch  along with flux rope eruption in the simulation are suspected to produce the apparent three ribbons later. Also, the footpoint evolutions in a slipping fashion follow the bright point as well as the ribbons as they expand later. Interestingly similar to our case, <cit.> in their work of a non-typical flare, have found a co-spatiality of flare ribbons and bald patch near the PIL from observations. In their numerical modeling, they found a remarkable agreement between the footpoints of bald patches to the ribbons. In our simulation, albeit the initial set-up of the evolution, the reconnection near the bald patch facilitates the formation of the multiple ribbons which eventually contributed to the jet base region alongside the untwining of flux rope. However, attributing the eruption of this flux rope from the jet base region to the consequent oscillation in the filament on the east side of the active region, as reported in <cit.>, is beyond the scope of the simulation presented here. Together, Figure <ref> and Figure <ref> show the onset of the jet eruption and the formation of the ribbons.





 §.§ More energetics during the evolution near all the ROIs


Topological evolution aside, the energetics near the transients are other hints of the energy release mechanism. Hence, in Figure <ref> in panel (a), we have plotted the time variation of volume averaged magnetic energy density (MED; (B^2/8π)) and kinetic energy density (KED; v^2/2,  also ρ =1)  throughout the simulation period for the whole of the computational box. Important is the trend of decrease in the magnetic energy while the kinetic energy is sharply increasing at the expense of it during the first ≈ 4 minutes of the evolution, decreases between ≈ 5-8 minutes, and then it increases till 10 minutes. Later, both are evolving toward an almost quasi-steady state. Also, notable is that the initial value of the KED is zero as suggested from the numerical setup of the model in the Section <ref>. In panel (b), we have focused on the base of the jet. The sharp change in the magnetic energy density can be noticed there, ≈ 13% decline to its original strength in the 15 minutes of evolution time scale. However, within the same period of evolution, kinetic energy density is increasing with a rate of ≈ 16% from a motionless state while showing a fall in between ≈ 4-8 minutes due to viscous dissipation. Further, in the panel (c), we noticed an interesting variation near the neighborhood of the bright point and observed remarkable changes in both MED and KED, where the transients are triggered. The overall decrease in MED during first ≈ 10 minutes is due to the continuous reconnections near the bald patch. This dissipation may be contributing to the increase in KED.  Noteworthy is the peak in the KED after ≈ 12 minutes of the start time. The early peak in KED might be responsible for enhancing the twist in the field lines near the bright point as well as in the neighborhood of the flux rope adding more flux to it, which can be seen from panel (b) of Figure <ref>. Subsequently, the increasing twist in the field lines assists in increasing the magnetic energy near those locations, which is discernible from the peak of magnetic energy t≈ 17 minutes in Figure <ref>(c) after ≈ 5 minutes to that of kinetic energy. Important is also the continuous reconnection near the bald patch which facilitates the dissipation of magnetic energy along with the triggering of the jet. The KED over every region is dissipated by the action of viscous relaxation. The feedback sharing between the MED and KED  is significant in the jet region and more in the bright point region. This suggests a higher rate of dissipation of magnetic energy and conversion to kinetic energy near the locality of the jet base region and the bright point location. Here, we also draw the attention of the reader toward the efficacy of the ILES scheme of EULAG-MHD model being majorly effective at the locations of primary reconnection sites, whereas not contributing to the other part of the computational box.




Next in Figure <ref>, we have plotted the profile of average of Lorentz force density (|𝐉× 𝐁|) for all ROIs. The profiles are denoted in solid line for the whole of AR, dotted line for the jet region, and the dashed line for the bright point. We observe an initial fall till the first couple of minutes as the |𝐉× 𝐁| is expended on the generation of the plasma flow which can be seen in the trend of KED in panel (a) of Figure <ref>. A strikingly rise |𝐉× 𝐁| is seen near the jet region. This might be due to the effect of prior increase of KED which pushes the field lines into the jet region while injecting more twist to flux rope as well as the surrounded arcades. But, it then decreases rapidly due to the unwinding of the flux rope and dissipation near the bright point. The trend of (|𝐉× 𝐁|) is however relatively similar to that of the complete AR till ≈ 17 minutes the evolution. The continuous reconnection near the bald patch may inhibit the sharp development of Lorentz force near the bright point. Furthermore, in Figure <ref>, we have plotted the variation of mean of twist parameter for the three ROIs. Similar to Figure <ref>, we have maintained the same legends for three ROIs. The twist in the whole region is decreasing throughout the simulation period whereas, in the case of the jet base region, it increases till ≈ 8 minutes of the total time, then decreases and maintains almost constancy till the end of the evolution.







Next, Figure <ref> shows the evolution of mean of free energy again for the same ROIs keeping the same legends as the above Figure. We have estimated the free energy in the following way;



    E_free =  1/V∫_V{(B^2/8π)_Sim|_t - (B^2/8π)_Pot|_t=0} dV,

where, V is the total volume of our simulation box. The first term on the R.H.S of Equation <ref> refers to the non-potential energy density at each time step. The second term refers to the potential energy density calculated at the initial time step, as the potential magnetic field is the minimum energy state <cit.>. The profiles for the whole AR and jet region are seen to decrease till the end of the evolution. However, for the bright point it shows an increement  between ≈ 7-12 minutes of the simulation period.





To investigate the energetics further near these locations during the onset and afterward, in Figure <ref> we have provided the evolution of volume averaged total current density (|𝐉|), z-component of 𝐉, |𝐉|/|𝐁|, and the J_trans or √(J_x^2+J_y^2) for all the locations stated in the above paragraphs. In panel (a), all the parameters show a monotonous decrease. Near the jet base region and the bright point location, in panel (b) and (c), J and J_z show a fall during the initial few minutes, then start picking up due to twisting of the field lines. Interestingly, J_z shows an increase in the bright point neighborhood during the peak of the eruption phase which was also the case for magnetic energy density in panel (c) of Figure <ref>. Interestingly, we observe similar growth in the average temperature over the complete activity period in the DEM averaged temperature map, in panel (c) of Figure <ref>. This survives till the start of the decay phase of the eruption. After ≈ 30 minutes of evolution, both are decreasing in a monotonous fashion. However, we fail to capture any abrupt changes in |𝐉|/|𝐁| during the evolution in any of the regions, particularly near the jet base region. The reason may be partially due to the lower spatial resolutions of computation.  However, locations of current sheet  are not found in the extrapolation resolution scales as well which was of nearly twice of the simulation. Another overt reason can be accredited to the boundary condition employed at the bottom boundary as line-tying.











§ SUMMARY



In this work, we investigate the triggering of an active region jet by magnetic field extrapolation, and its subsequent magnetohydrodynamics evolution through numerical simulation. We have also analyzed the temperature evolution based on the emission measurement from observation near the transients. Then we compare the temporal evolution of temperature (derived from observation) with the heating implication based on MHD simulation. Besides, we have explored the formation of three ribbons on the base of the jet which are atypical.

First, we used the non-force-free-field magnetic field extrapolation to obtain the initial magnetic field topology of the active region. The extrapolated field provides a non-zero Lorentz force in the region at least to a certain height in the lower atmosphere, afterward achieving a nearly force-free state in the upper atmosphere of the computational box. The top boundary of the computational box reaches ≈ 184 Mm height of the atmosphere. The overall magnetic topology agrees well with the loop morphology from the observation. Near the jet base region, we found a bald patch triggering the jetting process and multiple ribbons. Additionally, a flux rope, though low-lying,  was found in the vicinity of the bright point.

To understand the onset of the jet and the formation of the three ribbons, we utilized the NFFF extrapolated field as an initial state to a data-constrained EULAG-MHD simulation with a line-tying boundary condition on the bottom boundary covering the period of jet eruption. The findings from the simulation are summarized below.

 
  * The study provides topological evolution of magnetic field in the neighborhood of the jet, particularly field line evolutions near the flux rope and the bald patch. The untwisting of field lines releases the jet materials along the direction of the plasma flows and ultimately contributes toward the eruption of the jet. The onset of the jet starts due to reconnection in the bald patch which is different from a traditional jet onset mechanism.  
    
  * The simulation simultaneously focuses on the formation of multiple flare ribbons where the same reconnection near the bald patch and the eruption of the flux rope generate multiple ribbons due to the transfer of energetic particles to the lower atmosphere. The footpoint evolution of these topologies are in good agreement with the observed ribbon locations.
    
  * Further discussions on the variation of the magnetic and kinetic energies on different parts of the jet not only elucidate the triggering of the jetting process but also validate the credibility of the simulation performed. A detailed focus on the variation of energy densities and the active driving of Lorentz force and modulation by the plasma flow over time near both the jet base and bright point regions explain the relaxation of magnetic energy and kinetic energy in facilitating the triggering and eruption process.  These parametric study highlight the importance of the onset region in order to understand magnetic reconnection. 
    
  * Additional highlights are on the profile of the mean of total current and current components which shed light on the current dissipation near those sites. A qualitative comparison with dynamics from the DEM analysis also falls in the line of argument of onset of the jet. We observe a hotter base and cooler spire in the DEM analysis for the jet. The temperature and emission profile near the bright point from the observation, the topological evolution, and the profiles of the simulated energetics comply with the onset process of the jet remarkably. This motivates further to focus on both the magnetic and thermal structure of a jet in order to understand the energy transfer through the jet and as well as in the atmosphere globally.  We noticed a key aspect from the evolution of the current profiles is that the onset point is a maximum current carrying site than the overall transient affected area. However, the observational and/or numerical modeling artifacts are not ignored. Thus, directly quantifying the amount of dissipation into the ohmic heating in these locations may not be so accurate.
	
Yet, the presented simulation could not (a) capture the formation of the current sheet as similar to the 2D model proposed by <cit.>, and as discussed in <cit.>. This  may be due to the resolution limit and the boundary conditions adapted in our simulation, and (b) delineate the complete thermodynamics of the jet region due to utilization of an incompressible and thermally inactive model. The jets are ubiquitous and impact not only the local atmosphere but also influence the solar wind generation based on their source of origins. Hence, in our future work, we aim to attempt a more realistic simulation covering all these shortcomings. We will also approach a data-driven boundary condition to the simulation to understand the impact on regularization of both the onset and eruption process of such events, particularly focusing on the genesis of these multiple ribbons.

1cm

AIA/SDO, HMI/SDO
Solarsoft <cit.>, VAPOR/NCAR

We thank the anonymous referee for the insightful comments and suggestions which improved the manuscript considerably. The extrapolation and simulations are performed in the Bladerunner cluster located in the Center for Space Plasma and Aeronomic Research department of the University of Alabama in Huntsville. We acknowledge the use of the visualization software VAPOR (www.vapor.ucar.edu) for generating relevant graphics. Data and images are courtesy of NASA/SDO and the HMI and AIA science teams. SDO/HMI is a joint effort of many teams and individuals to whom we are greatly indebted for providing the data. We thank Ranadeep Sarkar for sharing the PIL calculation code. S.S.N. acknowledges NSF-AGS-1954503, NASA-LWS-80NSSC21K0003, and 80NSSC21K1671 grants. S.S. acknowledges support by the European Research Council through the Synergy Grant #810218 (“The Whole Sun”, ERC-2018-SyG). A.K.S. is supported by funds of the Council of Scientific & Industrial Research (CSIR), India, under file no. 09/079(2872)/2021-EMR-I. We are grateful to Qiang Hu for his valuable discussions on development of the manuscript. 





