Fragmentation of Filamentary Cloud Permeated by Perpendicular Magnetic Field

# Fragmentation of Filamentary Cloud Permeated by Perpendicular Magnetic Field

Tomoyuki Hanawa Center for Frontier Science, Chiba University
1-33 Yayoi-cho, Inage-ku, Chiba
Chiba 263-8522, Japan
Takahiro Kudoh Faculty of Education, Nagasaki University
1-14 Bonkyo-machi, Nagasaki
Nagasaki 852-8521, Japan
Kohji Tomisaka Division of Theoretical Astronomy, National Astronomical Observatory of Japan
Mitaka, Tokyo 181-8588, Japan, Department of Astronomical Science, School of Physical Sciences, SOKENDAI (The Graduate University for Advanced Studies), Mitaka, Tokyo 181-8588, Japan
###### Abstract

We examine the linear stability of an isothermal filamentary cloud permeated by a perpendicular magnetic field. Our model cloud is assumed to be supported by gas pressure against the self-gravity in the unperturbed state. For simplicity, the density distribution is assumed to be symmetric around the axis. Also for simplicity, the initial magnetic field is assumed to be uniform and turbulence is not taken into account. The perturbation equation is formulated to be an eigenvalue problem. The growth rate is obtained as a function of the wavenumber for fragmentation along the axis and the magnetic field strength. The growth rate depends critically on the outer boundary. If the displacement vanishes in the region very far from the cloud axis (fixed boundary), cloud fragmentation is suppressed by a moderate magnetic field, which means the plasma beta is below 1.67 on the cloud axis. If the displacement is constant along the magnetic field in the region very far from the cloud, the cloud is unstable even when the magnetic field is infinitely strong. The cloud is deformed by circulation in the plane perpendicular to the magnetic field. The unstable mode is not likely to induce dynamical collapse, since it is excited even when the whole cloud is magnetically subcritical. For both the boundary conditions the magnetic field increases the wavelength of the most unstable mode. We find that the magnetic force suppresses compression perpendicular to the magnetic field especially in the region of low density.

MHD — ISM clouds — ISM: magnetic fields
\AuthorCallLimit

=3

## 1 Introduction

It is widely accepted that stars form from the fragmentation of filamentary clouds (see, e.g., André et al., 2014, and the references therein). This idea is supported by observations showing that the filamentary structure is universal in interstellar clouds and that prestellar cores and newly formed stars are associated with the dense parts of filamentary clouds. It is also well known that filamentary clouds are unstable against fragmentation if they are gravitationally bound (see, e.g., Stodólkiewicz, 1963; Larson, 2003, and the references therein).

Fragmentation of a magnetized filamentary cloud has been studied extensively (see, e.g., Stodólkiewicz, 1963; Nakamura, Hanawa & Nakano, 1993; Hanawa et al., 1993; Fiege & Pudritz, 2000), particularly in 1990s. We note that the magnetic field was generally assumed to be either parallel or helical to the cloud axis in these studies. However, fragmentation depends on the magnetic field direction, since the magnetic force is perpendicular to the local magnetic field. We know that massive main filamentary clouds are often associated with a magnetic field perpendicular to the axis, although less massive sub-filaments are parallel to the magnetic field (see, e.g., Sugitani et al., 2011; Palmeirim et al., 2013; André et al., 2014; Kusune et al., 2016, and the references therein). Prestellar cores are often associated with the main filaments. Thus, we should examine the effects of the magnetic field perpendicular to the cloud axis on fragmentation.

It is not easy to construct a self-consistent model of a filamentary molecular cloud associated with a magnetic field perpendicular to the cloud axis. Consequently, it is not easy to analyze the stability against fragmentation. As shown in Tomisaka (2014), the cloud is likely to be flattened, like fettuccine pasta. His equilibrium model shows that the cloud is supported in part by the magnetic force against gravity. The numerically obtained density distribution is difficult to apply in a stability analysis. Hanawa & Tomisaka (2015) examined fragmentation of a filamentary cloud permeated by a perpendicular magnetic field under the approximation that the filamentary cloud is highly flattened in the direction perpendicular to the magnetic field. They neglected the density structure along the magnetic field.

With the abovementioned difficulties in mind, we examine the stability of a simple self-consistent model for a magnetized filamentary cloud. The model cloud is assumed to be isothermal and supported against the self-gravity solely by the gas pressure in equilibrium. The initial magnetic field is assumed to be uniform and perpendicular to the cloud axis so that the magnetic force vanishes in the initial state. Although turbulence is dominant over the thermal pressure, especially in the low density region, we ignore interstellar turbulence for simplicity. Nevertheless, some effects of turbulence may be taken into account effectively by replacing the sound speed with the typical turbulent velocity. The model cloud is also assumed to be infinitely long and the density to be a function of the radius from the axis. The magnetic field works against fragmentation through tension. We evaluate the effects of the magnetic field on fragmentation by conducting a linear stability analysis. The growth rate of the instability is obtained as a function of the wavelength and the magnetic field strength, i.e., the initial plasma beta at the cloud center.

