# Slightly embarassing questions about Rhie-Chow interpolation

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

 June 4, 2013, 02:10 Slightly embarassing questions about Rhie-Chow interpolation #1 New Member   Martin Join Date: Jun 2013 Posts: 5 Rep Power: 5 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

 June 4, 2013, 02:21 #2 Senior Member     Paolo Lampitella Join Date: Mar 2009 Location: Italy Posts: 613 Blog Entries: 17 Rep Power: 20 This is how the momentum interpolation is done in OpenFOAM: http://dx.doi.org/10.1016/j.compfluid.2011.11.003 RodriguezFatz likes this.

 June 4, 2013, 03:44 #3 Senior Member   Join Date: Dec 2011 Location: Madrid, Spain Posts: 134 Rep Power: 7 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.

June 4, 2013, 07:55
#4
New Member

Martin
Join Date: Jun 2013
Posts: 5
Rep Power: 5
Quote:
 Originally Posted by sbaffini 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?

June 4, 2013, 07:59
#5
New Member

Martin
Join Date: Jun 2013
Posts: 5
Rep Power: 5
Quote:
 Originally Posted by michujo 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.

 June 4, 2013, 08:56 #6 Senior Member   Join Date: Dec 2011 Location: Madrid, Spain Posts: 134 Rep Power: 7 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.

June 4, 2013, 09:40
#7
New Member

Martin
Join Date: Jun 2013
Posts: 5
Rep Power: 5
Quote:
 Originally Posted by michujo 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.

June 5, 2013, 02:33
#8
Senior Member

Philipp
Join Date: Jun 2011
Location: Germany
Posts: 1,297
Rep Power: 19
Quote:
 Originally Posted by sbaffini 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.

June 5, 2013, 06:49
#9
Senior Member

Join Date: Dec 2011
Posts: 134
Rep Power: 7
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.
Attached Images
 Pfield_cavityflow_Re20_LinearMomInterp.png (4.7 KB, 22 views) Pfield_cavityflow_Re20_RhieChowMomInterp.png (4.6 KB, 24 views) Velfield_cavityflow_Re20_RhieChowInterp.png (5.3 KB, 20 views)

 June 5, 2013, 07:12 #10 Senior Member     Philipp Join Date: Jun 2011 Location: Germany Posts: 1,297 Rep Power: 19 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.

June 5, 2013, 07:21
#11
New Member

Martin
Join Date: Jun 2013
Posts: 5
Rep Power: 5
Quote:
 Originally Posted by michujo 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

 June 5, 2013, 08:52 #12 Senior Member   Join Date: Dec 2011 Location: Madrid, Spain Posts: 134 Rep Power: 7 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.

 June 10, 2013, 09:03 #13 Senior Member     Philipp Join Date: Jun 2011 Location: Germany Posts: 1,297 Rep Power: 19 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.

 June 10, 2013, 11:51 #14 Senior Member   Join Date: Dec 2011 Location: Madrid, Spain Posts: 134 Rep Power: 7 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.

 June 11, 2013, 03:07 #15 Senior Member     Philipp Join Date: Jun 2011 Location: Germany Posts: 1,297 Rep Power: 19 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.

 June 11, 2013, 04:24 #16 Senior Member   Join Date: Dec 2011 Location: Madrid, Spain Posts: 134 Rep Power: 7 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.

 June 11, 2013, 04:35 #17 Senior Member     Philipp Join Date: Jun 2011 Location: Germany Posts: 1,297 Rep Power: 19 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.

 Tags rhie-chow simple

 Thread Tools Display Modes Linear Mode

 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 OffTrackbacks are On Pingbacks are On Refbacks are On Forum Rules

 Similar Threads Thread Thread Starter Forum Replies Last Post gerardosrez Main CFD Forum 9 March 12, 2014 17:35 Rachit Main CFD Forum 13 March 7, 2014 02:52 Javier Goicochea Main CFD Forum 0 November 27, 2003 19:01 Goicox Main CFD Forum 0 November 22, 2003 01:00 Barry Main CFD Forum 3 January 20, 2003 10:11

All times are GMT -4. The time now is 00:08.