Hi,

I'm using a second order colocated FV method w. either triangular or quadrilateral control volumes to solve flow over a laminar backward facing step (Re 100 only). Pressure-velocity coupling is handled via the Rhie-Chow interpolation method.

Using meshes of similar resolution, I am getting MUCH more accurate results with the quad meshes, and I cannot seem to figure out why. Has anyone ever seen/heard of this before?

Furthermore, while mesh anisotropy improves solution accuracy in the case of quad meshes, it actually adversely affects the numerical solution with tris.

I am using two measures to evaluate solution accuracy: a comparison of solutions obtained using moderate resolution (about 20 CVs in the transverse direction) to a mesh independent (MI) solution, and a comparison of the flow re-attachement points to the MI solution.

Any ideas and/or references would be appreciated.

Thanks in advance, Steve

What kind of finite volume method are you using (cell centered, or vertex centered)? What are your approximation functions (constant, linear...)

When you answer my questions, I might be able to give you a complete answer.

Hello,

Thanks for your reply.

I am using a cell centered method, with linearly exact (second order) approximations for everything. Linear corrections are deferred, and expressed in terms of cell gradients that are reconstructed using the least-squares approach by T. Barth. Viscous/diffusion terms are corrected for mesh non-orthogonality, but the correction is deferred.

Thanks in advance for any help you can offer. Please let me know if there is any additional information that may help your response.

Cheers, Steve