This paper is organized as follows. We describe our equilibrium model and methods of stability analysis in §2. We employ two types of boundary conditions, fixed and free, for the magnetic field in the region very far from the cloud center. The result of the stability analysis is presented in §3. The growth rate is shown to depend on the choice of boundary condition. We discuss the implication of our stability analysis in §4. In Appendix A, we prove that an unstable mode has a real growth rate. In Appendix B, we describe the elements of the matrixes used for our stability analysis. In Appendix C, we calculate the growth rate for an unmagnetized filamentary cloud in order to evaluate the accuracy of our numerical analysis.

## 2 Methods

### 2.1 Basic Equations

For our stability analysis, we employ the ideal magnetohydrodymamic (MHD) equations,

 ∂ρ∂t+\boldmath∇⋅(ρ\boldmathv)=0, (1) ρd\boldmathvdt=−\boldmath∇P+\boldmathj×\boldmathB+ρ\boldmath∇Φ, (2) \boldmathj=\boldmath∇×% \boldmathB4π, (3)

where , , denote the density, pressure, gravitational potential, velocity, magnetic field, and electric current density, respectively. The gas is assumed to be isothermal and so the equation of state is expressed as

 P = ρc2s, (4)

where denotes the sound speed. The self-gravity of the gas is taken into account through the Poisson equation

 ΔΦ=4πGρ, (5)

where denotes the gravitational constant. We ignore ambipolar diffusion and turbulence for simplicity.

### 2.2 Equilibrium Model

We consider an isothermal filamentary cloud in hydrostatic equilibrium. The cloud axis is assumed to be located on in Cartesian coordinates. Then the density distribution is expressed as

 ρ0(x,y,z) = ρc(1+x2+y28H2)−2, (6) H2 = c2s4πGρc, (7)

where and denote the density distribution in equilibrium and the value on the cloud axis, respectively (Stodólkiewicz, 1963; Ostriker, 1964).

We assume that the cloud in equilibrium is permeated by a uniform magnetic field in the -direction,

 \boldmath{B}0=B0\boldmathex, (8)

where denotes the unit vector in the -direction. This uniform magnetic field does affect the cloud stability but not the equilibrium.

To specify the initial magnetic field strength in the analysis, we use the plasma beta at the cloud center,

 β = 8πρcc2sB20. (9)

The plasma beta is related to the mass to flux ratio,

 f = ∫ρ0(x′,y,z)dx′B0 (10) = √2πρcHB0(1+y28H2)−3/2. (11)

The mass to flux ratio can be rewritten as

 f = 14√βG(1+y28H2)−3/2, (12)

by using Equation (7). The critical mass to flux ratio is . Thus, the cloud is subcritical when

 β < 0.405(1+y28H2)3. (13)

### 2.3 Perturbation Equation

We consider a small perturbation around the equilibrium in order to search for an unstable mode. The perturbation is described by the displacement defined by

 ξ = ξx(x,y)coskz\boldmathex+ξy(x,y)coskz\boldmathey (14) +ξz(x,y)sinkz\boldmathez,

where the perturbation is assumed to be sinusoidal in the -direction with the wavenumber . The density perturbation is derived from the equation of continuity,

 δρ=−\boldmath∇⋅(ρ0\boldmathξ). (15)

Substituting Equation (14) into Equation (15), we obtain

 δρ(x,y,z) = δϱ(x,y)coskz, (16) δϱ = −∂∂x(ρ0ξx)−∂∂y(ρ0ξy)−kρ0ξz. (17)

Similarly, we obtain the perturbation in the magnetic field from the induction equation,

 δ\boldmathB=\boldmath∇×(\boldmathξ×\boldmathB0). (18)

The induction equation is further expressed as

 δ\boldmathB(x,y,z) = bx(x,y)coskz\boldmathex+by(x,y)coskz% \boldmathey (19) +bz(x,y)sinkz\boldmathez, bx = −B0[∂∂yξy(x,y)+kξz], (20) by = B0∂ξy∂x, (21) bz = B0∂ξz∂x. (22)

We evaluate the change in the current density to be

 δ\boldmathJ=14π\boldmath∇×δ\boldmathB, (23)

using Equation (3). Each component of the current density is expressed as

 δ\boldmathJ(x,y,z) = jx(x,y)sinkz\boldmathex+jy(x,y)sinkz% \boldmathex (24) +jz(x,y)coskz\boldmathez, jx = 14π(∂bz∂y+kby), (25) jy = −14π(kbx+∂bz∂x), (26) jz = 14π(∂by∂x−∂bx∂y). (27)

