# Finite difference scheme for 2D sound propagation

I am simulating the sound wave propagation in non-rectangular and asymetric spaces using finite-difference method. I presume linear acoustics equations to be enough (i.e. $\Box p = 0$, $\Box \vec{v} = 0$) and I need to have open and closed boundaries of the system (e.g. node or antinode of acoustic pressure on the wall and the same for velocity).

The problem occurs when I have to deal with the Neumann boundary conditions. The method is unstable even for dense and computationally inconvenient grids.

The physics behind the maths allows 3+1 more or less straightforward possible ways how to do the simulation:

• Model the acoustic pressure (then the Neumann condition is needed at the open/loose boundaries: $\partial_x p = 0$ on the boundary, closed boundaries are OK)
• Model the velocity potential $\phi$ (not very convenient cause $\partial_x \phi = 0$ on the open boundary and $\partial_t \phi = 0$ on the closed boundary)
• Model the acoustic velocity (i.e. model d'Alembert wave equation for the vector quantity - not very elegant, that's why we have potential)
• Reformulate the physics to be e.g. hamiltonian (i.e. have only the first order derivatives), lattice Boltzmann etc. etc.

Any suggestions or experiences please?

Edit:

CFL condition is satisfied.

The problem actually might be the numerical dispersion. The picture shows a system with the source on the left oscillating from -1 to 1 sinusoidally. The values on the picture outside the source area vary between -0.1 and 0.1. On the top and bottom there are Dirichlet boundaries, there is Neumann boundary at the right.

FDM with centered time and space was used. In matlab (dx=dy):

A(ii,jj,n+1)=2*A(ii,jj,n)-A(ii,jj,n-1)+c1^2*dT^2/dX^2*(A(ii-1,jj,n)+A(ii+1,jj,n)+A(ii,jj-1,n)+A(ii,jj+1,n)-4*A(ii,jj,n));

• could you please elaborate a little more? what do you mean by "the method is unstable even for dense and computationally inconvenient grids"? When I did acoustic wave modelling using finite differences there was a well known formula to calculate the timestep and the grid discretization space so that the system would be stable. This was unrelated to the boundary condition. Are you sure your system is unstable? or could you have a bug in your code? I just saw the "fluid-dynamics" tag, so this may be a totally different problem to the one I was solving (earth-sciences related) – jbcolmenares Oct 9 '15 at 9:03
• @Victor Pira, What kind of time integrator are you using? I'm not used to deal with sound/wave propagation, but what I remenber is that FTCS is not that great for it, but Crank-Nicholson is fine. – Hydro Guy Oct 12 '15 at 0:02
• Edits have been made. I haven't find anything satisfactory with Crank-Nicholson. It's more suited for parabolic problems. – Victor Pira Oct 12 '15 at 14:05
• If you are concerned that the errors in your approximation are due to the dissipative and dispersive properties of your finite difference scheme you need to look in Fourier space so that you can quantify the errors and know what to expect when you run your simulations. The book by Strikwerda (Chap 5) gives a nice outline of how to perform such analysis. – Andrew Winters Oct 13 '15 at 6:58