[maka]

To reproduce the behavior:
openFOAM version: 1.3.
channel395 tutorial.

If we change fvSchemes dictionary to the description at the end of this message to achieve 4th order accuracy and rerun channel395 tutorial, will finish the simulation with max Courant number of: 0.99 while using 2nd order schemes it finishes with 0.35.

Q1) Is there something wrong in the selection of schemes made below?

Q2) Is the three fold increase in Courant number a bug or such increase is the normal price in terms of stability that we pay for higher order schemes?

Q3) Henry before in one thread mentioned that 2nd order schemes should be fine "particularly when the mesh is non-orthogonal". Can any body highlight the relation between the choice of scheme and the non-orthogonal correction in OpenFOAM.

------------------------------------
fvSchemes:

ddtSchemes
{
default backward;


gradSchemes
{
default Gauss cubic; //linear;
grad(p) Gauss cubic; //linear;
grad(U) Gauss cubic; //linear;

divSchemes
{
default none;
div(phi,U) Gauss cubic; //linear;
div(phi,k) Gauss limitedCubic 1; //limitedLinear 1;
div(phi,B) Gauss limitedCubic 1; limitedLinear 1;
div(B) Gauss cubic; //linear;
div(phi,nuTilda) Gauss limitedCubic 1; //limitedLinear 1;
div((nuEff*dev(grad(U).T()))) Gauss cubic; //linear;
}

laplacianSchemes
{
default none;
laplacian(nuEff,U) Gauss cubic corrected; //linear corrected;
laplacian(1|A(U),p) Gauss cubic corrected; //linear corrected;
laplacian(DkEff,k) Gauss cubic corrected; //linear corrected;
laplacian(DBEff,B) Gauss cubic corrected; //linear corrected;
laplacian(DnuTildaEff,nuTilda) Gauss cubic corrected; //linear corrected;
}

interpolationSchemes
{
default cubic; //linear;
interpolate(U) cubic; //linear;
}

snGradSchemes
{
default corrected;
}

fluxRequired
{
default no;
p;
}


}}

Best regards,
Maka.

[henry]

I would guess that the increase in Courant number is a consequence of instabilities in the velocity field. Have you looked at the results in detail? Are there unphysical spikes in the velocity when running cubic differencing?

Note that running such scheme does not provide overall 4th-order accuracy because the Gauss integrals are 2nd-order accurate.

H

[maka]

(1) Yes, there are localized increase of velocity near the core of the channel but due to diffusion it is hard to say if they are physical or not because of the turbulent nature of the flow. I have seen the same behavior in other problems using LES and tried to reproduce it with the the channel tutorial. I agree that in some case the spikes are very clear and are unphysical. If cubic scheme is unbounded then that explains such behavior.

(2) regarding the overall 4th-order accuracy (actually that was my objective): I used to think that when we replace volume integral with surface integral in FVM using Gauss theorem, the resulting equation is exact till we approximate the surface terms. I used to think that Gauss that we write in front of the scheme (which is the only option for divergence and Laplacian schemes) refers to such process of replacing the volume integrals but surface integral. I would be grateful if you explain where such limitation of overall 2nd order accuracy in Gauss integrals come from? Many thanks Henry for your help.

Best regards,
Maka.

[henry]

(1) All schemes other than upwind are potentially unbounded and the higher the order the more likely the unboundedness and the larger number of computational modes supported.

(2) Gauss theorem allows us to convert volume integrals into surface integrals which if evaluated exactly would allow us to evaluate the volume integrals exactly. However, currently we use interpolation to obtain values at the face-centres (on regular meshes) which are the Gauss-points for evaluating the surface integrals to second-order simply by multiplying the face-areas. Hence irrespective of the order of the interpolation the order of the volume-integrals is second on a regular mesh. Nevertheless it is still often beneficial to use higher-order schemes for interpolation, particularly for convection terms.

Note also that irregular meshes the above procedure is not usually even second-order because the simple interpolation does not provide values at the face-centres but this can be corrected for and we are developing a new set of polynomial-fit interpolation schemes which are formally second and higher order on arbitrary meshes.

H

[maka]

Many thanks for your explanations Henry.

Best regards,
Maka.

[ivan_cozza]

Prof. Weller,
one of my objective is to implement a linearized euler solver for acoustic propagation in OpenFOAM.
Now, as you said, it's not possible to achieve an integration order higher than 2 because of the Gauss integration. For my purposes, I need 4th order non dispersive schemes, so, with your polynomial fitting it's possible to achieve it in a FVM framework?
Thanks, Ivan

[henry]

It is possible to implement the FVM higher than second-order in general and in OpenFOAM by interpolating to more than one Gauss-point on the cell faces in order to obtain a higher-order surface integral. However, for our current purposes this has not proved necessary and we can achieve the accuracy we need with high-order interpolation to the cell face centres only. Do you think that this might be suitable for acoustic propagation or do you think that formally higher-order integration is also necessary?

Another option would be to use the polynomial fitting within a finite-difference framework to provide a formally higher-order scheme if local conservation is less important than formal accuracy.

H

[ivan_cozza]

Dear Prof Weller,
I have allready tried to implement an aeroacoustic solver (linearized euler) with the allready existing schemes (for istance, I tried linear and cubic), but I can't get good results, because linear is too dissipative and cubic is unstable and not so good with the dispersion error, so I gave up and wait for new evolutions.
The state of art in aeroacoustic is DRP (Dispersion Relation Preserving) 4th order finite difference schemes with filters (huge stencils!) or Discontinuous Galerkin schemes, there are people (Davidson at Chalmer, as I know) that work on finite difference - like FVM schemes, but they need structured cells.
So, conservation is not as fundamental as for istance in LES.
Can you give me some reference on polynomial fitting schemes?
Do you manage to put it in some new release of OF?
Thank you

[henry]

Do you need high-order centred or upwinded schemes? My guess is that cubic is unstable because it is a centred scheme and some kind of high frequency filtering would have to be applied.

We are developing both high-order centred and upwinded schemes and these will be included in the next release of OpenFOAM.

Some information on these schemes is contained in:

"Voronoi, Delaunay and Block Structured Mesh Refinement for Solution of the Shallow Water
Equations on the Sphere"
Hilary Weller, Henry Weller and Aime' Fournier
Submitted to Monthly Weather Review, 2008

http://www.met.reading.ac.uk/~sws02h...rnier_2008.pdf

and more details can be found in the references.

H

[ivan_cozza]

Someone uses centered schemes with selective filters, someone uses specific low dissipation upwind schemes. Now, me and my research group have a Fortran code based on Discontinuous Galerkin that is 4th order with an upwind treatment of the interfaces fluxes (Lax-Friedrichs fluxes), and with the acoustics test cases it works quite well.
It would be great to have High-order schemes in the next release, I hope it will be out as soon as possible!