Then the changes in the density and current density are expressed as an explicit function of .

The change in the gravitational potential is given as the solution of the Poisson equation

 \boldmath∇2δψ=4πGδρ. (28)

Thus, it can be regarded as an implicit function of .

We derive the equation of motion for the perturbation by taking account of the force balance,

 \boldmath∇(c2sρ0)+ρ0\boldmath% ∇ψ0=0, (29)

with no electric current density, , in equilibrium. Then the equation of motion is expressed as

 σ2ρ0\boldmathξ=−c2sρ0\boldmath∇(δρρ0)−ρ0∇δψ+δ\boldmathJ×\boldmathB0, (30)

where the last term represents the magnetic force. The linear growth rate, , is obtained as the eigenvalue of the differential equation (30), since the right-hand side is proportional to . As shown in Appendix A, the growth rate should be either real or pure imaginary.

The equilibrium model is symmetric with respect to the - and -axes. Thus, all eigenmodes should be either symmetric or anti-symmetric with respect to these axes. We restrict ourselves to the eigenmodes symmetric to both - and -axes, since the unstable mode has the same symmetry in the case of no magnetic field (Nakamura, Hanawa & Nakano, 1993). The choice of this symmetry is justified since we are interested only in the unstable mode. Using this symmetry, we can reduce the region of computation to the first quadrant, and . The variables describing the perturbation and their symmetries are summarized in Table 1.

We consider two types of the boundary conditions. The first one assumes that the displacement should vanish in the region very far from the filament center. We call this the fixed boundary since the magnetic field lines are fixed on the boundary. The second one allows the magnetic field lines to move while remaining straight and normal to the boundary. This restriction is expressed as

 (\boldmathB0⋅\boldmath∇)\boldmathξ=0. (31)

Thus, we assume on the boundary in the -direction and in the -direction. We refer to this as the free boundary condition. In both types of boundary conditions, we use the symmetries given in Table 1 to set the boundary conditions for and .

### 2.4 Numerical Methods

We solve the eigenvalue problem numerically by a finite difference approach. The differential equations are evaluated on the rectangular grid in the plane. We evaluate , , , , and at the points

 (xi,yj) = (iΔx,jΔy), (32)

where and specify the grid points, while and denote the grid spacing in the - and -directions, respectively (see Table 1). These variables are symmetric with respect to both the - and -axes. Using this symmetry, we consider the range and , where and specify the number of grid points in each direction. When or , the displacement is assumed to vanish for the fixed boundary and to have the same values at neighboring points in the computation domain for the free boundary condition. We use the indexes, and , to specify the position where the variables are evaluated, such as .

The variables and are evaluated at the points

 (xi−1/2,yj) = [(i−12)Δx,jΔy]. (33)

These variables are anti-symmetric with respect to and symmetric with respect to . Similarly, the variables and are evaluated at the points

 (xi,yj−1/2) = [iΔx,(j−12)Δy], (34)

since they are symmetric with respect to and anti-symmetric with respect to . Given that is anti-symmetric with respect to both and , it is evaluated at the points

 (xi−1/2,yj−1/2) = [(i−12)Δx,(j−12)Δy]. (35)

All these variables are evaluated in the region and . In other words, we use staggered grids to achieve second-order accuracy in space.

Using the variables defined on the grids, we rewrite the perturbation equations. Equation (17) is rewritten as

 δϱi,j = −ρ0,i+1/2,jξx,i+1/2,j−ρ0,i−1/2,jξx,i−1/2,jΔx (36) − ρ0,i,j+1/2ξy,i,j+1/2−ρ0,i,j−1/2ξy,i,j−1/2Δy − kρ0,i,jξz,i,j.

Equation (28), the Poisson equation, is expressed as

 δψi+1,j+δψi−1,jΔx2+δψi,j+1+δψi,j−1Δy2 −(2Δx2+2Δy2+k2)δψj,k=4πGδϱi,j. (37)

The solution of Equation (37) is expressed as

 δψi,j = ∑i′∑j′Gi,j,i′,j′δϱi′,j′, (38)

where denotes the Green’s function and the value is obtained by solving Equation (37) numerically. The change in the magnetic field is evaluated as

 bx,i,j = −B0(ξy,i,j+1/2−ξy,i,j−1/2Δy+kξz,i,j), (39)
 by,i−1/2,j−1/2 = B0(ξy,i,j−1/2−ξy,i−1,j−1/2Δx), (40) bz,i−1/2,j = B0(ξz,i,j−ξz,i−1,jΔx), (41)

