A&A 374, 337-347 (2001)
DOI: 10.1051/0004-6361:20010724

Spatial aspect of wave transformations in astrophysical flows

G. Bodo1 - S. Poedts2,[*] - A. Rogava3,4,[*] - P. Rossi1

1 - Osservatorio Astronomico di Torino, Strada dell'Osservatorio 20, 10025 Pino Torinese, Italy
2 - Centre for Plasma Astrophysics, K.U. Leuven, Celestijnenlaan 200B, 3001 Heverlee, Belgium
3 - Dipartimento di Fisica Generale, Universitá degli Studi di Torino, Via Pietro Giuria 1, 10125 Torino, Italy
4 - Abdus Salam International Centre for Theoretical Physics, 34014 Trieste, Italy

Received 28 March 2001 / Accepted 17 May 2001

The phenomenon of Shear Induced wave Transformations (SITs) (Chagelishvili et al. 1996), is a common feature of flows, sustaining n>1 mode of wave motion. Until now this "nonmodal'' phenomenon was described only in terms of the time evolution of individual Fourier harmonics of perturbations in the space of wave numbers (k-space). In this paper we present the results of the first, direct numerical simulations of SITs, giving compelling evidence in favor of the robust and recognizable character of SITs. The importance of these results for an understanding of physical processes in astrophysical shear flows is pointed out. Namely, concrete astrophysical situations (solar plasma flows, galactic gaseous disks, accretion disks and accretion-ejection flows), where this approach may help to shed some light on observational appearances of related objects, are indicated and discussed.

Key words: magnetohydrodynamics (MHD) - waves

1 Introduction

In astrophysics we are frequently challenged to consider physical processes in shear flows: flows of astrophysical plasmas with spatially inhomogeneous rates. Recent studies of the physics of shear flows have revealed a number of peculiar, shear-induced processes, that might be observationally important in the astrophysical context (see Rogava et al. 2000 and references therein). Therefore, the study of astrophysical shear flows (ASF) is important both from physical and astronomical points of view. For the physical theory, these studies could provide ample range of observational tests, since ASF are so diverse in terms of physical conditions and a spectrum of physical parameters. For astronomy, these studies could possibly provide a clue to the understanding of observational appearances of related astronomical objects.

In many cases of practical importance ASF maintain different modes of collective behavior: vortices, waves and instabilities. When a flow sustains more than one mode of collective behavior, for instance two or more different kinds of wave modes, then the velocity shear may couple these waves with each other and may lead to their mutual transformation. Originally, this effect was discovered for a magnetohydrodynamic (MHD), two-dimensional (2-D) shear flow: it was shown that the shear couples the slow magnetosonic waves (SMWs) and the fast magnetosonic waves (FMWs) and that under favorable conditions, namely, when the speed of sound $C_{\rm S}$ and the Alfvén speed $C_{\rm A}$ are of the same order of magnitude, the coupling leads to the efficient reciprocal transformations of these waves (Chagelishvili et al. 1996).

It was apparent from the very beginning that the SITs may easily occur in a number at astrophysical situations and may have diverse and characteristic observational appearances. One of the first attempts of astrophysical application was related to pulsar physics. Namely, it was shown that SITs couple Langmuir oscillations with electromagnetic waves in relativistic electron-positron plasmas and the process was suggested (Mahajan et al. 1997) as an alternative mechanism for the generation of pulsar radio emission. Later it was shown that in the solar wind there exist favorable conditions for the mutual transformation of MHD waves and it was argued that the observed large-amplitude Alfvén waves (AWs) in the solar wind may be generated as a result of the SITs of photospheric SMWs and FMWs (Poedts et al. 1998). It was surmised that the SITs may play an important role in the generation of MHD waves, coronal heating and the acceleration of the solar wind (Rogava et al. 2000). Another prominent astrophysical example of a magnetized ASF is a galactic gaseous disk that sustains both slow ( ${\sim}0.2\;$ kms-1) and fast ( ${\sim}20\;$ kms-1) MHD density waves (Fan & Lou 1996). Theoretical studies imply (Fan & Lou 1997; Rogava et al. 1999) that these waves may be coupled by the shear and may transform into each other.

The mathematical aspect of the SITs is nontrivial - it belongs to the so called "nonmodal processes''. These phenomena are not accessible by means of the traditional normal modes approach because of a fundamental quality of shear flow systems: non-self-adjointness, or "non-normality" (Trefethen et al. 1993) of their linear dynamics. These systems exhibit some interesting physics on finite time scales associated with existing modes of collective behavior and/or appearance of completely new processes (see, e.g., Rogava et al. 2000 and references therein). The so called nonmodal approach became a convenient instrument for the disclosure of these phenomena (Mahajan & Rogava 1999; Rogava et al. 2000). This method allows one to reduce the original system of partial differential equations for physical perturbation variables F(x,y,z;t) to the set of non-autonomous ordinary differential equations (ODEs) for the functions $\hat F[{\vec k}(t),t]$ - amplitudes of individual Fourier harmonics corresponding to the Fourier decomposition of perturbations.

