
[Sponsors] 
November 3, 2012, 18:13 
Do we really have a fourth order scheme?

#1 
Member
ville vuorinen
Join Date: Mar 2009
Posts: 67
Rep Power: 9 
Hi,
I'd like to find out an answer to some doubts regarding the cubic scheme. The aim is to evaluate in 1D fvc::div(phi, T) with phi=constant flux . This should give the gradient in 1D. For T I'm using a Gaussian centered at 0.5 and 0<=x<=1. The Gaussian is narrow enough to fit to the domain. OpenFOAM manual claims that the cubic scheme is fourth order. I'm not getting this decay rate but more likely second order decay with slightly lower error (also dispersion in advection problem is a bit reduced in my experience). In 1D Cartesian grid (uniform) a finite volume scheme can be viewed as a finite difference scheme and thereby the fourth order decay could be obtained here. I have some doubts that the cubic is not really fourth order or then I'm doing something wrong. My settings in fvSchemes are: divSchemes { default Gauss cubic; } gradSchemes { default Gauss cubic; } interpolationSchemes { default cubic; } snGradSchemes { default cubic corrected; } fluxRequired { default no; p ; } Any comments? I know that it is possible to construct a "good" fourth order scheme "easily" on unstructured grids by using face neighbor cell gradients i.e. essentially linear interpolation + correction based on the gradients. Is this what the cubic is aimed at? I'd thereby expect that constructing a fourth order scheme in OF is very possible and I'm wondering why we are not seeing the expected behavior. Best, Ville 

November 4, 2012, 05:20 

#2 
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 2,518
Rep Power: 31 
Hi,
your question has a trick ... FD and FV are not the same, also on uniform grids, if you go to accuracy greater than 2. The difference is in the resolved variable. For example: FD > d /dt (phi) + Div (f) = 0 FV > d /dt Int [V] (phi) dx + Int [SV] n.(f) dS = 0 You can see that in FV you solve for: phi_bar = (1/V) Int [V] (phi) dx which is analytically a second order approximation of phi (for a point centered around V). In conclusion, while in FD if you discretize space derivative to fourth order accuracy then is correct to see a fourth order slope in the error decay, in FV when you use a cubic reconstruction for the flux (fourth order accurate FV) is correct to see a second order slope in the error because you are comparing phi_bar to phi. If you want to see in FV a real fourth order slope you have to average the exact solution. This is because you compute the discretization error, not directly the local truncation error of the scheme Note that some FVbased ENO scheme are built in such a way to be equivalent to a fourth order FD derivative and you see an effective fourth order convergence. This issue is quite misleaded in the literature, but for example the book of Ferziger and Peric said clearly the difference between FD and FV, as well as some infos are in the book of LeVeque. I have also more references on this issue. I hope to be of some help .... 

November 4, 2012, 06:06 

#3 
Super Moderator

The OPs question is correctly posed. He is putting a smooth function and computing its derivative. He is not solving a PDE at this stage.
He is doing this in 1d on uniform grids. He should see the theoretical convergence rates and there is no issue of FD versus FV in this case. Actually it is not easy to construct stable and efficient higher than second order schemes on unstructured grids with FVM. 

November 4, 2012, 06:07 

#4 
Member
ville vuorinen
Join Date: Mar 2009
Posts: 67
Rep Power: 9 
Hi,
that is interesting. So what you say is that a "more right" way of formulating the problem would be something like this: FVmethod does not give us an approximation for dg/dx (g=Gaussian here) but more likely to dg_ave/dx . And to compare correctly means having to smooth the Gaussian a posteriori  assuming that the order is higher than 2. I'm not completely eating this yet but let me try it out Ville 

November 4, 2012, 06:32 

#5  
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 2,518
Rep Power: 31 
Quote:
Yes, you got it if you have some like this [f(i+1/2)  f(i1/2)]/dx (1) in a FV formulation this means that you reconstruct the fluxes in i+1/2 and i1/2 also with high order reconstruction (e.g. a cubic interpolation) but (1) is a discrete approximation of the surface integration that balances d/dt(phi_bar). If want, you can check the convergence error on the exact flux and you will recover a fourth order accuracy. Contact me if you need more referrences. Conversely 

November 4, 2012, 08:53 
Another try

#6 
Member
ville vuorinen
Join Date: Mar 2009
Posts: 67
Rep Power: 9 
It is unbelievable but it indeed seems that it might be important to really
explicitly consider the filtering. I repeated the test and stated that we are not seeking for approximation for dg/dx but instead for volume averaged quantity dg_ave/dx = (1/dx) * \int dg/dx which has the explicit "solution" for cell i i.e. analytically dg_ave/dx = 1/dx * (g(i+1/2)  g(i1/2) ) g=Gaussian function. What we note is that by comparing the finite volume approximation to this analytical expression we get quite a different outcome and we note the slope 4 for the cubic scheme. Although at the finest resolution the decay rate settles down a bit deviating from the 4th order decay. Quite a miracle and I'm still not eating this completely Btw, I'd really appreciate if you could send me some more references on the topic. Here or to my email ville.vuorinen at aalto.fi This really makes me think a bit further on the FVmethod. The volume averaging as such does not filter enough in highReynolds number flow when central differences are applied because the diagonal is missing from the stencil (well, on uniform grids). Best, Ville 

