Open and Closed Magnetic Configurations of Twisted Flux Tubes

The small-scale magnetic structure in the solar photosphere plays an important role in several phenomena, such as the evolution of active regions (Muller & Mena 1987; Aschwanden et al. 2000; Centeno et al. 2007), heating of the corona through the dissipation of waves (Ruzmaikin & Berger 1998; Srivastava et al. 2017), and reconnection between flux tubes (van Ballegooijen 1986; Muller et al. 1994). 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 magnetic reconnection, leading to energy release in the solar corona (Parker 1988; Peter et al. 2005; Thalmann et al. 2013). Therefore, modeling of the proverbial flux tube is a key aspect of understanding the various phenomena on the solar surface and the Sun's outer atmosphere.

There have been several previous attempts to construct a model of flux tubes for both twisted and untwisted magnetic fields. Schlüter & Temesváry (1958) 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 with respect to 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 in which the magnetic lines of force rise from a horizontal plane and do not return in the model domain. Yun (1971) 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}_{\phi }(r,z)$ was taken from the data obtained from observations (Stepanov 1965). 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 Schlüter & Temesváry (1958), Osherovitch (1982) 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.

Steiner et al. (1986) numerically studied a 2D model of an open single flux tube with a twist using the standard boundary conditions (BCs) 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 that span from the photosphere to the corona have been studied for the case of an untwisted magnetic field (Gent et al. 2013, 2014), in which an empirical form of the magnetic field components is motivated by a self-similar construction. A numerical model of flux tubes was studied by Murawski et al. (2015), where an empirical form of magnetic flux function was assumed; this was followed by a model to study the propagation of the MHD waves through 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, Vigeesh et al. (2009) 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. Fedun et al. (2009) studied the propagation of the acoustic wave through the solar atmosphere due to the periodic drivers at the photosphere, and Shelyag et al. (2010) 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 the Grad–Shafranov equation (GSE; Grad & Rubin 1958; Shafranov 1958). We have assumed a quadratic form of the flux function for the gas pressure and poloidal current, which has previously been used to study the equilibrium solution of terrestrial plasma (Atanasiu et al. 2004), 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 main goal of this paper is to show that the closed field model, reported in Sen & Mangalam (2018; hereafter SM18), 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 the 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, was reported by Solov'ev (1968). SM18 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 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 the transition region. After building the solutions semianalytically and applying appropriate BCs, we calculate the magnetic field structure and thermodynamic quantities inside the flux tube. As MBPs observed in the photosphere (Muller & Mena 1987; Centeno et al. 2007; Lagg et al. 2010; Shelyag et al. 2010) 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 Section 2, we apply the GSE to the cylindrical flux tube case and describe the common BCs that are physically realistic and used in the modeling of our flux tubes. In Section 3, 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 Section 4. In Section 5, 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 Section 6, the results obtained from the models are applied to the existing observations of MBPs and the simulations for other solar flux tubes. In Section 7, we compare the Coulomb function and self-similar models and find the regime of validity; we also discuss the advancements we have made and how the models for open and closed field flux tubes are useful for building realistic structures. Finally, in Section 8, we summarize and highlight the major points of the work and conclude how our findings may be useful for future numerical studies. A glossary of all the symbols used throughout the paper is provided in Table 1.

Table 1.  Glossary of Symbols Used in the Different Flux Tube Models

Symbols common to all the models
r	Radial coordinate in cylindrical geometry	ϕ	Azimuthal coordinate in cylindrical geometry
z	Vertical coordinate in cylindrical geometry	Br	Radial component of the magnetic field strength
Bϕ	Azimuthal component of the magnetic field strength	Bz	Vertical component of the magnetic field strength
B0	Magnetic field strength at the center of the flux tube	Ip	Poloidal current
Ψ	Total flux function	Ψb	Flux at the boundary
k	Scale factor of the pressure-scale height	p	Gas pressure
pT	Total pressure inside flux tube	p0	Gas pressure at the photosphere
pe	Pressure outside the flux tube in the solar atmosphere	pt	Pressure at the transition region
zt	Height of the transition region from the photosphere	ρ	Gas density inside the flux tube
T	Temperature inside the flux tube	g	Acceleration due to gravity at the photosphere
Rg	Universal gas constant	μe	Effective molar mass
$\bar{\mu }$	Mean effective molar mass
Coulomb function model
a, α, b, β, κ	Parameters of the GSE	ψh	Homogeneous solution of GSE, same as ${\psi }_{C}^{C}$
s	Radial part of the homogeneous solution	Z	Vertical part of homogeneous solution
ψp	Particular solution of GSE	${\psi }_{C}^{X}$	General Coulomb solution of the GSE
${\psi }_{C}^{O}$	Coulomb function open field solution	${\psi }_{C}^{C}$	Coulomb function closed field solution
jϕ	Sheet current along azimuthal direction at flux tube boundary	jz	Sheet current along vertical direction at flux tube boundary
p20	Gas pressure at the center of the flux tube	R	Radius of the flux tube
rb	Cutoff radius inside the flux tube
Self-similar model
ξ	Self-similar parameter	DX	Shape function of the flux function
${B}_{z0}^{{\prime} }$	Vertical gradient of magnetic field at the center of the flux tube	$\bar{\chi }$	Dimensionless twist parameter of the field lines
$\bar{f}$	Dimensionless shape function parameter	pc	Gas pressure at the center of the flux tube
${\psi }_{S}^{O}$	Self-similar open field flux	R90	Radius of the flux tube where 90% of the flux is enclosed
Generalized Gaussian flux tube model
ψG	Flux function	DG	Shape function of the flux function
nG	Power index of shape function profile	RG	Radius of the flux tube
Power-law flux tube model
${\psi }_{P}$	Flux function	DP	Shape function of the flux function
nP	Power index of the shape function profile	RP	Radius of the flux tube

Download table as:  ASCIITypeset image

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

$$\begin{eqnarray}&&-{\rm{\nabla }}p+\displaystyle \frac{1}{4\pi }({\rm{\nabla }}\times {\boldsymbol{B}})\times {\boldsymbol{B}}+\rho {\boldsymbol{g}}=0,\end{eqnarray} \tag{ 1 }$$

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, ${\rm{\Psi }}(r,z)={\int }_{0}^{r}r^{\prime} {B}_{z}(r^{\prime} ,z){dr}^{\prime}$, in the following way:

$$\begin{eqnarray}&&{B}_{r}=-\displaystyle \frac{1}{r}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial z};\quad {B}_{z}=\displaystyle \frac{1}{r}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial r};\quad {B}_{\phi }=\displaystyle \frac{{I}_{p}}{r},\end{eqnarray} \tag{ 2 }$$

where I_p represents the poloidal current. These forms of ${B}_{r},{B}_{\phi }$ and B_z automatically ensure the solenoidal condition for B. Using the axisymmetric condition we split the MHS equilibrium Equation (1) into the r- and z-directions and plug in the forms of magnetic field components from Equation (2), to find two scalar partial differential equations

$$\begin{eqnarray}&&\displaystyle \frac{\partial {\rm{\Psi }}}{\partial r}\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {z}^{2}}+\displaystyle \frac{\partial {\rm{\Psi }}}{\partial r}\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {r}^{2}}-\displaystyle \frac{1}{r}{\left(\displaystyle \frac{\partial {\rm{\Psi }}}{\partial r}\right)}^{2}+\displaystyle \frac{1}{2}\displaystyle \frac{\partial {I}_{p}^{2}}{\partial r}=-4\pi {r}^{2}\displaystyle \frac{\partial p}{\partial r},\end{eqnarray} \tag{ 3a }$$

$$\begin{eqnarray}&&\begin{array}{l}-\displaystyle \frac{\partial p}{\partial z}+\displaystyle \frac{1}{4\pi }\left[\displaystyle \frac{1}{r}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial z}\left(\displaystyle \frac{1}{{r}^{2}}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial r}-\displaystyle \frac{1}{r}\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {r}^{2}}\right)\right.\\ \quad -\left.\,\displaystyle \frac{1}{{r}^{2}}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial z}\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {z}^{2}}-\displaystyle \frac{1}{2{r}^{2}}\displaystyle \frac{\partial {I}_{p}^{2}}{\partial z}\right]-\rho g=0,\end{array}\end{eqnarray} \tag{ 3b }$$

where we assume the form of the gas pressure to be

$$\begin{eqnarray}&&p(r,z)={p}_{1}({\rm{\Psi }})+{p}_{2}(z);\end{eqnarray} \tag{ 4 }$$

this form is required in order to have a non-zero density (see SM18). The ϕ part of Equation (1) gives ∇Ψ × ∇I_p = 0, which implies ${I}_{p}={I}_{p}({\rm{\Psi }})$. We have the following form from Equations 3(a), (4) for the GSE to be given by:

$$\begin{eqnarray}&&\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {r}^{2}}-\displaystyle \frac{1}{r}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial r}+\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {z}^{2}}=-\displaystyle \frac{1}{2}\displaystyle \frac{\partial {I}_{p}^{2}({\rm{\Psi }})}{\partial {\rm{\Psi }}}-4\pi {r}^{2}\displaystyle \frac{\partial {p}_{1}({\rm{\Psi }})}{\partial {\rm{\Psi }}}.\end{eqnarray} \tag{ 5 }$$

From Equations 3(b), (4) we find

$$\begin{eqnarray}&&\begin{array}{l}-\displaystyle \frac{\partial {p}_{2}}{\partial z}-\displaystyle \frac{\partial {p}_{1}({\rm{\Psi }})}{\partial z}+\displaystyle \frac{1}{4\pi }\left[\displaystyle \frac{1}{{r}^{2}}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial z}\left(\displaystyle \frac{1}{r}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial r}-\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {r}^{2}}\right)\right.\\ \quad -\left.\,\displaystyle \frac{1}{{r}^{2}}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial z}\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {z}^{2}}-\displaystyle \frac{1}{2{r}^{2}}\displaystyle \frac{\partial {I}_{p}^{2}({\rm{\Psi }})}{\partial z}\right]-\rho g=0.\end{array}\end{eqnarray} \tag{ 6 }$$

Following SM18, by multiplying $4\pi {r}^{2}\tfrac{\partial z}{\partial {\rm{\Psi }}}$ on both sides of Equation (6) and using Equation (5), we obtain

$$\begin{eqnarray}&&\rho (z)=-\displaystyle \frac{1}{g}\displaystyle \frac{{{dp}}_{2}(z)}{{dz}}.\end{eqnarray} \tag{ 7 }$$

We will see later that the prescription of p₂(z) will lead ρ to be a positive quantity, 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:

$$\begin{eqnarray}&&T(r,z)=\displaystyle \frac{\bar{\mu }\ p(r,z)}{{R}_{g}\ \rho (z)},\end{eqnarray} \tag{ 8 }$$

where, ${R}_{g}=8.314$ J mol⁻¹ K⁻¹ represents the gas constant and

$$\begin{eqnarray}&&\bar{\mu }=\displaystyle \frac{1}{{z}_{t}}{\int }_{0}^{{z}_{t}}{\mu }_{e}(z){dz}=1.12\end{eqnarray} \tag{ 9 }$$

