Edoardo Gabrielli – Home page

Menu

Alfvén waves in one dimension

Fortran code that simulates how nonlinear interaction deforms Alfvén waves in a plasma, with a LaTeX report. I did it on my own in one month, between August and September 2022.

Stack: Fortran (gfortran), Python and Matplotlib, LaTeX.

The problem

Alfvén waves are magnetohydrodynamic waves. In the Elsässer variables, Z+ and Z-, they are two waves travelling in opposite directions: if the equation were linear they would cross without disturbing each other; with the nonlinear terms they deform each other. I wanted to see how, by solving the equations myself.

Evolution of Z+ and Z- at five times, with and without nonlinearity

How it works

  1. 01

    The scheme

    Compact finite differences with a Lele filter; the derivatives require solving tridiagonal and pentadiagonal systems. About 525 lines of Fortran in separate modules.

  2. 02

    The report

    It derives the equations in the Elsässer variables and discusses the results (title: "Nonlinear effects on Alfvén waves in one dimension").

  3. 03

    Building and plots

    Compiled with gfortran and a Makefile; the plots are made in Python with Matplotlib.

Evolution of Z+ and Z- at five times, with and without nonlinearity
Z+ (top) and Z- (bottom) at five times; dashed, the linear solution. At t = 3.48 the profile steepens and steps appear: the filter is needed.
Animation of current density and energy spectrum in the Orszag-Tang vortex
Orszag-Tang at 512×512: current density j (left) and E(k) (right), up to t = 4.

2D extension (2026). In 2026 I extended the code to two dimensions with the help of Claude Code. I kept the 1D method: compact derivatives and Lele filter, applied one direction at a time on a periodic grid, with the same fourth-order Runge-Kutta. The variables are the Elsässer vorticities of Z+ and Z-. The only non-compact part is the Poisson solver, which is spectral. 512×512 grid, viscosity 5·10-4. Same code, two initial conditions: Orszag-Tang forms current sheets, a shear layer (weak magnetic field B0 = 0.1 along the flow) gives Kelvin-Helmholtz vortices. At this resolution the spectrum follows roughly k-3/2 only over a short range: it is not a fully developed turbulent cascade.

Animation of the Kelvin-Helmholtz instability with magnetic field lines
Shear layer with a weak magnetic field along the flow: vorticity and field lines, up to t = 6. The domain is shown twice along x.

  • What I did

    Everything: the code, the report and the plots. 21 commits.

  • Limits

    It is study code, not an application. I did not compare the results with other codes.

Code on GitHub
← All projects