Here is an example of my data... This is del(l) = total displacement of the system with time. >> attached u_t.png .... The small oscillations,.. Are they numerical artefacts?

No, I do not have an analytic solution to compare with... I will tell you my equations,

but what I was actually looking for is like if there is any way to be sure about numerical dispersion in the solution...

will post the eq s in a few hours...

Yes, it would be better to check your equation in order to say something. But as Mao-Lin say, commonly you compare you're discrete scheme with an analytical solution.

Comparing with the analytical solution is good for classroom situations.

If no analytic solution is available, I would try and compare different discretisation methods and/or different discretisation parameters (such as meshwidth).

The solutions would reflect the same solution behaviour with different amounts of numerical dispersion.

Central Difference in space DOES NOT work if your Peclet no.is larger than 2.

You have to refine the Grid swo that the cell peclet no is less than 2 OR preferably use upwind scheme

Mr. Hemker,

Analytical solutions doesn't exist for general situations, but it actually does for simple geometries and linear differential equations. The idea of study the dispersion and dissipatiob of a numerical scheme is to control the behavior of the wave in order to make a comparison.

See here , these pictures shows the early time developments in the variables >>

Another variable >>

Dear Deb,

What equation are you solving for?

I have two coupled pde s- described in the attached doc >>

Here

I see, your equation is more like a diffusion problem, so the term dispersion is more like a mathematical one. In diffusion problem the "dispersion relation" is not a straight line (cone) and then one always expects that a "signal" distorts while traveling. Have you try your numerical scheme in a linearized version of your PDE? Is there the dispersion fine?

In these two references the author touch this subject:

http://pauli.uni-muenster.de/tp/fileadmin/lehre/NumMethoden/WS0910/ScriptPDE/Heat.pdf

http://people.maths.ox.ac.uk/trefethen/5all.pdf

If these things are going fine... I would do something like Mr. Henker suggests, but instead of changing the mesh I would keep the mesh fix and change the signal (bandwidth and central frequency), then would try to characterize the numerical dispersion from these simulations.

@ Nicolas >> Thanks for the ref. I will go through them.. will take some time...

I computed some solutions for your system of equations (with a straightforward Method Of Lines) but I don't observe any of these wiggles you show. Can you specify a set of parameters (A,D1,D2,D3,ku0,ε) for which you find these (spurious) rapid oscillations?

@Mr. Hemker, by method of lines do you mean Euler's method?

I don't remember exactly parameter values of those plots I uploaded...

but I will generate some more results with dispersion and then let you know of the parameter values and the results. Thank you very much....

@Mr. Hemker here is the data...

Thanks for the data. Solving the problem for your parameters (and bell-shaped initial conditions for Y1[0] and Y2[0]=0), I used a variable-order variable-stepsize BDF-formula for time-integration. I do not observe the spurious wiggles, but it appears that the Jacobian becomes singular at about t=0.496. (see picture attached). This appearance of the singularity becomes earlier for a finer space-discretisation. So now I'm questioning the stability of the (analytic form of the) system of differential equations.

If we take A=-1 (instead of A=-2) there is no problem for 0

Yes... the parameter space for this system has three phases.... for A ~ small the solution will converge to the fixed point (0,0)... I think for the set of values I gave you

will show problem in numerical soln at A < -1.5 or so... that value (and the set I gave you), should give a advection dominated region.... and if you decrease A like A

ok, backward difference, i got it, I used central difference up to 6th order(fixed)...and adaptive time stepping like rkf45... but looked like did not work... to be sure once I have to check with matlab ODE45 result..

If you can use Wolfram's Mathematica, I can send you the program (34 lines, including graphics)

oh yes, I know some of mathematica. It would be a great help.

For a while the program will be downloadable from:

https://dl.dropboxusercontent.com/u/72331160/Deb_Banerjee_20130918-c2.nb

Let me know (by email) when you got a copy.

Central differences will oscillate whenever the Peclet no. is outside the range -2 to 2

Peclet No is the ratio of convection to diffusion. Therefore central differences are useless for highly convective flows. Either use upwind differences or refine your grid to reduce the Peclet no.

@ Mr. Markatos

Can you give me hints how should I be calculating Peclet number for a system of two coupled PDEs? Any trick to even estimate it. Any reference is also welcome.

Thank you,

Regards,

Deb.

Density times velocity times the length of your finite difference cell divided by your diffusion coefficient.Reference for example the book by Malalasekera

In the case of these (nonlinear) equations you should not care so much about the Peclet number. That isn't the problem here. The problem is the (inherent) stability of the system of equations. With A=-1 the stability is no problem (over the range of interest): we can solve it by discretisation with even a large h=0.5 (and smaller, of course). If we replace the constant of convection, A, by A=-1.5, the instability appears for t>1, even for a meshwidth h=0.1 (and smaller).

The fact that, if we increase the convection constant by a factor 1.5 and decrease the meshwidth by a factor 5, the problem changes from stable to instable, makes clear that the Peclet number is not the problem.

Moreover, if we compute the cell-Peclet number (which is (x,t)-dependent in this nonlinear case) we see that its absolute value remains less than 0.015, for A=-1.5 and h=0.1, the other constants as above in this discussion).

I'm attaching a picture of the cell-Pecletnumber in this case, where you also see the instability appear after some time.

Conclusion: it is the instability of the system of equations itself and (most probably) not of the numerical method.

Here is the attachment which shows the cell-Pecletnumber.

I have no idea of the numbers, but if Peclet=0.015 obviously it is Not the problem

Dear Sankar Banerjee

I have read (Num Meth. for Scientist and Engineers, J.D. Hoffman 2° ed.) that numerical dispersion can be spotted if the modified finite difference equation (MDE), which is used to assess consistency of the finite difference scheme you are using, has a factor (either delta t, delta x or a combination of these) multiplying an odd-order derivative and not going to zero as both delta t and delta x do. Should this happen with even-order derivatives, we would have the so-called numerical diffusion (also dissipation or damping).

Hope this helps.

Regards

Dear Mauricio Barrera ,

Yes, convection terms (1st order derivative in space) will always have some numerical dispersion error in many of the discretization  schemes. I finally used a Finite Volume method with van leer / min-mod flux limiter. The gist of this method is that it uses 2nd order accurate difference form for less steep parts of the density profile (the variable we usually solve for) and 1st order accurate diffusive difference form in parts where the density has high gradient values.. this diffusive form helps preventing dispersive and spurious peaks that result in non physical negative density values.

I would like to thank you all for helping me in the problem. The work related to this problem can be found here :

https://www.nature.com/articles/s41467-017-01130-1

Finally the equations were solved using Finite volume method with flux limiters.