CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > Main CFD Forum

Slightly embarassing questions about Rhie-Chow interpolation

Register Blogs Members List Search Today's Posts Mark Forums Read

Like Tree1Likes
  • 1 Post By sbaffini

Reply
 
LinkBack Thread Tools Display Modes
Old   June 4, 2013, 02:10
Lightbulb Slightly embarassing questions about Rhie-Chow interpolation
  #1
New Member
 
Martin
Join Date: Jun 2013
Posts: 5
Rep Power: 4
mercator is on a distinguished road
Hi community,

as the title says, I am trying to understand Rhie-Chow interpolation, and I am unfortunately quite confused and in need of some advice.

(I have read Rhie Chow interpolation, but I am still not sure how it is actually implemented.)

I realize these are very basic questions, sorry for this...

1. If using SIMPLE, which step is actually modified by RC? Normally I do:
a) Solve convection-diffusion to obtain u^*.
b) Calculate divergence of u^*, and solve pressure Poisson eqn.
c) Calculate u^n+1 by adding pressure gradients.
My understanding is that RC *only* affects (b). Is that correct or not?

2. To derive the actual matrix used in (b), I use the same approach as for
SIMPLE (assume rho = delta_t = 1 here):
- start with u^n+1 = u^* - nabla P
- use continuity eqn. on both sides: div u^n+1 = div u^* - div nabla P
- require div u^n+1 = 0 => div nabla P = div u^*.
Now, in the book of Versteeg, formula (11.83) details the face velocity
for RC interpolation as (e.g.)
u^e = (u_P + u_E) / 2 + 1/2 (1/a_P + 1/a_E) (p_P - p_E)
- 1/4 1/a_P (p_W - p_E) - 1/4 1/a_E (p_P - p_EE)
So my understanding, again, is that we actually use this to discretize
the div u^* term on the rhs of the pressure equation. Is that so?

