Numerical MHD modelling of waves in solar coronal

  • Slides: 40
Download presentation
Numerical MHD modelling of waves in solar coronal loops Petr Jelínek 1 and Marian

Numerical MHD modelling of waves in solar coronal loops Petr Jelínek 1 and Marian Karlický 2 1 University 2 Astronomical of South Bohemia, Department of Physics České Budějovice Institute, Academy of Sciences of the Czech Republic Ondřejov

Outline n Motivation of numerical studies n Equations of magnetohydrodynamics (MHD) n Numerical methods

Outline n Motivation of numerical studies n Equations of magnetohydrodynamics (MHD) n Numerical methods & solutions n Results 1 D model – impulsively generated acoustic waves ¨ 1 D model – gravitational stratification ¨ 2 D model – impulsively generated acoustic waves ¨ 2 D model – modelling of wavetrains ¨ n Conclusions

Motivation of numerical studies – I. n Oscillations in solar coronal loops have been

Motivation of numerical studies – I. n Oscillations in solar coronal loops have been observed for a few decades n The importance of such oscillations lies in their potential for the diagnostics of solar coronal structure (magnetic field, gas density, etc. ) n The various oscillation modes in coronal loops were observed with highly sensitive instruments such as SUMER (So. HO), TRACE n The observed oscillations include propagating and slow magnetosonic waves. There also observations of fast magnetosonic waves, kink and sausage modes of waves

Motivation of numerical studies – II. n Coronal loop oscillations were studied analytically but

Motivation of numerical studies – II. n Coronal loop oscillations were studied analytically but these studies are unfortunately applicable only onto highly idealised situations n The numerical simulations are often used for solutions of more complex problems – these studies are based on numerical solution of the full set of MHD equations n Mentioned studies of coronal loop oscillations are very important in connection with the problem of coronal heating, solar wind acceleration and many unsolved problems in solar physics n Magnetohydrodynamic coronal seismology is one of the main reasons for studying waves in solar corona

MHD equations n In our models we describe plasma dynamics in a coronal loop

MHD equations n In our models we describe plasma dynamics in a coronal loop by the ideal magnetohydrodynamic equations n The plasma energy density n The flux vector

Numerical solution of MHD equations n n The MHD equations (1) – (4) are

Numerical solution of MHD equations n n The MHD equations (1) – (4) are transformed into a conservation form For the solution of the equations in conservation form exist many numerical algorithms including professional software such as NIRVANA, ATHENA, FLASH, . . (www. astro-sim. org)

Numerical methods – I. n There exist a lot of numerical methods used for

Numerical methods – I. n There exist a lot of numerical methods used for the solution of equations in conservation form in numerical mathematics n Generally we can use the two types of numerical methods ¨ explicit methods – calculate the state of a system at a later time from the state of the system at the current time easy to programming unstable in many cases ¨ implicit methods – find the solution by solving an equation involving both the current state of the system and the later one unconditionally stable difficult to programming (tridiagonal matrix; solution by Thomas algorithm)

Numerical methods – II. n We use only explicit methods in our calculations for

Numerical methods – II. n We use only explicit methods in our calculations for this reason we must use the artificial smoothing for the stabilisation of the numerical scheme n Some mathematical definitions of numerical methods for PDEs ¨ Consistency – the numerical scheme is called consistent if ¨ Convergence – the numerical method is called convergent if

Numerical methods – III. n For the solution of the MHD equations in a

Numerical methods – III. n For the solution of the MHD equations in a conservation form the methods of so-called flux limiters are used n These numerical methods are able to jump down the oscillations near sharp discontinuities and jumps n Generally, for the solution of PDE in conservation form in 1 D we can write

Numerical methods – IV. n Many authors often use the linear methods ¨ upwind

Numerical methods – IV. n Many authors often use the linear methods ¨ upwind scheme ¨ Lax-Wendroff scheme (downwind slope) ¨ Beam-Warming scheme (upwind slope) ¨ Fromm scheme (centered slope)

Numerical methods – V.

Numerical methods – V.

Numerical methods – VI. n n To avoid the “overshoots” we limit the slope

Numerical methods – VI. n n To avoid the “overshoots” we limit the slope by flux limiter methods ¨ minmod ¨ superbee ¨ MC ¨ van Leer And many others – van Albada, OSPRE, UMIST, MUSCL schemes

Numerical methods – VII.

Numerical methods – VII.

1 D model of acoustic standing waves n There exists a lot of types

1 D model of acoustic standing waves n There exists a lot of types of oscillations in solar coronal loop acoustic oscillations ¨ kink and sausage oscillations ¨ fast and slow propagating waves, . . . ¨ n Acoustic oscillations are easy to simulate, they can be modelled in 1 D, without magnetic field, etc. n Kink and sausage oscillations were directly observed (SOHO, TRACE) and there are many unanswered questions – excitation and damping mechanisms, etc. n We focused on the impulsively generated acoustic standing waves in coronal loops

