[FMDenaro]

Originally Posted by TurbJet

I didn't get to test Lax-Wendroff. But for the upwind scheme, the convergence order for L_inf is almost zero (with fixed dt); The reason I think is due to the dispersion at these non-smooth regions (highlighted by red circle).

However, lowering dx/dt at the same time (where CFL = 0.9) can give the 1st order convergence in L_inf.

Attachment 65764

You are using a wrong periodic condition! if you see the profile it does not repeat periodically

[TurbJet]

Originally Posted by FMDenaro

You are using a wrong periodic condition! if you see the profile it does not repeat periodically

You are saying using the entire sine wave on the interval, instead of single period of it?

[FMDenaro]

Originally Posted by TurbJet

You are saying using the entire sine wave on the interval, instead of single period of it?

If you use periodic BC, the initial signal must be periodic, the function and all derivatives must be the same.
Therefore, given L the lenght of periodicity you can choose a wavenumber k=1,2... and setting the sine at 2*pi*k/L.

[TurbJet]

Originally Posted by FMDenaro

If you use periodic BC, the initial signal must be periodic, the function and all derivatives must be the same.
Therefore, given L the lenght of periodicity you can choose a wavenumber k=1,2... and setting the sine at 2*pi*k/L.

Ok, I fixed the problem. Now it gives 1st order convergence. But how does it relate to my original 4th-order scheme?

[FMDenaro]

Originally Posted by TurbJet

Ok, I fixed the problem. Now it gives 1st order convergence. But how does it relate to my original 4th-order scheme?

You create a discontinuity in the derivatives therefore the local truncation error cannot be defined and the accuracy order cannot be estimated.

[TurbJet]

Originally Posted by FMDenaro

You create a discontinuity in the derivatives therefore the local truncation error cannot be defined and the accuracy order cannot be estimated.

Ok, I'll try this sine wave for my original scheme and see if it can give me 4th order convergence. And I'll get back to you.

[TurbJet]

Originally Posted by FMDenaro

You are using a wrong periodic condition! if you see the profile it does not repeat periodically

I tried with the new sine wave. It gives me 1st order convergence. I looked into the solution in detail and it makes me feel that my periodic BC has some issues. I attach the schematic for my mesh below (sorry for the poor image quality).

As in the figure below (I use collocated mesh), I have N+1 interior points and the boundary points, 6 ghost points in total (3 on each sides).
As can see, point 4 and N+1 are boundary points, points 1~3 and N+5~N+7 are ghost points. I only evaluate the derivative on the interior points 5 ~ N+3.

As for the periodic BC, I do (assuming f is the original function I am solving)

f(N+4 : N+7) = f(5 : 8);
f(1 : 4)     = f(N : N+3);

Would you please take a look where it goes wrong? I just can't see it. Thanks a lot!

Screen Shot 2018-09-26 at 10.28.01 PM.jpg

[FMDenaro]

Assume that the periodicity lenght L goes from node i=1 at x=0 and i=N at x=L. Ath the boundaries, then you solve the equation only at the node i=1 (or i=N) and set the peridocity prescribing the same value at i=N. Of course you get


