Open and closed magnetic configurations of twisted flux tubes
Samrat Sen A. Mangalam





A. Mangalam
mangalam@iiap.res.in

0000-0003-1546-381X]Samrat Sen
Indian Institute of Astrophysics, Sarjapur road, Koramangala 2nd Block, Bangalore-560034, India
samrat@iiap.res.in

0000-0001-9282-0011]A. Mangalam
Indian Institute of Astrophysics, Sarjapur road, Koramangala 2nd Block, Bangalore-560034, India
aasjournal


We construct two classes of magnetohydrostatic (MHS) equilibria for an axisymmetric vertical flux tube spanning from the photosphere to the lower part of the transition region within a realistic stratified solar atmosphere subject to solar gravity. We assume a general quadratic expression  of the magnetic flux function for the gas pressure and poloidal current and solve the Grad-Shafranov equation analytically. The solution is a combination of a homogeneous and a particular part where the former is separable by a Coulomb function in r and exponential in z, while the particular part is an open configuration that has no z dependence.  We also present another open field solution by using a self-similar formulation with two different profile functions and incorporating stratified solar gravity to maintain the magnetohydrostatic equilibria, which is a modification of earlier self-similar models with a twist. We study the admitted parameter space that is consistent with the conditions in the solar atmosphere and derive magnetic and the thermodynamic structures inside the flux tube that are reasonably consistent with the photospheric magnetic bright points (MBPs) for both open and closed field Coulomb function and self-similar models as estimated from observations and simulations. The obtained open and closed field flux tube solutions can be used as the background conditions for the numerical simulations for the study of the wave propagation through the flux tubes. The solutions can also be used to construct realistic magnetic canopies. 







§ INTRODUCTION

The small scale magnetic structure in the solar photosphere plays an important role in several phenomena like the evolution of active regions <cit.>, heating of corona through the dissipation of waves <cit.> and reconnection between flux tubes <cit.>. Magnetic flux tubes span from the  photosphere to the higher atmosphere and are observed in the form of small scale magnetic structures. The topological rearrangement of these flux tubes due to the motion of the photospheric footpoints gives rise to the magnetic reconnection leading to the energy release in the solar corona <cit.>. Therefore, the modeling of the proverbial flux tube is one of the key aspects to understand various phenomena on the solar surface and its outer atmosphere.

Several attempts have been made earlier to construct the model of flux tubes for both twisted and untwisted magnetic fields. <cit.> studied a two dimensional (2D) axisymmetric flux tube model without twist for sunspots  using self–similar structure, where a self-similar parameter was defined as a combination of r and z, and the relative vertical magnetic field strength at any arbitrary point w.r.t. the magnetic field strength at the axis is scaled with a Gaussian profile function of the self-similar parameter. This model is valid for open field lines where the magnetic lines of force rise from a horizontal plane and do not return in the model domain. <cit.> implemented a twist in the self–similar structure to model the  sunspots. In this model, an empirical form of the azimuthal magnetic field strength B_ϕ(r,z) was taken from the data obtained from observations <cit.>. By solving for the  variation of the pitch angle and  gradient of the pitch angle, the thermodynamic quantities with the depth were calculated. Motivated by the model and the self–similar structure proposed by <cit.>, <cit.> assumed a quadratic form of the flux function for the gas pressure to model a closed field flux tube, where the magnetic lines of force rise and return to the same horizontal plane.

<cit.> have numerically studied a 2D model of open single flux tube with a twist using the standard boundary conditions including sheet current to study the magnetic field line structure within and outside the flux tube. The magnetic and thermodynamic structure for both single and multiple flux tubes which span from  photosphere to corona have been studied for the case of untwisted magnetic field <cit.>,  where an empirical form of the magnetic field components is motivated by a  self-similar construction. A numerical model of flux tubes has been studied by <cit.>, where an empirical form of magnetic flux function has been assumed; this was followed by a model to study the propagation of the MHD waves through the flux tubes with an azimuthal velocity perturbation. The steady structure of the 2D flux tube was used as a background initial condition to study the propagation of the MHD waves. For example, <cit.> assumed an empirical form of gas pressure for investigating the wave propagation and energy transport through the flux tube. Other interesting results of wave behavior in the solar atmosphere have been presented by several authors. <cit.> have studied the propagation of the acoustic wave through the solar atmosphere due to the periodic drivers at the photosphere, and <cit.> have modeled the wave propagation through the photospheric magnetic bright points (MBPs).

In this work, we have constructed two different models of flux tubes with twisted magnetic field for open and closed field lines by solving Grad–Shafranov equation (GSE) <cit.>. Here, we have assumed a quadratic form of the flux function for the gas pressure and poloidal current which has been used to study the equilibrium solution of terrestrial plasma <cit.>, and we extend it to solar flux tubes. As the MHD waves follow the magnetic field lines, it is important to model flux tubes with open field lines, so that MHD waves propagate through the flux tube and dissipate in the upper atmosphere, which is a key aspect of the coronal heating. A key aspect of this paper is to show that the closed field model, reported in <cit.> (SM18 hereafter), is a special case of the open field model with a twisted field line. The flux tube we build is axisymmetric in structure and spans vertically upward from photosphere to the transition region. The case of a linear form of the flux function for the gas pressure and poloidal current, an equilibrium solution near the magnetic axis of a plasma torus has been reported by <cit.>. SM18 have studied the homogeneous solution of GSE which is a special case of the general solution of the quadratic case to model a flux tube with closed field lines with a twist. Here, we present the full solution of the GSE including both homogeneous and the particular parts to model a twisted open field flux tube. The other model we have built is a self–similar magnetic structure with twist, with a generalized Gaussian (or power law) incorporated into the magnetic shape functions; the gas pressure and poloidal current  are taken to be quadratic functions of the flux function. The self–similar flux tube model expands with height which spans from the photosphere to transition region. After building the solutions semi–analytically and applying appropriate boundary conditions (BCs), we calculate the magnetic field structure and thermodynamic quantities inside the flux tube. As magnetic bright points (MBPs) observed in the photosphere <cit.> are likely to be flux tubes, we compare our model with the existing observations and simulations of MBPs.

The paper is organized as follows. In <ref>, we apply the GSE to the cylindrical flux tube case and describe the common BCs which are physically realistic and used in modeling of our flux tubes. In <ref>, we present the Coulomb function model for open and closed fields, the appropriate BCs, and show how the open field Coulomb model generalizes the Coulomb field closed model. The solution of the self-similar model and the appropriate BCs are presented in <ref>. In <ref>, the results of our simulations and the variation of the magnetic and thermodynamic profile functions are presented for Coulomb function and self-similar models, and in <ref>, the results obtained from the models are applied to the existing observations of MBPs and the simulations for other solar flux tubes. In <ref>, we have compared between the Coulomb function and self-similar models and find the regime of the validity; we have also discussed the advancement made and how the models for open and closed field flux tubes are useful for building realistic structures. Finally, in <ref>, we summarize and highlight the major points of the work and conclude how the work may be useful for future numerical studies. A glossary of all the symbols used throughout the paper is provided in Table <ref>.  
      










§ GRAD–SHAFRANOV EQUATION FOR THE CYLINDRICAL FLUX TUBE 


In a magnetic medium of field strength B, with gas (or plasma) pressure p and mass density ρ, the magnetohydrostatic (MHS) pressure balance equation is given by

    -∇ p+1/4π(∇× B)× B+ρ g=0,
      
where g denotes the acceleration due to gravity at the solar surface. The individual components of B can be expressed in terms of the poloidal flux function, Ψ(r,z)=∫_0^r r' B_z(r',z)  dr', in the following way

    B_r=-1/r∂Ψ/∂ z ;   B_z=1/r∂Ψ/∂ r ;    B_ϕ=I_p/r,

where I_p represents the poloidal current. These forms of B_r,B_ϕ and B_z automatically ensure the solenoidal condition for B. Using the axisymmetric condition we split the MHS equilibrium eqn (<ref>) into r and z direction and plug in the forms of magnetic field components from eqn (<ref>), to find two scalar partial differential equations


    ∂Ψ/∂ r∂^2Ψ/∂ z^2+∂Ψ/∂ r∂^2Ψ/∂ r^2-1/r(∂Ψ/∂ r)^2+1/2∂ I_p^2/∂ r=-4π r^2∂ p/∂ r
    
    -∂ p/∂ z+1/4π[1/r∂Ψ/∂ z(1/r^2∂Ψ/∂ r-1/r∂^2Ψ/∂ r^2)-1/r^2∂Ψ/∂ z∂^2Ψ/∂ z^2-1/2r^2∂ I_p^2/∂ z]-ρ g=0,

  