However, an important disadvantage of all "nonmodal" schemes is that although they help to disclose new phenomena, providing us with basic understanding of a temporal span of involved processes; they do not give us any information about the spatial appearance of these phenomena. We do know that SITs are taking place in shear flows, but we still know next to nothing about the appearance of these phenomena in real physical space. This is a non-trivial step ahead, because, when we talk about real space, we have to consider the superposition of many Fourier components that undergo the process of mode conversion at different moments in time and we have also to take into account the effects of wave propagation. Therefore, before addressing in detail the concrete situations pertaining to the different astrophysical contexts described above, we think that an essential step is to consider an ideal case in which we can examine whether mode transformation can indeed be observed also in real space and for which we can develop the tools for analysing the process. This is the aim of the present paper and in the next section we describe a simple 2-D MHD setup, corresponding to Chagelishvili et al. (1996) where SITs were originally identified, which allows as to study mutual transformations of SMWs and FMWs. We write down the system of partial differential equations for the perturbations of physical variables and describe briefly the "nonmodal" evolution of the amplitudes of individual Fourier harmonics. The third section contains a description of the numerical simulation of the initial set of partial differential equations. Detailed analysis of the results shows that the initial compact "package'' of SMWs undergoes an efficient transformation into FMWs. The latter waves manage to persist and even augment their strength despite the presence of nonzero viscous dissipation in the system. The latter fact clearly indicates that the transformed FMWs, in full accordance with the nonmodal behaviour of Fourier amplitudes, efficiently extract energy from the mean flow.

In the final section we discuss the importance of these results that give the first direct and firm evidence, obtained by means of numerical simulations, in favor of the efficiency of SITs. Possible applications of these results in astrophysics are outlined and critically discussed. We think that solar plasma flows, galactic gaseous disks, accretion-ejection flows and pulsar magnetospheric plasmas are astrophysical instances where one may look for SITs and their observational appearances.

2 General formalism

We consider a standard MHD system:

\begin{displaymath}D_{\rm t}{\rho}+{\rho}{\nabla}{\cdot}\vec{V}=0,
\end{displaymath} (1)

\begin{displaymath}{\rho}D_{\rm t}\vec{V}=-{\nabla}{\left(P+{{B^2}\over{8{\pi}}}\right)}
+{1 \over {4{\pi}}}(\vec{B}{\cdot} {\nabla}){\vec
\end{displaymath} (2)

\begin{displaymath}D_{\rm t}\vec{B}=(\vec{B}{\cdot}{\nabla})\vec{V}-\vec{B}{\nabla}{\cdot} \vec{V},
\end{displaymath} (3)

\end{displaymath} (4)

where $D_{\rm t} \equiv
\partial_{\rm t}+\vec{V}{\cdot}{\nabla}$ is the convective (material) derivative. The above set of equations needs to be complemented by "an equation of state" relating the thermal pressure to the density of the plasma. In this study we adopt an adiabatic MHD fluid approximation and employ $D_{\rm t}(P{\rho}^{-{\gamma}})=0$ as a closure equation. This equation of state implies that the pressure and the density perturbations are linked through $p=C_{\rm S}^2{\rho}$ where $C_{\rm S}$ is the homogeneous speed of sound.

We study the evolution of 2D waves in compressible, unbounded parallel flow with uniform velocity shear contained in an external regular magnetic field $\vec{B}_0{\parallel} \vec{U}_0$. The mean velocity field is given by:

\begin{displaymath}\vec{U}_0 (Ay;~~0),
\end{displaymath} (5)

where the constant A (without loss of generality, we adopt A>0) parameterizes the linear inhomogeneity of the background flow. The linear phase of the temporal evolution of the system is governed by the following set of partial differential equations (Chagelishvili et al. 1996):

\begin{displaymath}{\cal D}_{\rm t}d + \partial_{x}u_{x} +
\partial_{y}u_{y} =0,
\end{displaymath} (6)

$\displaystyle %
{\cal D}_{\rm t}u_{x} + Au_{y} =
-C_{\rm s}^{2}\partial_{x}d,$     (7a)

$\displaystyle %
{\cal D}_{\rm t}u_{y}=
-C_{\rm s}^{2}\partial_{y}d +C_{\rm A}^2{\left ( \partial_{x}b_{y}
-\partial_{y}b_{x} \right )},$     (7b)

\begin{displaymath}{\cal D}_{\rm t}b_{y} =
\end{displaymath} (8)

\begin{displaymath}\partial_{x}b_{x} + \partial_{y}b_{y} = 0,
\end{displaymath} (9)

where ${\cal D}_{\rm t} \equiv \partial_{\rm t}+Ay\partial_x$, $d \equiv
{\rho}'/{\rho}_{0}$ and $\vec{b} \equiv \vec{B}{'}/\vert \vec{B}_0 \vert$, while $C_{\rm s}$ and $C_{\rm A}$ are the speed of sound and the Alfvén speed, respectively. In (7a-7b) the polytropic equation of state is used to express the pressure perturbation by means of the density perturbation: ${\nabla}p^{'}=C_{\rm s}^{2}{\nabla}{\rho}^{'}$.

The set of Eqs. (6)-(9) forms our basic system of equations. Before employing nonmodal techniques, we can derive from (6)-(9), still staying in the "physical space'', a pair of second order partial differential equations, which describes the coupling between the SMWs and the FMWs. Taking into account (6)-(9) we easily get:

$\displaystyle {\cal D}_{\rm t}^2 b_{y}-C_{\rm A}^2[\partial_x^2+\partial_y^2]
b_y=-C_{\rm s}^2\partial_x\partial_yd,$     (10)

