Quick von Neumann stability analysis

The von Neumann stability analysis is a great tool to assess stabilty of discretizing schemes for PDEs. But too often, imho, the discussion gets blurred with trigonometric expressions. Here, I try to provide a shortcut.

BTW: it’s “fon no ee man”

Continue reading


Discrete Fourier transform and Fourier series

This is quite silly, but the relationship between the discrete Fourier transform (DFT) and the Fourier series (FS) is clouded by annoying factors. I will try to connect both in this article. The motivation is to employ DFT techniques in a computer simulation. In the latter, one usually has a finite simulation box, which makes Fourier series the most interesting (a connection to the Fourier transform may also be made, see below).

Continue reading

Sound waves with attenuation

Just a simple derivation of the role of attenuation in the standard sound wave equation. Original work: Stokes, 1845.

Starting with the Navier-Stokes momentum equation

\frac{\partial }{\partial t} \mathbf{u} + \mathbf{u} \nabla \mathbf{u} = - \frac{1}{\rho} \nabla p + \frac{\mu}{\rho} \nabla^2 \mathbf{u} + \left(\frac{\lambda+\mu}{\rho}\right)\nabla (\nabla\cdot\mathbf{u}) ,

where \lambda is a Lamé viscosity coefficient. The bulk viscosity coeficient  is defined as  \zeta = \lambda + (2/3) \mu. The last term  is often neglected, even in compressible flow, but sound attenuation is one of the few cases where it may have some influence. All viscosities are assumed to be constant, but in this case this is a safe assumption, since we are going to assume small departures about equilibrium values.

Continue reading

Sparse Poisson problem in eigen

OwlgenBack to scientific computing. Lately, I have been using the Eigen libraries for linear algebra. They were recommended by the CGAL project, and they indeed share a common philosophy. Thanks to the rather high C++ level they can accomplish this sort of magic:

  int n = 100;
  VectorXd x(n), b(n);
  SpMat A(n,n);



  // solve Ax = b
  Eigen::BiCGSTAB<SpMat> solver;
  //ConjugateGradient<SpMat> solver;


  x = solver.solveWithGuess(b,x0);

Notice that A is a sparse matrix! I am next describing how to use this in order to solve the 1D Poisson equation.

Continue reading

The weak form


One can write a general differential equation as:


An example, which we will often consider is

ex. A(u)=-\frac{d}{dx}\left( a(x)\frac{d u }{dx}\right),

to model e.g. heat conduction with a space-dependent heat conduction coefficient a(x) and a heat production term f(x). (Notice the “ex” when we are talking about a particular example.)

This is pretty technical, so please stop reading here if you are already lost!

Continue reading