November 4, 2012, 09:44 

#7 
Super Moderator

I stand corrected. The averaged derivative can indeed be shown to be fourth order accurate by doing a truncation error analysis, at least in the interior. I dont know what openfoam would do at boundaries. It may not be fourth order accurate at the boundaries.


November 4, 2012, 10:11 

#8  
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 2,518
Rep Power: 31 
Quote:
I am happy you found the trick FV and FD are actually different things! If you are interested:  G. De Stefano, F.M. Denaro, G. Riccardi High Order Filtering Formulation for Control Volumes Simulation of Turbulent Flows, Int. J. Num. Methods in Fluids, 37, 7, 2001.  P.Iannelli, F.M. Denaro, G. De Stefano,A Deconvolutionbased Fourth Order Finite Volume Method for Incompressible Flows On NonUniform Grids Int. J. Num. Methods in Fluids, 43, 4, 2003, 431462.  F.M. Denaro, What does Finite Volumebased implicit filtering really resolve in largeeddy simulations? J. Comp.Phys., 230, 10, May 2011, 38493883, DOI 10.1016/j.jcp.2011.02.011 

November 4, 2012, 15:49 

#9  
Member
Francesco Capuano
Join Date: May 2010
Posts: 81
Rep Power: 8 
Quote:


November 4, 2012, 16:22 

#10 
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 2,518
Rep Power: 31 
Note that the accuracy order of a FV scheme is defined on the basis of the local truncation error of the flux reconstruction, which differs much from the accuracy order of a FD scheme, based on the local truncation error of the derivatives 

November 4, 2012, 17:40 

#11 
Member
ville vuorinen
Join Date: Mar 2009
Posts: 67
Rep Power: 9 
Lots of reading for tomorrow.
Does anyone have further info on what the cubic scheme in OF is actually is based on? I assume that what OF is doing for cubic in the background is to use the explicitly evaluated (second order) gradients and the flux would become something like 0.5*( T_i + T_i+1) + 0.5*(grad(T_i) dot (Rface  Ri) + grad(T_{i+1}) dot (Rface R_{i+1} ) i.e. uses the gradient/Jacobian to correct the flux on the cell face center. In 1D central (finite) difference framework (2nd order) this strategy would lead eventually to fourth order approximate for the derivative which is seen easily. The source for cubic.H could be something like that (since gaussGrad is used in a mystical way) but the correction is encrypted into the syntax. V 

November 4, 2012, 18:16 

#12  
Member
Francesco Capuano
Join Date: May 2010
Posts: 81
Rep Power: 8 
Quote:
Quote:
http://www.tandfonline.com/doi/abs/1...07790500314947 in which the vairous schemes are catalogued with the same OpenFOAM's logic. I'm investigating OF numerical schemes too, so let's keep each other updated in case we find out something Francesco 

November 5, 2012, 05:53 
Just an update on convergence rate

#13 
Member
ville vuorinen
Join Date: Mar 2009
Posts: 67
Rep Power: 9 
To avoid any confusion I attach a final picture of L2error. The error cleanly
goes down with slope=4 (earlier it didn't, purely a matter of boundary effects which are avoided by making the Gaussian more narrow so that it really fits into the domain). V 

November 5, 2012, 07:31 

#14 
Senior Member
Filippo Maria Denaro
Join Date: Jul 2010
Posts: 2,518
Rep Power: 31 
good to see that OF works!
As Francesco said, the matter of distinguish FV from FD is quite unexplored, I tried in some papers to give inlights since also in the book of Leveque is necessary to "extrapolate" the concept. This matter become relevant both for peaople working with Euler flows and thos working on LES (FV are a way to filter) Regards Filippo 

Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Second order scheme with buoyantSimpleFoam  mighelone  OpenFOAM Running, Solving & CFD  2  September 26, 2012 10:28 
How to compute UDS fourth order gradient  Emma66  FLUENT  1  April 27, 2010 05:44 
How to find order of the scheme  madhunitturi  OpenFOAM Running, Solving & CFD  3  April 12, 2007 05:53 
FVMhigher order scheme  K S Chang  Main CFD Forum  6  December 9, 2003 11:26 
second order FD upwind scheme  Heinz Wilkening  Main CFD Forum  2  November 3, 1998 15:33 