and defining an auxiliary variable $ {\Psi} \equiv
\partial_xd+\partial_yb_y$ we find that the problem reduces to the solution of the following pair of second-order equations:

D}_{\rm t}^2b_{y}-[C_{\rm A}^2\partial_x^2+
(C_{\rm s}^2+C_{\rm A}^2)\partial_y^2]b_y=-C_{\rm s}^2\partial_y{\Psi},
\end{displaymath} (11)

\begin{displaymath}{\cal D}_{\rm t}^2{\Psi}-C_{\rm s}^2\partial_x^2{\Psi}=-C_{\rm s}^2
\end{displaymath} (12)

3 Transformations as seen in ${\vec k}$-space

At this stage we can employ standard "nonmodal" approach (see e.g., Chagelishvili et al. 1996). For the flow described by Eq. (5) the only nonzero component of the shear matrix (Mahajan & Rogava 1999) in this case is ${\cal S}_{12} =A$and, therefore, kx stays constant while the transverse component of the wave vector varies linearly in time:

Ky(t)=ky - Akxt. (13)

The amplitude of an individual Fourier harmonic of our auxiliary function ${\Psi}$ is related to the amplitudes of d and by via the relation $
{\hat {\Psi}}=i(k_x{\hat d}+K_y(t){\hat b}_y)$, while for the amplitudes of ${\Psi}$ and by from (11-12) we get the following pair of second-order ordinary differential equations:
$\displaystyle \partial_{\rm t}^2{\hat b}_y+[C_{\rm A}^2k_x^2+(C_{\rm s}^2+C_{\rm A}^2)K_y^2(t)]{\hat
b}_y= -iC_{\rm s}^2K_y(t){\hat{\Psi}},$     (14a)

$\displaystyle \partial_{\rm t}^2{\hat{\Psi}}+C_{\rm s}^2k_x^2{\hat{\Psi}}=iC_{\rm s}^2k_x^2
K_y(t){\hat b}_y,$     (14b)

which we can rewrite, switching to new variables $X \equiv {\hat b}_y,$and $Y \equiv -i{\hat{\Psi}}/k_x$, in a more symmetric form (and without imaginary coefficients):
$\displaystyle \partial_{\rm t}^2X+[C_{\rm A}^2k_x^2+(C_{\rm s}^2+C_{\rm A}^2)K_y^2(t)]X=C_{\rm s}^2k_xK_y(t)Y,$     (15a)

$\displaystyle \partial_{\rm t}^2Y+C_{\rm s}^2k_x^2Y=C_{\rm s}^2k_xK_y(t)X.$     (15b)

In this form this system resembles a standard pair of equations describing coupled oscillations (Magnus 1976; Morse 1981):
$\displaystyle \partial_{\rm t}^2X+W_X^2X+{\cal C}Y=0,$     (16a)

$\displaystyle \partial_{\rm t}^2Y+W_Y^2Y+{\cal C}X=0,$     (16b)

with WX and WY being the eigenfrequencies of the X and the Y oscillators, respectively; while ${\cal C}$ is the coupling coefficient. Obviously in our case the roles of effective eigenfrequencies and the coupling coefficient are played by the following quantities:
$\displaystyle W_X(t)=[C_{\rm A}^2k_x^2+(C_{\rm s}^2+C_{\rm A}^2)K_y^2(t)]^{1/2},$     (17a)

$\displaystyle W_Y=C_{\rm s}k_x,$     (17b)

$\displaystyle {\cal C}(t)=-C_{\rm s}^2k_xK_y(t).$     (17c)

Certainly this is not a standard coupled system with two degrees of freedom, because one of our eigenfrequencies and the coupling coefficient vary in time. However, the analogy with mechanics may be used for the preparation of initially "pure" normal modes (for the mathematical details of this procedure see Poedts et al. (1998) or the Appendix of this paper).

\par\resizebox{8.4cm}{!}{\includegraphics[angle=90]{MS1313f1.eps}} \end{figure} Figure 1: The evolution of the density (solid line) and the longitudinal (bx) magnetic field (dashed line) perturbations. The values of parameters are: ky/kx=3, $C_{\rm A}/C_{\rm s}=1$, and $R \equiv A/C_{\rm s}k_x=0.1$.
Open with DEXTER

\par\resizebox{8.1cm}{!}{\includegraphics[angle=90]{MS1313f2.eps}} \end{figure} Figure 2: The evolution of the total energy of perturbations (solid line) plotted together with the graphs for the time-dependent dispersion curves of the SMW (dashed line) and the FMW (dashed-dotted line) effective normal frequencies. The values of the parameters are the same as for Fig. 1.
Open with DEXTER

In order to illustrate how this method of the initial separation of wave modes works we present here an example of the solution when ky/kx=3, $C_{\rm A}/C_{\rm s}=1$, and $R \equiv A/C_{\rm s}k_x=0.1$. Figure 1 displays the evolution of the density (solid line) and the longitudinal (bx) magnetic field (dashed line) perturbations. It is clear that the initial perturbation is a pure SMW with phase difference being almost ${\pi}/2$ between these two oscillations. This "anticorrelation" of the density and the magnetic field perturbations is a "bona fide'' signature of the SMW (Fan & Lou 1997) and helps to distinguish the SMW from the FMW. The latter MHD mode, on the contrary, is characterized by an exact correlation between these oscillations. Obviously Fig. 1 clearly shows that around the time $t_{*} \equiv C_{\rm s}k_y/A=30$ a drastic change of the oscillation mode occurs and the initial pure SMW transforms into the FMW.