where we assume the form of the gas pressure to be

    p(r,z)=p_1(Ψ)+p_2(z);

this form is required in order to have a non-zero density (see SM18). The ϕ part of eqn (<ref>) gives ∇Ψ×∇ I_p =0, which implies I_p=I_p(Ψ). We have the following form from eqns (<ref>, <ref>) for the GSE to be given by

    ∂^2Ψ/∂ r^2-1/r∂Ψ/∂ r+∂^2Ψ/∂ z^2=-1/2∂ I_p^2(Ψ)/∂Ψ-4π r^2∂ p_1(Ψ)/∂Ψ.

From eqns (<ref>, <ref>) we find

    -∂ p_2/∂ z-∂ p_1(Ψ)/∂ z+1/4π[1/r^2∂Ψ/∂ z(1/r∂Ψ/∂ r-∂^2Ψ/∂ r^2)-1/r^2∂Ψ/∂ z∂^2Ψ/∂ z^2-1/2r^2∂ I_p^2(Ψ)/∂ z]-ρ g=0.

Following SM18, by multiplying 4 π r^2∂ z/∂Ψ on both sides of eqn (<ref>) and using eqn (<ref>), we obtain

    ρ(z)=-1/g d p_2(z)/ d z.

We will see later that the prescription of p_2(z) will lead ρ to be a positive quantity, and hence the density within the flux tube is independent of the radial distance r but varies with height z. The temperature, T, inside the flux tube is calculated by the ideal gas law according to the following form

    T(r,z)=μ̅ p(r,z)/R_g ρ(z),

where, R_g=8.314 J mol^-1 K^-1 represents the gas constant and 

    μ̅=1/z_t∫_0^z_tμ_e(z) dz=1.12

is the mean effective molar mass from photosphere to transition region given by the empirical relation, μ_e(z)=1.288 [ 1-0.535 ( z/2.152)^3 ] <cit.> in the domain of 0<z<2.152 Mm. The formulary of the derived functions for the Coulomb function helical flux tube model are summarized in the Table <ref>. In <ref> and <ref>, we reduce the GSE for different models of flux tubes having open or closed field line structures. A flowchart of the solutions of the two  different flux tube models obtained is shown in Fig. <ref>.

Before we solve for the various cylindrical structures, we discuss the boundary conditions below that are crucial to the models, applicable to both open and closed field flux tubes. The magnetic field lines that rise from a horizontal plane and do not return to the same plane within the domain of interest, are called open field lines (see Figs. [<ref>, <ref>, <ref>]). On the other hand, the field lines that rise and return to the same horizontal plane are called the closed field lines (see Fig. <ref>). We take an idealized case in which the flux tube is embedded in a magnetic field free region where there is no current outside the flux tube. We apply the following standard BCs which are used by several authors [e.g. <cit.>], that [B_r(r=0,z)=0, B_ϕ(r=0,z)=0] which implies that the magnetic field line is vertical at the axis of the flux tube. At the boundary, the radial component vanishes i.e. B_r(r=R,z)=0. We also use the BCs that the total pressure at the boundary of the flux tube matches with the external pressure and the radial average of the internal gas pressure at the transition region (z=z_t) is equal to p_t, where the pressure at the photosphere (z=0) outside the flux tube is taken to be p_0=1.228× 10^5 dyne cm^-2 and at the transition region (z_t=2 Mm), it is p_t=0.1488 dyne cm^-2; these are taken from Avrett-Loeser model <cit.>. We specify the appropriate BCs to model both open and closed field flux tubes below: 




    BC 1:    B_r(r=0,z) =0 
       BC 2:    B_ϕ(r=0,z)=0 
       BC 3:    B_r(R,z)=0 
       BC 4:    p_T(R,z) =p_e(z) 
       BC 5:  1/R∫_0^R p(r,z_t)  dr =p_t.




The BCs that distinguishes between the closed and open field flux tubes is the following

    B_ϕ(R,z)
    
    =0;   closed field
    ≠ 0;   open field,

which reduces to the condition,

    Ψ(R,z)=Ψ_b
    
    =0;   closed field
    ≠ 0;   open field,

which is derived in <ref>. The open (general) solution is obtained in <ref> and it is reduced to the special case of the closed solution by taking Ψ_b=0 is presented in <ref>.




[hbt!]
0.8



decision = [diamond, draw,  
    text width=6em, text badly centered, node distance=5cm, inner sep=0pt]
block = [rectangle, draw, 
    text width=6.5em, text centered, rounded corners, minimum height=2em]
line = [draw, -latex']
cloud = [draw, ellipse,node distance=3cm,
    minimum height=2em]
    

[node distance = 2cm, auto]

block = [rectangle, draw, 
    text width=9.2em, text centered,rounded corners, minimum height=1em]
[block] (c1) <ref>
 BCs:
 B_r(0,z)=0

B_ϕ(0,z)=0

B_r(R,z)=0

p_T(R,z)=p_e(z)

1/R∫_0^Rp(r,z_t) dr=p_t;

block = [rectangle, draw, 
    text width=7em, text centered,rounded corners, minimum height=1em]
[block, left of=c1, node distance=5cm] (c2) <ref>
 Cylindrical flux tube general solutions to GSE with B_ϕ≠ 0;

block = [rectangle, draw, 
    text width=5em, text centered,rounded corners, minimum height=3em]
[block, left of=c2, below of=c2, node distance=5.5cm](l1)<ref>
 Coulomb function model;

block = [rectangle, draw, 
    text width=11em, text centered,rounded corners, minimum height=1em]
[block, right of=l1,node distance=12cm](r1)<ref>
 Self-similar model
 (Open field solution)
 ψ_b ≠ 0
 B_ϕ(r →∞)=0;

block = [rectangle, draw, 
    text width=7em, text centered,rounded corners, minimum height=1em]
[block, below of=l1,node distance=3cm](l1b1)ψ=ψ_h+ψ_p
j_ϕ s≠ 0;

block = [rectangle, draw, 
    text width=5.5em, text centered,rounded corners, minimum height=1em]
[block, left of=l1b1,node distance=4cm](l1b1l1)<ref>
 Closed Field solution
 ψ_p=0
 ψ_b=0
(SM18);

block = [rectangle, draw, 
    text width=7em, text centered,rounded corners, minimum height=1em]
[block, below of=l1b1l1,node distance=3.9cm](l1b2l1)ψ(ϖ=1)=0
 B_ϕ(ϖ=1)=0
 j_zs=0;



block = [rectangle, draw, 
    text width=5em, text centered,rounded corners, minimum height=1em]
[block, right of=l1b1,node distance=3.5cm](l1b1r1)<ref>
 Open Field solution
 ψ_p ≠ 0
 ψ_b ≠ 0;

block = [rectangle, draw, 
    text width=7em, text centered,rounded corners, minimum height=1em]
[block, below of=l1b1r1,node distance=3.9cm](l1b2r1)ψ(ϖ=1)=ψ_b
 B_ϕ(ϖ=1) ≠ 0
 j_zs≠ 0;

block = [rectangle, draw, 
    text width=17em, text centered,rounded corners, minimum height=1em]
[block, left of=r1, below of=r1, node distance=3.3cm](l1b1r2)<ref>
 Generalized Gaussian shape function 
D_G(ξ)=D_G0exp(-ξ^n_G)
(n_G>0)

<cit.> for n_G=2;

block = [rectangle, draw, 
    text width=12em, text centered,rounded corners, minimum height=1em]
[block, right of=l1b1r2,node distance=6.8cm](l1b1r3)<ref>
 Power law shape function
D_P(ξ)=D_P0 (1+ξ)^-n_P
(n_P>2)

(current paper);

[->,thick,line width=0.5mm,draw=blue!100](c2.south)to[out=-90,in=90]node [text width=1.9cm,above=0.17cm]R= finite(l1.north);

[->,thick,line width=0.5mm,draw=blue!100](c2.south)to[out=-90,in=90]node [text width=1.9cm,above=0.17cm]R →∞(r1.north);
[line] (c1)– (c2);
[line] (l1b1)– (l1b1l1);
[line] (l1b1)– (l1b1r1);
[line] (l1b1l1)– (l1b2l1);
[line] (l1b1r1)– (l1b2r1);
[line] (l1)– (l1b1);
[->,thick,line width=0.5mm,draw=blue!88](r1.south)to[out=-90,in=90](l1b1r2.north);
[->,thick,line width=0.5mm,draw=blue!88](r1.south)to[out=-90,in=90](l1b1r3.north);


The two families of solutions presented in various sections are indicated above: the Coulomb function model and the self-similar model are shown along with the applicable boundary conditions.







§ COULOMB FUNCTION SOLUTION OF HELICAL FLUX TUBE MODEL
 
For the magnetohydrostatic equilibria with uniform solar gravity and axisymmetric condition, we have split p(r,z)=p_1(Ψ(r,z))+p_2(z) in order to have a non-zero density (SM18). If we take the form of p_1(Ψ) as a linear function of Ψ, we found that the BCs 1–5 [eqn (<ref>)–(<ref>)], which are crucial for our model, will not be satisfied for arbitrary R values. On the other hand, the quadratic function of Ψ, which is more general than the linear form, is the simplest allowed form for p_1 and I_p, satisfies all the BCs [eqns (<ref>)–(<ref>)], where R becomes a free parameter, and can be chosen any value within the domain of our interest. Therefore, we assume p_1(Ψ) and I_p^2(Ψ) to be polynomials of Ψ upto second order <cit.>,


    p_1(Ψ)=1/4 π(a'/2Ψ^2+b' Ψ), 
        I_p^2(Ψ)=α' Ψ^2+ 2β' Ψ+I_0^2,

  
where the parameters a', b',α', β', I_0 are to be determined by appropriate boundary conditions (BCs) and the function p_2(z) will be evaluated later. Plugging eqns (<ref>, <ref>) into eqn (<ref>) we obtain a second order scalar linear partial differential equation

    ∂^2Ψ/∂ r^2-1/r∂Ψ/∂ r+∂^2Ψ/∂ z^2= -(a'r^2+α')Ψ-(b'r^2+β').
 
We define the dimensionless parameters (in the LHS) by introducing the scaling relations,

    ϖ=r/R,    τ=R/z_0,    ψ_b= Ψ_b √(a')/B_0,   a=a' R^4/4, 
       α=α'/4√(a'),    b=b'/B_0 √(a'),    β=R^2 β'/Ψ_b,