3. When doing so, the "rhs" suddenly contains many new pressure coefficients
in the linear system. I have now two choices: Use the old pressure values
(p^* in the SIMPLE lingo), which is essentially deferred correction, or
use the new values (p' = p as above). When using deferred correction,
I have the problem that the resulting velocity field after the first iteration
won't actually be divergence-free. When coupling p into the system, I
have a wide stencil which I don't really like when using multigrid methods.
So how is this done in practice?

4. And finally, to complete my confusion, I have tried to understand how
OpenFOAM does it from
http://www.tfd.chalmers.se/~hani/kur...7/rhiechow.pdf
They claim that OpenFOAM does not "actually" do RC (which I get since
no wide p stencil is used anywhere). But the paper is quite technical,
and the author never actually points out *where* OpenFOAM differs
from RC. It appears to me that they temporarily use face-based velocities
and then interpolate back, which would make it actually a MAC approach.
But is that true?

Again sorry for the lengthy post. Any insight from the experts would be greatly appreciated!

Cheers
mercator is offline   Reply With Quote

Old   June 4, 2013, 02:21
Default
  #2
Senior Member
 
sbaffini's Avatar
 
Paolo Lampitella
Join Date: Mar 2009
Location: Italy
Posts: 519
Blog Entries: 14
Rep Power: 17
sbaffini will become famous soon enough
This is how the momentum interpolation is done in OpenFOAM:

http://dx.doi.org/10.1016/j.compfluid.2011.11.003
RodriguezFatz likes this.
sbaffini is offline   Reply With Quote

Old   June 4, 2013, 03:44
Default
  #3
Senior Member
 
Join Date: Dec 2011
Location: Madrid, Spain
Posts: 133
Rep Power: 6
michujo is on a distinguished road
Hi, I'm facing a similar problem when trying to implement Rhie-Chow in my little code here. I also don't like at all the wide stencil, especially when you approach the walls and you have to code variations of that expression depending on the cell type.

As I understand it, the interpolated face velocity will be used to calculate the convective fluxes across the east, west, north and south cell faces (Fe,Fw,Fn and Fs following notation from Versteeg&Malalasekera).

The convective fluxes will be used to compute the discretization coefficients aP, aE, ... based on your differentiation scheme (upwind, QUICK, ...), both for solving the convection-diffusion equation for u* and v*, and the corresponding aP_u and aP_v. I do not think you have to include the whole expression for the calculation of the pressure correction, you just use the old values of pressure in the calculation of ue, uw, un, us. Anyway you will repeat the calculation in the SIMPLE loop so you will eventually converge to the right values. What do you think?

Cheers,
Michujo.
michujo is offline   Reply With Quote

Old   June 4, 2013, 07:55
Default
  #4
New Member
 
Martin
Join Date: Jun 2013
Posts: 5
Rep Power: 4
mercator is on a distinguished road
Quote:
Originally Posted by sbaffini View Post
This is how the momentum interpolation is done in OpenFOAM:

http://dx.doi.org/10.1016/j.compfluid.2011.11.003
Thanks for the link. Unfortunately our uni does not have access to this journal, but I will try to use the library service. Amyway, this could take a few weeks.

Any insights on the other questions?
mercator is offline   Reply With Quote

Old   June 4, 2013, 07:59
Default
  #5
New Member
 
Martin
Join Date: Jun 2013
Posts: 5
Rep Power: 4
mercator is on a distinguished road
Quote:
Originally Posted by michujo View Post
Hi, I'm facing a similar problem when trying to implement Rhie-Chow in my little code here. I also don't like at all the wide stencil, especially when you approach the walls and you have to code variations of that expression depending on the cell type.

As I understand it, the interpolated face velocity will be used to calculate the convective fluxes across the east, west, north and south cell faces (Fe,Fw,Fn and Fs following notation from Versteeg&Malalasekera).

The convective fluxes will be used to compute the discretization coefficients aP, aE, ... based on your differentiation scheme (upwind, QUICK, ...), both for solving the convection-diffusion equation for u* and v*, and the corresponding aP_u and aP_v. I do not think you have to include the whole expression for the calculation of the pressure correction, you just use the old values of pressure in the calculation of ue, uw, un, us. Anyway you will repeat the calculation in the SIMPLE loop so you will eventually converge to the right values. What do you think?

Cheers,
Michujo.
Hi Michujo,

you're right, for ultimately convergence it shouldn't matter whether we use old or new values. Maybe it affects the speed of convergence though..

You could also be right that convective fluxes could be affected, but it seems to me that this is not necessary (?). After all, RC is needed to avoid checkerboard problems when projecting, so I guess one could get away with only doing it during projection. I am not sure though.
mercator is offline   Reply With Quote

Old   June 4, 2013, 08:56
Default
  #6
Senior Member
 
Join Date: Dec 2011
Location: Madrid, Spain
Posts: 133
Rep Power: 6
michujo is on a distinguished road
Hi. Yes I agree.

However I am not sure how to calculate the convective fluxes for the calculation of the discretization coefficients in the solution of u* and v*. If I do not need to use Rhie-Chow interpolation here, how shall I calculate the face velocities then? Do you think a linear interpolation between adjacent cells will suffice?

So far I have been trying to apply Rhie-Chow in the calculation of the convective fluxes and subsequent solution of u* and v*, but I could be completely wrong. I'll keep reading about it.

Cheers,
Michujo.
michujo is offline   Reply With Quote

Old   June 4, 2013, 09:40
Default
  #7
New Member
 
Martin
Join Date: Jun 2013
Posts: 5
Rep Power: 4
mercator is on a distinguished road
Quote:
Originally Posted by michujo View Post
Hi. Yes I agree.

However I am not sure how to calculate the convective fluxes for the calculation of the discretization coefficients in the solution of u* and v*. If I do not need to use Rhie-Chow interpolation here, how shall I calculate the face velocities then? Do you think a linear interpolation between adjacent cells will suffice?
That's what I do right now. Apart from the checkerboard problem (which I think has to do with the projection) that should work fine.
mercator is offline   Reply With Quote

Old   June 5, 2013, 02:33
Default
  #8
Senior Member
 
RodriguezFatz's Avatar
 
Philipp
Join Date: Jun 2011
Location: Germany
Posts: 1,097
Rep Power: 16
RodriguezFatz will become famous soon enough
Quote:
Originally Posted by sbaffini View Post
This is how the momentum interpolation is done in OpenFOAM:

http://dx.doi.org/10.1016/j.compfluid.2011.11.003
Thank you for the link, I just read it. I wish all papers were as good as this one. If you know the authors tell them they are great!
__________________
The skeleton ran out of shampoo in the shower.
RodriguezFatz is offline   Reply With Quote

Old   June 5, 2013, 06:49
Default
  #9
Senior Member
 
Join Date: Dec 2011
Location: Madrid, Spain
Posts: 133
Rep Power: 6
michujo is on a distinguished road
Hi, I finally got my little code to work.

1) Lid-driven cavity flow, 2D, incompressible, Reynolds=20 based on the top wall velocity.