One of the notable features of the SITs is that in the case of sufficiently low shearing rates ( $R\ {\ll}\ 1$) the temporal evolution of the perturbations stays adiabatic both before and after the transformations. This is related to the fact that the shear-induced change of the oscillation frequencies is a much slower process than the oscillation time scales associated with both modes. This is illustrated by Fig. 2. The $E(\tau)$ graph is plotted together with the graphs for the time-dependent dispersion curves of the SMW (dashed line) and the FMW (dashed-dotted line) effective normal frequencies. The energy evolution is adiabatic, following the SMW/FMW dispersion curve before/after the transformation.

4 Transformations as seen in "real" space: Methods

The theory of SITs in Fourier space, summarized in the previous Sections, gives us important information on the process. However, a complete and vivid picture of the phenomenon can be obtained by examining the behavior of the solutions in real space for a wave-packet - superposition of different Fourier harmonics.

In this section we adopt this perspective and we approach the problem with two different methods: by a numerical solution of the full set of MHD equations and by a semi-analytical approach based on asymptotic expansions (Whitham 1974). These two methods will be described in the next two subsections.

4.1 Numerical solutions

We have solved numerically the fully nonlinear viscous and resistive MHD equations by using a central differencing scheme with fourth order spatial accuracy.

As an initial condition, we superimpose on the basic shear flow (5) a wave packet, that we define in k-space centered around a value (kx(0) = kx0, ky(0) = ky0). The form of the perturbations for the different physical variables in real space can be then obtained through the inverse Fourier transform and expressed in the following form:

\begin{displaymath}d(x,y) = \epsilon \int D(k_x, k_y)
F(k_x, k_y) \exp(ik_x x + i k_y y) {\rm d}k_x {\rm d}k_y
\end{displaymath} (18a)

\begin{displaymath}u_x(x,y) =
\epsilon \int U(k_x, k_y) F(k_x, k_y) \exp(ik_x x + i k_y y) {\rm d}k_x
{\rm d}k_y,
\end{displaymath} (18b)

\begin{displaymath}u_y(x,y)\, {=}\, \epsilon \int V(k_x, k_y) F(k_x, k_y)
\exp(ik_x x + i k_y y) {\rm d}k_x {\rm d}k_y
\end{displaymath} (18c)

\begin{displaymath}b_x(x,y)\, {=}\, \epsilon
\int B_1(k_x, k_y) F(k_x, k_y) \exp(ik_x x + i k_y y) {\rm d}k_x {\rm d}k_y
\end{displaymath} (18d)

\begin{displaymath}b_y(x,y)\, {=}\, \epsilon \int B_2(k_x, k_y) F(k_x, k_y) \exp(ik_x x + i k_y y)
{\rm d}k_x {\rm d}k_y,
\end{displaymath} (18e)

where F(kx, ky) is the shape of the wave packet in k-space and the functions $D, \, U, \, V$ and B give the relative amplitudes of the other physical variables, for the different Fourier components, and determine the relative strength of SMW and FMW for each of the components of the wave packet. In particular, if we want to excite initially a packet of "pure'' slow or fast magnetosonic waves, we have to use the expressions given in the Appendix. The explicit form for F(kx, ky) is also given in the Appendix. Finally $\epsilon$ gives the scaling of the perturbation, which will be taken to be sufficiently small to examine the behavior of linear waves.

The system we are considering does not have any characteristic scale, we have only the constraint that the characteristic length scale of perturbations $l \sim \vert k_0\vert^{-1}$ has to be much less than the size of the computational domain L. This is because we want to avoid the influence of domain boundaries on the shear-induced processes, which we are intending to simulate. We have therefore chosen an arbitrary unit of length and in this unit we have fixed kx0 = 2 and ky0 = 6, so that the central wavenumber of the wave-packet corresponds to the case we have presented in the previous section. Our unit of time will then be the sound crossing time over the unit of length. The Reynolds and magnetic Reynolds numbers are both equal to 500. The computational domain, defined by (-L < x < L, -L < y < L), with L = 40, is divided in $1024 \times 1024$ grid points and the boundary conditions are periodic along the direction of the basic flow and outflow in the transverse direction.

4.2 Semi-analytical method

We can follow the propagation of a wave packet using the theory of ray tracing described by Whitham (1974) and Lighthill (1978), under the assumption that the wavelengths are much smaller than the scale of variation of the basic flow. In this case we can describe the solution in terms of a slowly varying wavetrain, writing each quantity in the form

\begin{displaymath}{\cal A}(\vec{r}, t) \exp \left[ i \theta (\vec{r}, t) \right],
\end{displaymath} (19)

where $\theta ({\vec r}, t)$ is the phase function and its derivatives represent the local frequency and wavenumber of the wave packet, i.e.

\begin{eqnarray*}\omega = \frac {\partial\theta}
{\partial t} \;, \quad \vec{k} = - \nabla \theta\;.

In a shear flow the frequency $\omega_{\rm r}$ in a local reference frame moving with the local fluid velocity satisfies locally the dispersion relation