where Ψ_b, R, B_0 are the boundary flux, radius and the magnetic field strength at the center of the flux tube respectively, and z̅=z/z_0, where z_0 is a constant.
To solve this equation, we split ψ=Ψ/Ψ_b into homogeneous ψ_h and particular part ψ_p, i.e. ψ=ψ_h+ψ_p. We plug ψ into eqn (<ref>) and separate  out the homogeneous and particular parts to obtain the following dimensionless equations

    ∂^2ψ_h/∂ϖ^2-1/ϖ∂ψ_h/∂ϖ+τ^2 ∂^2ψ_h/∂z̅^2= -4√(a) (√(a)ϖ^2+2α)ψ_h,  
       ∂^2ψ_p/∂ϖ^2-1/ϖ∂ψ_p/∂ϖ= -(4aϖ^2+8√(a)α)ψ_p-(4ab/ψ_bϖ^2+β).
 
The solution of eqn (<ref>) is separable and given by ψ_h(ϖ,z̅)=s(ϖ) Z(z̅) which has been shown in SM18 to be given by 

    s(ϖ)= c  F_0(-α-κ^2,√(a)ϖ^2)       (with a>0),

where F_0(-α-κ^2,√(a)ϖ^2) represents the Coulomb function <cit.> and κ=k R/2√(2)a^1/4, where the value of k is evaluated later. The z–part solution of eqn (<ref>) is given by  

    Z(z̅)=exp(-2√(2)κ a^1/4z̅/τ).

The homogeneous solution takes the following form

    ψ_h(ϖ,z̅)=c exp(-2√(2)κ a^1/4z̅/τ)F_0(-α-κ^2,√(a)ϖ^2).

The solution of the inhomogeneous eqn (<ref>) is given by a power series solution

    ψ_p(ϖ)=-β/8α+i √(a)ϖ^2/2 e^-i √(a)ϖ^2(β/8α-b/ψ_b)∑_n=0^∞F_2^1(n+2,1;n+2+iα;1/2)(i√(a)ϖ^2)^n/(n+1+iα)n!.

A similar but different homogeneous solution which is oscillatory along z-direction has been used for the cases of both D-shaped and toroidally diverted laboratory plasma <cit.>. The general solution of the GS eqn (<ref>) is given by ψ=ψ_h+ψ_p. Since ψ(ϖ,z̅) and ψ^*(ϖ,z̅), its complex conjugate, are the valid solutions of eqn (<ref>), we construct a real solution of eqn (<ref>) by redefining ψ_h(ϖ,z̅)+ψ_h^*(ϖ,z̅)/2→ψ_h(ϖ,z̅) and ψ_p(ϖ)+ψ_p^*(ϖ)/2→ψ_p(ϖ) which leads to ψ(ϖ,z̅)+ψ^*(ϖ,z̅)/2 →ψ(ϖ,z̅)=ψ^X_C= ψ_C^C=ψ_h   (closed field)
 ψ_C^O=ψ_h+ψ_p   (open field).  The solution ψ_h alone gives the closed field structure of flux tube (SM18), which we denote as ψ^C_C; the general solution is a combination of ψ_h and ψ_p and we denote the open field flux tube structure as ψ^O_C.

The total flux function ψ^X_C(ϖ,z̅) is given by

    ψ^X_C(ϖ,z̅)=s(ϖ) Z(z̅)+ψ_p(ϖ),
 
where s(ϖ) and ψ_p(ϖ) are given by eqns (<ref>) and (<ref>) respectively. Now, ψ^X_C(ϖ,z̅) has to be zero at the axis (i.e. ϖ=0) for all z̅, to keep the field finite at the origin. Since s(0)=0, which satisfies the BC 1 [eqn (<ref>)], we obtain from eqn (<ref>), ψ_p(ϖ=0)=0. From eqn (<ref>), we have ψ_p(ϖ=0)=-β/8α. Therefore, we obtain β=0. From eqn (<ref>), ψ_p(ϖ) reduces to

    ψ_p(ϖ)=i√(a)bϖ^2/4ψ_b[  e^i √(a)ϖ^2∑_n=0^∞F_2^1(n+2,1;n+2-iα;1/2)(-i√(a)ϖ^2)^n/(n+1-iα)n!
     -e^-i √(a)ϖ^2∑_n=0^∞F_2^1(n+2,1;n+2+iα;1/2)(i√(a)ϖ^2)^n/(n+1+iα)n!].

From eqns (<ref>) and (<ref>) we have 

    B_z(ϖ,z̅)=B_0ψ_b/2√(a)ϖ(∂ψ_h/∂ϖ+∂ψ_p/∂ϖ),

and we obtain the explicit form for

    B_z(ϖ,z̅)=   B_0/8ϖ[{1-ib/2(F_2^1(1,-iα,2-iα,-1)/1-iα-F_2^1(1,iα,2+iα,-1)/1+iα)}
       ·exp(-2√(2)κ a^1/4z̅/τ)  d/ dϖ[F_0(-α-κ^2,√(a)ϖ^2)+F^*_0(-α-κ^2,√(a)ϖ^2)]
       +ib  d/ dϖ{ϖ^2 e^i √(a)ϖ^2∑_n=0^∞F_2^1(n+2,1,n+2-iα,1/2)(-i√(a)ϖ^2)^n/(n+1-iα)n!
       -ϖ^2 e^-i √(a)ϖ^2∑_n=0^∞F_2^1(n+2,1,n+2+iα,1/2)(-i√(a)ϖ^2)^n/(n+1+iα)n!}].

The steps to obtain eqn (<ref>) are given in Appendix <ref>. From eqns (<ref>, <ref>) we have

    B_r(ϖ,z̅)=-B_0ψ_b τ/2√(a) s(ϖ) Z'(z̅),

whose explicit form is given by

    B_r(ϖ,z̅)=   B_0 κ/2√(2)a^1/4ϖexp(-2√(2)κ a^1/4z̅/τ)[F_0(-α-κ^2,√(a)ϖ^2)+F^*_0(-α-κ^2,√(a)ϖ^2)] 
       ·[1-ib/2(F_2^1(1,-iα,2-iα,-1)/1-i α- F_2^1(1,iα,2+iα,-1)/1+i α)].

