1D wave equation - bizarre problem

[smolloy]

1D wave equation -- bizarre problem!

I am trying to write a solver for a 1D wave equation in python, and I have run into a bizarre problem that I just can't find a way out of.

I start with the wave equation, and then discretise it, to arrive at the following,

phi(i,j+1) = delta_t²/delta_x² * [phi(i+1,j) - 2phi(i,j) + phi(i-1,j)] + 2phi(i,j) - phi(i,j-1)

I'm pretty sure this is correct.

I set boundary conditions such that the string is fixed at both ends (i.e. phi(0,t) = phi(L,t) = 0). As initial conditions I use phi(x,0)=0, and phi(x,delta_t) = dirac_delta (i.e. all zeros, except for a spike in the middle).

With delta_t and delta_x set such that the solver should be stable, I then run my code.

The problem is that the central point (the point at which I place the Dirac delta) moves upwards (as it should), before slowing down (as it should), stopping (as it should), and then remains there indefinitely (as it most definitely should *not*!).

Thinking a little bit more about the structure of the wave equation, this kinda makes sense. It is an equality between the net force due to the spring-constant tension on the curved string, and the F=ma like momentum of the string. Thus, a negative curvature should accelerate downwards, thus generating a positive curvature, which accelerates upwards.

In other words, the string never comes back to zero until the reflections drag it that way!

This, surely, isn't the way the string should behave -- surely it should break down into two separate waves moving in opposite directions, with the initial disturbance moving back to zero.

What have I done wrong, and where should I go from here?? :D

PS: Thanks for your attention to this long rambling post!

 

[Geomancer]

I suspect your intuition is more at fault than you numerics. Consider the function f=1 for positive x, -1 for negative x. Then the sum of a left moving copy of f (or a close smooth approximation) and a right moving copy of -f will give the behavior you describe.

 

[smolloy]

Thanks for your reply.

I think the problem is that I didn't really use the right initial conditions. When I tested the code with a sinusoidal initial condition, I got a standing wave as the result -- which means that the code is correct.

Then, when I expanded it to two spatial dimensions, I got something that looked much more as I imagined.

Thanks for your time.