\begin{eqnarray*}\omega_{\rm r} = \omega_{\rm r} \left( - \frac { \partial\theta}
{\partial x}, - \frac {\partial\theta} {\partial y} \right),

with the two solutions of FMW and SMW, while

\begin{eqnarray*}\omega = k_x
U_0(y) + \omega_{\rm r},


\begin{eqnarray*}\frac {\partial\theta}
{\partial t} = \omega_{\rm r} \left( - \...
...artial y} \right) - \frac {
\partial\theta} {\partial x} U_0(y),

differentiating this equation with respect to x and y, and using the identity $\frac{\partial k_x}{\partial y} = \frac{\partial k_y}{\partial x}$, we get

\begin{displaymath}\frac {\partial
k_x} {\partial t} + \left( U_0 + \frac {\part...
...{\rm r}} {\partial k_y} \frac {\partial k_x} {\partial
y}= 0,
\end{displaymath} (20)

\begin{displaymath}\frac {\partial k_y} {\partial t} + \left(
U_0 + \frac {\part...
...} {\partial y} = - k_x \frac {\partial
U_0} {\partial y}\cdot
\end{displaymath} (21)

Considering the initial conditions described in the previous subsection, we have at t = 0 a set of wavenumbers concentrated near the origin. We can then follow their propagation using Eqs. (20)-(21), which can be written in characteristic form

\begin{displaymath}\frac {{\rm d}k_x} {{\rm d}t} = 0 ,
\end{displaymath} (22)

\begin{displaymath}\frac {{\rm d}k_y} {{\rm d}t} =
-Ak_x ,
\end{displaymath} (23)

on the rays defined by the equations

\begin{displaymath}\frac {{\rm d}x} {{\rm d}t} = U_0 + V_{{\rm g}x} ,
\end{displaymath} (24)

\begin{displaymath}\frac {{\rm d}y}
{{\rm d}t} = - V_{{\rm g}y} ,
\end{displaymath} (25)

where $V_{{\rm g}x}$ and $V_{{\rm g}y}$ are the x and the y components of the group velocity, respectively.

This method allows us to capture the propagation properties of FMWs and SMWs but not their transformation, since the solutions of the previous section show that the transformation process occurs quite abruptly when ky = 0, thus violating the adiabatic assumption made in the derivation above. Following the wave packet we have, therefore, assumed that each wave number initially starts as a FMW or a SMW (appropriate for the considered case) and then, during the propagation, when ky = 0 it switches to the other kind of wave.

5 Transformations as seen in "real'' space: Results

The study of the dynamics of a single Fourier harmonic, described in Sect. 2, allows us to analyze the process of mode conversion, but cannot capture the effects related to the wave propagation and their interplay with the process of conversion. This can only be achieved by studying the behavior of solutions in real space and in this Section we will indeed follow, as an example, the evolution of a wave packet, analyzing propagation and conversion effects. As discussed in Sect. 3, this has been done through a numerical solution of the full system of MHD equations. The initial conditions for the wave packet are described in Sect. 3. Initially, as we have discussed, the packet is composed mainly of slow waves and the amplitude of the perturbation of the x component of velocity is $2 \times
10^{-2}$ (small enough to consider it linear). Saying "mainly" we want to emphasize that the method of initial preparation of "pure" normal modes is valid only approximately. It is exact only for oscillatory systems with constant characteristics, while the parameters of our system vary adiabatically. Therefore, the initial wave package, "baked" according to the approximate "recipe", will be mainly a package of SMWs with a small dash of FMWs. Since group velocities of these two wave modes are quite different, this initial FMW impurity will drift much faster than the bulk of the package, composed of SMWs, and the two of them will be spatially separated.

\par {\includegraphics[angle=90,width=16.4cm]{MS1313f3.eps} }
\end{figure} Figure 3: Real-space evolution of initial (predominantly SMW) perturbation package, given by eight panels for eight different moments of time: t = 0, 8, 16, 24, 32, 40, 48 and 56, respectively.
Open with DEXTER

\par {\includegraphics[angle=90,width=16.7cm]{MS1313f4.eps} }
\end{figure} Figure 4: The shapes of the wave-packets, calculated by means of analytical approximation, at the same moments of time, for which the numerical solution is displayed in Fig. 3. Different tonalities of grey are used to distinguishing between: SMW transformed into FMW (the darkest), SMW that stays SMW (lighter), FMW transformed into SMW (still lighter), and FMW staying FMW (the lightest).
Open with DEXTER

The general evolution can be seen in Fig. 3, where we show images of the distribution of the y-component of velocity at eight different times. From the figure, we can see that initially the wave packet splits in four parts. The two central parts slowly propagate along the positive and negative x-direction. These two packets, as they propagate, also rotate in the clockwise direction. The rotation of the packet is associated with the drift of ky induced by the shear and discussed in the previous section (see Eqs. (22)-(23)). The group velocity for slow waves with wavenumber equal to the central value of the packet is mainly directed along the x-direction and this explains the propagation of the packet in this direction. The Fourier harmonics composing the packet evolve separately, each one follows its ky-drift described by Eq. (23) and, at different times, they will reach the value ky = 0 around which transformation to fast waves occurs. Fast waves have a large component of the group velocity along the y-direction and, therefore, as soon as transformation begins to occur, around t = 30, the packet begins to get an elongated component in the y-direction. The fast waves propagating along the y-direction are also dragged forward for positive y (or backward for negative y) by the increasing x-velocity component of the shear flow and thus the packet takes a curved shape (see panels 5-6-7-8 in Fig. 3).

In panel 2-3, alongside the slow packets that propagate along the x-direction, we observe also the presence of two other packets with a higher y-component of their propagation velocity. Their propagation can also be followed in the successive panels. The presence of these additional packets is due to the fact that, as we have discussed, the initial conditions we have imposed do not correspond exactly to pure slow waves but contain a small mixture of fast waves.

The changes in the shape of the wave-packet can be also captured using the ray-tracing approach described in the previous section. We start assuming that initially, at t = 0, all the wave number components of the packet are concentrated at the origin. Then we follow the propagation of each wave number according to Eqs. (22)-(25), and the shape of the packet at a given time is reconstructed from the position of all wave numbers at that time. We assume that initially in the wave-packet both slow and fast waves are present and we follow their propagation and their mutual transformations. Initially, the propagation velocities for each wave number of the slow and fast components are given by Eqs. (24)-(25) where $v_{\rm g}$ is taken to be the slow or the fast group velocity for that wave number, respectively.

As time progresses, Eq. (23) shows that ky changes its value and, when ky = 0, we assume that the transformation of slow to fast or of fast to slow waves may occur (as discussed above). Accordingly, hereafter, the group velocity is switched from slow to fast or from fast to slow. The transformation, however, is not complete. Slow waves are not totally transformed to fast waves or, conversely, fast waves are not totally transformed to slow waves. Therefore, we have followed also the propagation of slow and fast waves that do not undergo the transformation process. In Fig. 4 we show the shapes of the wave-packets at the same moments of time, for which the numerical solution is displayed in Fig. 3, so that a comparison between the exact numerical solution and the analytical approximation can be performed. Originally, when transformation has not yet occurred, we observe only the presence of two pairs of wave-packets (fast and slow waves), but, from the time when transformation begins to occur, we observe the presence of four pairs of wave packets (slow transformed into fast, slow that stays slow, fast transformed into slow and fast that stays fast). To distinguish between the different wave packets we have used for each of them a different tonality of grey: the darkest for slow transformed into fast, lighter for slow that stays slow, still lighter for fast transformed into slow and the lightest for fast that stays fast.

In comparing the numerical and analytical solutions, we have to keep in mind that we have somewhat different initial conditions, since for the analytical solution we assume essentially a point source at t = 0. However, we can see that the characteristics of the solution are captured quite well by the analytical approach. The comparison also allows a better interpretation of the numerical results since it enables the identification of the different components. For example, we can notice that in the numerical solution it is also possible to observe the transformation of the fast wave-packets in slow wave-packets.

\par {\includegraphics[angle=90,width=17.1cm]{MS1313f5a.eps} }
{\includegraphics[angle=90,width=17.1cm]{MS1313f5b.eps} }
\end{figure} Figure 5: The distribution of the y-component of the perturbed velocity for four different times, t = 5, 20, 35, 50, respectively. The distribution at each time (shown in the first column) is decomposed in its slow (second column) and fast (third column) components. The arrows and letters show the points where we have taken the cuts that are shown in Fig. 6.
Open with DEXTER

\par {\includegraphics[angle=90,width=18cm]{MS1313f6.eps} }
\end{figure} Figure 6: Cuts of the density and Bx perturbations, in the positions and along the directions shown by the arrows marked in Figure 5. In the positions of slow waves there is a substantial phase difference between density and magnetic field perturbations (see panels A, B, C, E), while at the location of fast waves the perturbations are almost in phase (see panels D, F).
Open with DEXTER

The firm evidence of the occurrence of the SITs is demonstrated by Fig. 5, where we show images of the distribution of the y-component of the perturbed velocity for four different times. The distribution at each time (shown in the first column) has been decomposed in to its slow and fast components and these are shown respectively in the second and third column of images of the figure. The decomposition has been obtained first projecting the Fourier components of the perturbation onto the eigenvectors for the slow and fast modes and then transforming back in real space the slow and fast components separately. We can see that initially the four wave-packets can be easily separated in their slow and fast components. At later times, we see that the elongated arcs extending in the y direction can be effectively identified as fast waves and a remnant of the slow mode is still found slowly "creeping" in the x direction.

Another strong signature of the SITs - drastic change of the phase characteristics of MHD waves (Fan & Lou 1997; Rogava et al. 1999) before and after the transformation - is illustrated in Fig. 6, where we plot cuts of the density and Bx perturbations, in the positions and along the directions shown by the arrows marked in Fig. 5. From the plots we can see that in the positions of slow waves we indeed find a substantial phase difference (opposition of phases) between density and magnetic field perturbations, in complete accordance with the theoretical expectations (see panels A, B, C, E). While at the location of fast waves, density and magnetic field perturbations are almost in phase (see panels D, F).

In Fig. 7 we show the behavior of the total perturbation energy defined as

\begin{displaymath}E = {1 \over 2} \int (d^2 + \vert\vec{u}\vert^2 + \vert\vec{B}\vert^2 ) {\rm d}x {\rm d}y,
\end{displaymath} (26)

as a function of time. We can compare the behavior shown in the figure with the behavior of single harmonics shown in Fig. 2. For single harmonics, energy stays constant before the transformation, and a strong increase follows after the transformation. In the case of the numerical calculations, we have to take into account that dissipation is present and that we are observing a packet of Fourier harmonics undergoing transformation at different times. Initially, we have a strong decrease of energy due to the action of dissipation. Around t = 30, we observe a flat minimum, when the different SMW harmonics are transforming into FMWs, and afterwards the energy keeps increasing at a lower rate than that observed in Fig. 2. In this last part we have again to take into account the action of dissipation that competes against the energy growth of FMWs, but this latter effect seems to prevail, implying persistence of the generated FMWs.

6 Discussion

Several years ago "nonmodal" study of MHD perturbations led to the discovery of velocity shear-induced wave transformations (SITs). This phenomenon, despite its linear nature, was suspected to be quite persistent and universal. However, direct numerical simulations that could confirm these expectations, could display SITs explicitly and could show their appearance within real physical space, were dramatically lacking.

In this article we report the first real-space simulations of the SITs. We found that the effect is indeed robust - it persists and it is strongly pronounced even in the presence of the nonzero dissipation. This allows us to assume that in real shear flows - laboratory MHD devices, terrestrial and astrophysical situations - this effect will not be obliterated by the presence of other physical phenomena.

Where, in astrophysics, could we expect to observe manifestations of the SITs? Below we discuss three classes of ASF - solar plasma flows, galactic gaseous disks and accretion-ejection flows - where theoretical modeling of these objects and the interpretation of their observational features might lead to the identification of the SITs.

The solar atmosphere is known to be particularly rich in diverse kinds of complicated oscillatory, transient or sustained motions throughout different kinds of flow patterns within the solar photosphere, chromosphere, transition region and corona (Pneuman & Orrall 1986; Thomas 1996). Solar plasma flows of different nature and geometry can sustain all three kinds of MHD waves, but only Alfvén waves seem to be capable of making the whole trip through the chromosphere, penetrate through the transition region and reach the low coronal layers. Recent observational results (Grall et al. 1996) indicate that the acceleration of the polar wind is almost completed by $10\,R_\odot$, implying that the acceleration of the solar wind and the heating of the solar corona occur in approximately the same region, where plasma motions, coupled with magnetic fields, are having significant impact on physical processes.

\par\resizebox{8.8cm}{!}{\includegraphics{MS1313f7.eps}} \end{figure} Figure 7: Temporal evolution of the total perturbation energy normalized to its initial value.
Open with DEXTER

Therefore, it is reasonable to admit that solar plasma flows might exhibit SITs (Poedts et al. 1998). Recently it was argued (Rogava et al. 2000) that FMWs, produced by high-frequency photospheric motions within isolated intense magnetic flux tubes, could become subject to shear-induced FMW-AW transformations in the chromosphere. The Alfvén waves, being weakly dissipative, and channelized by spicules, could penetrate through the transition region and reach those "heights" of the solar corona, where the main bulk of the coronal heating and the solar wind acceleration processes take place. The pressure of these waves, as well as compressional damping of "secondary" fast waves (the outcome of the opposite (AW-FMW) transformation processes) might contribute to the acceleration of the solar wind.

Now, having first encouraging results of direct MHD simulations, showing robust character of SITs, we could develop a detailed, 3-D numerical MHD model and check whether the regime of shear-induced MHD wave oscillations can indeed be realized in the solar wind (Poedts et al. 1998), in the inner corona and in various kinds of solar flux tube flows (Rogava et al. 2000). Can SITs help to generate large amplitude Alfvén waves penetrating through the transition region? The analysis of this problem involves a quite broad program, because it requires one to perform high-resolution simulations of the SITs with a properly chosen equilibrium setup, boundary conditions and realistic values of the solar plasma physical parameters.

Another challenging task would be to see how SITs behave? in galactic gaseous disks. It is believed that shear-induced processes play a significant role in the generation and maintenance of the galactic MHD density waves for gas-rich spiral galaxies (Fan & Lou 1996, 1997), where physical conditions seem to be favorable for the SITs. Thermal and magnetic energy densities are of comparable magnitude: for typical nearby galaxies the Alfvén ($C_{\rm A}$) and sound ($C_{\rm s}$) speeds are comparable and fall in the range $\sim $10-20 km s-1. This supports the existence of both fast ($\sim $20 km s-1) and slow ($\sim $0.2 km s-1) MHD density waves in such a disk system. Phase characteristics of mass density and azimuthal magnetic field perturbations help to distinguish galactic MHD density waves from one another (Fan & Lou 1997). In particular, for a FMW the enhancement of the surface density and parallel magnetic field are more or less in phase, while for a SMW there is a significant ( ${\ge}{\pi}/2$) phase difference between these perturbations. The "Whirlpool galaxy" M 51, for instance, has been identified as a clear case of a galactic disk with FMWs, while in NGC 6946 coherent large-scale magnetic spiral arms are located between the optical spiral arms leading us to believe that here we observe SMWs. Moreover, there is at least one example of a spiral galaxy where optical and radio continuum observations indicate the presence of both kinds of MHD density waves! This is the galaxy IC 342, where observations show the presence of SMWs in the inner (rigidly rotating) region of the galactic gaseous disk and FMWs in the outer (differentially rotating) region. Recently it was argued (Poedts & Rogava 2000) that this galaxy features the very first observational evidence of the SITs. Real-space simulations, with due consideration of the rotation curve and all its physical characteristics, taking into account the influence of Coriolis forces and the self-gravity, might result in the numerical spatial portrait of the wave evolution, which could be compared with existing observational data about magnetic field structure of this galaxy.

Yet another broad class of ASF, where SITs can be modeled numerically, are accretion-ejection flows. The possible role of SITs in the generation of ``mixed'' MHD turbulence in accretion disks has already been indicated (Rogava et al. 2000). Numerical models may help to understand what role is played by the SITs in the generation of the disk wind and the appearance of ejectional outflows in different kinds of accretion-ejection systems (Balbus & Hawley 1998; Ferrari 1998).

Therefore, apart from the general importance of the obtained results, which constitute the first, real-space simulations and verification of the SITs, they may serve as a steppingstone for the further, comprehensive study of shear-induced phenomena in different astrophysical situations. The attractive advantage of the real-space simulations in the astrophysical context is that they are much closer to the real observational picture than the analysis of the same phenomena in the space of wave numbers. In this way we may come much closer to the interpretation and understanding of some puzzling observational features and appearances of related astronomical objects.


Our work was supported, in part, by the INTAS grant No. 97-0504. Andria Rogava wishes to thank Abdus Salam International Centre for Theoretical Physics and Universitá degli Studi di Torino for supporting him, in part, through the Regular Associate Membership Award and the Assegno di Ricerca A02.070./99, respectively. The calculations have been performed on the Cray T3E at CINECA in Bologna, Italy, thanks to the support of CNAA.

Appendix A: Preparation of "pure" initial modes

In order to excite "pure'' modes in a coupled linear oscillatory system with two degrees of freedom, described by (16), one needs to calculate fundamental (normal) frequencies for t=0, given by the standard expressions (Magnus 1976):

\begin{displaymath}{\Omega}_{\mp}={1 \over
\sqrt{2}}{\sqrt{W_X^2+W_Y^2 {\mp}{\sqrt{{(W_X^2-W_Y^2)^2+4{\cal
\end{displaymath} (A.1)

It is easy to see that in this particular case ${\Omega}_{+}$ and ${\Omega}_{-}$ correspond to the standard frequencies of FMWs and SMWs, respectively.

When coefficients are constant, the general solution of (16) is a certain mixture of the normal modes:

$\displaystyle X=X_+\cos(\Omega_+t-\phi_+)+ X_-\cos(\Omega_-t-\phi_-),$     (A.2a)

$\displaystyle Y=\sigma_+X_+\cos(\Omega_+t-\phi_+)+\sigma_-X_-\cos(\Omega_-t-\phi_-),$     (A.2b)

where the auxiliary quantities $\sigma_\pm$ relate the oscillation amplitudes of the two normal modes to one other and are defined as:
$\displaystyle \sigma_\pm\equiv{\Omega_\pm^2-W_X^2\over{\cal C}}= {{\cal
C}\over\Omega_\pm^2-W_Y^2},$     (A.3)

while the $\phi_\pm$ are the initial phases of the coupled oscillators.

From these expressions we see that the system exhibits simple harmonic motion with one of the fundamental (normal) frequencies (${\Omega}_{+}$ or ${\Omega}_{-}$) when either X+ or X- is equal to zero. It immediately follows that:

$X_+\ne 0$ and X-=0, when $Y=\sigma_+X$, and $\partial_{\rm t}
Y=\sigma_+\partial_{\rm t} X$

$X_-\ne0$ and X+=0, when $Y=\sigma_-X$, and $\partial_{\rm t}
Y=\sigma_-\partial_{\rm t} X$.

Thus in order to excite initially a "pure" oscillation with one of two normal frequencies we should make an arbitrary choice for X(0) and $\partial_{\rm t}X(0)$ and take for Y(0) and $\partial_{\rm t}Y(0)$:

\begin{displaymath}Y(0)=\sigma_{-}X(0),~~~~\partial_{\rm t}Y(0)=\sigma_{-}
\partial_{\rm t}X(0).
\end{displaymath} (A.4)

This set of initial conditions will give us pure SMWs, while for the excitation of an initial pure FMWs we have to repeat the same procedure, but replace $\sigma_{-}$ by $\sigma_{+}$ in (A.4). If we express these conditions in terms of the physical variables, we obtain $d = (\sigma_- - k_y/k_x) b_y$, ux= iAby/kx, uy = 0, bx = -kyby/kx. While for the functions $D, \, U, \, V, \, B_1,$ and B2 we get:

\begin{eqnarray*}D(k_x, k_y) = \sigma_- - \frac{k_y}{k_x},
V(kx, ky) = 0,

\begin{eqnarray*}B_1(k_x, k_y) = -\frac{k_y}{k_x},

B2(kx, ky) = 1.

Finally the explicit form that we have chosen for the wave packet is

\begin{displaymath}F(k_x, k_y) =
... [4pt]
0 &\quad \hbox{if} \;\;\xi < 0
\end{array} \right. ,
\end{displaymath} (A.5)


\begin{displaymath}\xi = \left( 1 - \left( k_x - k_{x0} \right)^2 -
\left( k_y - k_{y0} \frac{k_x}{k_{x0}} \right)^2 \right).
\end{displaymath} (A.6)



Copyright ESO 2001