From eqn (<ref>) we obtain the toroidal component 

    B_ϕ(ϖ,z̅)=√(2)B_0α^1/2ψ_b/a^1/4ϖ(ψ_h+ψ_p),

whose explicit form is given by

    B_ϕ(ϖ,z̅) = B_0 α^1/2a^-1/4/4√(2)[1/ϖ{ 1-ib/2(F_2^1(1,-iα,2-iα,-1)/1-iα-F_2^1(1,iα,2+iα,-1)/1+iα)}
       ·exp(--2√(2)κz̅/τ) [F_0(-α-κ^2,√(a)ϖ^2)+F^*_0(-α-κ^2,√(a)ϖ^2)] 
        +ib ϖ{ e^i √(a)ϖ^2∑_n=0^∞F_2^1(n+2,1,n+2-iα,1/2)(-i√(a)ϖ^2)^n/(n+1-iα)n!
        -e^-i √(a)ϖ^2∑_n=0^∞F_2^1(n+2,1,n+2+iα,1/2)(-i√(a)ϖ^2)^n/(n+1+iα)n!}].

Applying the BC 3 [eqn (<ref>)] and using eqn (<ref>), we find s(ϖ=1)=0. From BC 2 [eqn (<ref>)] and eqn (<ref>) we find I_0=0. 
We assume that the external pressure from photosphere to transition region decreases exponentially as

    p_e(z)=p_0 exp(-2kz),

where k is the pressure scale height, which is determined by the relation, k=1/2  z_tln(p_0/p_t) =3.405 Mm^-1, where p_0=1.22 × 10^5 dyne cm^-2, p_t=0.148 dyne cm^-2, and z_t=2 Mm. By matching the pressure scale heights inside and outside the flux tube, we see that p_2(z) follows 

    p_2(z)=p_20exp(-2kz),
    
where p_20 is evaluated later. Finally, we have the expression of p(r,z) from eqns (<ref>, <ref>) to be given by

    p(r,z)=1/4π(a'/2Ψ^2+b'Ψ)+p_20exp(-2kz),

whose explicit form is given by

    p(ϖ,z̅)=  B_0^2[(ψ_b^2 s^2(ϖ)/8π+p̅_20)Z^2(z̅)+(ψ_b^2 s(ϖ)ψ_p(ϖ)/4π+bψ_b s(ϖ)/2√(2a))Z(z̅)+(ψ_b^2 ψ_p^2/8π+b ψ_b ψ_p/2√(2a))],

where p̅_20=p_20/B_0^2, and s(ϖ), Z(z̅), ψ_p(ϖ) are given by eqns (<ref>, <ref>, <ref>) respectively. We now calculate the total pressure at the boundary of the flux tube that includes the  contribution due to gas pressure and the magnetic forces due to the presence of the sheet currents j_ϕ (SM18) and j_z. The pressure and radial component of the MHS force balance eqn (<ref>) yields

    -∂ p/∂ r|_r=R+1/4π(B_r∂ B_r/∂ r+B_z∂ B_r/∂ z)|_r=R -∂/∂ r(B^2/8π)|_r=R
     +j_ϕ(R)B_z(R)-j_z(R) B_ϕ(R)=0.
  
The sheet currents j_ϕ and j_z take the forms


    j_ϕ(r)=j_ϕ s δ(r-R),
    
    j_z(r)=j_z s δ(r-R).


Integrating  eqn (<ref>) w.r.t. r from r=R-ϵ to r=R+ϵ where ϵ is an infinitesimal positive quantity we obtain

    -∫_R-ϵ^R+ϵ∂ p/∂ rdr+1/4π(∫_R-ϵ^R+ϵB_r∂ B_r/∂ rdr+∫_R-ϵ^R+ϵB_z∂ B_r/∂ zdr)
     -∫_R-ϵ^R+ϵ∂/∂ r(B^2/8π)dr+ ∫_R-ϵ^R+ϵj_ϕ(r)B_z(r) dr-∫_R-ϵ^R+ϵj_z(r)B_ϕ(r) dr=0,

which leads to

    p_i(R,z)-p_e(z)+j_ϕ s B_z(R)-j_z s B_ϕ(R)
     +1/4π[ B_r ∂ B_r/∂ r+B_z∂ B_r/∂ z]_R+ B_i^2(R,z)-B_e^2(R,z)/8π=0,

where [...]_R denotes the jump condition at the boundary and {  B_i, p_i } and  {
B_e, p_e } are the internal and external magnetic fields and gas pressures in the flux tube respectively. To calculate j_ϕ s and j_zs, we assume an infinitesimal current loop of vertical height L and radial extent R-ϵ to R+ϵ placed at the boundary of the flux tube (see Fig. <ref>) and by applying the line integral along the loop, we obtain 



    B_z(R) L=4π L ∫_R-ϵ^R+ϵj_ϕ sδ(r-R) dr
    
    -B_ϕ(R) L=4π L ∫_R-ϵ^R+ϵj_z sδ(r-R) dr,
   

which implies


    j_ϕ s=B_z(R)/4π
        j_z s=-B_ϕ(R)/4π.


The total internal magnetic field is given by, B_i^2=B_r^2+B_ϕ^2+B_z^2. Applying the BC 3 [eqn (<ref>)] and B_e=0, we have from eqn (<ref>),   

    p_i(R,z)-p_e(z)+3 B_z^2(R,z)/8π+3 B_ϕ^2(R,z)/8π=0.

By expanding eqn (<ref>) we obtain

    ψ_b^2/3(1+6α/√(a))+2/3b ψ_b +ψ_b^2/4a[1/ϖ^2(s'(ϖ)exp(-2√(2)κ a^1/4z̅/τ)+ψ'_p(ϖ))^2]_ϖ=1
       =p̅exp(-4√(2)κ a^1/4z̅/τ),