Ok, Steve. First, It is manditory (As you do), to correct the diffusion flux for a cell centered FVM when the mesh is not orthogonal (Sorry but I don't know the Barth, approach). But, your results suggest a wrong approximation of those diffusion flux. Have you tried running the same case with an triangular othogonal mesh? Or tried solving a simple diffusion problem with a known analytical solution on an isotropic triangular mesh?

Second of all, Are you using some king of upwinding or artificial viscosity to avoid unwanted oscillations (I guess yes)?

Good luck.

Hi Steve,

That's intresting. I observed the same. I was using SIMPLE on Unstructured grids (Hybrid Scheme). The results on quad mesh were much better even if there was some skewness. The triangle mesh gave results of comparable accuracy only when the mesh was refined. I guess the diffusion term is the problem on triangular mesh. I evaluate mu*du_i/dx_j at face of triangle by greens theorem by constructing a auxilary control voulume having edge as one daigonal & line joining centers of left & right cells as other daigonal. Do u do in in same way? If u do it this way u can split ther diff terms in to orthogonal & non-ortho. parts. I printed values of non-ortho terms they are significant. so reduce accuracy. Let me know your experience.

Thanks

Abhijit Tilak

Hello Sebastian,

The stress/diffusion terms are roughly implemented as an active node-to-node derivative + the correct lagged derivative - a lagged node-to-node derivative. An alternate implementation, the use of auxiliary nodes that are orthogonal to the face has not been used. For clarification, T. Barth's work is only used to reconstruct the nodal solution derivatives used in the lagged part of the approximation. The implemented stress/diffusion term comes from Ferziger and Peric's text.

I have run the same problem with 2:1 aspect ratio (in SW direction) quads, and isotropic and anisotropic triangular meshes. All meshes have approximately the same number of mesh points in the transverse direction. In the order of best to worste results (by all measures) are quads, isotropic tris, and anisotropic tris.

For what its worth, the anisotropic tris are stretched in the streamwise direction where d2u/dy2 is large (jets above and below the shear layer), and are stretched in the transverse direction through the shear layer (where d2u/dy2 is small). The greatest solution errors, and the greatest differences between the quad and tri generated results occur in the shear layer - this is surprising because the linearly exact discrete approximations should very accurately resolve this region of the flow where the solution varies nearly linearly.

Thanks again in advance for any thoughts. Steve

Hello Abhijit,

I feel bad to say, but I am happy to hear you have a similar problem.

I think we have different discretizations. I treat the du_i/dx_j term similar to how scalar diffusion would be handled. As noted in my reply to Seb, this term is approximated using a node-to-node derivative that assumes orthogonality (for active coefficients), plus a lagged non-orthogonal correction. The correction essentially adds the correct normal derivative, and removes the node-to-node derivative.

I suspect that the magnitude of the non-orth correction is not so important, provided that you can get the system to converge.

However, just today I realized that this discretization does not account for the error introduced if the node-to-node vector (I call 's') does not cross the face at the face's mid point. This can adversely affect the order accuracy of the mid-point-rule quadrature approximation... and the overall force approximation.

An auxiliary node approach will account for this error. In this approach, the true nodes are projected onto the face normal vector to obtain auxiliary nodes. Nodal solution values (u,v,p, etc) are extapolated from the true nodes to the aux. nodes. These aux nodes are both orthogonal to the mesh face AND the resulting aux. 's' vector passes through the face at the integration point.

I will implement this aux point approximation, and first output the differences between the original discretization and aux node discretization. These may be signficant. If they are, I will use the aux node implementation to solve the flow and see what happens.

I will post my observations since you are have a similar problem.

Cheers, Steve

Dear Steve,

As part of my experience, the collocated mesh scheme is not fully conservative. There are actually two terms responsible for spoiling the kinetic energy conservation. A term du to the pressure (since the cartesian velocities at the center do not eunsure mass exactly), and an another term du to the interpolation method. So if you are using a second order centered operateur, the error is present and will significantly degrade the result the results if you have more than one non-uniform(and/or curved) direction. The easiest way of getting rid of the interpolation error is to use a First order centered interpolation (which mean that you still the points on both side of the face, but you tune the weighted average to be mesh independant: U(i+1/2)= 0.5*( U(i) + U(i+1) ) ). If you do this, you will drop the interpolation error. In addition, the scheme is 1st order locally, but still second order globally. I actually wrote a paper about these issues. you can download the paper at the following address: http://utacfdi.uta.edu/~felten/resume/taicdl_art.ps

So just try the simple mesh independant average, and you should see an improvement in the solution. Sincerely,

Frederic Felten.

Hello Frederic,

Thank you for your reply. Using the terminology presented in your paper, I am using a first order accurate (mesh independent) upwind approximation plus a deferred linear correction from the upwind node to the integration point (the mesh dependent part).

I re-ran the test case with an isotropic triangular mesh, using the mean mesh independent *advected* velocity approximation (i.e. the integration point value) sugested in your paper. The mesh dependent part was removed. For this Reynolds number (100), my mesh independent flow re-attachment point is approximately 3.22 step heights downstream of the step. Using the original (second order) implementation, the predicted value was 3.18, and using the modified (first order) implementation, the predicted value became 3.14 (i.e. became worse).

I suspected that this would occur because of the diffusive nature of the first order scheme.

Please let me know if I have missed something, or if I have incorrectly interpreted your ideas/paper.

Thanks in advance, Steve

hi there,

Instead of using an upwind interpolation still use the two points on both side of the face to compute the average but just use 0.5 of one + 0.5 of the other. F(i+1/2)= 0.5*( u(i,j) + u(i+1,j) ). So you are not upwind since you are grabing info on both side of the face, but you are not second order exactly since the weighted average is independant of the actual mesh distances. This procedure is 1st order locally, but second order globally. (Check Veldman and Verstappen, 1992, 1998).

If you look at what we did in the paper, we compared this "1st order Centered" interpolation to the "2nd order Centered" and the "Second order Upwind". The result was, that you don't need upwinding. I hope this helps, Sincerely,

Frederic.

Dear all,

I am trying to write the numerical solution for this step problem using C language.For that i am using HYBRID method and SIMPLER algorithem using finite volume method. But i could not succeed because i cannot go beyond 100 simultenious equations as it is the max size for a two dimentional array in C.

If anybody know how to do it please mail me at

rambabu_sirisetty@yahoo.com

This sounds like your C compiler is only using 16 bits for array indices. If you are using an old DOS C compiler it will probably have a switch to use more bits. A better alternative would proably be to install and use gcc:

http://gcc.gnu.org/

Hello,

Perhaps other programming languages may be more appropriate (C++ or Fortran). The limit of 100x100 may be a number that you can change in the compiler that you are using - this seems quite limiting, and will provide only a very inaccurate solution for this problem; I have realized a mesh independent solution for this problem with about 200x80 (i.e. 16000) CVs, and thus 48000 equations.

Another alternative is to take advantage of sparse matrix storage techniques - The discrete equation sets are quite sparse, and the resulting coefficients are easily stored in 3 vectors (rather than one array), for which you may not have the noted size constraints. One vector contains the coefficients, one vector contains the column information, and one vector contains row starting indices - a search on sparse matrix storage techniques will provide additional information. Perhaps someone else has a specific reference?

Good luck to you, Steve R