1 D model – initial conditions n Initial condtitions n The length of the

1 D model – initial conditions n Initial condtitions n The length of the coronal loop was L = 50 Mm which corresponds to loop radius about 16 Mm. n The loop footpoints were settled at positions x = 0 and x = L

1 D model – perturbations n Perturbations n In the view of our interest

1 D model – perturbations n Perturbations n In the view of our interest to study impulsively generated waves in the solar coronal loops, we have launched a pulse in the pressure and mass density n The pulse had the following form

1 D model – numerical solution n The numerical region was covered by a

1 D model – numerical solution n The numerical region was covered by a uniform grid with 2 500 cells and open boundary conditions that allow a wave signal freely leave the region were applied n The time step used in our calculations satisfied the Courant. Friedrichs-Levy stability condition in the form n In order to stabilize of numerical methods we have used the artificial smoothing as the replacing all the variables at each grid point and after each full time step as

Results – 1 D model Time evolution of velocity v(x = L/4, t), mass

Results – 1 D model Time evolution of velocity v(x = L/4, t), mass density r(x = L/4, t) (top panels) and spatial profiles of velocity v(x, Dt), v(x, 7. 12 T 1) (bottom panels); all for mass density contrast d = 108, pulse width w = L/40, and initial pulse position x 0 = L/4.

Results – 1 D model Time evolution of velocity v(x = L/4, t), mass

Results – 1 D model Time evolution of velocity v(x = L/4, t), mass density r(x = L/4, t) (top panels) and spatial profiles of velocity v(x, Dt), v(x, 7. 89 T 2) (bottom panels); all for mass density contrast d = 108, pulse width w = L/40, and initial pulse position x 0=L/2.

Results – 1 D model Time evolution of velocity v(x = L/4, t), mass

Results – 1 D model Time evolution of velocity v(x = L/4, t), mass density r(x = L/4, t) (top panels) and spatial profiles of velocity v(x, Dt), v(x, 11. 00 T 1) (bottom panels); all for mass density contrast d = 108, pulse width w = L/40, and initial pulse position x 0 = L/50.

Results – 1 D model Time evolution of average pressure, increased by the factor

Results – 1 D model Time evolution of average pressure, increased by the factor 103, initial pulse position x 0 = L/4 (left top panel), x 0 = L/2 (right top panel) and x 0 = L/50 (bottom panel), mass density contrast d = 108 and pulse width w = L/40; note that x-axis is in the logarithmic scale.

Results – 1 D model Fourier power spectra of velocities v for initial pulse

Results – 1 D model Fourier power spectra of velocities v for initial pulse position x 0 = L/2 (left) and x 0 = L/4 (right), mass density contrast d = 105 (top panels) and d = 108 (bottom panels) and pulse width w = L/40. The amplitude of the power spectrum A(P) is normalized to 1.

Results – 1 D model Time evolution of total (red), pressure (blue) and kinetic

Results – 1 D model Time evolution of total (red), pressure (blue) and kinetic (green) energies for various positions in numerical box. Left upper panel – whole simulation region, left upper panel – “transition region”, bottom panel – “coronal region”. The initial pulse position x 0 = L/4, d = 108 and pulse width w = L/40.

Results – 1 D model Time evolution of total (red), pressure (blue) and kinetic

Results – 1 D model Time evolution of total (red), pressure (blue) and kinetic (green) energies for various positions in numerical box. Left upper panel – whole simulation region, left upper panel – “transition region”, bottom panel – “coronal region”. The initial pulse position x 0 = L/2, d = 108 and pulse width w = L/40.

1 D – gravitational stratification n To create more realistic model the gravitational stratification

1 D – gravitational stratification n To create more realistic model the gravitational stratification was added n We consider a semi-circular loop with the curvature radius RL, in this model we incorporate the effect of loop plane inclination the shift of circular loop centre from the baseline was omitted

1 D – gravitational stratification – I. n The MHD equation of motion has

1 D – gravitational stratification – I. n The MHD equation of motion has the following form n The gravitational acceleration at a distance s measured from the footpoint along the loop, is n For the plasma pressure in the loop we can write

1 D – gravitational stratification – II. n The temperature profile was calculated by

1 D – gravitational stratification – II. n The temperature profile was calculated by means of this formula n The mass density was calculated from n The length of the coronal loop was L = 100 Mm in this case which corresponds to loop radius about 32 Mm.