is the mean effective molar mass from the photosphere to the transition region given by the empirical relation, ${\mu }_{e}(z)=1.288\left[1-0.535{\left(\tfrac{z}{2.152}\right)}^{3}\right]$ (Solov'ev & Kirichek 2015), in the domain of 0 < z < 2.152 Mm. The formularies of the derived functions for the Coulomb function helical flux tube model are summarized in the Table 6. In Sections 3 and 4, we reduce the GSE for different models of flux tubes with open or closed field line structures. A flowchart of the obtained solutions of the two different flux tube models is shown in Figure 1.

Zoom In Zoom Out Reset image size

Before we solve for the various cylindrical structures, we discuss the BCs 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 Figures 6, 12, 13). The field lines that rise and return to the same horizontal plane are called the closed field lines (see Figure 7). 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 that have been used by several authors (e.g., Mangalam & Krishan 2000; Solov'ev & Kirichek 2015; Sen & Mangalam 2018): $[{B}_{r}(r=0,z)=0$, ${B}_{\phi }(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 BCs for which 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\times {10}^{5}$ dyne cm⁻² and at the transition region (${z}_{t}=2$ Mm), it is p_t = 0.1488 dyne cm⁻²; these are taken from the Avrett–Loeser model (Avrett & Loeser 2008). We specify the appropriate BCs to model both the open and closed field flux tubes below:

$$\begin{eqnarray}&&\mathrm{BC}1:\quad {B}_{r}(r=0,z)=0,\end{eqnarray} \tag{ 10a }$$

$$\begin{eqnarray}&&\mathrm{BC}2:\quad {B}_{\phi }(r=0,z)=0,\end{eqnarray} \tag{ 10b }$$

$$\begin{eqnarray}&&\mathrm{BC}3:\quad {B}_{r}(R,z)=0,\end{eqnarray} \tag{ 10c }$$

$$\begin{eqnarray}&&\mathrm{BC}4:\quad {p}_{T}(R,z)={p}_{e}(z),\end{eqnarray} \tag{ 10d }$$

$$\begin{eqnarray}&&\mathrm{BC}5:\quad \displaystyle \frac{1}{R}{\int }_{0}^{R}p(r,{z}_{t}){dr}={p}_{t}.\end{eqnarray} \tag{ 10e }$$

The BCs that distinguish between the closed and open field flux tubes are as follows:

$$\begin{eqnarray}\begin{array}{c}{B}_{\phi }(R,z)\left\{\begin{array}{ll}=0; & {\rm{closedfield}}\\ \ne 0; & {\rm{openfield}},\end{array}\right.\end{array}\end{eqnarray} \tag{ 11 }$$

which reduces to the condition,

$$\begin{eqnarray}{\rm{\Psi }}(R,z)={{\rm{\Psi }}}_{b}\left\{\begin{array}{ll}=0; & \mathrm{closedfield}\\ \ne 0; & \mathrm{openfield},\end{array}\right.\end{eqnarray} \tag{ 12 }$$

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

For magnetohydrostatic equilibria with uniform solar gravity and axisymmetric conditions, we have split $p(r,z)\,={p}_{1}({\rm{\Psi }}(r,z))+{p}_{2}(z)$ in order to have a non-zero density (SM18). If we take the form of p₁(Ψ) as a linear function of Ψ, we find that BCs 1–5 (Equations 10(a)–(e)), 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₁ and I_p, satisfying all the BCs (Equations 10(a)–(e)), in which R becomes a free parameter, and can be can be any chosen value within the domain of our interest. Therefore, we assume p₁(Ψ) and ${I}_{p}^{2}({\rm{\Psi }})$ to be polynomials of Ψ up to second order (Atanasiu et al. 2004),

$$\begin{eqnarray}&&{p}_{1}({\rm{\Psi }})=\displaystyle \frac{1}{4\pi }\left(\displaystyle \frac{a^{\prime} }{2}{{\rm{\Psi }}}^{2}+b^{\prime} {\rm{\Psi }}\right),\end{eqnarray} \tag{ 13a }$$

$$\begin{eqnarray}&&{I}_{p}^{2}({\rm{\Psi }})=\alpha ^{\prime} {{\rm{\Psi }}}^{2}+2\beta ^{\prime} {\rm{\Psi }}+{I}_{0}^{2},\end{eqnarray} \tag{ 13b }$$

where the parameters $a^{\prime} ,b^{\prime} ,\alpha ^{\prime} ,\beta ^{\prime} ,{I}_{0}$ are to be determined by appropriate BCs and the function ${p}_{2}(z)$ will be evaluated later. Plugging Equations 13(a), (b) into Equation (5) we obtain a second order scalar linear partial differential equation:

$$\begin{eqnarray}&&\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {r}^{2}}-\displaystyle \frac{1}{r}\displaystyle \frac{\partial {\rm{\Psi }}}{\partial r}+\displaystyle \frac{{\partial }^{2}{\rm{\Psi }}}{\partial {z}^{2}}=-(a^{\prime} {r}^{2}+\alpha ^{\prime} ){\rm{\Psi }}-(b^{\prime} {r}^{2}+\beta ^{\prime} ).\end{eqnarray} \tag{ 14 }$$

We define the dimensionless parameters (on the LHS) by introducing the scaling relations,

$$\begin{eqnarray}\begin{array}{rcl}\varpi & = & r/R,\ \tau =R/{z}_{0},\ {\psi }_{b}=\displaystyle \frac{{{\rm{\Psi }}}_{b}\sqrt{a^{\prime} }}{{B}_{0}},\ a=\displaystyle \frac{a^{\prime} {R}^{4}}{4},\\ \alpha & = & \displaystyle \frac{\alpha ^{\prime} }{4\sqrt{a^{\prime} }},\ b=\displaystyle \frac{b^{\prime} }{{B}_{0}\sqrt{a^{\prime} }},\ \beta =\displaystyle \frac{{R}^{2}\beta ^{\prime} }{{{\rm{\Psi }}}_{b}},\end{array}\end{eqnarray} \tag{ 15 }$$

where Ψ_b, R, B₀ are the boundary flux, radius and the magnetic field strength at the center of the flux tube, respectively, and $\bar{z}=z/{z}_{0}$, where z₀ is a constant. To solve this equation, we split $\psi ={\rm{\Psi }}/{{\rm{\Psi }}}_{b}$ into homogeneous ψ_h and particular part ψ_p, i.e., $\psi ={\psi }_{h}+{\psi }_{p}$. We plug ψ into Equation (14) and separate out the homogeneous and particular parts to obtain the following dimensionless equations:

$$\begin{eqnarray}&&\displaystyle \frac{{\partial }^{2}{\psi }_{h}}{\partial {\varpi }^{2}}-\displaystyle \frac{1}{\varpi }\displaystyle \frac{\partial {\psi }_{h}}{\partial \varpi }+{\tau }^{2}\displaystyle \frac{{\partial }^{2}{\psi }_{h}}{\partial {\bar{z}}^{2}}=-4\sqrt{a}(\sqrt{a}{\varpi }^{2}+2\alpha ){\psi }_{h},\end{eqnarray} \tag{ 16 }$$

$$\begin{eqnarray}&&\displaystyle \frac{{\partial }^{2}{\psi }_{p}}{\partial {\varpi }^{2}}-\displaystyle \frac{1}{\varpi }\displaystyle \frac{\partial {\psi }_{p}}{\partial \varpi }=-(4a{\varpi }^{2}+8\sqrt{a}\alpha ){\psi }_{p}-\left(\displaystyle \frac{4{ab}}{{\psi }_{b}}{\varpi }^{2}+\beta \right).\end{eqnarray} \tag{ 17 }$$

The solution of Equation (16) is separable and given by ${\psi }_{h}(\varpi ,\bar{z})=s(\varpi )Z(\bar{z})$, which has been shown in SM18 to be given by

$$\begin{eqnarray}&&s(\varpi )={{cF}}_{0}(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})\qquad (\mathrm{with}\ a\gt 0),\end{eqnarray} \tag{ 18 }$$

where ${F}_{0}(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})$ represents the Coulomb function (Abramowitz & Stegun 1972) and $\kappa =\tfrac{{kR}}{2\sqrt{2}{a}^{1/4}},$ where the value of k is evaluated later. The z-part solution of Equation (16) is given by