f(0)=f(N-1)
f(-1)=f(N-2)
....
f(N+1)=f(2)
f((N+2)=f(3)


the periodicity must work this way. Otherwise the bug is in your 4th order discretization.

[TurbJet]

Originally Posted by FMDenaro

Assume that the periodicity lenght L goes from node i=1 at x=0 and i=N at x=L. Ath the boundaries, then you solve the equation only at the node i=1 (or i=N) and set the peridocity prescribing the same value at i=N. Of course you get


f(0)=f(N-1)
f(-1)=f(N-2)
....
f(N+1)=f(2)
f((N+2)=f(3)


the periodicity must work this way. Otherwise the bug is in your 4th order discretization.

What do you mean by

solve the equation only at the node i=1 (or i=N)

My understanding is that I have to determine which way the wave goes: if the wave goes to the right, need to evaluate the derivatives on the right boundary point i=N, then I do

f(1) = f(N)

if the wave goes to the left, then I need to do

f(N) = f(1)

Am I correct?

[FMDenaro]

Originally Posted by TurbJet

What do you mean by


My understanding is that I have to determine which way the wave goes: if the wave goes to the right, need to evaluate the derivatives on the right boundary point i=N, then I do

f(1) = f(N)

if the wave goes to the left, then I need to do

f(N) = f(1)

Am I correct?

No. You need to solve the wave equation into the domain using the upwinding and that depends on the sign of the velocity in deciding the stencil for the derivative.
That has nothing to do with the setting of periodic BC.s that is totally independent from the scheme that is adopted and does not depend on the velocity.


Setting the correct periodicity means that if you solve the equation at i=1 and the velocity is, for example, positive, the derivative is (f(1)-f(0))/dx and the correct setting of the periodiciy make this exactly equal to (f(1)-f(N-1))/dx. Once computed the updated at i=1, just set the update at i=N without solving the equations.

[TurbJet]

Originally Posted by FMDenaro

No. You need to solve the wave equation into the domain using the upwinding and that depends on the sign of the velocity in deciding the stencil for the derivative.
That has nothing to do with the setting of periodic BC.s that is totally independent from the scheme that is adopted and does not depend on the velocity.


Setting the correct periodicity means that if you solve the equation at i=1 and the velocity is, for example, positive, the derivative is (f(1)-f(0))/dx and the correct setting of the periodiciy make this exactly equal to (f(1)-f(N-1))/dx. Once computed the updated at i=1, just set the update at i=N without solving the equations.

mmm, I see your point. I'll fix it and get back to you. Thanks a lot!

[TurbJet]

Originally Posted by FMDenaro

No. You need to solve the wave equation into the domain using the upwinding and that depends on the sign of the velocity in deciding the stencil for the derivative.
That has nothing to do with the setting of periodic BC.s that is totally independent from the scheme that is adopted and does not depend on the velocity.


Setting the correct periodicity means that if you solve the equation at i=1 and the velocity is, for example, positive, the derivative is (f(1)-f(0))/dx and the correct setting of the periodiciy make this exactly equal to (f(1)-f(N-1))/dx. Once computed the updated at i=1, just set the update at i=N without solving the equations.

Finally! I figured out where I went wrong! I forgot to apply BCs after each sub-time step of the RK scheme. Now I have 4th order convergence (with dx and dt lowering at the same time, CFL=1).
And I can get 4th order in space with fixed small time step.

But for the temporal convergence, at the beginning, the error is already somewhere 1E-8 (with fixed small grid size). And getting smaller gives the same magnitude of error, i.e., the order is O(0). How can I fix this?

[TurbJet]

Originally Posted by FMDenaro

So, you have a linear wave equation? You can compare the exact solution at any time as f(x,t)=f(x-ut,0). However, if you want to see the slope of the error curve ONLY for the time accuracy, you have to use a very refined spatial grid to eliminate the effect of the spatial term in the local truncation error.

Even though I have a 4th order convergence in space, it still tricky to get the temporal convergence.

I set grid size to a small value to avoid the error swamp the temporal error. However, since dt is bounded by CFL and dx. If set dx small, even with CFL=1, the error is already push to 1E-12. Lowering dt doesn't give a much better result.

If I am getting 4th order convergence with lowering dt dx at the same time, is it safe to assume that my time integral is 4th order as well? Since it's not compromising the total order.

[FMDenaro]

Originally Posted by TurbJet

Even though I have a 4th order convergence in space, it still tricky to get the temporal convergence.

I set grid size to a small value to avoid the error swamp the temporal error. However, since dt is bounded by CFL and dx. If set dx small, even with CFL=1, the error is already push to 1E-12. Lowering dt doesn't give a much better result.

If I am getting 4th order convergence with lowering dt dx at the same time, is it safe to assume that my time integral is 4th order as well? Since it's not compromising the total order.

Consistence requires that the LTE vanishes as the time and space steps go to zero therefore you are assessing that there is no bug in the scheme.
Start again the time accuracy check using the one-step error, this way you can use any CFL since there is no stability problem. You should be able to see the error slope.

[TurbJet]

Originally Posted by FMDenaro

Consistence requires that the LTE vanishes as the time and space steps go to zero therefore you are assessing that there is no bug in the scheme.
Start again the time accuracy check using the one-step error, this way you can use any CFL since there is no stability problem. You should be able to see the error slope.

I set (with convective velocity = 1)

CFL = 1000, error_max ~ 1E-4;
CFL = 100, error_max ~ 1E-8;
CFL = 10, error_max ~ 1E-12.

It does has 4th order in time now.

It just feel wired that CFL went to 1000....

[FMDenaro]

Originally Posted by TurbJet

I set (with convective velocity = 1)

CFL = 1000, error_max ~ 1E-4;
CFL = 100, error_max ~ 1E-8;
CFL = 10, error_max ~ 1E-12.

It does has 4th order in time now.

It just feel wired that CFL went to 1000....

How about the error going up to CFL=1 ?

[TurbJet]

Originally Posted by FMDenaro

How about the error going up to CFL=1 ?

When CFL = 1, error_max ~1E-16

Is it safe to say that my time integral does have 4th order and it's all good now??

[TurbJet]

Originally Posted by FMDenaro

How about the error going up to CFL=1 ?

Another issue, I can't really get 1E-16 for small grid size. The error hanging around 1E-12 ~ 1E-13. Maybe the machine precision already taken over?

[FMDenaro]

Originally Posted by TurbJet

Another issue, I can't really get 1E-16 for small grid size. The error hanging around 1E-12 ~ 1E-13. Maybe the machine precision already taken over?

yes, you are at the limit of double precision... It seems ok now the accuracy

[TurbJet]

Originally Posted by FMDenaro

yes, you are at the limit of double precision... It seems ok now the accuracy

Great! And thanks very much for all your help. Really appreciate it!