from Equations (20) through (22). The current density is evaluated as

 jy,i,j = −14π(kbx,i,j+bz,i+1/2,j−bz,i−1/2,jΔx), (42)
 jz,i,j−1/2 = 14π(by,i+1/2,j−1/2−by,i−1/2,j−1/2Δx (43) −bx,i,j−bx,i,j−1Δy).

The -component of the current density, , is not evaluated, since it does not appear in the equation of motion. The fixed boundary conditions are expressed as

 ξx,nx+1/2,j = 0, (44) ξy,nx+1,j−1/2 = 0, (45) ξz,nx+1,j = 0, (46) ξx,i−1/2,ny+1 = 0, (47) ξy,i,ny+1/2 = 0, (48) ξz,i,ny+1 = 0. (49)

When the free boundary is applied, the conditions are replaced with

 ξx,nx+1/2,j = ξx,nx−1/2,j, (50) ξy,nx+1,j−1/2 = ξy,nx,j−1/2, (51) ξz,nx+1,j = ξz,nx,j, (52) ξx,i−1/2,ny+1 = ξx,i−1/2,ny, (53) ξy,i,ny+1/2 = ξy,i,ny−1/2, (54) ξz,i,ny+1 = ξz,i,ny. (55)

The equation of motion (30) is expressed as

 σ2ρ0,i−1/2,jξx,i−1/2,j = −c2sρ0,i−1/2,jΔx(δϱi,jρ0,i,j−δϱi−1,jρ0,i−1,j) (56) −ρ0,i−1/2,jΔx(δψi,j−δψi−1,j). σ2ρ0,i,j−1/2ξy,i,j−1/2 = −c2sρ0,i,j−1/2Δy(δϱi,jρ0,i,j−δϱi,j−1ρ0,i,j−1) (57) −ρ0,i,j−1/2Δy(δψi,j−δψi,j−1) +B0jz,i,j−1/2. σ2ρ0,i,jξz,i,j = −kc2sδϱi,j−kρ0,i,jδψi,j (58) −B0jy,i,j.

Equations (56) through (58) are summarized in the form

 σ2\boldmathB\boldmathζ=(% \boldmathA+B20\boldmathC)\boldmathζ (59)

by using Equations (36), and (38) through (43). Here, denotes an array of components, , , and for all the combinations of and . The matrix elements of , , and are evaluated numerically as a function of . See Appendix B for further details. Then the growth rate is given as the solution of

 det(σ2\boldmathB−\boldmathA−B20\boldmathC)=0. (60)

We use the subroutine DGGEVX of LAPACK (see, Anderson et al., 1999, for the software) to solve Equation (60). The subroutine returns all the eigenvalues .

The matrixes , , and have dimension . Thus, we obtain eigenmodes. However, we select only one unstable mode () for a given and . The remaining eigenmodes denote oscillation of the filamentary cloud. In the following, we restrict ourselves to the unstable mode.

## 3 Results

### 3.1 Growth Rate

Figure 1 shows the eigenvalue for the fixed boundary condition as a function of . Each curve denotes in units of for a given . The grid spacing is set as and the computation region is specified as . Thus, the outer boundary is set as and . The growth rate is evaluated over the interval . Figure 1 shows the eigenvalue only when , i.e., . The imaginary part of is omitted since it is negligibly small. We find at most one growing mode for a given pair of and .

The growth rate has its maximum around in the absence of the magnetic field (). The growth rate is reduced greatly by the magnetic field having its root fixed at infinity. We find an unstable mode only in the range for . The unstable mode disappears when , i.e., when the plasma beta is smaller than 1.67 ().

The stabilization due to the magnetic field is shown clearly in Figure 2, where each curve is as a function of for a given . The square of the growth rate, , decreases in proportion to , i.e., in the range . The proportional constant is larger for a larger . The dispersion relation is similar to that for the MHD fast wave,

 ω2=(c2s+\boldmathB204πρ0)k2, (61)

where the wave is assumed to propagate normal to the magnetic field .

The growth rate shows a different dependence on the magnetic field strength when we apply the free boundary condition. Figure 3 is the same as Figure 1 except for the free boundary. Each curve in Figure 3 denotes the growth rate as a function of for a given whose value is designated by the same color. The dashed line denotes the growth rate when , and the value is obtained by extrapolation. As the magnetic field strength increases, the growth rate decreases but remains positive for , and therefore the cloud is unstable even when the magnetic field is infinitely strong.

The result shown in Figure 3 seems to contradict the naive expectation that a subcritical cloud is stable (see, e.g., Nakano & Nakamura, 1978). We explain the apparent discrepancy in the next subsection.

Figure 4 is the same as Figure 2 except for the free boundary. The growth rate approaches a certain value in the limit of large . The growth rate is well approximated by

 σ2(k,B0) = [σ∞(k)]2+[σ′(k)]2β. (62)