$$\begin{eqnarray}&&Z(\bar{z})=\exp \left(-\displaystyle \frac{2\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right).\end{eqnarray} \tag{ 19 }$$

The homogeneous solution takes the following form:

$$\begin{eqnarray}&&{\psi }_{h}(\varpi ,\bar{z})=c\exp \left(-\displaystyle \frac{2\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right){F}_{0}(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2}).\end{eqnarray} \tag{ 20 }$$

The solution of the inhomogeneous Equation (17) is given by a power=series solution

$$\begin{eqnarray}\begin{array}{rcl}{\psi }_{p}(\varpi ) & = & -\displaystyle \frac{\beta }{8\alpha }+\displaystyle \frac{i\sqrt{a}{\varpi }^{2}}{2}{e}^{-i\sqrt{a}{\varpi }^{2}}\left(\displaystyle \frac{\beta }{8\alpha }-\displaystyle \frac{b}{{\psi }_{b}}\right)\\ & & \times \,\displaystyle \sum _{n=0}^{\infty }\displaystyle \frac{{F}_{2}^{1}\left(n+2,1;n+2+i\alpha ;\tfrac{1}{2}\right){\left(i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1+i\alpha )n!}.\end{array}\end{eqnarray} \tag{ 21 }$$

A similar but different homogeneous solution that is oscillatory along the z-direction has been used for the cases of both D-shaped and toroidally diverted laboratory plasma (Atanasiu et al. 2004). The general solution of the GS Equation (14) is given by ψ = ψ_h + ψ_p. Since $\psi (\varpi ,\bar{z})$ and ${\psi }^{* }(\varpi ,\bar{z})$, its complex conjugate, are the valid solutions of Equation (14), we construct a real solution of Equation (14) by redefining $\tfrac{{\psi }_{h}(\varpi ,\bar{z})+{\psi }_{h}^{* }(\varpi ,\bar{z})}{2}\to {\psi }_{h}(\varpi ,\bar{z})$ and $\tfrac{{\psi }_{p}(\varpi )+{\psi }_{p}^{* }(\varpi )}{2}\to {\psi }_{p}(\varpi )$, which leads to

$$\begin{eqnarray*}&&\begin{array}{l}\displaystyle \frac{\psi (\varpi ,\bar{z})+{\psi }^{* }(\varpi ,\bar{z})}{2}\to \psi (\varpi ,\bar{z})={\psi }_{C}^{X}\\ \ =\,\left\{\begin{array}{l}{\psi }_{C}^{C}={\psi }_{h}\quad (\mathrm{closed}\ \mathrm{field})\\ {\psi }_{C}^{O}={\psi }_{h}+{\psi }_{p}\quad (\mathrm{open}\ \mathrm{field})\end{array}\right..\end{array}\end{eqnarray*}$$

The solution ψ_h alone gives the closed field structure of the flux tube (SM18), which we denote as ${\psi }_{C}^{C};$ the general solution is a combination of ψ_h and ψ_p and we denote the open field flux tube structure as ${\psi }_{C}^{O}$.

The total flux function ${\psi }_{C}^{X}(\varpi ,\bar{z})$ is given by

$$\begin{eqnarray}&&{\psi }_{C}^{X}(\varpi ,\bar{z})=s(\varpi )Z(\bar{z})+{\psi }_{p}(\varpi ),\end{eqnarray} \tag{ 22 }$$

where s(ϖ) and ${\psi }_{p}(\varpi )$ are given by Equations (18) and (21), respectively. Now, ${\psi }_{C}^{X}(\varpi ,\bar{z})$ has to be zero at the axis (i.e., ϖ = 0) for all $\bar{z}$, to keep the field finite at the origin. Since s(0) = 0, which satisfies the BC 1 (Equation 10(a)), we obtain from Equation (22), ${\psi }_{p}(\varpi =0)=0$. From Equation (21), we have ${\psi }_{p}(\varpi =0)=-\tfrac{\beta }{8\alpha }$. Therefore, we obtain β = 0. From Equation (21), ${\psi }_{p}(\varpi )$ reduces to

$$\begin{eqnarray}&&\begin{array}{l}{\psi }_{p}(\varpi )=\displaystyle \frac{i\sqrt{a}b{\varpi }^{2}}{4{\psi }_{b}}\\ \ \times \,\left[{e}^{i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\displaystyle \frac{{F}_{2}^{1}\left(n+2,1;n+2-i\alpha ;\tfrac{1}{2}\right){\left(-i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1-i\alpha )n!}\right.\\ \ -\,\left.{e}^{-i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\displaystyle \frac{{F}_{2}^{1}\left(n+2,1;n+2+i\alpha ;\tfrac{1}{2}\right){\left(i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1+i\alpha )n!}\right].\end{array}\end{eqnarray} \tag{ 23 }$$

From Equations (2) and (22) we have

$$\begin{eqnarray}&&{B}_{z}(\varpi ,\bar{z})=\displaystyle \frac{{B}_{0}{\psi }_{b}}{2\sqrt{a}\varpi }\left(\displaystyle \frac{\partial {\psi }_{h}}{\partial \varpi }+\displaystyle \frac{\partial {\psi }_{p}}{\partial \varpi }\right),\end{eqnarray} \tag{ 24 }$$

and we obtain the explicit form for

$$\begin{eqnarray}&&\begin{array}{l}{B}_{z}(\varpi ,\bar{z})=\displaystyle \frac{{B}_{0}}{8\varpi }\\ \quad \times \,\left[\left\{1-\displaystyle \frac{{ib}}{2}\left(\displaystyle \frac{{F}_{2}^{1}(1,-i\alpha ,2-i\alpha ,-1)}{1-i\alpha }\right.\right.\right.\\ \quad -\left.\left.\,\displaystyle \frac{{F}_{2}^{1}(1,i\alpha ,2+i\alpha ,-1)}{1+i\alpha }\right)\right\}\\ \quad \cdot \,\exp \left(-\displaystyle \frac{2\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right)\displaystyle \frac{d}{d\varpi }\\ \quad [{F}_{0}(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})+{F}_{0}^{* }(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})]\\ \quad +\,{ib}\displaystyle \frac{d}{d\varpi }\left\{{\varpi }^{2}{e}^{i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\right.\\ \quad \displaystyle \frac{{F}_{2}^{1}(n+2,1,n+2-i\alpha ,1/2){\left(-i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1-i\alpha )n!}\\ \quad -\,{\varpi }^{2}{e}^{-i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\\ \quad \left.\left.\displaystyle \frac{{F}_{2}^{1}(n+2,1,n+2+i\alpha ,1/2){\left(-i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1+i\alpha )n!}\right\}\right].\end{array}\end{eqnarray} \tag{ 25 }$$

The steps to obtain Equation (25) are given in Appendix A. From Equations (2), (22) we have

$$\begin{eqnarray}&&{B}_{r}(\varpi ,\bar{z})=-\displaystyle \frac{{B}_{0}{\psi }_{b}\tau }{2\sqrt{a}}s(\varpi )Z^{\prime} (\bar{z}),\end{eqnarray} \tag{ 26 }$$

whose explicit form is given by

$$\begin{eqnarray}&&\begin{array}{l}{B}_{r}(\varpi ,\bar{z})=\displaystyle \frac{{B}_{0}\kappa }{2\sqrt{2}{a}^{1/4}\varpi }\exp \left(\displaystyle \frac{-2\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right)\\ \quad \times \,[{F}_{0}(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})+{F}_{0}^{* }(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})]\\ \quad \cdot \,\left[1-\displaystyle \frac{{ib}}{2}\left(\displaystyle \frac{{F}_{2}^{1}(1,-i\alpha ,2-i\alpha ,-1)}{1-i\alpha }\right.\right.\\ \quad -\,\left.\left.\displaystyle \frac{{F}_{2}^{1}(1,i\alpha ,2+i\alpha ,-1)}{1+i\alpha }\right)\right].\end{array}\end{eqnarray} \tag{ 27 }$$

From Equation (2) we obtain the toroidal component

$$\begin{eqnarray}&&{B}_{\phi }(\varpi ,\bar{z})=\displaystyle \frac{\sqrt{2}{B}_{0}{\alpha }^{1/2}{\psi }_{b}}{{a}^{1/4}\varpi }({\psi }_{h}+{\psi }_{p}),\end{eqnarray} \tag{ 28 }$$

whose explicit form is given by

$$\begin{eqnarray}&&\begin{array}{l}{B}_{\phi }(\varpi ,\bar{z})=\displaystyle \frac{{B}_{0}{\alpha }^{1/2}{a}^{-1/4}}{4\sqrt{2}}\\ \quad \times \,\left[\displaystyle \frac{1}{\varpi }\left\{1-\displaystyle \frac{{ib}}{2}\left(\displaystyle \frac{{F}_{2}^{1}(1,-i\alpha ,2-i\alpha ,-1)}{1-i\alpha }\right.\right.\right.\\ \quad -\,\left.\left.\displaystyle \frac{{F}_{2}^{1}(1,i\alpha ,2+i\alpha ,-1)}{1+i\alpha }\right)\right\}\\ \quad \cdot \,\exp \left(-\displaystyle \frac{-2\sqrt{2}\kappa \bar{z}}{\tau }\right)\\ \quad \times \,[{F}_{0}(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})+{F}_{0}^{* }(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})]\\ \quad +\,{ib}\varpi \left\{{e}^{i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\right.\\ \quad \displaystyle \frac{{F}_{2}^{1}(n+2,1,n+2-i\alpha ,1/2){\left(-i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1-i\alpha )n!}\\ \quad -\,{e}^{-i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\\ \quad \left.\left.\displaystyle \frac{{F}_{2}^{1}(n+2,1,n+2+i\alpha ,1/2){\left(-i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1+i\alpha )n!}\right\}\right].\end{array}\end{eqnarray} \tag{ 29 }$$

Applying BC 3 (Equation 10(c)) and using Equation (26), we find $s(\varpi =1)=0$. From BC 2 (Equations 10(b)) and 13(b) we find I₀ = 0. We assume that the external pressure from photosphere to transition region decreases exponentially as

$$\begin{eqnarray}&&{p}_{e}(z)={p}_{0}\exp (-2{kz}),\end{eqnarray} \tag{ 30 }$$

where k is the pressure-scale height, which is determined by the relation, $k=\tfrac{1}{2{z}_{t}}\mathrm{ln}\left(\tfrac{{p}_{0}}{{p}_{t}}\right)=3.405$ Mm⁻¹, where p₀ = 1.22 × 10⁵ dyne cm⁻², p_t = 0.148 dyne cm⁻², and z_t = 2 Mm. By matching the pressure-scale heights inside and outside the flux tube, we see that p₂(z) follows

$$\begin{eqnarray}&&{p}_{2}(z)={p}_{20}\exp (-2{kz}),\end{eqnarray} \tag{ 31 }$$

where p₂₀ is evaluated later. Finally, the expression of $p(r,z)$ from Equations (4), 13(a) is given by

$$\begin{eqnarray}&&p(r,z)=\displaystyle \frac{1}{4\pi }\left(\displaystyle \frac{a^{\prime} }{2}{{\rm{\Psi }}}^{2}+b^{\prime} {\rm{\Psi }}\right)+{p}_{20}\exp (-2{kz}),\end{eqnarray} \tag{ 32 }$$

whose explicit form is given by

$$\begin{eqnarray}\begin{array}{rcl}p(\varpi ,\bar{z}) & = & {B}_{0}^{2}\left[\left(\displaystyle \frac{{\psi }_{b}^{2}{s}^{2}(\varpi )}{8\pi }+{\bar{p}}_{20}\right){Z}^{2}(\bar{z})\right.\\ & & +\,\left(\displaystyle \frac{{\psi }_{b}^{2}s(\varpi ){\psi }_{p}(\varpi )}{4\pi }+\displaystyle \frac{b{\psi }_{b}s(\varpi )}{2\sqrt{2a}}\right)Z(\bar{z})\\ & & +\,\left.\left(\displaystyle \frac{{\psi }_{b}^{2}{\psi }_{p}^{2}}{8\pi }+\displaystyle \frac{b{\psi }_{b}{\psi }_{p}}{2\sqrt{2a}}\right)\right],\end{array}\end{eqnarray} \tag{ 33 }$$

where ${\bar{p}}_{20}={p}_{20}/{B}_{0}^{2}$, and $s(\varpi ),Z(\bar{z}),{\psi }_{p}(\varpi )$ are given by Equations (18), (19), (23), 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 Equation (1) yields

$$\begin{eqnarray}&&\begin{array}{l}{\left.-\displaystyle \frac{\partial p}{\partial r}\right|}_{r=R}{\left.+\displaystyle \frac{1}{4\pi }\left({B}_{r}\displaystyle \frac{\partial {B}_{r}}{\partial r}+{B}_{z}\displaystyle \frac{\partial {B}_{r}}{\partial z}\right)\right|}_{r=R}{\left.-\displaystyle \frac{\partial }{\partial r}\left(\displaystyle \frac{{B}^{2}}{8\pi }\right)\right|}_{r=R}\\ \quad +\,{j}_{\phi }(R){B}_{z}(R)-{j}_{z}(R){B}_{\phi }(R)=0.\end{array}\end{eqnarray} \tag{ 34 }$$

The sheet currents j_ϕ and j_z take the forms

$$\begin{eqnarray}&&{j}_{\phi }(r)={j}_{\phi s}\ \delta (r-R),\end{eqnarray} \tag{ 35a }$$

$$\begin{eqnarray}&&{j}_{z}(r)={j}_{{zs}}\ \delta (r-R).\end{eqnarray} \tag{ 35b }$$

Integrating Equation (34) with respect to r from $r=R-\epsilon$ to $r=R+\epsilon$ where is an infinitesimal positive quantity we obtain

$$\begin{eqnarray}&&\begin{array}{l}-{\displaystyle \int }_{R-\epsilon }^{R+\epsilon }\displaystyle \frac{\partial p}{\partial r}{dr}+\displaystyle \frac{1}{4\pi }\left({\displaystyle \int }_{R-\epsilon }^{R+\epsilon }{B}_{r}\displaystyle \frac{\partial {B}_{r}}{\partial r}{dr}+{\displaystyle \int }_{R-\epsilon }^{R+\epsilon }{B}_{z}\displaystyle \frac{\partial {B}_{r}}{\partial z}{dr}\right)\\ \quad -\,{\displaystyle \int }_{R-\epsilon }^{R+\epsilon }\displaystyle \frac{\partial }{\partial r}\left(\displaystyle \frac{{B}^{2}}{8\pi }\right){dr}+{\displaystyle \int }_{R-\epsilon }^{R+\epsilon }{j}_{\phi }(r){B}_{z}(r){dr}\\ \quad -\,{\displaystyle \int }_{R-\epsilon }^{R+\epsilon }{j}_{z}(r){B}_{\phi }(r){dr}=0,\end{array}\end{eqnarray} \tag{ 36 }$$

which leads to

$$\begin{eqnarray}&&\begin{array}{l}{p}_{i}(R,z)-{p}_{e}(z)+{j}_{\phi s}\ {B}_{z}(R)-{j}_{{zs}}\ {B}_{\phi }(R)\\ +\,\displaystyle \frac{1}{4\pi }{\left[{B}_{r}\displaystyle \frac{\partial {B}_{r}}{\partial r}+{B}_{z}\displaystyle \frac{\partial {B}_{r}}{\partial z}\right]}_{R}+\displaystyle \frac{{B}_{i}^{2}(R,z)-{B}_{e}^{2}(R,z)}{8\pi }=0,\end{array}\end{eqnarray} \tag{ 37 }$$

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}_{\phi s}$ and j_zs, we assume an infinitesimal current loop of vertical height L and radial extent $R-\epsilon$ to $R+\epsilon$ placed at the boundary of the flux tube (see Figure 2), and by applying the line integral along the loop, we obtain

$$\begin{eqnarray}&&{B}_{z}(R)L=4\pi L{\int }_{R-\epsilon }^{R+\epsilon }{j}_{\phi s}\delta (r-R){dr},\end{eqnarray} \tag{ 38a }$$

$$\begin{eqnarray}&&-{B}_{\phi }(R)L=4\pi L{\int }_{R-\epsilon }^{R+\epsilon }{j}_{{zs}}\delta (r-R){dr},\end{eqnarray} \tag{ 38b }$$

which implies

$$\begin{eqnarray}&&{j}_{\phi s}=\displaystyle \frac{{B}_{z}(R)}{4\pi }\end{eqnarray} \tag{ 39a }$$

$$\begin{eqnarray}&&{j}_{{zs}}=-\displaystyle \frac{{B}_{\phi }(R)}{4\pi }.\end{eqnarray} \tag{ 39b }$$

The total internal magnetic field is given by ${B}_{i}^{2}={B}_{r}^{2}\,+{B}_{\phi }^{2}+{B}_{z}^{2}$. Applying the BC 3 (Equation 10(c)) and B_e = 0, from Equation (37) we obtain

$$\begin{eqnarray}&&{p}_{i}(R,z)-{p}_{e}(z)+\displaystyle \frac{3{B}_{z}^{2}(R,z)}{8\pi }+\displaystyle \frac{3{B}_{\phi }^{2}(R,z)}{8\pi }=0.\end{eqnarray} \tag{ 40 }$$

By expanding Equation (40) we obtain

$$\begin{eqnarray}&&\begin{array}{l}\displaystyle \frac{{\psi }_{b}^{2}}{3}\left(1+\displaystyle \frac{6\alpha }{\sqrt{a}}\right)+\displaystyle \frac{2}{3}b{\psi }_{b}+\displaystyle \frac{{\psi }_{b}^{2}}{4a}\\ \ \times \,{\left[\displaystyle \frac{1}{{\varpi }^{2}}{\left(s^{\prime} (\varpi )\exp \left(-\displaystyle \frac{2\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right)+{\psi }_{p}^{{\prime} }(\varpi )\right)}^{2}\right]}_{\varpi =1}\\ =\,\bar{p}\exp \left(-\displaystyle \frac{4\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right),\end{array}\end{eqnarray} \tag{ 41 }$$

where ${\psi }_{b}={\psi }_{p}(\varpi =1)$ and $\bar{p}=\tfrac{8\pi ({p}_{0}-{p}_{20})}{3{B}_{0}^{2}}$. By equating the coefficients of $\exp \left(-\tfrac{2\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right)$, $\exp \left(-\tfrac{4\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right)$ and the constant quantity between both sides of Equation (41), we obtain

$$\begin{eqnarray}&&{\left(\displaystyle \frac{1}{{\varpi }^{2}}s^{\prime} (\varpi ){\psi }_{p}^{{\prime} }(\varpi )\right)}_{\varpi =1}=0,\end{eqnarray} \tag{ 42a }$$

$$\begin{eqnarray}&&\displaystyle \frac{{\psi }_{b}^{2}}{4a}{\left(\displaystyle \frac{s{{\prime} }^{2}(\varpi )}{{\varpi }^{2}}\right)}_{\varpi =1}=\bar{p},\end{eqnarray} \tag{ 42b }$$

$$\begin{eqnarray}&&\displaystyle \frac{{\psi }_{b}^{2}}{3}\left(1+\displaystyle \frac{6\alpha }{\sqrt{a}}\right)+\displaystyle \frac{2b{\psi }_{b}}{3}+\displaystyle \frac{{\psi }_{b}^{2}}{4a}{\left[\displaystyle \frac{{\psi }_{p}^{{\prime} 2}(\varpi )}{{\varpi }^{2}}\right]}_{\varpi =1}=0.\end{eqnarray} \tag{ 42c }$$

From Equations 42(a), 42(b) we find that ${\left(\tfrac{{\psi }_{p}^{{\prime} }(\varpi )}{\varpi }\right)}_{\varpi =1}=0$, as $\bar{p}\ne 0$, hence, from Equation 42(c),

$$\begin{eqnarray}&&b=-\displaystyle \frac{{\psi }_{b}}{2}\left(1+\displaystyle \frac{6\alpha }{\sqrt{a}}\right).\end{eqnarray} \tag{ 43 }$$

We summarize the set of equations we need to solve numerically for the open field model as follows:

$$\begin{eqnarray}&&s(\varpi =1)=0,\end{eqnarray} \tag{ 44a }$$

$$\begin{eqnarray}&&{\left(\displaystyle \frac{{\psi }_{p}^{{\prime} }(\varpi )}{\varpi }\right)}_{\varpi =1}=0,\end{eqnarray} \tag{ 44b }$$

$$\begin{eqnarray}&&{\left(\displaystyle \frac{s{{\prime} }^{2}(\varpi )}{{\varpi }^{2}}\right)}_{\varpi =1}=\displaystyle \frac{4a\bar{p}}{{\psi }_{b}^{2}},\end{eqnarray} \tag{ 44c }$$

$$\begin{eqnarray}&&b=-\displaystyle \frac{{\psi }_{b}}{2}\left(1+\displaystyle \frac{6\alpha }{\sqrt{a}}\right),\end{eqnarray} \tag{ 44d }$$

$$\begin{eqnarray}&&{\int }_{0}^{1}p(\varpi ,{\bar{z}}_{t})d\varpi ={p}_{t}.\end{eqnarray} \tag{ 44e }$$

The explicit forms of Equations 44(a)–44(e) are given by Equations (79)–(83) in Appendix B. The five Equations (79)–(83) consist of seven unknown variables, $\{a,\alpha ,b,{\psi }_{b},\bar{p},R,{B}_{0}\};$ so there is a unique solution to the Coulomb function open field model for a given pair of unknown variables. Equations (79) and (80) contain three variables, a, α, and R, and we use these two equations to obtain a(R) and α(R). From Equation (83), we calculate ${\psi }_{b}(R,{B}_{0},b)$ and then find $b(R,{B}_{0})$ from Equation (82), and hence ${\psi }_{b}(R,{B}_{0})$, and then evaluate $\bar{p}(R,{B}_{0})$ from Equation (81). As a result, the complete solution for the open field flux tube depends only on R and B₀, which are the free parameters of the model. Per BC1–BC5 (Equations 10(a)–10(e)), which are used for the open field flux tube model, the magnetic field component at the boundary of the flux tube is given by

$$\begin{eqnarray}&&{B}_{\phi }(\varpi =1,\bar{z})=\displaystyle \frac{\sqrt{2\alpha }{B}_{0}{\psi }_{b}^{2}}{{a}^{1/4}}.\end{eqnarray} \tag{ 45 }$$

If we demand additionally that ${B}_{\phi }(\varpi =1,\bar{z})=0$, then from Equation 44(d), ψ_b = 0; also, Equation 44(d) gives b = 0. Therefore, from Equation (23) we obtain

$$\begin{eqnarray}&&{\psi }_{p}(\varpi )=0.\end{eqnarray} \tag{ 46 }$$

This represents the solution of the homogeneous part ψ_h of GSE, which has been discussed in SM18, and is applicable for the closed field flux tube model. Therefore, we need to solve Equations 44(a), (c), (e), which are given in their explicit forms by Equations (79), (81), (83) (with b = ψ_b = 0), numerically to find the parameters a, α, and $\bar{p}$ in terms of $\{R,{B}_{0}\}$, which are the free parameters of the closed field model. The formularies of the derived functions for the Coulomb function helical flux tube model are summarized in Table 6. We discuss the various configurations of the Coulomb function open and closed field structure of flux tubes in Section 5.1.

Zoom In Zoom Out Reset image size

The basic formulation of a self-similar model of a flux tube is based on Schlüter & Temesváry (1958; hereafter ST58). For an axially symmetric cylindrical geometry (r, ϕ, z), where ϕ is ignorable, the magnetic field components are given by Equation (2). 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., ${\rm{\Psi }}(r,z)={\rm{\Psi }}(\xi )$ (ST58). We define the dimensionless parameters (in the LHS) by introducing the scaling relations

$$\begin{eqnarray}\begin{array}{rcl}\varpi & = & r/R,\ \ \tau =R{/z}_{0},\ \ \psi ={{\rm{\Psi }}/{\rm{\Psi }}}_{b},\ \ {\psi }_{b}=\displaystyle \frac{{{\rm{\Psi }}}_{b}}{{B}_{0}{z}_{0}^{2}},\\ {\bar{p}}_{1} & = & {p}_{1}{/B}_{0}^{2},\ \ {\bar{I}}_{p}=\displaystyle \frac{{I}_{p}}{{B}_{0}R},\ \ \bar{\chi }=\chi {z}_{0}^{2},\end{array}\end{eqnarray} \tag{ 47 }$$

where, Ψ_b, R, B₀ are the boundary flux, radius, and magnetic field strength at the center of the flux tube, respectively, and $\bar{z}=z/{z}_{0}$, where z₀ is a constant length. From ST58, the self-similar parameter ξ is defined by

$$\begin{eqnarray}&&\xi =\zeta (\bar{z})\varpi ,\end{eqnarray} \tag{ 48 }$$

which describes the radial size distribution of the flux tube with height from the base. Plugging in Equation (48), we can rewrite the GS Equation (5) in the form

$$\begin{eqnarray}&&\begin{array}{l}\displaystyle \frac{{\psi }_{b}^{2}}{2{\tau }^{2}}\displaystyle \frac{d}{d\xi }{\left(\displaystyle \frac{d\psi }{d\xi }\right)}^{2}\zeta {{\prime} }^{2}(\bar{z})+\displaystyle \frac{{\psi }_{b}^{2}}{{\tau }^{2}}\displaystyle \frac{1}{\xi }{\left(\displaystyle \frac{d\psi }{d\xi }\right)}^{2}\zeta ^{\prime} (\bar{z})\zeta ^{\prime\prime} (\bar{z})\\ \ +\,\displaystyle \frac{{\psi }_{b}^{2}}{2{\tau }^{4}}\displaystyle \frac{d}{d\xi }{\left(\displaystyle \frac{1}{\xi }\displaystyle \frac{d\psi }{d\xi }\right)}^{2}{\zeta }^{4}(\bar{z})+\displaystyle \frac{1}{2{\xi }^{2}}\displaystyle \frac{d{\bar{I}}_{p}^{2}}{d\xi }{\zeta }^{2}=-4\pi \displaystyle \frac{\partial {\bar{p}}_{1}}{\partial \xi },\end{array}\end{eqnarray} \tag{ 49 }$$

and the z-part of GS Equation 3(b) gives the expression of ρ, Equation (7), which is self-consistent for both the Coulomb function and self-similar models. We define a quantity that is called the magnetic shape function given by

$$\begin{eqnarray}&&{D}_{X}(\xi )=\displaystyle \frac{1}{\xi }\displaystyle \frac{d\psi }{d\xi }.\end{eqnarray} \tag{ 50 }$$

Plugging Equation (50) into Equation (49) we obtain

$$\begin{eqnarray}&&\begin{array}{l}\displaystyle \frac{{\psi }_{b}^{2}}{{\tau }^{2}}\xi {D}_{X}^{2}(\xi )\zeta \zeta ^{\prime\prime} +\displaystyle \frac{{\psi }_{b}^{2}}{2{\tau }^{2}}\displaystyle \frac{d}{d\xi }[{\xi }^{2}{D}_{X}^{2}(\xi )]\zeta {{\prime} }^{2}\\ \quad +\,\displaystyle \frac{{\psi }_{b}^{2}}{2{\tau }^{4}}\displaystyle \frac{d}{d\xi }\left({D}_{X}^{2}(\xi )\right){\zeta }^{4}+\displaystyle \frac{1}{2{\xi }^{2}}\displaystyle \frac{d{\bar{I}}_{p}^{2}}{d\xi }{\zeta }^{2}=-4\pi \displaystyle \frac{\partial {\bar{p}}_{1}}{\partial \xi },\end{array}\end{eqnarray} \tag{ 51 }$$

and integrating Equation (51) with respect to ξ from 0 to $\infty$ we write

$$\begin{eqnarray}&&\begin{array}{l}\displaystyle \frac{{\psi }_{b}^{2}}{{\tau }^{2}}\zeta \zeta ^{\prime\prime} {\displaystyle \int }_{0}^{\infty }\xi {D}_{X}^{2}(\xi )d\xi +\displaystyle \frac{{\psi }_{b}^{2}\zeta {{\prime} }^{2}}{2{\tau }^{2}}{\left[{\xi }^{2}{D}_{X}^{2}(\xi )\right]}_{\xi =0}^{\infty }\\ \quad +\,\displaystyle \frac{{\psi }_{b}^{2}{\zeta }^{4}}{2{\tau }^{4}}[{D}_{X}^{2}(\infty )-{D}_{X}^{2}(0)]\\ \quad +\,\displaystyle \frac{{\zeta }^{2}}{2}{\displaystyle \int }_{0}^{\infty }\displaystyle \frac{1}{{\xi }^{2}}\displaystyle \frac{d{\bar{I}}_{p}^{2}}{d\xi }d\xi =-4\pi {\displaystyle \int }_{0}^{\infty }\displaystyle \frac{\partial {\bar{p}}_{1}}{\partial \xi }d\xi .\end{array}\end{eqnarray} \tag{ 52 }$$

Following ST58, we define

$$\begin{eqnarray}&&{y}^{2}(\bar{z})=\displaystyle \frac{{\psi }_{b}{D}_{0}}{\tau }{\zeta }^{2}(\bar{z}),\end{eqnarray} \tag{ 53 }$$

where, $y(\bar{z})\equiv {\left(\tfrac{{B}_{z}(0,\bar{z})}{{B}_{0}}\right)}^{1/2}$, ${B}_{0}\equiv {B}_{z}(0,0)$, and ${D}_{0}\equiv {D}_{X}(\xi =0)$. Next, using Equations (52) and (53), we obtain

$$\begin{eqnarray}&&\begin{array}{l}\displaystyle \frac{{\psi }_{b}}{\tau }\displaystyle \frac{{yy}^{\prime\prime} }{{D}_{0}}{\displaystyle \int }_{0}^{\infty }\xi {D}_{X}^{2}(\xi )d\xi +\displaystyle \frac{{\psi }_{b}}{\tau }\displaystyle \frac{y{{\prime} }^{2}}{2{D}_{0}}{\left[{\xi }^{2}{D}_{X}^{2}(\xi )\right]}_{\xi =0}^{\infty }\\ \quad +\,\displaystyle \frac{{y}^{4}}{2{D}_{0}^{2}{\tau }^{2}}[{D}_{X}^{2}(\infty )-{D}_{0}^{2}]\\ \quad +\,\displaystyle \frac{\tau }{{\psi }_{b}}\displaystyle \frac{{y}^{2}}{2{D}_{0}}{\displaystyle \int }_{0}^{\infty }\displaystyle \frac{1}{{\xi }^{2}}\displaystyle \frac{d{\bar{I}}_{p}^{2}}{d\xi }d\xi =-4\pi {\displaystyle \int }_{0}^{\infty }\displaystyle \frac{\partial {\bar{p}}_{1}}{\partial \xi }d\xi .\end{array}\end{eqnarray} \tag{ 54 }$$

To solve Equation (54), we need to specify the functional form of p₁, I_p, and ${D}_{X}(\xi )$ to study the flux tube model with twisted magnetic field under the similarity assumption. The functional form of ${p}_{1}=\tfrac{f}{2}{{\rm{\Psi }}}^{2}$ is taken from Osherovitch (1982), where f is the shape function parameter, and the poloidal current, I_p, defined by Yun (1971) and Osherovitch (1979), and motivated from the observations of Stepanov (1965). Hence, the form of gas pressure p and poloidal current I_p are taken to be

$$\begin{eqnarray}&&p={p}_{c}\exp (-2{kz})+\displaystyle \frac{f}{2}{{\rm{\Psi }}}^{2},\end{eqnarray} \tag{ 55 }$$

$$\begin{eqnarray}&&{I}_{p}^{2}={{\rm{\Psi }}}_{b}^{2}\chi {\xi }^{4}{D}_{X}^{2}(\xi ),\end{eqnarray} \tag{ 56 }$$

for the positivity of ρ(z) at all z, which is given by Equation (7). Here, ${p}_{2}(z)={p}_{c}\exp (-2{kz})$ denotes the gas pressure at the flux tube axis, with p_c being the pressure at the center of the flux tube, and $\chi ={\left(\tfrac{{B}_{\phi }}{{{rB}}_{z}}\right)}^{2}$ being a constant pitch angle parameter. We deviate from (Yun 1971; Osherovitch 1982) by employing the extra term, ${p}_{c}\exp (-2{kz})$, with p₁ in Equation (55), to maintain the hydrostatic vertical pressure balance condition under the influence of solar gravity, with two options for the shape function ${D}_{X}(\xi )$ specified by

$$\begin{eqnarray}&&\begin{array}{l}{D}_{X}(\xi )\\ =\,\left\{\begin{array}{l}{D}_{G}(\xi )={D}_{G0}\exp (-{\xi }^{{n}_{G}});\quad ({n}_{G}\gt 0):\quad \mathrm{Generalized}\ \mathrm{Gaussian}\\ {D}_{P}(\xi )={D}_{P0}{\left(1+\xi \right)}^{-{n}_{P}};\quad ({n}_{P}\gt 2):\quad \mathrm{Power}\ \mathrm{law}\end{array}\right.\end{array}\end{eqnarray} \tag{ 57 }$$

where

$$\begin{eqnarray}&&{D}_{G0}=\displaystyle \frac{{n}_{G}}{{\rm{\Gamma }}(2/{n}_{G})},\end{eqnarray} \tag{ 58a }$$

$$\begin{eqnarray}&&{D}_{P0}=({n}_{P}-1)({n}_{P}-2).\end{eqnarray} \tag{ 58b }$$

We see that both the shape functions (Equation (57)) vanish asymptotically at $\xi \to \infty ;$ hence from Equation (54) we obtain

$$\begin{eqnarray}&&\begin{array}{l}\displaystyle \frac{{\psi }_{b}}{\tau }\displaystyle \frac{{yy}^{\prime\prime} }{{D}_{0}}{\displaystyle \int }_{0}^{\infty }\xi {D}_{X}^{2}(\xi )d\xi -\displaystyle \frac{{y}^{4}}{2{\tau }^{2}}\\ \quad +\,\displaystyle \frac{\tau }{{\psi }_{b}}\displaystyle \frac{{y}^{2}}{{D}_{0}}{\displaystyle \int }_{0}^{\infty }\displaystyle \frac{1}{2{\xi }^{2}}\displaystyle \frac{d{\bar{I}}_{p}^{2}}{d\xi }d\xi =-4\pi {\displaystyle \int }_{0}^{\infty }\displaystyle \frac{\partial {\bar{p}}_{1}}{\partial \xi }d\xi .\end{array}\end{eqnarray} \tag{ 59 }$$

Next, we evaluate the integrals of Equation (59) for both generalized Gaussian and power-law shape functions. We will see later that, from Equation (66), the flux function ${\psi }_{P}(\xi )$ varies as ${\xi }^{2-{n}_{P}}$, in the domain $0\lt \xi \lt \infty ;$ therefore ψ_P(ξ) will converge to a finite value at $\xi \to \infty$, if n_p > 2. The results of the integrals are provided in Table 2.

Table 2.  Expression of the Integrals in Equation (59) and λ_X in Equation (60), for Generalized Gaussian (X = G), and Power-law (X = P) Shape Function, Equation (57)

Functions	Shape Function ${D}_{G}(\xi )\quad ({n}_{G}\gt 0)$	Shape Function ${D}_{P}(\xi )\quad ({n}_{P}\gt 2)$
${\int }_{0}^{\infty }\xi {D}_{X}^{2}(\xi )d\xi$	$\tfrac{{n}_{G}}{{\rm{\Gamma }}(2/{n}_{G}){2}^{2/{n}_{G}}}$	$\tfrac{({n}_{P}-1){\left({n}_{P}-2\right)}^{2}}{2(2{n}_{P}-1)}$
${\int }_{0}^{\infty }\tfrac{1}{2{\xi }^{2}}\tfrac{d{\bar{I}}_{p}^{2}}{d\xi }d\xi$	$\tfrac{\bar{\chi }{\psi }_{b}^{2}{n}_{G}}{{2}^{2/{n}_{G}}{\tau }^{2}{\rm{\Gamma }}(2/{n}_{G})}$	$\tfrac{\bar{\chi }{\psi }_{b}^{2}{({n}_{P}-2)}^{2}({n}_{P}-1)}{2{\tau }^{2}(2{n}_{P}-1)}$
${\int }_{0}^{\infty }\tfrac{\partial {\bar{p}}_{1}}{\partial \xi }d\xi$	$\bar{f}{\psi }_{b}^{2}/2$	$\bar{f}{\psi }_{b}^{2}/2$
${\lambda }_{X}(n)$	${2}^{2/{n}_{G}-1}$	$\tfrac{2{n}_{P}-1}{{n}_{P}-2}$

Download table as:  ASCIITypeset image

Using the values of the integrals from Table 2 and redefining $y^{\prime} (\bar{z}=0)\equiv {y}_{0}^{{\prime} }$, we reduce Equation (59) to the form

$$\begin{eqnarray}\begin{array}{rcl}\displaystyle \frac{{dy}}{d\bar{z}} & = & \left[\displaystyle \frac{{\lambda }_{X}(n)}{2{\psi }_{b}\tau }({y}^{4}-1)-2\bar{\chi }({y}^{2}-1)\right.\\ & & -\,{\left.8\pi {\psi }_{b}\bar{f}{\lambda }_{X}(n)\tau \mathrm{ln}(y)+\displaystyle \frac{{\bar{B^{\prime} }}_{z0}^{2}}{4}\right]}^{1/2},\end{array}\end{eqnarray} \tag{ 60 }$$

in which ${\bar{B^{\prime} }}_{z0}=\tfrac{{B}_{z0}^{{\prime} }{z}_{0}}{{B}_{0}}$, $\bar{\chi }=\chi {z}_{0}^{2}$, and $\bar{f}={{fz}}_{0}^{4}$. The form of Equation (60) 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 2. Here we have used the notation ${B}_{z0}^{{\prime} }={B}_{z}^{{\prime} }(0,0)$, which represents the vertical gradient of ${B}_{z}(0,z)$ at the center. From Equation (60), we have the following integral relation:

$$\begin{eqnarray}&&\bar{z}(y)={\int }_{1}^{y}\displaystyle \frac{{dy}^{\prime} }{G(y^{\prime} )},\end{eqnarray} \tag{ 61 }$$

where the function G(y) is given by

$$\begin{eqnarray}\begin{array}{rcl}G(y) & = & \left[\displaystyle \frac{{\lambda }_{X}(n)}{2{\psi }_{b}\tau }({y}^{4}-1)-2\bar{\chi }({y}^{2}-1)\right.\\ & & -\,{\left.8\pi {\psi }_{b}\bar{f}{\lambda }_{X}(n)\tau \mathrm{ln}(y)+\displaystyle \frac{{\bar{B^{\prime} }}_{z0}^{2}}{4}\right]}^{1/2}.\end{array}\end{eqnarray} \tag{ 62 }$$

We evaluate the integral (61) numerically, which gives $\bar{z}=\bar{z}(y)$. Thereafter, inverting the function between $\bar{z}$ and y, we evaluate $y=y(\bar{z})$. From Equations (48), (53) we obtain

$$\begin{eqnarray}&&\xi =\sqrt{\displaystyle \frac{\tau }{{\psi }_{b}{D}_{0}}}\varpi y(\bar{z}).\end{eqnarray} \tag{ 63 }$$

Using the similarity assumption, ${B}_{z}(r,z)=\tfrac{{B}_{0}{\psi }_{b}}{\tau }{\zeta }^{2}(z){D}_{X}(\xi )$ (ST58), and Equations (2), (48), (63) we calculate the magnetic field components, representing the most general self-similar solution, to be,

$$\begin{eqnarray}&&{B}_{z}(\varpi ,\bar{z})=\displaystyle \frac{{B}_{0}}{{D}_{0}}{y}^{2}(\bar{z}){D}_{X}(\xi ),\end{eqnarray} \tag{ 64a }$$

$$\begin{eqnarray}&&{B}_{r}(\varpi ,\bar{z})=-\displaystyle \frac{{B}_{0}\varpi }{{D}_{0}}y(\bar{z})y^{\prime} (\bar{z}){D}_{X}(\xi ),\end{eqnarray} \tag{ 64b }$$

$$\begin{eqnarray}&&{B}_{\phi }(\varpi ,\bar{z})=\displaystyle \frac{\sqrt{\bar{\chi }}{B}_{0}}{{D}_{0}}\varpi {y}^{2}(\bar{z}){D}_{X}(\xi ).\end{eqnarray} \tag{ 64c }$$

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

$$\begin{eqnarray}&&{\psi }_{S}^{O}(\xi )={\int }_{0}^{\xi }\xi ^{\prime} {D}_{X}(\xi ^{\prime} )d\xi ^{\prime} .\end{eqnarray} \tag{ 65 }$$

Employing Equations (50), (57), we obtain the open flux function for generalized Gaussian, ψ_G, and power-law, ψ_P, models to be given by

$$\begin{eqnarray}&&{\psi }_{S}^{O}(\xi )=\left\{\begin{array}{l}{\psi }_{G}(\xi )=1-\displaystyle \frac{{\rm{\Gamma }}(2/{n}_{G},{\xi }^{{n}_{G}})}{{\rm{\Gamma }}(2/{n}_{G})};\quad ({n}_{G}\gt 0):\quad \mathrm{Generalized}\ \mathrm{Gaussian}\\ {\psi }_{P}(\xi )=1-{\left(1+\xi \right)}^{1-{n}_{P}}\left(1+\xi ({n}_{P}-1)\right);\quad ({n}_{P}\gt 2):\quad \mathrm{Power}\ \mathrm{law}.\end{array}\right.\end{eqnarray} \tag{ 66 }$$

From Equation (66) it is seen that ${\psi }_{G}(\xi )$ and ${\psi }_{P}(\xi )$ converge to unity for $\xi \to \infty$. 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 that can make a partition with the external solar atmosphere. In other words, the self-similar flux tube is embedded in a continuous magnetic medium that has a maximum field strength at the base of the axis of the flux tube and the 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 refer to this radius as R₉₀. The total flux is zero at the axis and it increases asymptotically with r. The explicit forms of the magnetic field components, obtained from Equations 64(a), (b), (c), using Equations (57), 58(a), (b), (63), are

$$\begin{eqnarray}&&\begin{array}{l}{B}_{z}(\varpi ,\bar{z})\,=\,\left\{\begin{array}{l}{B}_{0}{y}^{2}(\bar{z})\exp \left[-{\left(\sqrt{\displaystyle \frac{\tau {\rm{\Gamma }}(2/{n}_{G})}{{n}_{G}{\psi }_{b}}}y(\bar{z})\varpi \right)}^{{n}_{G}}\right],\quad ({n}_{G}\gt 0):\quad \mathrm{Generalized}\ \mathrm{Gaussian}\\ {B}_{0}{y}^{2}(\bar{z}){\left[1+\sqrt{\displaystyle \frac{\tau }{({n}_{P}-1)({n}_{P}-2){\psi }_{b}}}y(\bar{z})\varpi \right]}^{-{n}_{P}},\quad ({n}_{P}\gt 2):\quad \mathrm{Power}\ \mathrm{law}\end{array}\right.\end{array}\end{eqnarray} \tag{ 67 }$$

$$\begin{eqnarray}&&{B}_{r}(\varpi ,\bar{z})=\left\{\begin{array}{l}-{B}_{0}y(\bar{z})y^{\prime} (\bar{z})\varpi \exp \left[-{\left(\sqrt{\displaystyle \frac{\tau {\rm{\Gamma }}(2/{n}_{G})}{{n}_{G}{\psi }_{b}}}y(\bar{z})\varpi \right)}^{{n}_{G}}\right],\quad ({n}_{G}\gt 0):\quad \mathrm{Generalized}\ \mathrm{Gaussian}\\ -{B}_{0}y(\bar{z})y^{\prime} (\bar{z})\varpi {\left[1+\sqrt{\displaystyle \frac{\tau }{({n}_{P}-1)({n}_{P}-2){\psi }_{b}}}y(\bar{z})\varpi \right]}^{-{n}_{P}},\quad ({n}_{P}\gt 2):\quad \mathrm{Power}\ \mathrm{law}\end{array}\right.\end{eqnarray} \tag{ 68 }$$

$$\begin{eqnarray}&&{B}_{\phi }(\varpi ,\bar{z})=\left\{\begin{array}{l}{B}_{0}\sqrt{\bar{\chi }}{y}^{2}(\bar{z})\varpi \exp \left[-{\left(\sqrt{\tfrac{\tau {\rm{\Gamma }}(2/{n}_{G})}{{n}_{G}{\psi }_{b}}}y(\bar{z})\varpi \right)}^{{n}_{G}}\right],\quad ({n}_{G}\gt 0):\quad \mathrm{Generalized}\ \mathrm{Gaussian}\\ {B}_{0}\sqrt{\bar{\chi }}{y}^{2}(\bar{z})\varpi {\left[1+\sqrt{\tfrac{\tau }{({n}_{P}-1)({n}_{P}-2){\psi }_{b}}}y(\bar{z})\varpi \right]}^{-{n}_{P}},\quad ({n}_{P}\gt 2):\quad \mathrm{Power}\ \mathrm{law}.\end{array}\right.\end{eqnarray} \tag{ 69 }$$

The magnetic field components ${B}_{r}(\varpi ,\bar{z})$ and ${B}_{\phi }(\varpi ,\bar{z})$ for the self-similar model follow the BCs (1, 2, 3; Equations 10(a), (b), (c)) for $R=\infty$. ${B}_{z}(\varpi ,\bar{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 (Equation 10(d)) at z = 0, for $r\to \infty$, so

$$\begin{eqnarray}&&p(r\to \infty ,0)={p}_{0}.\end{eqnarray} \tag{ 70 }$$

From Equation (66), we see that the flux function for both a generalized Gaussian and power law converges to unity at $\varpi \to \infty$, i.e.,

$$\begin{eqnarray}&&{\psi }_{S}^{O}(\varpi \to \infty )=1.\end{eqnarray} \tag{ 71 }$$

Using Equations (55), (70), we obtain

$$\begin{eqnarray}&&\bar{f}=2({\bar{p}}_{0}-{\bar{p}}_{c}),\end{eqnarray} \tag{ 72 }$$

and the explicit form of $p(\varpi ,\bar{z})$ for both the generalized Gaussian and the power-law models is given by

$$\begin{eqnarray}&&p(\varpi ,\bar{z})={B}_{0}^{2}\left(\displaystyle \frac{\bar{f}}{2}{\psi }^{2}+{\bar{p}}_{c}{e}^{-2\bar{k}\bar{z}}\right),\end{eqnarray} \tag{ 73 }$$

where ${\bar{p}}_{0}={p}_{0}/{B}_{0}^{2}$ and ${\bar{p}}_{c}={p}_{c}/{B}_{0}^{2}$. The formularies of the derived functions for the self-similar flux tube model are summarized in Table 6. A flowchart of the solutions to the Coulomb function and self-similar models is shown in Figure 1.

5.1. Coulomb Function Helical Flux Tube Model

This magnetohydrostatic Coulomb function helical flux tube model consists of the free parameters R and B₀ and its functional dependence through a(R), α(R), κ(R), $b(R,{B}_{0})$, ${\psi }_{b}(R,{B}_{0})$, and $\bar{p}(R,{B}_{0})$. We choose the parameter ranges, 1 kG ≤ B₀ ≤ 1.5 kG and 100 km ≤ R ≤ 180 km, consistent with the observations of MBP size and field strength distributions (Utz et al. 2009, 2013). In Table 3, we show the solutions for combinations of the free parameters $\{R,{B}_{0}\}$, where we notice the following trends:

• 1.  
  The boundary flux ψ_b decreases with R for the same B₀, and with B₀ for the same R within the parameter space of runs $C1-C21$.
• 2.  
  Due to the pressure balance at the boundary of the flux tube, $\bar{p}$ increases with R for the same B₀, but there is no fixed trend with B₀ for the same R within the parameter space of runs $C1-C21$.

As an example, we show the solution of ${\psi }_{C}^{O}$, and the magnetic and thermodynamic structure of the flux tube for run C4. The radial variation of the solution of ${\psi }_{C}^{O}$, magnetic components, and pressure inside the flux tube are shown in Figures 3–5, respectively. Examples of 3D the configuration of the magnetic field lines for open and closed field are shown in Figures 6 and 7 for runs C4 and C10. A 2D vertical projection of the magnetic field lines for ${\psi }_{C}^{O}$ inside the flux tube along r − z plane is shown in Figure 8. The density inside the flux tube is constant along the radial direction, but it decreases along z, whereas the temperature varies along the r-direction and is nearly constant along the z-direction at the axis. The vertical variations of B_z, p, and ρ are shown in Figure 9. Conclusions from figures and tables are discussed in Section 7.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 3.  Numerical Values of the Different Parameters Obtained from the Coulomb Function Open Field Flux Tube Model for Different Combinations of R and B₀ Are Shown; the Units of the Various Quantities Are in the Square Brackets at the Top

Run #	B0 [kG]	R [km]	${\psi }_{b}$ $[{10}^{-3}]$	a	α $[{10}^{-2}]$	κ [106]	b $[{10}^{-3}]$	$\bar{p}$
C1	1	100	2.57	9.390	2.74	4.85	1.350	0.109
C2	1.2	100	2.22	9.390	2.74	4.85	1.170	0.105
C3	1.5	100	1.80	9.390	2.74	4.85	0.949	0.104
C4	1	120	1.92	9.388	2.54	6.99	1.008	0.159
C5	1.2	120	1.62	9.388	2.54	6.99	0.849	0.165
C6	1.5	120	1.31	9.388	2.54	6.99	0.692	0.163
C7	1	130	1.69	9.389	2.43	8.21	0.880	0.184
C8	1.2	130	1.42	9.389	2.43	8.21	0.744	0.181
C9	1.5	130	1.15	9.389	2.43	8.21	0.603	0.182
C10	1	140	1.50	9.383	2.31	9.52	0.783	0.205
C11	1.2	140	1.25	9.383	2.31	9.52	0.661	0.204
C12	1.5	140	1.00	9.383	2.31	9.52	0.535	0.208
C13	1	150	1.38	9.378	2.18	10.93	0.723	0.233
C14	1.2	150	1.16	9.378	2.18	10.93	0.606	0.235
C15	1.5	150	0.94	9.378	2.18	10.93	0.492	0.237
C16	1	160	1.31	9.388	1.98	12.43	0.665	0.276
C17	1.2	160	1.11	9.388	1.98	12.43	0.577	0.276
C18	1.5	160	0.89	9.388	1.98	12.43	0.465	0.279
C19	1	180	1.14	9.395	1.72	15.73	0.587	0.402
C20	1.2	180	0.96	9.395	1.72	15.73	0.497	0.407
C21	1.5	180	0.78	9.395	1.72	15.73	0.405	0.409

Download table as:  ASCIITypeset image

5.2. Self-similar Model

The self-similar model we developed consists of the dimensionless parameters ${\psi }_{b},{\bar{B^{\prime} }}_{z0}$, $\bar{f}$, and $\bar{\chi }$, which are the functions of the input parameter set $\{{{\rm{\Psi }}}_{b},{B}_{0},{B}_{z0}^{{\prime} },{p}_{c}$, $\chi \}$. 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, ${{\rm{\Psi }}}_{b}={10}^{17}\mbox{--}{10}^{18}$ Mx (Zhang et al. 1998; Hagenaar et al. 1999; Guglielmino et al. 2011), B₀ = 1–2 kG (Zhang et al. 1998), and ${B}_{z0}^{{\prime} }$ in the range 1–2 G km⁻¹ (Wittmann 1974; Pahlke & Wiehr 1990; Balthasar & Schmidt 1993), ${p}_{c}\lt {p}_{0}$ (Shelyag et al. 2010 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 Section 4 that for the power-law profile, the flux function converges to a 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 $\bar{\chi }=0.01\mbox{--}100$ for different combinations of the other parameter sets $\{{\bar{{\rm{\Psi }}}}_{b},{\bar{B^{\prime} }}_{z0}$, $\bar{f}\}$, which are shown in runs S1–S19 of Table 4. For the parameter set of runs S1–S19, we find the following results:

• 1.  
  For the same ${{\rm{\Psi }}}_{b}$ and B₀, with the increase of n_G and n_P, R_G and R_P decrease, respectively.
• 2.  
  For the same n_G, n_P, and Ψ_b, with the increase of B₀, R_G and R_P decrease; whereas for the same n_G, n_P, and B₀, R_G and R_P increase with increasing Ψ_b.
• 3.  
  For a fixed parameter set $\{{\psi }_{b}$, ${\bar{B^{\prime} }}_{z0}$ $\bar{f}$, $\bar{\chi }\}$, we see that ${R}_{G}\gt {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 fall off more quickly than those of the generalized Gaussian profiles for higher values of n_G and n_P.

Table 4.  Different Combinations of the Input Parameters and the Dimensionless Parameters for the Self-similar Model where R_G and R_P Represent the Radii of the Flux Tubes for the Generalized Gaussian and Power-law Profiles, Respectively

Run #	${{\rm{\Psi }}}_{b}$	B0	pc	${B}_{z0}^{{\prime} }$	χ	f	nG	nP	${\psi }_{b}$	$\bar{f}$	${\bar{B}}_{z0}^{{\prime} }$	$\bar{\chi }$	RG	RP
	$({10}^{17}$ Mx)	(kG)	(105 dyne cm−2)	(G km−1)	(cm−2)	$({10}^{-30}$ cm−4)							(km)	(km)
S1	1	1	1	1	10−16	4.56	2	3	0.01	456	1	1	214	261
S2	1	1	1	1	10−14	4.56	2	3	0.01	456	1	100	214	261
S3	1	1	0.8	1.5	10−18	8.56	2.5	3.5	0.01	856	1.5	0.01	196	138
S4	1	1	0.5	2	10−16	14.56	3	4	1	1456	2	1	186	100
S5	1	2	1	1	10−14	4.56	2	3	0.005	456	0.5	100	151	184
S6	1	2	0.8	1.5	10−18	8.56	2.5	3.5	0.005	856	0.75	0.01	139	98
S7	1	2	0.5	2	10−16	14.56	3	4	0.005	1456	0.25	1	131	71
S8	5	1	1	1	10−14	0.182	2	3	0.005	18.2	1	100	479	584
S9	5	1	0.8	1.5	10−18	0.342	2.5	3.5	0.005	34.2	1.5	0.01	439	308
S10	5	1	0.5	2	10−16	0.582	3	4	0.005	58.2	2	1	416	225
S11	5	2	1	1	10−14	0.182	2	3	0.025	18.2	0.5	100	339	413
S12	5	2	0.8	1.5	10−18	0.342	2.5	3.5	0.025	34.2	0.75	0.01	310	218
S13	5	2	0.5	2	10−16	0.582	3	4	0.025	58.2	0.25	1	294	159
S14	10	1	1	1	10−14	0.0456	2	3	0.1	4.56	1	100	678	826
S15	10	1	0.8	1.5	10−18	0.0856	2.5	3.5	0.1	8.56	1.5	0.01	621	436
S16	10	1	0.5	2	10−16	0.145	3	4	0.1	14.56	2	1	589	318
S17	10	2	1	1	10−14	0.0456	2	3	0.05	4.56	0.5	100	479	584
S18	10	2	0.8	1.5	10−18	0.0856	2.5	3.5	0.05	8.56	0.75	0.01	439	308
S19	10	2	2	0.5	10−16	0.145	3	4	0.05	14.56	0.25	1	416	225

Download table as:  ASCIITypeset image

As an example, we show the solution of ${\psi }_{S}^{O}$ and the magnetic and thermodynamic structures for run S1 in Table 4. The values of the magnetic and thermodynamic quantities obtained from the self-similar model are reported in Table 5, for both the Gaussian and power-law shape function profiles. The radial variations of the generalized Gaussian and power-law flux function are shown in Figure 10 for different values of n_G and n_P, and the variation along the r − z plane is shown in Figure 11. The 3D configurations of the field lines for the generalized Gaussian and power-law self-similar models are shown in the Figures 12 and 13 for the parameter sets for runs S1 and S2 given in Table 4. The radial and vertical distributions of the magnetic field components are shown in Figures 14 and 15, respectively, for both the Gaussian and power-law models, whereas the density inside the flux tube does not vary along the r-direction but decreases along the z-direction which is shown in the Figure 16. The variations of p and T in the r − z plane obtained from the self-similar model are shown in Figures 19 and 20 for Gaussian and power-law shape function profiles. Conclusions drawn from Figures 6, 7, 9, 12, 13, 16, 18, 19, 20, 21, and Tables 3–5 are discussed in Section 7.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Table 5.  The Values of the Magnetic Field Strength and Thermodynamic Quantities Obtained from the Coulomb Function Open Field Flux Tube Model for the Parameter Set of Run C4 in Table 3, where r_b = 84 km; and the Self-similar Model with Gaussian Profile, where R_G = 214 km, and Power-law Profile, with n_P = 3, where R_P = 261 km, for the Parameter Set of Run S1 Corresponding to Table 4, Are Shown

Models	r	z (Mm)	Bz(G)	p (dyne cm−2)	ρ (g cm−3)	T (K)
	0	0	1000	$1.030\times {10}^{5}$	$2.44\times {10}^{-7}$	5656
Coulomb function open field	0	2	2.61	0.234	$5.56\times {10}^{-13}$	5656
	rb	0	0	$1.04\times {10}^{5}$	$2.44\times {10}^{-7}$	5690
	rb	2	0	0.2445	$5.56\times {10}^{-13}$	5890
	0	0	103	$1.0\times {10}^{5}$	$2.37\times {10}^{-7}$	5630
Generalized Gaussian	0	2	3.44	0.227	$5.44\times {10}^{-13}$	5630
	RG	0	6.73	$1.17\times {10}^{5}$	$2.37\times {10}^{-7}$	6620
	RG	2	2.2	1.54	$5.44\times {10}^{-13}$	38000
	0	0	103	$1.0\times {10}^{5}$	$2.37\times {10}^{-7}$	5630
Power law	0	2	75	0.227	$5.44\times {10}^{-13}$	5630
	RP	0	50	$1.17\times {10}^{5}$	$2.37\times {10}^{-7}$	6620
	RP	2	19	1.75	$5.44\times {10}^{-13}$	43000

Download table as:  ASCIITypeset image

We compare our models with the observations reported by high-resolution and high-cadence instruments. The small-scale magnetic structures in the solar photosphere are often found in the forms of MBPs, which are small-scale magnetic flux tubes with open field lines (Berger et al. 1995; Centeno et al. 2007; Lagg et al. 2010). Therefore, the MBPs are the best candidates to compare our open field flux tube models with the observations. MBPs can be identified by spectropolarimetric measurements or they can be seen by the G-band filtergrams (Utz et al. 2009, 2013; Yang et al. 2016). Next, we compare the observed magnetic field strength, size, and thermodynamic quantities of the MBPs with those 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 cutoff 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 cutoff radius is always positive. From the recent observations by Utz et al. (2009, 2013), 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 Figure 4 we see that B_z vanishes at r_b = 84 km, where R = 120 km is the entire radial simulation domain. We choose the parameter ranges, 1 kG ≤ B₀ ≤ 1.5 kG and 100 km ≤ R ≤ 180 km, for which the magnetic and thermodynamic quantities obtained from our model are in reasonable agreement with the solar atmosphere (Vernazza et al. 1981), and the selection of the $\{R,{B}_{0}\}$ parameter space is also consistent with the observations of MBP size and field strength distributions (Utz et al. 2009, 2013). The values of the magnetic and thermodynamic quantities obtained from the Coulomb function open field model are reported in Table 5. 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¹⁷–10¹⁸ Mx which is a typical flux value for MBPs(Zhang et al. 1998; Hagenaar et al. 1999; Guglielmino et al. 2011). According to previous studies by Shelyag et al. (2010) and SM18, the gas pressure at the axis of MBP is less than its boundary gas pressure, so we have chosen the parameter p_c < p₀. The field strengths of the magnetic footpoints observed in the photosphere for MBPs are ∼1 kG, with a distribution peak at 1.3 kG (Utz et al. 2013). Thus, we use the value of B₀ in the typical range of 1–2 kG (Zhang et al. 1998) in our model. The vertical gradient of the magnetic field strength at the photosphere is ∼1 G km⁻¹ (Wittmann 1974; Pahlke & Wiehr 1990; Balthasar & Schmidt 1993). Hence, we use the value of ${B}_{z0}^{\prime}$ in the range of 1–2 G km⁻¹ in our model. We have reported combinations of the free parameters and the corresponding input parameters in Table 4. Within the parameter sets of runs S1–S19 in Table 4, we note that the minimum and maximum radii of the flux tubes are 151 and 678 km, respectively, for the Gaussian model, and 71 and 826 km, respectively, for the power-law model, which are in reasonable agreement with the observations of MBP size distributions (Utz et al. 2009). The values of the magnetic and thermodynamic quantities obtained from a self-similar model are reported in Table 5, which is also in reasonable agreement with the solar atmosphere reported by Vernazza et al. (1981).

Zoom In Zoom Out Reset image size

The main findings of our simulations are as follows:

• 1.  
  The Coulomb function model is easier to implement for the numerical studies, as it consists of two free parameters (R, B₀); on the other hand, the self-similar model consists of five free parameters (Ψ_b, B₀, p_c, ${B}_{z0}^{{\prime} }$, and χ). From Table 5, we see that the rise of the gas pressure along the radial direction from the 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 the 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 as where B_z vanishes; on the other hand, for the self-similar model, B_z(R) $\ne$ 0, whereas B_z reduces along the radial direction from the axis to the boundary for the Gaussian model faster than it does for the power-law model.
• 2.  
  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 note that, in the domains of 100 km ≤ R ≤ 180 km, and 1 kG ≤ B₀ ≤ 1.5 kG, the magnetic and thermodynamic quantities are in reasonable agreement with the solar atmosphere (Vernazza et al. 1981), which also validates the MBP size and magnetic field strength distribution (Utz et al. 2009, 2013). For the self-similar model, the radial sizes of the flux tubes depend on the choice of the dimensionless input parameters $\{{\psi }_{b},\bar{f},\bar{{B}_{z0}^{{\prime} }},\bar{\chi }\}$. In the domain of the selected parameter space (see Table 4), the maximum and minimum radii of the flux tubes are 678 and 151 km as 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 observations of MBP size distribution by Utz et al. (2009).
• 3.  
  For the Coulomb function model, we note that the value of α decreases with R (see Table 3), which lowers the poloidal current I_p and the twist of the field lines. The 3D geometries of the field lines for different twists are shown in Figures 6 and 7 for open and closed field Coulomb function models, respectively. In the self-similar model, the twist of the field lines increases with $\bar{\chi }$, as shown in Figures 12 and 13 for the Gaussian and power-law profiles, respectively, which follows from Equation (66).
• 4.  
  The gas pressure for both the Coulomb function and self-similar models increases along the radial direction from the axis to the boundary (see Figures 14 and 17), whereas it decreases along the vertical height from the photosphere to the transition region (see Figures 18 and 19), which is similar to the result obtained by Shelyag et al. (2010) for MBPs, where the gas pressure inside the MBPs increases radially though the change is not significant, and decreases vertically. Gent et al. (2013, 2014) studied 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 Figures 9 and 16). Our model predicts that the atmosphere inside the flux tube is nearly plane-parallel, which is comparable to the model obtained for MBPs by Shelyag et al. (2010).
• 5.  
  In the solar atmosphere, the temperature in the transition region rises perhaps because the shock dissipation of waves plays a dominant role, which is not included in our model. We have also not implemented the temperature profile by Vernazza et al. (1981) (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 in which our input external atmosphere model is valid. Both the flux tube models we built are non-isothermal, and the temperature increases along the radial direction for both Coulomb function open field and self-similar models (see Figures 17, 18, and 20). 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 Figures 18 and 20).
• 6.  
  Hewitt et al. (2014), Uitenbroek & Criscuoli (2013), and Riethmüller & Solanki (2017) reported the simulation results of MBPs using MuRAM and Copenhagen-Stagger code in which the obtained values of magnetic field strength, pressure, density, and temperature inside the flux tube are in reasonable agreement with our predictions.
• 7.  
  The 2D simulations of the propagation of linear and non-linear magnetoacoustic waves through an open magnetic flux tube, embedded in the solar atmosphere from the photosphere to the corona, were carried out by Fedun et al. (2011). We can incorporate our solutions as the background condition for such numerical studies of waves and their kinematic properties by taking realistic inputs of the field strength and pressure distribution observed in the solar atmosphere.
• 8.  
  The Coulomb function model gives both open and closed field flux tube solutions, which can be co-added to build a canopy structure. A cartoon diagram of the magnetic canopy is shown in Figure 21, where the closed field lines (red), ${{\rm{\Psi }}}_{C}^{C}$, are present between the open field flux tubes and obtained from the Coulomb function, where the open field lines (blue), ${{\rm{\Psi }}}_{C}^{O}$ and ${{\rm{\Psi }}}_{S}^{O}$, of the neighboring flux tubes merge with each other to form a canopy structure. This is similar to the structures assumed in the numerical simulations by Gent et al. (2014), 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 that is embedded in a continuous magnetic medium and spans up to 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.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

More advanced observations of the magnetic and thermodynamic structures of the MBPs will provide a better selection of parameter inputs to discriminate between our models.

In this work, we have constructed two different models of flux tubes with twisted magnetic fields, the Coulomb function helical flux tube model and the self-similar model, by solving GSE semianalytically. We tabulate the expressions of the magnetic and thermodynamic functions for the Coulomb and self-similar models in Table 6, and highlight the novel features of this work as follows:

• 1.  
  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 that is independent of z.
• 2.  
  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₀), ${\psi }_{b}(R,{B}_{0})$ and $\bar{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 the photosphere and transition region obtained from Avrett & Loeser (2008) model. The values of the parameters for the Coulomb function model are listed in Table 3.
• 3.  
  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 5. 3D visualizations of both the open and closed field lines are shown in Figures 6 and 7 for the parameter sets of runs C4 and C10 corresponding to Table 3.
• 4.  
  In the self-similar model, we have employed an extra term, ${p}_{c}\exp (-2{kz})$, with p₁ in Equation (55), to maintain the hydrostatic pressure balance under the influence of stratified solar gravity, and taken two options for the shape functions: ${D}_{X}(\xi )$ from Equation (57), which is the extension of previous models by ST58, Yun (1971), and Osherovitch (1979, 1982). 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₀, p_c, ${B}_{z0}^{{\prime} }$, and χ (see Table 4), that is 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. These are given in Table 5, for the parameter set of run S1 corresponding to Table 4.
• 5.  
  Preliminary calculations using the constraint of relative helicity based on the formulations given in Prasad et al. (2014) and Prasad & Mangalam (2016), and applying the constrained energy minimization principle (Taylor 1974; Finn & Antonsen 1983; Mangalam & Krishan 2000), 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, that can be verified by 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 coronal heating by waves.

Table 6.  Formularies of Different Functions Obtained for the Coulomb Function Helical Flux Tube and Self-similar Model

Functions	Formulae for the Coulomb Function Model
$s(\varpi )$	$c({F}_{0}\left(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2}\right)+{F}_{0}^{* }\left(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2}\right)$
$Z(\bar{z})$	$\exp \left(-\tfrac{2\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right)$
${\psi }_{p}(\varpi )$	$\tfrac{i\sqrt{a}b{\varpi }^{2}}{4{\psi }_{b}}\left[{e}^{i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\tfrac{{F}_{2}^{1}\left(n+2,1;n+2-i\alpha ;\tfrac{1}{2}\right){\left(-i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1-i\alpha )n!}-{e}^{-i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\tfrac{{F}_{2}^{1}\left(n+2,1;n+2+i\alpha ;\tfrac{1}{2}\right){\left(i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1+i\alpha )n!}\right]$
${\psi }_{C}(\varpi ,\bar{z})$	$s(\varpi )Z(\bar{z})+{\psi }_{p}(\varpi )$
${B}_{r}(\varpi ,\bar{z})$	$\tfrac{{B}_{0}\sqrt{2}{\psi }_{b}\kappa }{{a}^{1/4}}s(\varpi )Z(\bar{z})$
${B}_{z}(\varpi ,\bar{z})$	$\tfrac{{B}_{0}{\psi }_{b}}{2\sqrt{a}\varpi }[s^{\prime} (\varpi )Z(\bar{z})+{\psi }_{p}^{{\prime} }(\varpi )]$
${B}_{\phi }(\varpi ,\bar{z})$	$\tfrac{{B}_{0}\sqrt{2}{\alpha }^{1/2}{\psi }_{b}}{{a}^{1/4}}[s(\varpi )Z(\bar{z})+{\psi }_{p}(\varpi )]$
$p(\varpi ,\bar{z})$	${B}_{0}^{2}\left[\left(\tfrac{{\psi }_{b}^{2}{s}^{2}(\varpi )}{8\pi }+{\bar{p}}_{20}\right){Z}^{2}(\bar{z})+\left(\tfrac{{\psi }_{b}^{2}s(\varpi ){\psi }_{p}(\varpi )}{4\pi }+\tfrac{b{\psi }_{b}s(\varpi )}{2\sqrt{2a}}\right)Z(\bar{z})+\left(\tfrac{{\psi }_{b}^{2}{\psi }_{p}^{2}}{8\pi }+\tfrac{b{\psi }_{b}{\psi }_{p}}{2\sqrt{2a}}\right)\right]$
$\rho (\bar{z})$	$\tfrac{4\sqrt{2}\kappa {a}^{1/4}{p}_{20}}{{gR}}{Z}^{2}(\bar{z})$
$T(\varpi ,\bar{z})$	$\tfrac{\bar{\mu }{B}_{0}^{2}{gR}}{4\sqrt{2}{R}_{g}\kappa {a}^{1/4}{p}_{20}}\left[\left(\tfrac{{\psi }_{b}^{2}{s}^{2}(\varpi )}{8\pi }+{\bar{p}}_{20}\right)+\left(\tfrac{{\psi }_{b}^{2}s(\varpi ){\psi }_{p}(\varpi )}{4\pi }+\tfrac{b{\psi }_{b}s(\varpi )}{2\sqrt{2a}}\right)\tfrac{1}{Z(\bar{z})}+\left(\tfrac{{\psi }_{b}^{2}{\psi }_{p}^{2}}{8\pi }+\tfrac{b{\psi }_{b}{\psi }_{p}}{2\sqrt{2a}}\right)\tfrac{1}{{Z}^{2}(\bar{z})}\right]$
Functions	Formulae for the Self-similar Model
$\xi (\varpi ,\bar{z})$	$\sqrt{\tfrac{\tau }{{\psi }_{b}{D}_{0}}}\varpi y(\bar{z})$
${\psi }_{G}(\varpi ,\bar{z})$	$1-\tfrac{{\rm{\Gamma }}(2/{n}_{G},{\xi }^{{n}_{G}})}{{\rm{\Gamma }}(2/{n}_{G})};\quad ({n}_{G}\gt 0)$
${\psi }_{P}(\varpi ,\bar{z})$	$1-{(1+\xi )}^{1-{n}_{P}}\left(1+\xi ({n}_{P}-1)\right);\quad ({n}_{P}\gt 2)$
${B}_{r}(\varpi ,\bar{z})$	$-\tfrac{{B}_{0}\varpi }{{D}_{0}}y(\bar{z})y^{\prime} (\bar{z}){D}_{X}(\xi )$
${B}_{z}(r,z)$	$\tfrac{{B}_{0}}{{D}_{0}}{y}^{2}(z){D}_{X}(\xi )$
${B}_{\phi }(\varpi ,\bar{z})$	$\tfrac{\sqrt{\bar{\chi }}{B}_{0}}{{D}_{0}}\varpi {y}^{2}(\bar{z}){D}_{X}(\xi )$
$p(\varpi ,\bar{z})$	${B}_{0}^{2}\left[\tfrac{\bar{f}}{2}{\psi }^{2}+{\bar{p}}_{c}{e}^{-2\bar{k}\bar{z}}\right]$
$\rho (\bar{z})$	$\tfrac{2\bar{k}\bar{z}{\bar{p}}_{c}{B}_{0}^{2}}{{{gz}}_{0}}{e}^{-2\bar{k}\bar{z}}$
$T(\varpi ,\bar{z})$	$\tfrac{\bar{\mu }{{gz}}_{0}}{2{R}_{g}\bar{k}\bar{z}{\bar{p}}_{c}}\left(\tfrac{\bar{f}}{2}{\psi }^{2}{e}^{2\bar{k}\bar{z}}+{\bar{p}}_{c}\right)$

Note. Here, $s(\varpi )$, $Z(\bar{z})$, and ${\psi }_{p}(\varpi )$ are given by Equations (18), (19), (23), and ${\bar{p}}_{20}={p}_{20}/{B}_{0}^{2}$. $\xi$ is a self-similar parameter in which $y(\bar{z})$ is obtained by solving Equation (60) and ${\bar{p}}_{0}={p}_{0}/{B}_{0}^{2}$, ${\bar{p}}_{c}={p}_{c}/{B}_{0}^{2}$, $\bar{k}={{kz}}_{0}$. The value of the constants are $\bar{\mu }=1.12$, $g=2.74\times {10}^{4}\,\mathrm{cm}$ s⁻², $k=3.4\times {10}^{-8}$ cm⁻¹, and ${z}_{0}={10}^{8}\,\mathrm{cm}$.

Download table as:  ASCIITypeset image

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 their hospitality during our visits, and Saikat Das for help with Figures 2 and 21.

The homogeneous solution $s(\varpi )$ which is given by the Equation (18), can be represented in terms of the Whittaker–M function (SM18), where the Whittaker–M function can be expressed in terms of the hypergeometric function by the standard relation (Dixit & Moll 2015)

$$\begin{eqnarray}&&{M}_{t,m}(\nu )={e}^{-\nu /2}{\nu }^{m+\tfrac{1}{2}}{F}_{1}^{1}\left(\displaystyle \frac{1}{2}+m-t,1+2\,m,\nu \right),\end{eqnarray} \tag{ 74 }$$

where ${F}_{1}^{1}$ represents the hypergeometric function with the arguments $t,m$, and ν. Taking the real part of ${\psi }_{h}(\varpi ,\bar{z})$ from Equation (20) and ${\psi }_{p}(\varpi )$ from (23), and using Equations ((24), (74)) we obtain

$$\begin{eqnarray}\begin{array}{rcl}{B}_{z}(\varpi ,\bar{z}) & = & {B}_{0}{\psi }_{b}c\exp \left(-\displaystyle \frac{2\sqrt{2}\kappa {a}^{1/4}\bar{z}}{\tau }\right)\\ & & \times \,\left[8(1+i\sqrt{a}{\varpi }^{2}){F}_{1}^{1}(1+i\alpha ,2,2i\sqrt{a}{\varpi }^{2})\right.\\ & & -\,\left.8\sqrt{a}{\varpi }^{2}{F}_{1}^{1}(1+i\alpha ,3,2i\sqrt{a}{\varpi }^{2})\right]\\ & & +\,\displaystyle \frac{{B}_{0}{\psi }_{b}}{2\sqrt{a}}\left(\displaystyle \frac{{\psi }_{p}^{{\prime} }(\varpi )+{\psi }_{p}^{{\prime} * }(\varpi )}{2\varpi }\right)\end{array}\end{eqnarray} \tag{ 75 }$$

where, ${B}_{0}\equiv {B}_{z}(0,0)$, and from Equation (75) we obtain,

$$\begin{eqnarray}&&1=8{\psi }_{b}c+\displaystyle \frac{{\psi }_{b}}{2\sqrt{a}}{\left[\displaystyle \frac{{\psi }_{p}^{{\prime} }(\varpi )+{\psi }_{p}^{{\prime} * }(\varpi )}{2\varpi }\right]}_{\varpi =0},\end{eqnarray} \tag{ 76 }$$

where the identity ${F}_{1}^{1}(1+i\alpha ,2,0)=1$ is applied. By expanding the last term on the RHS of Equation (76), we obtain

$$\begin{eqnarray}\begin{array}{rcl}1 & = & 8{\psi }_{b}c+\displaystyle \frac{{ib}}{2}\left[\displaystyle \frac{{F}_{2}^{1}(1,-i\alpha ,2-i\alpha ,-1)}{1-i\alpha }\right.\\ & & -\,\left.\displaystyle \frac{{F}_{2}^{1}(1,i\alpha ,2+i\alpha ,-1)}{1+i\alpha }\right]\end{array}\end{eqnarray} \tag{ 77 }$$

and the expression for c is given by

$$\begin{eqnarray}\begin{array}{rcl}c & = & \displaystyle \frac{1}{8{\psi }_{b}}\left[1-\displaystyle \frac{{ib}}{2}\left(\displaystyle \frac{{F}_{2}^{1}(1,-i\alpha ,2-i\alpha ,-1)}{1-i\alpha }\right.\right.\\ & & -\,\left.\left.\displaystyle \frac{{F}_{2}^{1}(1,i\alpha ,2+i\alpha ,-1)}{1+i\alpha }\right)\right].\end{array}\end{eqnarray} \tag{ 78 }$$

Hence, from Equation (75), we obtain the explicit form for ${B}_{z}(\varpi ,\bar{z})$ that is given in Equation (25).

The explicit forms of Equations 44(a)–44(e) are

$$\begin{eqnarray}&&{F}_{0}(-\alpha -{\kappa }^{2},\sqrt{a})+{F}_{0}^{* }(-\alpha -{\kappa }^{2},\sqrt{a})=0,\end{eqnarray} \tag{ 79 }$$

$$\begin{eqnarray}&&\begin{array}{l}\left[\displaystyle \frac{1}{{\varpi }^{2}}\displaystyle \frac{d}{d\varpi }\left({\varpi }^{2}{e}^{i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\right.\right.\\ \quad \displaystyle \frac{{F}_{2}^{1}(n+2,1,n+2-i\alpha ,1/2){\left(-i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1-i\alpha )n!}\\ \quad -\,{\varpi }^{2}{e}^{-i\sqrt{a}{\varpi }^{2}}\displaystyle \sum _{n=0}^{\infty }\\ \quad {\left.\left.\displaystyle \frac{{F}_{2}^{1}(n+\mathrm{2,1},n+2+i\alpha ,1/2){\left(i\sqrt{a}{\varpi }^{2}\right)}^{n}}{(n+1+i\alpha )n!}\right)\right]}_{\varpi =1}=0,\end{array}\end{eqnarray} \tag{ 80 }$$

$$\begin{eqnarray}&&\begin{array}{l}\left[1-\displaystyle \frac{{ib}}{2}\left(\displaystyle \frac{{F}_{2}^{1}(1,-i\alpha ,2-i\alpha ,-1)}{1-i\alpha }\right.\right.\\ \quad -\,{\left.\left.\displaystyle \frac{{F}_{2}^{1}(1,i\alpha ,2+i\alpha ,-1)}{1=i\alpha }\right)\right]}^{2}\\ \quad \cdot \,{\left(\displaystyle \frac{d}{d\varpi }{\left[{F}_{0}(\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})+{F}_{0}^{* }(\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})\right]}_{\varpi =1}\right)}^{2}\\ =\,\bar{p}a,\end{array}\end{eqnarray} \tag{ 81 }$$

$$\begin{eqnarray}&&b=-\displaystyle \frac{{\psi }_{b}}{2}\left(1+\displaystyle \frac{6\alpha }{\sqrt{a}}\right),\end{eqnarray} \tag{ 82 }$$

$$\begin{eqnarray}&&\begin{array}{l}{\displaystyle \int }_{0}^{1}\left[{e}^{-2{\kappa }_{t}}\left(1-\displaystyle \frac{{ib}}{2}\left[\displaystyle \frac{{F}_{2}^{1}(1,-i\alpha ,2-i\alpha ,-1)}{1-i\alpha }\right.\right.\right.\\ \quad -\,{\left.\left.\displaystyle \frac{{F}_{2}^{1}(1,i\alpha ,2+i\alpha ,-1)}{1+i\alpha }\right]\right)}^{2}\\ \cdot \,{\left({F}_{0}(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})+{F}_{0}^{* }(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})\right)}^{2}\\ \quad +\,8{e}^{-{\kappa }_{t}}{\psi }_{b}\left({\psi }_{p}+{\psi }_{p}^{* }+\displaystyle \frac{2b}{{\psi }_{b}}\right)\\ \quad \cdot \,\left(1-\displaystyle \frac{{ib}}{2}\left[\displaystyle \frac{{F}_{2}^{1}(1,-i\alpha ,2-i\alpha ,-1)}{1-\alpha }\right.\right.\\ \quad -\,{\left.\left.\displaystyle \frac{{F}_{2}^{1}(1,i\alpha ,2+i\alpha ,-1)}{1+\alpha }\right]\right)}^{2}\\ \cdot \,\left({F}_{0}(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})+{F}_{0}^{* }(-\alpha -{\kappa }^{2},\sqrt{a}{\varpi }^{2})\right)\\ \quad +\,\left.\displaystyle \frac{{\psi }_{b}^{2}}{4}{\left({\psi }_{p}+{\psi }_{p}^{* }\right)}^{2}+b{\psi }_{b}({\psi }_{p}+{\psi }_{p}^{* })\right]d\varpi =\bar{p}\end{array}\end{eqnarray} \tag{ 83 }$$

respectively, where, ${\psi }_{p}(\varpi )$ is given by Equation (23), ${\kappa }_{t}={{kz}}_{t}$, and ${\psi }_{p}(\varpi =1)={\psi }_{b}$.