2) Collocated grid, 15x15, uniform.

3) Upwind discretization scheme.

I tried two different methods to interpolate the cell face velocity: linear interpolation between adjacent cells and Rhie&Chow.
I attach two pictures showing the pressure fields. It can be seen that the pressure field obtained using Rhie&Chow is much smoother than the one obtained using linear interpolation.

Coming back to the original question:

1) I applied Rhie&Chow to calculate the cell face velocities and therefore the convective fluxes, which I used to solve the discretized u* and v* equations.

2) I did not explicitly introduce the expression for the cell face velocities in the pressure correction equation, which would have led to a wider stencil. The additional pressure terms appear implicitly in the mass imbalance term through the interpolated cell face velocities. I guess that a explicit substitution of the cell face velocities expressions in the pressure correction equation would probably improve convergence though, but so far this worked fine.

I also noticed that Rhie&Chow requires less than 1/3 of the SIMPLE iterations needed using linear interpolation for convergence. Anyone noticed something similar?

Cheers,
Michujo.
michujo is offline   Reply With Quote

Old   June 5, 2013, 07:12
Default
  #10
Senior Member
 
RodriguezFatz's Avatar
 
Philipp
Join Date: Jun 2011
Location: Germany
Posts: 1,097
Rep Power: 16
RodriguezFatz will become famous soon enough
What criterion did take for "convergence"? I am asking because the absolute value of residuals can change when you change the equations. Thus, the residual can appear lower without having a better solution.
__________________
The skeleton ran out of shampoo in the shower.
RodriguezFatz is offline   Reply With Quote

Old   June 5, 2013, 07:21
Default
  #11
New Member
 
Martin
Join Date: Jun 2013
Posts: 5
Rep Power: 4
mercator is on a distinguished road
Quote:
Originally Posted by michujo View Post
Hi, I finally got my little code to work.

/snip

Cheers,
Michujo.
Very interesting! I think I will try something similar then. Maybe I was confused because the textbook claimed that the continuity eqn was modified, which should not affect convection (pure momentum eqn). Otoh maybe you could do both as you indicated, at the expense of having wider stencils for the pressure eqn.

Just to be sure, so I understand: in the incompressible case the convective flux at the east face is simply F_e = u_e. And that is the velocity you replace by the RC term, using the guess pressure p^* as a corrector. Right?

Cheers
mercator is offline   Reply With Quote

Old   June 5, 2013, 08:52
Default
  #12
Senior Member
 
Join Date: Dec 2011
Location: Madrid, Spain
Posts: 133
Rep Power: 6
michujo is on a distinguished road
Hi.

1) For convergence I just track the residuals calculated as the difference between the old and new values of u, v and p within the SIMPLE loop. I take the maximum of those three and accept convergence when this is lower than my tolerance (I took 1e-5). Do you think I could implement a better procedure?

2) Hi, yes I agree with you. However I tried to avoid the wider stencil in the pressure correction equation. I'm assuming that the influence of the pressure field is implicitly accounted for in the mass imbalance term calculated with Rhie&Chow.

Yes, I first calculate the cell face velocities with Rhie&Chow (u_east, u_west, v_north and v_south) and with them the convective fluxes as you wrote. Yes, I use p* as the corrector term in the R&C expression.

Maybe it is better to use R&C expression explicitly in the pressure correction equation as well, but so far it seems to be working...

Thanks for your comments.

Cheers,
Michujo.
michujo is offline   Reply With Quote

Old   June 10, 2013, 09:03
Default
  #13
Senior Member
 
RodriguezFatz's Avatar
 
Philipp
Join Date: Jun 2011
Location: Germany
Posts: 1,097
Rep Power: 16
RodriguezFatz will become famous soon enough
Dear all,
I also busied myself with these implementation details. Unfortunately all explanations