The growth rate is evaluated by a spline fit to the growth rate in the range . The value of is shown by the dashed curve in Figure LABEL:growthM.

When is smaller, the growth rate is close to the asymptotic value for a weaker magnetic field. This result implies that the magnetic field is bent or compressed little by the unstable perturbation when the initial magnetic field is strong. Otherwise, the distorted magnetic field would induce a strong force acting on the gas. At the same time, the displacement should be appreciable even in the region very far from the cloud center. Otherwise, the growth rate would be independent of the boundary condition in the -direction. We confirm this expectation in the next subsection.

The wavenumber of the most unstable mode is smaller for a larger . The wavenumber is for and for a small . The magnetic field decreases the growth rate and increases the wavelength of a typical perturbation.

The filamentary cloud is stable for any perturbation having a wavenumber larger than when . When the free boundary is applied, the critical wavenumber is slightly reduced, i.e., for a large . The critical wavenumber is evaluated from the high resolution computation of and .

### 3.2 Eigenmode

First, based on our stability analysis, we review the unstable perturbation in the absence of the magnetic field. When , both the growth rate and the eigenmode depend little on the boundary condition. Figure 5 shows the unstable mode for the fixed boundary when and . The values are normalized so that at the origin. The arrows indicate the displacement , while the contours are of . The density perturbation is representing according to the color scale given in the right panel. The displacement is very small in the region when the fixed boundary is set as and . Thus, the growth rate depends little on the outer boundary condition. The gas concentrates towards the origin and the density decreases in the surrounding gas.

Low-density gas is affected by a magnetic field even when the magnetic field is very weak. Figures 6 and 7 are the same as Figure 5, except the initial plasma betas are and 16.67 at the cloud center, respectively. The displacement in the -direction is greatly suppressed, especially in the region of large , since it is perpendicular to the magnetic field. Also, the displacement along the cloud axis (: contour lines) is restricted in the region when 16.67 . It should be noted that is large in the more extended region in the -direction (along the initial magnetic field) than in the absence of the magnetic field. The gas in the region far from the cloud axis is anchored through the magnetic field by the gas around the cloud center. The magnetic fields in Figures 6 and 7 are loosely bent because the roots are fixed on the outer boundary. The magnetic tension works against fragmentation as expected from the variational principle shown in Appendix A. When the initial magnetic field is stronger, the density perturbation is smaller. We again note that the unstable mode is normalized so that the displacement in the -direction is .

Figure 8 is the same as Figure 5 except for the equilibrium model, . The growth rate is much smaller than that for . The displacement is negligibly small in the region , where the magnetic flux tube is subcritical. Note that the displacement is toward the origin along the -axis, while it is away from the origin in the -direction. The region is concentrated around the -axis. This density enhancement is mainly due to the - and -components of the displacement. The -component elongates the density enhancement along the -axis.

Not only the growth rate but also the eigenfunction depends on the boundary condition. Figure 9 is the same as Figure 6 except the boundary is free and . The -component of the displacement has a large value even in the region of . The -component of the displacement changes its sign around .

When the free boundary is applied, the magnetic field lines can rotate freely around the -axis, as shown in Figure 10. Figure 10 shows the perturbation in the plane for the eigenmode depicted in Figure 9. The arrows indicate the displacement, while the color indicates the relative change in the density, . The contours depict the change in the -component of the magnetic field in intervals of . The displacement is essentially incompressible circulation, which deforms the shape of the filamentary cloud. The cloud diameter increases in the -direction at the points where the density increases. In other words, each dense fragment expands in the -direction. It should be noted that the cloud changes its form in the opposite way when the magnetic field is absent or longitudinal. In that case, the cloud diameter decreases at the points where the density increases (see Figure 5). The density enhancement is mainly due to the displacement along the density gradient . The compression by the displacement perpendicular to the magnetic field is very weak. Compared with the relative change in the density, the change in the magnetic field is much smaller, i.e., .

Figures 11 is the same as Figure 9 except = 2. The -component of displacement, , changes sign around . In contrast to the case of the fixed boundary, the -component of displacement is appreciably large in the low-density region. This is because magnetic tension does not arise by displacement when the free boundary is applied.

Figure 12 is the same as Figure 9 except . Although the model cloud is subcritical [see Equation (13)], it is unstable. In other words, it is unstable though the magnetic energy dominates over the gravitational energy.

The free boundary allows rearrangement of the magnetic flux tubes in the plane while keeping them straight and therefore without induction of magnetic force. The rearrangement reduces the gravitational energy if it gathers more massive tubes, i.e., those tubes initially located near the cloud axis. The gravitational energy release by the rearrangement depends little on the magnetic field. Thus, the growth rate also depends little on . This mode is not taken into account in Hanawa & Tomisaka (2015), where the mass-to-flux ratio is assumed to be constant in the stability analysis for simplicity.