Gravitational stratification – first results in 1 D Time evolution of velocity v(x =

Gravitational stratification – first results in 1 D Time evolution of velocity v(x = L/4, t), mass density contrast d = 102, pulse width w = L/80, and initial pulse position x 0 = L/4 and x 0 = L/2, inclination angle a = 0° (blue line) and a = 45° (red line).

2 D modelling of magnetoacoustic standing waves n n We consider a coronal slab

2 D modelling of magnetoacoustic standing waves n n We consider a coronal slab with a width w = 1 Mm and mass density ri, embedded in a environment of mass density re The pressure, mass density, temperature and initial pulses in pressure and mass density are calculated similarly as in 1 D model

Numerical solution in 2 D n For the solution of 2 D MHD equations

Numerical solution in 2 D n For the solution of 2 D MHD equations the Lax-Wendroff numerical scheme was used, this method is often used for the solutions of MHD by many authors n Step 1 n Step 2 n The stability condition

Results – 2 D model Time evolution of velocity v(x = L/4, y =

Results – 2 D model Time evolution of velocity v(x = L/4, y = 0, t) (left top panel). Spatial profile of x-component of velocity – vx at time t = 8. 17 T 1 (right top panel) and the corresponding slices of vx along y = H/2 (x = L/2) – bottom left (right) panel; all for mass density contrast d = 108, pulse width w = L/40, and initial pulse position x 0 = L/2.

Results – 2 D model Time evolution of velocity v(x = L/4, y =

Results – 2 D model Time evolution of velocity v(x = L/4, y = 0, t) (left top panel). Spatial profile of x-component of velocity – vx at time t = 6. 15 T 2 (right top panel) and the corresponding slices of vx along y = H/2 (x = L/4) – bottom left (right) panel; all for mass density contrast d = 108, pulse width w = L/40, and initial pulse position x 0 = L/4.

Modelling of wave trains in 2 D n The wave trains were directly observed

Modelling of wave trains in 2 D n The wave trains were directly observed and discovered by SECIS (Solar Eclipse Coronal Imaging System) n Observed in Ondrejov in radio waves n The theoretical description is needed – the comparison of observed and modelled tadpoles → what type of waves are present

Modelling of wave trains in 2 D n We study impulsively generated propagating along

Modelling of wave trains in 2 D n We study impulsively generated propagating along a coronal loop magnetoacoustic wave trains n The problem is modelled by means of 2 D model presented before, but magnetic field is parallel to the y axis n The equilibrium is perturbed by a pulse in velocity, situated at L/4 of the numerical domain

Wave trains – first results – I. The spatial profile of the velocity vx

Wave trains – first results – I. The spatial profile of the velocity vx at time t = 30 s from initial pulse (left upper panel), and corresponding slices of vx along x axis (y = H/2) (right upper panel) and along y axis (x = L/4). Initial pulse position x 0 = L/4, mass density contrast d = 108, pulse width w = L/40

Wave trains – first results – II. Time evolution of mass density r(x =

Wave trains – first results – II. Time evolution of mass density r(x = L/2, t), (top panel) and corresponding wavelet analysis (bottom panel); all for mass density contrast d = 108, pulse width w = L/40, and initial pulse position x 0 = L/4.

Conclusions – I. n Computer modelling seems to be very useful tool for the

Conclusions – I. n Computer modelling seems to be very useful tool for the understanding of processes in solar coronal loops n The next step in our research will be the extension of current model to three dimensions (by means of mentioned software – Athena, Nirvana, FLASH. . . ), including the source terms such as cooling term, heating term, gravitational stratification, etc. n By means of this model we could investigate effects like attenuation of waves in coronal loops, plasma energy leakage by the dissipation into solar atmosphere and more very interesting problems in solar coronal physics. . .

Conclusions – II. n More informations about 1 D or 2 D models can

Conclusions – II. n More informations about 1 D or 2 D models can be found in ¨ Jelínek P. , Karlický, M. : Numerical Modelling of Slow Standing Waves in a Solar Coronal Loop, Proc. 12 th ESPM, Freiburg, Germany, 2008 ¨ Jelínek, P. , Karlický, M. : Computational Study of Implusively Generated Standing Slow Acoustic Waves in a Solar Coronal Loop, Eur. Phys. J. D, after revisions.

References n [1] M. Aschwanden, Physics of the Solar Corona (Springer, Praxis Publ. ,

References n [1] M. Aschwanden, Physics of the Solar Corona (Springer, Praxis Publ. , Chichester UK 2004). n [2] T. J. Chung, Computational Fluid Dynamics (Cambridge University Press, New York USA 2002). n [3] E. R. Priest, Solar Magnetohydrodynamics (D. Reidel Publishing Company, London England 1982). n [4] M. Selwa, K. Murawski, S. K. Solanki, A&A 436, 701 (2005). n [5] Tsiklauri, D. , Nakariakov, V. M. , A&A, 379, 1106 (2001). n [6] Nakariakov, V. M. et al. : Mon. Not. R. Astron. Soc. , 349, 705 (2004).