I found (here and in books) weren't satisfying and I hope someone here can shed light on this:

1) In "Ferziger / Peric" Computational Methods for Fluid Dynamics (3rd Edition), the authors explain how modification of the Poisson-equation for the pressure prevents pressure oszillations in colocated arrangements. (page 196, et sqq.). As I understand it, the "outer" derivative (which stems from continity eq.) of the Poisson equation is evaluated using face values of pressure derivatives of the current cell. Thus, only "dx" and "dy" derivatives occour, instead of "2dx", "2dy". Now they argue, that they "introduced an inconsistency in the treatment of the pressure gradient in the momentum and pressure equations" (page 200). I don't get it, can someone tell me, what the inconsinstency is about?

Philipp.
__________________
The skeleton ran out of shampoo in the shower.
RodriguezFatz is offline   Reply With Quote

Old   June 10, 2013, 11:51
Default
  #14
Senior Member
 
Join Date: Dec 2011
Location: Madrid, Spain
Posts: 133
Rep Power: 6
michujo is on a distinguished road
Hi, I think the inconsistency he is talking about is the fourth order pressure term that is added to the momentum equations when interpolating values at the cell faces (where variables are not stored in a collocated arrangement), which is essentially the Rhie&Chow interpolation. Being this extra term fourth order, its effect on the accuracy is negligible for second or third order discretization schemes, yet it damps out pressure oscillations with a 2*deltax pattern (checkerboarding).

Anyone else has a better explanation? Thanks.

Cheers,
Michujo.
michujo is offline   Reply With Quote

Old   June 11, 2013, 03:07
Default
  #15
Senior Member
 
RodriguezFatz's Avatar
 
Philipp
Join Date: Jun 2011
Location: Germany
Posts: 1,097
Rep Power: 16
RodriguezFatz will become famous soon enough
Thank you for your answer, michujo! I don't see how this can be called "inconsistent". For me, "inconsistent" in this context would mean to treat the same term differently in the two equations. However, I think this isn't done here, right?
__________________
The skeleton ran out of shampoo in the shower.
RodriguezFatz is offline   Reply With Quote

Old   June 11, 2013, 04:24
Default
  #16
Senior Member
 
Join Date: Dec 2011
Location: Madrid, Spain
Posts: 133
Rep Power: 6
michujo is on a distinguished road
Hi, I just read Ferziger again, more carefully this time. As I understood it, the inconsistency consists in using a central difference approximation using points 2*deltax apart in Eq. 7.120 and, on the other hand, using a central difference approximation using points separated deltax in Eq.7.126 (and then interpolating the values of the fluxes at the cell faces). The difference between both approaches is the fourth order term.

I personally find it hard understanding Ferziger's book compared to other textbooks, it's so dense...

Cheers,
Michujo.
michujo is offline   Reply With Quote

Old   June 11, 2013, 04:35
Default
  #17
Senior Member
 
RodriguezFatz's Avatar
 
Philipp
Join Date: Jun 2011
Location: Germany
Posts: 1,097
Rep Power: 16
RodriguezFatz will become famous soon enough
Yes, that this is what he states. But these are just different approximations of the same equation, of course there is a difference. I am just wondering why he calls that "inconsistent".

Additionally a second thing I am thinking about:
After deriving the "dx"-based pressure equation, he states that
"Use of this eliminates the oscillation in the pressure field".
If this is already sufficient, why all the "Rhie Chow" stuff for the face velocities? Is this (Rhie Chow) made for a second kind of oscillations?
__________________
The skeleton ran out of shampoo in the shower.

Last edited by RodriguezFatz; June 13, 2013 at 01:57.
RodriguezFatz is offline   Reply With Quote

Reply

Tags
rhie-chow simple

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
rhie chow interpolation gerardosrez Main CFD Forum 9 March 12, 2014 17:35
Rhie Chow interpolation Rachit Main CFD Forum 13 March 7, 2014 02:52
Rhie - Chow Interpolation ? Javier Goicochea Main CFD Forum 0 November 27, 2003 19:01
Rhie - Chow speed up factors Goicox Main CFD Forum 0 November 22, 2003 01:00
Rhie & Chow interpolation Barry Main CFD Forum 3 January 20, 2003 10:11


All times are GMT -4. The time now is 19:20.