Although it is due to the self-gravity of the gas, the instability may not result in dynamical collapse, since the displacement is dominated by circulation. The change in gravity is due to the change in the cloud shape. This instability is similar to that found by Nagai, Inutsuka & Miyama (1998). They found that a magnetized sheet-like cloud is unstable due to the self-gravity even when the cloud magnetic field is very strong. We discuss the choice of boundary conditions in the next section.

## 4 Discussions

We used a very simplified model for the filamentary cloud in our stability analysis. The medium is approximated by an isothermal gas and turbulence is not taken into account. These simplifications make the MHD stability analysis feasible. When the magnetic field is perpendicular to the cloud axis, we need to distinguish the direction parallel to the field line from that perpendicular to both the cloud axis and the magnetic field. Thus, the normal mode analysis cannot be reduced to a one-dimensional problem, unlike the previous MHD stability analyses for a filamentary cloud by Hanawa et al. (1993); Fiege & Pudritz (2000). Our analysis is formulated as a two-dimensional problem and its computational cost is much higher than that for a one-dimensional problem.

Our analysis has shown that the stability depends critically on the boundary assumed for . The dependence comes from the fact that the filamentary clouds are connected with the outer boundary through the magnetic field line. The cloud is stable if the magnetic field is substantial and the field lines are fixed in a region far from it. On the other hand, the cloud is unstable even when the magnetic field is extremely strong if the field lines are free to move in a region far from it. This implies that the magnetic field outside the cloud is important for the stability. Such a dependence on the outer boundary was not found in the previous MHD analysis (e.g., Hanawa et al., 1993), since the magnetic field lines were confined within the cloud.

In this section, we discuss numerical accuracy in §4.1 and the effects of magnetic field orientation in §4.2.

### 4.1 Numerical Accuracy

First we examine the accuracy of the growth rate obtained in our numerical analysis. The numerically obtained values should contain some numerical errors due to the finite spatial resolution, . The spatial resolution is only moderately high because the full width at half maximum (FWHM) of the model cloud is . Despite the moderate spatial resolution, matrix has dimension

 nA = 3nxny+2nx+2ny+1=4961, (63)

so solving the eigenvalue problem takes considerable computation time. The computation time is roughly proportional to since the computation time increases in proportion to . Hence, it is not easy to achieve higher spatial resolution to reduce the numerical error in our 2D analysis.

We examine the case of to assess the error due to the finite spatial resolution. When , the equilibrium model is symmetric around the axis and therefore the perturbation is also symmetric. The eigenvalue problem reduces to a 1D problem, as shown in Appendix B. We can achieve very high spatial resolution to eliminate the numerical error due to the finite spatial resolution.

As shown in Appendix C, our numerical method has second-order accuracy in space. The numerical error is proportional to the square of the spatial resolution, , in our 1D analysis and is evaluated as

 Δ(σ2)4πGρc≃3.4×10−4(Δr0.1H)2, (64)

from Figure 2 for the mode of .

Table 13 compares the growth rate squared, , obtained in our 1D and 2D analyses. The spatial resolutions are and in the 1D and 2D analyses, respectively. Despite the large difference in the formal spatial resolutions, the difference in the obtained growth rates is as small as

 −1.1×10−3<(σ2)2D−(σ2)1D4πGρc<−3×10−4. (65)

This difference is comparable to the expected numerical error for the growth rate obtained in our 1D analysis. Thus, we estimate that the numerical error of our 2D analysis is about . Although the formal spatial resolution is , the effective spatial resolution is likely to be higher since the displacement is evaluated not only on the - and -axes but at various points in the plane.

The accuracy mentioned above is consistent with the spatial extent of the displacement. As shown in Figure 5, the displacement is large in the region , which is covered by grid points. The error due to the finite difference is estimated to be , since the difference equations are of second-order accuracy. The factor 2 is introduced to take account of the region of negative and .

The presence of the magnetic field may increase the numerical error since the magnetic field confines the displacement in an especially narrow region when the fixed boundary condition is applied. Nevertheless, the error is likely to be moderate since the eigenmode is well resolved, as shown in the previous section. We expect that the numerical error is smaller than and does not have serious effects on our analysis.

### 4.2 Effects of Magnetic Field Perpendicular to the Cloud

Next we discuss the dependence of the instability on the direction of the magnetic field. As shown in the previous section, the magnetic field works against fragmentation and reduces the growth rate of instability when it is perpendicular to the cloud axis. On the other hand, the magnetic field does not suppress the instability when it is parallel to the cloud axis (see, e.g., Stodólkiewicz, 1963; Hanawa et al., 1993). A toroidal or helical magnetic field stabilizes fragmentation of the cloud and reduces the growth rate of the sausage (axisymmetric) mode (Fiege & Pudritz, 2000), although it induces the kink (non-axisymmetric) mode. These modes are characterized by the geometry of the magnetic field, namely, whether the magnetic field is bent or compressed by the fragmentation.