where  ψ_b=ψ_p(ϖ=1) and p̅=8π (p_0-p_20)/3 B_0^2. By equating the coefficients of exp(-2√(2)κ a^1/4z̅/τ), exp(-4√(2)κ a^1/4z̅/τ) and the constant quantity between both sides of eqn (<ref>), we obtain 


    (1/ϖ^2s'(ϖ)ψ_p'(ϖ))_ϖ=1=0,
       ψ_b^2/4a(s'^2(ϖ)/ϖ^2)_ϖ=1=p̅,   
       ψ_b^2/3(1+6α/√(a))+2b ψ_b/3+ψ_b^2/4a[ψ^'2_p(ϖ)/ϖ^2]_ϖ=1=0.


From eqns (<ref>, <ref>) we find that (ψ'_p(ϖ)/ϖ)_ϖ=1=0, as p̅≠ 0 and hence from eqn (<ref>)

    b=-ψ_b/2(1+6α/√(a)).
  
We now summarize the set of equations we need to solve numerically for the open field model to be given by


    s(ϖ=1)=0,
       (ψ'_p(ϖ)/ϖ)_ϖ=1=0,
       (s'^2(ϖ)/ϖ^2)_ϖ=1=4ap̅/ψ_b^2,
        b=-ψ_b/2(1+6α/√(a)),
       ∫_0^1 p(ϖ,z̅_t)  dϖ=p_t.

 
 The explicit forms of the eqns (<ref>–<ref>) are given by eqns (<ref>–<ref>) in Appendix <ref>. The five eqns (<ref>–<ref>) consists of seven unknown variables, { a, α, b, ψ_b, p̅, R, B_0 }; so there is a unique solution to the Coulomb function open field model for a given pair of the unknown variables. The eqns (<ref>) and (<ref>) contain three variables a, α and R, and we use these two equations to obtain a(R) and α(R). From eqn (<ref>), we calculate ψ_b(R,B_0,b) and then find b(R,B_0) from eqn (<ref>), and hence ψ_b(R,B_0) and then evaluate p̅(R,B_0) from eqn (<ref>). As a result, the complete solution for the open field flux tube depends only on R and B_0 which are the free parameters of the model. As per BC1–BC5 [eqns (<ref>–<ref>)], which are used for the open field flux tube model, the magnetic field component at the boundary of the flux tube is given by

    B_ϕ(ϖ=1,z̅)=√(2 α)B_0 ψ_b^2/a^1/4.

If we demand additionally, that B_ϕ(ϖ=1,z̅)=0, then, from eqn (<ref>), ψ_b=0; also eqn (<ref>) gives b=0. Therefore from eqn (<ref>), we obtain 

    ψ_p(ϖ)=0.

This represents the solution of the homogeneous part ψ_h of GSE, which has been discussed in SM18, that is applicable for closed field flux tube model. Therefore we need to solve eqns (<ref>, <ref>, <ref>), which are given by the explicit forms by eqns (<ref>, <ref>, <ref>) [with b=ψ_b=0], numerically to find the parameters a, α and p̅ in terms of { R, B_0}, which are  the free parameters of the closed field model. The formulary of the derived functions for the Coulomb function helical flux tube model are summarized in the Table <ref>. We discuss the various configurations of Coulomb function open and closed field structure of flux tubes in <ref>.






§ SELF-SIMILAR MODEL


The basic formulation of self-similar model of a flux tube is based on <cit.> (ST58 hereafter). For an axially symmetric cylindrical geometry (r, ϕ, z), where ϕ is ignorable, the magnetic field components are given by eqn (<ref>). The coordinates r and z are combined together into a new dimensionless variable ξ which is called the self-similar parameter and as a consequence, the flux function Ψ can be expressed only as a function of ξ i.e., Ψ(r,z)=Ψ(ξ) (ST58). We define the dimensionless parameters (in the LHS) by introducing the scaling relations,

    ϖ   =r/R,  τ =R/z_0,  ψ =Ψ/Ψ_b,  ψ_b =Ψ_b/B_0 z_0^2,
    p̅_1    =p_1/B_0^2,  I̅_p =I_p/B_0 R,  χ̅ =χ z_0^2,
 
where, Ψ_b, R, B_0 are the boundary flux, radius and the magnetic field strength at the center of the flux tube respectively, and z̅=z/z_0, where z_0 is a constant length. From ST58, the self-similar parameter ξ is defined by

    ξ= ζ(z̅) ϖ,
  
which describes the radial size distribution of the flux tube with height from the base. Plugging in eqn (<ref>), we can rewrite the GS eqn (<ref>) in the following form 

    ψ_b^2/2 τ^2 d/ dξ( dψ/ dξ)^2 ζ'^2(z̅)+ψ_b^2/τ^21/ξ( dψ/ dξ)^2 ζ'(z̅)ζ”(z̅)+ψ_b^2/2 τ^4 d/ dξ(1/ξ dψ/ dξ)^2 ζ^4(z̅)+1/2 ξ^2 dI̅_p^2/ dξζ^2=-4π∂p̅_1/∂ξ,

and the z-part of GS eqn (<ref>), gives the expression of ρ, eqn (<ref>), which is self consistent for both the Coulomb function and self-similar models. We define a quantity which is called the magnetic shape function given by

    D_X(ξ)=1/ξ dψ/ dξ.

Plugging eqn (<ref>) into eqn (<ref>) we obtain 

    ψ_b^2/τ^2ξ D_X^2(ξ) ζζ”+ψ_b^2/2τ^2 d/ dξ[ξ^2 D_X^2(ξ)]ζ'^2+ψ_b^2/2τ^4 d/ dξ(D_X^2(ξ))ζ^4+1/2ξ^2 dI̅_p^2/ dξζ^2=-4π∂p̅_1/∂ξ,

and integrating eqn (<ref>) w.r.t. ξ from 0 to ∞ we write

    ψ_b^2/τ^2ζζ”∫_0^∞ξ D_X^2(ξ)  dξ+ψ_b^2 ζ'^2/2 τ^2[ ξ^2 D_X^2(ξ)]_ξ=0^∞+ψ_b^2 ζ^4/2 τ^4[D_X^2(∞)-D_X^2(0)]
        +ζ^2/2∫_0^∞1/ξ^2 dI̅_p^2/ dξ dξ=-4π∫_0^∞∂p̅_1/∂ξ dξ.
 
Following ST58, we define

    y^2(z̅)=ψ_b D_0/τζ^2(z̅),

where, y(z̅) ≡(B_z(0,z̅)/B_0)^1/2, B_0 ≡ B_z(0,0) and D_0≡ D_X(ξ=0). Next, using eqns (<ref>) and (<ref>) we obtain 

    ψ_b/τy y”/D_0∫_0^∞ξ D_X^2(ξ)  dξ+ψ_b/τy'^2/2D_0[ ξ^2 D_X^2(ξ)]_ξ=0^∞+y^4/2D_0^2 τ^2[D_X^2(∞)-D_0^2]
        + τ/ψ_by^2/2 D_0∫_0^∞1/ξ^2 dI̅_p^2/ dξ dξ=-4π∫_0^∞∂p̅_1/∂ξ dξ.
  
To solve the eqn (<ref>), we need to specify the functional form of p_1, I_p and D_X(ξ) to study the flux tube model with twisted magnetic field under the similarity assumption. The functional form of  p_1=f/2Ψ^2 is taken from <cit.>, where f is the shape function parameter, and the poloidal current, I_p defined by <cit.> and <cit.>, and motivated from the observations of <cit.>. Hence the form of gas pressure p and poloidal current I_p are taken to be 

    p    =p_c exp(-2kz)+f/2Ψ^2
    
    I_p^2   =Ψ^2_b χξ^4 D_X^2(ξ),

for the positivity of ρ(z) at all z which is given by eqn (<ref>). Here, p_2(z)=p_c exp(-2kz) denotes the gas pressure at the flux tube axis, with p_c is the pressure at the center of the flux tube, and χ=(B_ϕ/r B_z)^2, is a constant pitch angle parameter. We deviate from [<cit.>; <cit.>] by employing the extra term, p_c exp(-2kz) with p_1 in eqn (<ref>), to maintain the hydrostatic vertical pressure balance condition under the influence of solar gravity, with two options for the shape function D_X(ξ) specified by

    D_X(ξ)=
      
      D_G(ξ)=D_G0exp(-ξ^n_G);    (n_G>0) :  Generalized Gaussian
     
     D_P(ξ)=D_P0 (1+ξ)^-n_P;    (n_P>2) :  Power law

where
 

    D_G0   =n_G/Γ (2/n_G), 
    
    D_P0   =(n_P-1)(n_P-2),


We see that both the shape functions, eqn (<ref>) vanish asymptotically  at ξ→∞; hence from eqn (<ref>) we obtain

    ψ_b/τy y”/D_0∫_0^∞ξ D_X^2(ξ)  dξ-  y^4/2 τ^2+τ/ψ_by^2/D_0∫_0^∞1/2ξ^2 dI̅_p^2/ dξ dξ=-4π∫_0^∞∂p̅_1/∂ξ dξ.

Next, we evaluate the integrals of the eqn (<ref>) for both generalized Gaussian and power law shape functions. We will see later that, from eqn (<ref>), that the flux function ψ_P(ξ) varies as ξ^2-n_P, in the domain 0< ξ< ∞; therefore ψ_P(ξ) will converge to a finite value at ξ→∞, if n_p>2. The results of the integrals are provided in Table  <ref>. 





Using the values of the integrals from Table <ref> and redefining, y'(z̅=0)≡ y'_0, we reduce eqn (<ref>) to the following form 

    dy/ dz̅=[λ_X(n)/2 ψ_b τ(y^4-1)-2χ̅(y^2-1)-8πψ_b f̅λ_X(n) τln(y)+B̅'̅^2_z0/4]^1/2,

where, B̅'̅_z0=B'_z0z_0/B_0, χ̅=χ z_0^2 and f̅=f z_0^4. The form of eqn (<ref>) is common for both generalized Gaussian (X=G) and power law (X=P) shape functions, where the functions λ_X(n) for both shape functions are defined in Table <ref>. Here we have used the notation B'_z0=B'_z(0,0), which represents the vertical gradient of B_z(0,z) at the center. From eqn (<ref>), we have the following integral relation

    z̅(y)=∫_1^y  dy'/G(y'),
 
where the function G(y) is given by 

    G(y)=[λ_X(n)/2 ψ_b τ(y^4-1)-2χ̅(y^2-1)-8πψ_b f̅λ_X(n) τln(y)+B̅'̅^2_z0/4]^1/2.

We evaluate the integral (<ref>) numerically which gives z̅=z̅(y). Thereafter, inverting the function between z̅ and y, we evaluate y=y(z̅). From eqns (<ref>, <ref>) we obtain

    ξ=√(τ/ψ_b D_0)ϖ y(z̅).

Using the similarity assumption, B_z(r,z)=B_0 ψ_b/τζ^2(z) D_X(ξ) (ST58), and eqns (<ref>, <ref>, <ref>) we calculate the magnetic field components, representing the most general self-similar solution, to be 


    B_z(ϖ,z̅)=B_0/D_0y^2(z̅) D_X(ξ)
        B_r(ϖ,z̅)=-B_0 ϖ/D_0y(z̅) y'(z̅)D_X(ξ)   
        B_ϕ(ϖ,z̅)=√(χ̅)B_0/D_0ϖ y^2(z̅) D_X(ξ).


The flux function for the self-similar model is obtained by integrating the shape function  

    ψ^O_S(ξ)=∫_0^ξξ' D_X(ξ')  dξ'.

Employing eqns (<ref>, <ref>), we obtain the open flux function for generalized Gaussian, ψ_G, and power law, ψ_P, models to be given by


    ψ_S^O(ξ)=
    ψ_G(ξ)= 1-Γ(2/n_G,ξ^n_G)/Γ(2/n_G)  ;    (n_G>0):   Generalized Gaussian
    ψ_P(ξ)=1-(1+ξ)^1-n_P(1+ξ(n_P-1)) ;    (n_P>2):   Power law.

From eqn (<ref>), it is seen that, ψ_G(ξ) and ψ_P(ξ) converges to unity for ξ→∞. The structure of the self-similar model of a flux tube is that the magnetic field decreases asymptotically in the radial direction to zero at infinity. The flux tube does not have any sharp boundary which can make a partition with the external solar atmosphere. In other words, the self-similar flux tube is embedded in a continuous magnetic medium which has the maximum field strength at the base of the axis of the flux tube and radius of the flux tube is infinity. We take the effective radius of the flux tube  as the distance from the axis on the z=0 plane, which makes a circular area where 90 % of the total flux is enclosed. We call this radius as R_90. The total flux is zero at the axis and it increases asymptotically with r. The explicit forms of the magnetic field components, obtained from the eqns (<ref>, <ref>, <ref>) by using eqns (<ref>, <ref>, <ref>, <ref>) are:


    B_z(ϖ,z̅)=
    
    B_0 y^2(z̅) exp[ -( √(τΓ(2/n_G)/n_G ψ_b) y(z̅) ϖ)^n_G],    (n_G>0):   Generalized Gaussian
     
    B_0 y^2(z̅) [1+√(τ/(n_P-1)(n_P-2)ψ_b) y(z̅) ϖ]^-n_P,    (n_P>2):   Power law



    B_r(ϖ,z̅)=
    
    -B_0 y(z̅)y'(z̅)ϖexp[ -( √(τΓ(2/n_G)/n_G ψ_b) y(z̅) ϖ)^n_G],    (n_G>0):   Generalized Gaussian
    
    -B_0 y(z̅)y'(z̅)ϖ[1+√(τ/(n_P-1)(n_P-2)ψ_b) y(z̅) ϖ]^-n_P,    (n_P>2):   Power law



    B_ϕ(ϖ,z̅)=
    
    B_0 √(χ̅)  y^2(z̅)ϖexp[ -( √(τΓ(2/n_G)/n_G ψ_b) y(z̅) ϖ)^n_G],    (n_G>0):   Generalized Gaussian
    
    B_0 √(χ̅) y^2(z̅) ϖ[1+√(τ/(n_P-1)(n_P-2)ψ_b) y(z̅) ϖ]^-n_P,    (n_P>2):   Power law.

The magnetic field components B_r(ϖ,z̅) and B_ϕ(ϖ,z̅) for the self-similar model follow the BCs (1, 2, 3) [eqns (<ref>, <ref>, <ref>)] for R = ∞. B_z(ϖ,z̅) decreases monotonically with ϖ and converges to zero at infinity. The total pressure far from the flux tube axis is only due to the gas pressure p. We use BC 4 [eqn (<ref>)] at z=0, for r →∞, so that

    p(r →∞,0)=p_0.

From eqn (<ref>), we see that the flux function for both generalized Gaussian and power law, converges to unity at ϖ→∞, i.e.

    ψ_S^O(ϖ→∞)= 1.

Using eqns (<ref>, <ref>), we obtain 

    f̅=2(p̅_0-p̅_c),

and the explicit form of p(ϖ, z̅) for both the generalized Gaussian and the power law models is given by 

    p(ϖ,z̅)=B_0^2(f̅/2ψ^2+p̅_c e^-2k̅z̅),

where p̅_0=p_0/B_0^2 and p̅_c=p_c/B_0^2. The formulary of the derived functions for the self-similar flux tube model are summarized in the Table <ref>. The flowchart of the solutions to the Coulomb function and self-similar models are shown in the Fig. <ref>.





§ RESULTS OBTAINED FROM THE MODELS



 §.§ Coulomb function helical flux tube model
 
This magnetohydrostatic Coulomb function helical flux tube model consists of the free parameters R and B_0 and its functional dependence through a(R), α(R), κ(R), b(R,B_0), ψ_b(R,B_0) and p̅(R,B_0). We choose the parameter range, 1 kG ≤ B_0 ≤ 1.5 kG and 100 km ≤ R ≤ 180 km, consistent with the observations of MBP size and field strength distributions <cit.>. In Table <ref>, we show the solutions for combinations of the free parameters { R, B_0}, where we notice the following trends:


  * The boundary flux ψ_b decreases with R for same B_0, and  with B_0 for same R within the parameter space of runs C1-C21.

  * Due to the pressure balance at the boundary of the flux tube, p̅ increases with R for same B_0, but there is no fixed trend  with B_0 for same R within the parameter space of runs C1-C21.

As example, we show the solution of ψ^O_C, and the magnetic and thermodynamic structure of the flux tube for run C4. The radial variation of the solution of ψ^O_C, magnetic components and pressure inside the flux tube are shown in the Figs. <ref>, <ref> and <ref> respectively. Examples of 3D configuration of the magnetic field lines for open and closed field are shown in the Figs. <ref> and <ref> for runs C4 and C10. 2D vertical projection of the magnetic field lines for ψ^O_C inside the flux tube along r-z plane is shown in the Fig. <ref>. The density inside the flux tube is constant along the radial direction but it decreases along z whereas the temperature varies along r direction and is nearly constant along z direction at the axis. The vertical variation of B_z, p and ρ are shown in the Fig. <ref>. Conclusions from figures and tables are discussed in <ref>.





















 §.§ Self similar model


The self-similar model we developed consists of the dimensionless parameters ψ_b, B̅'̅_z 0, f̅ and χ̅ which are the functions of the input parameter set {Ψ_b, B_0, B'_z0, p_c, χ}. The self-similar flux tube solutions are spanned by these parameters but the structures remain similar. We use the values of these input parameters in the range, Ψ_b=10^17–10^18 Mx <cit.>, B_0=1–2 kG <cit.>, B'_z0 in the range 1–2 G-km^-1 <cit.>, p_c<p_0 [<cit.> and SM18], which are observed for small scale magnetic structures in the photosphere. The generalized Gaussian profile reduces to the Gaussian profile for n_G=2, and it has been shown in <ref> that, for the power law profile, the flux function converges to finite value, at infinite radius, only for n_P>2. We study the different cases for n_G=2–3, n_P=3–4 and χ̅=0.01–100 for different combinations of the other parameter sets {Ψ̅_b, B̅'̅_z 0, f̅},  which are shown in runs S1–S19 of Table <ref>. For the parameter set of runs S1–S19, we find the following results:







  * For same Ψ_b and B_0, with the increase of n_G and n_P, R_G and R_P decrease respectively. 


  * For same n_G, n_P and Ψ_b, with the increase of B_0, R_G and R_P decrease; whereas for same n_G, n_P and B_0; R_G and R_P increase with the increasing of Ψ_b.    


  * For a fixed parameter set {ψ_b, B̅'̅_z 0 f̅, χ̅}, we notice that R_G>R_P for n_G=2 and n_P=3, but for values n_G ≥ 2.5 and n_P ≥ 3.5, R_G<R_P; this means that the radii of the flux tubes for the power law profiles falls off more quickly than those of the generalized Gaussian profiles for higher values of n_G and n_P. 


As an example, we show the solution of ψ^O_S and the magnetic and thermodynamic structures for run S1 of Table <ref>. The values of the magnetic and thermodynamic quantities obtained from the self-similar model are reported in Table <ref>, for both the Gaussian and power law shape function profiles.  The radial variation of the generalized Gaussian and power law flux function are shown in the Fig. <ref> for different values of n_G and n_P, and the variation along the r-z plane is shown in Fig. <ref>. The 3D configuration of the field lines for the generalized Gaussian and power law self-similar models are shown in the Figs. <ref> and <ref> for the parameter sets for runs S1 and S2 given in Table <ref>. The radial and vertical distribution of the magnetic field components are shown in Figs. <ref> and <ref> respectively for both the Gaussian and power law models, whereas the density inside the flux tube does not vary along r-direction but decreases along the z-direction which is shown in the Fig. <ref>. The variation of p and T in the r-z plane obtained from the self-similar model are shown in Figs. <ref> and <ref> for Gaussian and power law shape function profiles. Conclusions drawn from Figures [<ref>, <ref>, <ref>, <ref>, <ref>, <ref>, <ref>, <ref>, <ref>, <ref>] and Tables [<ref>, <ref>, <ref>] are discussed in <ref>. 




























§ COMPARING OUR MODELS WITH OBSERVATIONS


We compare our models with the observations reported by the high resolution and high cadence instruments. The small scale magnetic structures in the solar photosphere are often found in the forms of the magnetic bright points (MBPs) which are small scale magnetic flux tubes with open field lines <cit.>. Therefore the MBPs are the best candidates to compare our open field flux tube models with the observations. MBPs can be identified by spectro polarimetric measurements or they can be seen by the G-band filtergrams <cit.>. Next, we compare the observed magnetic field strength, size and the thermodynamic quantities of the MBPs with that obtained from our models. The MBPs are observed as a region of the unipolar flux concentration, therefore, in the Coulomb function model, we construct a cylindrical boundary of cut-off radius r_b inside the total simulation domain, where the line of sight magnetic field B_z vanishes. The magnetic field strength inside the cylinder of the cut-off radius is always positive. From the recent observations by <cit.>, it has been reported that the MBPs number distribution for the size, peaks in the range 160-200 km and the magnetic field strength is at ∼ 1.4 kG. From the Fig. <ref> we see that the B_z vanishes at r_b=84 km, where R=120 km is the entire radial simulation domain. We choose the parameter range, 1 kG ≤ B_0≤ 1.5 kG and 100 km ≤ R≤ 180 km, for which the magnetic and thermodynamic quantities obtained from our model is in reasonable agreement with the solar atmosphere <cit.>, and the selection of the { R, B_0 } parameter space is also consistent with the observations of MBP size and field strength distributions <cit.>. The values of the magnetic and thermodynamic quantities obtained from the Coulomb function open field model is reported in the Table <ref>. For the self-similar model, the choice of the parameter space is consistent with the MBPs. We take the flux value in the range of 10^17–10^18 Mx which is the typical flux value for MBPs <cit.>. According to the previous studies by <cit.> and SM18, the gas pressure at the axis of MBP is lesser than its boundary gas pressure, so we have chosen the parameter p_c< p_0. The field strength of the magnetic footpoints observed in the photosphere for MBPs are ∼1 kG with a distribution peak at 1.3 kG <cit.>. Thus, we use the value of B_0 in the typical range of 1–2 kG <cit.> in our model. The vertical gradient of the magnetic field strength at the photosphere is ∼ 1 G-km^-1 <cit.>. Hence we use the value of B_z0' in the range of 1–2 G-km^-1 in our model. We have reported the combinations of the free parameters and the corresponding input parameters in Table <ref>. Within the parameter sets of runs S1–S19 in Table <ref>, we notice that the minimum and maximum radii of the flux tubes  are 151 and 678 km respectively for Gaussian model, and 71 and 826 km respectively for the power law model, which are in the reasonable agreement  with the observations of MBP size distributions <cit.>. The values of the magnetic and thermodynamic quantities obtained from the self-similar model are reported in the Table <ref>, which is also in a reasonable agreement with the solar atmosphere reported by <cit.>.





  


 







§ DISCUSSION OF THE MODELS


We discuss the findings of our simulations below:


  * The Coulomb function model is easier to implement for the numerical studies as it consists of two free parameters (R, B_0); on the other hand, the self-similar model consists of five free parameters (Ψ_b, B_0,  p_c, B'_z0 and χ). From Table <ref>, we see that the rise of the gas pressure along the radial direction from axis to the boundary is higher for the self-similar model than the Coulomb function model at higher z. The density within the flux tube does not vary with r; hence, the rise of the temperature from axis to the boundary at higher z is also higher for the self-similar model relative to the Coulomb function model. For the Coulomb function model, the radial boundary of the flux tube is defined where B_z vanishes; on the other hand, for the self-similar model, B_z(R) ≠ 0, whereas B_z reduces along the radial direction from axis to the boundary for the Gaussian model faster than the Power law model. 


  * The radial size and the magnetic field strength at the center of the flux tube are the free parameters in the Coulomb function model. The magnetic and thermodynamic structure of the flux tube remain similar for different values of the free parameters, whereas the magnitude of the magnetic and thermodynamic quantities vary. We have explored the parameter space and notice that, in the domain of 100 km ≤ R ≤ 180 km, and 1 kG ≤ B_0 ≤ 1.5 kG, the magnetic and thermodynamic quantities are in reasonable agreement with the solar atmosphere <cit.>, which also validate the MBP size and magnetic field strength distribution <cit.>. For the self-similar model, the radial sizes of the flux tubes depend on the choice of the dimensionless input parameters {ψ_b, f̅, B̅'̅_̅z̅0̅, χ̅}. In the domain of the selected parameter space (see Table <ref>), the maximum and minimum radii of the flux tubes are 678 and 151 km obtained from the Gaussian model, whereas for the power law model with n_P=3, the maximum and minimum values of the radii are calculated to be 826 and 184 km respectively which are also in reasonable agreement with the observation of MBP size distribution by <cit.>.


  * For the Coulomb function model, we notice that the value of α decreases with R (see Table <ref>), which lowers the poloidal current I_p and the twist of the field lines. The 3D geometry of the field lines for different twists are shown in the Figs. <ref> and <ref> for open and closed field Coulomb function models respectively. In the self-similar model, the twist of the field lines increases with χ̅ and are shown in Figs. <ref> and <ref> for Gaussian and power law profiles respectively, which follows from eqn (<ref>).


  * The gas pressure for both Coulomb function and self-similar models increases along radial direction from axis to the boundary, whereas it decreases along the vertical height from photosphere to the transition region (see Figs. <ref> and <ref>) which is similar to the result obtained by <cit.> for MBPs, where the gas pressure inside the MBPs increases radially though the change is not significant, and decreases vertically. <cit.> have studied for the cases of single and multiple flux tubes, where the internal gas pressure is nearly same along the radial distance but decreases with height. The density within the flux tube does not change radially but it decreases along z, for both Coulomb function and self-similar models (see Figs. <ref> and <ref>).  Our model predicts that the atmosphere inside the flux tube is nearly plane parallel which is comparable to the model obtained for MBPs by <cit.>. 



  * In the solar atmosphere, the temperature in the transition region rises perhaps because the shock dissipation of waves play a dominant role, which is not included in our model. We have also not implemented the temperature profile by <cit.> (VAL model); however our model is self consistent, obtained by solving the GSE without shock heating. Therefore, we do not see the drastic rise of the temperature with height. Our vertical simulation domain is restricted from the photosphere to the transition region where our input external atmosphere model is valid. Both the flux tube models we built are non-isothermal where the temperature increases along the radial direction for both Coulomb function open field and self-similar models. The vertical variation of the temperature is constant at the axis but it increases with height away from the axis for the Coulomb function open field and self-similar models (see Figs. <ref> and <ref>).  



  * <cit.>, <cit.>, <cit.> have reported the simulation results of MBPs by using MuRAM and Copenhagen-Stagger code where the obtained values of magnetic field strength, pressure, density and temperature inside the flux tube are in reasonable agreement with our predictions. 




  * The 2D simulations of the propagation of linear and non-linear magneto acoustic wave through an open magnetic flux tube, embedded in the solar atmosphere from photosphere to corona were carried out by <cit.>. We can incorporate our solutions as the background condition for such numerical studies of waves and their kinematic properties taking realistic inputs of field strength and pressure distribution observed in the solar atmosphere.



  * The Coulomb function model gives both open and closed field flux tube solutions, which can be co-added to build the canopy structure. A cartoon diagram of the magnetic canopy is shown in Fig. <ref>, where the closed field lines (red), Ψ^C_C, are present between the open field flux tubes and obtained from the Coulomb function, where the open field lines (blue), Ψ^O_C and Ψ^O_S, of the neighboring flux tubes merge to each other to form a canopy structure. This is similar to structures assumed in the numerical simulations by <cit.>, constructed by a different self-similar flux tube solution. We can use our solutions for inputs to simulations to build such canopy structures. The self-similar flux tube model gives an open field structure of the flux tube which is embedded in a continuous magnetic medium and span upto  infinity in the radial direction. The magnetic and thermodynamic quantities we estimated from both Coulomb function and self-similar models are nearly similar, whereas there are some differences in the structures of the magnetic and thermodynamic profiles.





Future advancement of the observations of magnetic and thermodynamic structures of the MBPs will provide better selection of the parameter inputs and discriminate between our models. 







§ SUMMARY AND CONCLUSIONS


In this work, we have constructed two different models of flux tubes with twisted magnetic fields which are the Coulomb function helical flux tube model and self-similar model by solving GSE semi-analytically. We tabulate the expressions of magnetic and thermodynamic functions for Coulomb and self-similar models in Table <ref>, and highlight the novel features of this work below.



  * By incorporating the form of gas pressure and poloidal current we have solved GSE to obtain the flux function for the Coulomb function model. The solution of the Coulomb function model is the combination of a homogeneous part and a particular part. The homogeneous part with closed geometry is separable with a Coulomb function in r whereas the z part decreases exponentially with height, and  the particular part with open geometry is a power series of r which is independent of z. 


  * Using appropriate BCs and employing the presence of the sheet current at the boundary of the flux tube, we have determined the parameters a(R), α(R), κ(R), b(R,B_0) , ψ_b(R,B_0) and p̅(R,B_0) in terms of the input parameters { R, B_0 }, which are the free parameters in the model, and k is calculated from the pressure values at photosphere and transition region obtained from <cit.> model. The values of the parameters for Coulomb function model are listed in Table <ref>.



  * In the Coulomb function model, the solution consisting of homogeneous and particular parts together represents an open field flux tube solution, where the field lines rise from the photosphere. The  homogeneous solution depicts a closed field flux tube model which is discussed in SM18. The values of the magnetic field strength and thermodynamic quantities inside the flux tube are calculated and are summarized in Table <ref>. The 3D visualization of both open and closed field lines are shown in the Figs. <ref> and <ref> for the parameter set of run no. C4 and C10 corresponding to Table <ref>.



  * In the self-similar model, we have employed an extra term p_c exp(-2kz) with p_1 in eqn (<ref>), to maintain the hydrostatic pressure balance under the influence of stratified solar gravity, and taken two options for the shape functions, D_X(ξ) from eqn (<ref>), which is the extension of previous models by ST58; <cit.>. We have incorporated the resulting two different shape functions, generalized Gaussian and power law profiles, to obtain open field flux tube solutions. We have taken a range of the parameters Ψ_b, B_0, p_c, B'_z0 and χ (see Table <ref>), that are consistent with the solar atmosphere to study the structure and the properties of the flux tubes. The size of the flux tubes and the magnitude of the thermodynamic and magnetic field strengths depend on the choice of the input parameters, but the magnetic and thermodynamic structures remain similar. We have calculated the magnetic field strength and the thermodynamic quantities inside the flux tube which are given in Table <ref>, for the parameter set of run no. S1 corresponding to Table <ref>. 


  * Preliminary calculations using the constraint of relative helicity based on the formulations given in <cit.> and applying the constrained energy minimization principle <cit.> indicate that stable configurations are possible for some regions in the parameter space of { B_0, R }. We plan a complete solution of this allowed region and test it with numerical simulations in a paper in preparation. 


The flux tube models presented here give useful estimates of the magnitude and the distribution of the magnetic field strength and thermodynamic quantities from the photosphere to the transition region which can be verified by the future observations. Work on self-similar closed and twisted field structure is in progress. The solutions we obtained for different flux tubes can be used for the dynamical simulation of wave propagation through the flux tubes, which is important for studying the coronal heating by waves.



We thank V. Fedun, R. Erdélyi and S. Shelyag for useful discussions. We thank the anonymous referee for insightful comments and helpful suggestions. We also thank the support staff of the IIA HPC facility, VBO (IIA, Kavalur) staff for hospitality during our visits, and Saikat Das for help with Figs. <ref> and <ref>.






§ DERIVATION OF THE EXPLICIT FORM OF B_Z FOR COULOMB FUNCTION MODEL


The homogeneous solution s(ϖ) which is given by the eqn (<ref>), can be represented in terms of the Whittaker–M function (SM18), where the Whittaker–M function can be expressed in terms of hypergeometric function by the standard relation <cit.>

    M_t,m(ν)=e^-ν/2ν^m+1/2 F_1^1(1/2+m-t,1+2m,ν),

where F_1^1 represents the hypergeometric function with the arguments t,m and ν. Taking the real part of ψ_h(ϖ,z̅) from eqn (<ref>) and ψ_p(ϖ) from (<ref>), and using eqns [(<ref>), (<ref>)] we obtain 

    B_z(ϖ,z̅)=    B_0 ψ_b c exp(-2√(2)κ a^1/4z̅/τ) [ 8(1+i√(a)ϖ^2)F_1^1(1+iα,2,2i√(a)ϖ^2)-
        8√(a)ϖ^2F_1^1(1+iα,3,2i√(a)ϖ^2)]+B_0ψ_b/2√(a)(ψ'_p(ϖ)+ψ^'*_p(ϖ)/2ϖ)

where, B_0 ≡ B_z(0,0), and  from eqn (<ref>) we obtain, 

    1=8ψ_b c+ψ_b/2√(a)[ψ_p'(ϖ)+ψ^'*_p(ϖ)/2ϖ]_ϖ=0,

where the identity F_1^1(1+iα,2,0)=1 is applied. By expanding the last term on the RHS of eqn (<ref>), we obtain

    1=8 ψ_b c+ib/2[ F_2^1(1,-iα,2-iα,-1)/1-i α - F_2^1(1,iα,2+iα,-1)/1+i α]
 
and the expression for c is given by

    c=1/8ψ_b[1-ib/2(F_2^1(1,-iα,2-iα,-1)/1-i α - F_2^1(1,iα,2+iα,-1)/1+i α)].

Hence from eqn (<ref>), we obtain the explicit form for B_z(ϖ,z̅), given in eqn (<ref>). 


§ EXPLICIT FORMS OF THE BCS FOR COULOMB FUNCTION MODEL

   
The explicit forms of the eqns (<ref>–<ref>) are 

    F_0(-α-κ^2,√(a))+F^*_0(-α-κ^2,√(a))=0, 
       [1/ϖ^2 d/ dϖ(ϖ^2e^i√(a)ϖ^2∑_n=0^∞F_2^1(n+2,1,n+2-iα,1/2)(-i√(a)ϖ^2)^n/(n+1-iα)n!
       -ϖ^2e^-i√(a)ϖ^2∑_n=0^∞F_2^1(n+2,1,n+2+iα,1/2)(i√(a)ϖ^2)^n/(n+1+iα)n!)]_ϖ=1=0, 
       [1-ib/2(F_2^1(1,-iα,2-iα,-1)/1-iα-F_2^1(1,iα,2+iα,-1)/1=iα)]^2 
       ·(  d/ dϖ[F_0(α-κ^2,√(a)ϖ^2)+F^*_0(α-κ^2,√(a)ϖ^2)]_ϖ=1)^2 =p̅ a, 
       b=-ψ_b/2(1+6α/√(a)),
       ∫_0^1[e^-2κ_t(1-i b/2[F_2^1(1,-iα,2-iα,-1)/1-iα-F_2^1(1,iα,2+iα,-1)/1+iα])^2
       ·(F_0(-α-κ^2,√(a)ϖ^2)+F^*_0(-α-κ^2,√(a)ϖ^2))^2+ 8e^-κ_tψ_b (ψ_p+ψ^*_p+2b/ψ_b)
       ·(1-ib/2[F_2^1(1,-iα,2-iα,-1)/1-α-F_2^1(1,iα,2+iα,-1)/1+α])^2
       ·(F_0(-α-κ^2,√(a)ϖ^2)+F^*_0(-α-κ^2,√(a)ϖ^2))+ψ_b^2/4(ψ_p+ψ^*_p)^2+bψ_b (ψ_p+ψ^*_p)] dϖ=p̅

respectively, where, ψ_p(ϖ) is given by eqn (<ref>), κ_t= k z_t, and ψ_p(ϖ=1)=ψ_b.
   

   
   
   
   
   

