Hello Folks, I'm currently trying to write a code to solve the one-way wave equation (1st order wave-equation):

u_x + u_t =0

with boundary condition: u(0,t)=sin(t) and initial condition: u(x,0)= sinx

I have to solve this on a very large domain (x=0 -> 500) with a small timestep (dt approx 0.001). I have a stretched grid with very small values of dx near the x=0 end of the domain and very large step sizes towards the x=500 end.

Could someone please suggest a scheme that would be capable of dealing with this problem? I've tried a couple of schemes but have found instabilities when I try to adapt them to work on a stretched grid. To clarify, my grid is cartesian and the x-spacing is linearly stretched such that dx(I+1)=k*dx(I) where K is approximately 1.05.

If anyone could advise then I'd be very grateful, F

Well first, check to to see what your largest delta x value is at the outflow boundary. Make sure that if you are solving the equation explicitly, that your CFL number is within the range that you expect.

Also, obviously make sure that you are not using a forward differencing scheme, this will cause instabilities instantly.

If for some reason, you have growing waves or anything like that, you could use a backward difference scheme: [F(i) - F(i-1)] / deltax will help the solution damp itself.

Other than that, exactly what type of instabilities are you getting?

Have a look at ICASE/LaRC workshop.It is available on NASA server.Many different numerical schemes have been discussed for application to the wave equation.

Could you write the schemes used? This may help more than a series of question.

As I understand, you want to solve a one_space-dimentional equation:

au/ax+au/at=0 (a means partial derivative)

this equation have a very simply analitical solution:

Any function moving with unitary speed.

Take any function u(x), and make it move in time composing with the function x-t:

Solution=u(x-t).

Since your initial condition is u(x-0)=sin(x), make this sin(x) function move: u(x,t)=sin(x-t).

You are in bad luck u(x=0,t)=sin(-t)=-sin(t). It does not match your boundary conditions. But because the information flows from left to right, it does not matter. take this function=

u=sin(-x-t) (for x<=0, t=0) u=sin(x-t) (for x=>0, t=0)

and move it to left with speed dx/dt=1

Hello folks, I appreciate the comments that people have left - it's a big help having a forum like this where we can discuss such things. I've not had computer access for the past few days so this is my first chance to reply to your posts.

I realise that there is an analytical solution to the problem that I posed in my original post. The problem I presented, however, is merely a model problem I'm trying to solve numerically with the upstream boundary condition u(0,t) = sin(-t). My hope is that as this solver will also work when the upstream is not prescribed but varies in a similar way to the sine curve. Guillermo was correct in stating that the upstream BC i stated previously should have been "u(0,t)=-sin(t)"; forgetting to include the -ve sign that was a typo on my part.

I have tried an explicit scheme (which blew up because of the large dx spacing downstream made everything unstable); a first order implicit scheme (backwards time, backwards space) as well as the crank nicolson scheme.

Plotting my numerical curves against the analytical solutions I find that near the upstream edge where the mesh is fine I have good agreement but as the mesh spacing increases downstream the numerical solution is not a very good match.

Does anyone know of a good scheme to solve this problem on a stretched cartesian grid or perhaps a way that I could adapt an existing scheme to make it suitable for this grid-stretching?

Many thanks, F

Do you want to do it as an exercise in finite difference, or really need to solve this problem, as part of some other problem?

Explicit schemes generally fail when the timestep At is bigger enough to allow the fluid to travel more distance than the space-step Ax (A meaning uppercase greek delta letter).

If your speed is v=1, try At<Ax.

If you need an efficient solver, on a constant Ax grid, for any boundary condition, I have developped an efficient discrete solver using the fast fourier transform, for the Muskingham method, that is very similar.

I have writted it in spanish, and plan to traslate it to english, and to wikipedia. Today is a ms word document (and an excel worksheet example).

If you are interested, just tell.

I want to generalize it to a Saint Vennant solver, but can't manage the nasty av²/ax term...

I wanted to write:

"If your speed is v=1, try At<Ax"

Also, you need to differentiate betwen "schemes with point values" (that need to be interpolated betwen points), and "averaged over each discrete step". They are similar algorithms, but need subtle adaptations.

I wanted to write:

"If your speed is v=1, try At less than Ax"

It looks like a wiki-bug does not accept the character "less than"