Before applying our analysis to observed filamentary molecular clouds, we need to discuss the boundary conditions employed in our analysis. As shown in the previous section, the boundary conditions greatly affect the instability. When the free boundary condition is applied, the filamentary cloud is unstable even when the magnetic field is extremely strong. In other words, no magnetic field can suppress the instability if the displacement does not vanish in the region very far from the cloud. The unstable mode exists because the density decreases in proportion to , where denotes the radial distance from the cloud axis. In our model, the displacement travels from the cloud center to infinity through the Alfvén wave in a finite time,

 τA = ∫∞0√4πρ0B0dx=2π√2πρcHB0. (66)

The Alfvén transit time may be much longer in reality because the density is higher or the magnetic field is weaker in the region far from the cloud. Palmeirim et al. (2013) reported that the mean radial column density fits a Plummer-like density distribution,

 ρP = ρc[1+(r/R2flat)]p/2, (67)

where denotes the radius of the inner flat region. In the best-fit model, and for the B211 filament. The density slope is significantly shallower than that for the isothermal equilibrium, . Our model cloud is highly idealized, especially in the region very far from the cloud axis. The shallower slope might be due to turbulence, which is likely to be larger in a lower density gas. We again note that the striations running perpendicular to the main filament (see, e.g., Sugitani et al., 2011; Palmeirim et al., 2013) are unlikely to move easily by magnetic force since the Alfvén velocity is lowered by the relatively high density. If the Alfvén transit timescale is longer than the free-fall timescale, the fixed boundary condition may be appropriate. We need to account for the low-density region surrounding the filamentary cloud when applying our analysis.

Tomisaka (2014) made a model for a filamentary cloud in which the magnetic field lines lie in the plane perpendicular to the cloud axis. In his equilibrium model, the cloud is supported in part by the magnetic tension. The magnetic field is weaker at a greater distance from the cloud axis. In other studies with 3D MHD simulations, the magnetic fields are not uniform around the filamentary clouds formed (see, e.g., Nakamura & Li, 2008; Klassen, Pudritz & Kirk, 2017). If the density is very low in the periphery of the filamentary cloud, the magnetic field should be quasi force-free and thus nearly straight. It is also well known that the magnetic field is stronger where the column density is higher (see, e.g., Li et al., 2014).

Our result for the free boundary is quite similar to those obtained by Nagai, Inutsuka & Miyama (1998) and Fiege & Pudritz (2000), who studied the stability of an isothermal gas layer and a filament confined by the outer gas pressure, respectively.

Nagai, Inutsuka & Miyama (1998) assumed that the magnetic field is uniform and parallel to the gas midplane. This gas layer suffers from two types of instability: fragmentation due to displacement along the magnetic field and distortion of the surface layer. When the gas layer is thick and consequently mostly gravitationally bound, the instability due to displacement dominates. When the gas layer is thin and consequently pressure bound, the distortion dominates and the gas layer fragments into filaments. The magnetic fields remain straight as the instability due to distortion develops. As a result, the magnetized gas sheet is unstable and the growth rate remains constant in the limit of a strong magnetic field, as in the case of the free boundary condition in our stability analysis. Although the latter instability is due to the self-gravity, the filaments formed do not collapse directly, since the filaments expand in the direction normal to the initial sheet by the instability. It should again be noted that this instability develops even when the sheet is extremely thin and the self-gravity is not important. It seems that the this instability is analogous to the instability for the free boundary.

Fiege & Pudritz (2000) showed that their truncated Ostriker model is unstable against fragmentation. The model filamentary cloud is assumed to be surrounded by hot gas of negligibly small density. Their initial model was similar to that of Nagai, Inutsuka & Miyama (1998) except for the geometry and the magnetic field. When the pressure of the hot gas is comparable to that at the filament center, circulation dominates in the unstable mode, as in our model for the free boundary (see their Figure 4). Each dense clump expands in the radial direction as a result of the instability. The change in the gravitational potential is mainly due to that in the boundary between hot and cold gases (i.e., the boundary of the dense gas), as in our free boundary model. It seems that a relatively strong magnetic field increases the effective temperature of the gas through magnetic pressure and therefore we obtain similar results.

As shown in §2, the model cloud is subcritical when the plasma beta is as low as 0.405. It is well known that molecular clouds collapse dynamically only when the cloud is supercritical. Thus, the subcritical cloud forms clumps supported in part by the magnetic force, not by dynamically collapsing cores. In other words, fragmentation does not result in direct star formation.

The above argument implies that subcritical clumps can be formed via fragmentation of a filamentary cloud if the free boundary is the case. A subcritical filamentary cloud may suffers from fragmentation instability. The instability is not likely to result in the dynamical collapse to form gravitationally bound clumps supported mainly by magnetic field. Such clumps may form stars through quasi-static contraction, i.e., after substantial amount of magnetic field is liberated through the ambipolar diffusion.

Our analysis suggests to study low density region surrounding a filamentary cloud. Motion in the high-density region may propagate into the outer region through the Alfvén wave, when the filamentary cloud fragments. It would be meaningful to observe the velocity field around the filamentary structure. The motion perpendicular to the magnetic field may contain information on the motion in the dense region. The Alfvén speed is proportional to the inverse square root of the density (). Thus the motion is faster when the density is lower. The density distribution may be controlled by turbulence. If it is the case, turbulence may affect the stability of a filamentary cloud indirectly through the boundary condition.

We thank an anonymous referee for useful comments to clarify the arguments. This work was supported by JSPS KAKENHI Grant Number JP15K05032. \softwareLAPACK, Linear Package Algebra (Anderson et al., 1999)

## Appendix A Variational Principle

In this appendix, we prove that the growth rate is either real or pure imaginary. For this purpose, we integrate the inner product of Equation (30) and over the volume, where the asterisk denotes the complex conjugate. The integral is expressed as

 σ2I = WT+WG+WM, (A1) I = ∫+∞−∞∫+∞−∞∫2π/k0ρ|\boldmathξ|2dxdydz, (A2) WT = −c2s∫+∞−∞∫+∞−∞∫2π/k0ρ0\boldmathξ∗⋅\boldmath∇(δρρ0)dxdydz, (A3) WG = −∫+∞−∞∫+∞−∞∫2π/k0ρ0\boldmathξ∗⋅\boldmath∇δψdxdydz, (A4) WM = ∫+∞−∞∫+∞−∞∫2π/k0\boldmathξ∗⋅(δ\boldmathJ×\boldmathB0)dxdydz, (A5)

The integral is real and positive for any . The integral is also proved to be real and negative, since

 WT = c2s∫+∞−∞∫+∞−∞∫2π/k0\boldmath∇⋅(ρ0% \boldmathξ∗)(δρρ0)dxdydz (A6) = −c2s∫+∞−∞∫+∞−∞∫2π/k0|δρ|2ρ0dxdydz. (A7)

Here we use the boundary condition at infinity and Equation (15) in the integration by parts. Similarly, the integral is proved to be real and positive for any perturbation, since

 WG = ∫+∞−∞∫+∞−∞∫2π/k0\boldmath∇⋅(ρ0\boldmathξ∗)δψdxdydz (A8) = −∫+∞−∞∫+∞−∞∫2π/k0δρ∗δψdxdydz (A9) = −14πG∫+∞−∞∫+∞−∞∫2π/k0Δδψ∗δψdxdydz (A10) = 14πG∫+∞−∞∫+∞−∞∫2π/k0∣∣\boldmath∇δψ∣∣2dxdydz. (A11)

The integral is proved to be real and negative, since

 WM = −∫+∞−∞∫+∞−∞∫2π/k0δ\boldmathJ⋅(δ\boldmathξ∗×\boldmathB0)dxdydz (A12) = −14π∫+∞−∞∫+∞−∞∫2π/k0\boldmath∇×δ% \boldmathB⋅(δ\boldmathξ∗×\boldmathB0)dxdydz (A13) = −14π∫+∞−∞∫+∞−∞∫2π/k0δ\boldmathB⋅\boldmath∇×(δ\boldmathξ∗×\boldmathB0)dxdydz (A14) = −14π∫+∞−∞∫+∞−∞∫2π/k0∣∣δ\boldmathB∣∣2dxdydz. (A15)

Either or is assumed to vanish at infinity. corresponds to the fixed boundary, whereas corresponds to the free boundary. Thus, the square of the growth rate, , should be real since

 σ2 = WT+WG+WMI. (A16)

It is also clear from Equation (A16) that the filament can be unstable only through the self-gravity, .

## Appendix B Matrix Elements

This appendix describes each component of the vector , and matrixes , , and .

The -th element of the vector is set so that

 ζm = ξx,i−1/2,j,(i=1,2,…,nx and j=0,1,…,ny) (B1) m = i+nxj, (B2)

for . Indexes and are derived from by

 i = [(m−1)(modnx)]+1, (B3) j = m−inx. (B4)

Other elements of the vector, , are set as

 ζm = ξy,i,j−1/2,(i=0,1,…,nx and j=1,2,…,ny) (B5) m = i+1+(nx+1)(j−1)+nx(ny+1), (B6) i = {[m−nx(ny+1)−1](% modnx+1)}+1 (B7) j = m−i−nx(ny+1)nx+1, (B8)

for