[odellar]

Hi,

For a wall bounded incompressible flow, with no slip condition velocity boundary conditions at the walls, what is the right pressure boundary condition to apply, and why?

I've seen people use zeroGradient implying that at, for example a 'bottom' wall, dp/dy=0. This is what I typically use, but I'm not completely sure why.

I've seen a derivation of why one uses this BC based on the Navier Stokes momentum equations posted in this forum (but cannot find a link now), but don't agree with it - I think there's a fundamental flaw in it. So I've not seen a good reason for it.

This BC is even used for pressure in the lid driven cavity flow in the OpenFOAM tutorials, but they don't specify why.

Any help would be much appreciated.

[Blanco]

We have to consider that the pressure gradient is already present in the N-S equations, therefore using a neumann condition for pressure is quite "natural" for every boundary where pressure is not imposed by something external to the computational domain, i.e. everywhere you don't have a fixed-pressure boundary.
Neumann condition is natural to be imposed on a wall boundary, just think about the physics: a fixed wall doesn't exert any force on a fluid in its normal direction, apart the reaction to the force exerted by the fluid to the wall in the same direction. Therefore using dp/dn = 0, where n is the direction normal to the wall, is a straightforward way to get wall reaction to the fluid pressure applied to the fluid.

If we consider, as an example, two infinite flat plate placed one over the other at a distance h(x)=const., and a fluid flow occuring between them, and we consider this flow to be steady-state, incompressible, no gravity, viscid, single phase, if the two plates are fixed in space we can assume that the flow will be planar (2D to simplify) and that the velocity vector will be a function of y axis only, because we have nothing varying along the x-axis to justify a velocity change in that direction and we have to keep volumetric flow rate constant, so u=u(y), v=v(y) (where x axis and u component are aligned with the plates, y axis and v component are perpendicular to those plates).
Anyway, incompressible flow means

du/dx+dv/dy=0

but u=u(y) therefore the first term is 0 and we get dv/dy=0 -> v is constant along the y axis. We know that v(0)=v(h)=0 therefore v(y)=0. If we know write the momentum equations for the y direction we have that:

dp(x,y)/dy=0

therefore pressure is constant along the y axis -> p=p(x). As you can see, the neumann conditions arise "naturally" in this case, because you have dp/dy=0 in the flow core and you would impose that conditions also on the wall boundary, in order to have the flow "free" to develop it's internal behavior in the y direction.

Another topic: we have to consider that as long as incompressible flows are concerned, pressure is NOT representative of the real thermodynamic pressure but it acts just as a Lagrange multiplier. So pressure is not to be considered representative of the "real fluid pressure" but it's just a scalar quantity used to close the problem. The velocity field is the physical variable in this case.

Hope this helps,
Best

[murx]

Originally Posted by Blanco

If we consider, as an example, two infinite flat plate placed one over the other at a distance h(x)=const., and a fluid flow occuring between them, and we consider this flow to be steady-state, incompressible, no gravity, viscid, single phase, if the two plates are fixed in space we can assume that the flow will be planar (2D to simplify) and that the velocity vector will be a function of y axis only, because we have nothing varying along the x-axis to justify a velocity change in that direction and we have to keep volumetric flow rate constant, so u=u(y), v=v(y) (where x axis and u component are aligned with the plates, y axis and v component are perpendicular to those plates).
Anyway, incompressible flow means

du/dx+dv/dy=0

but u=u(y) therefore the first term is 0 and we get dv/dy=0 -> v is constant along the y axis. We know that v(0)=v(h)=0 therefore v(y)=0. If we know write the momentum equations for the y direction we have that:

dp(x,y)/dy=0

therefore pressure is constant along the y axis -> p=p(x). As you can see, the neumann conditions arise "naturally" in this case, because you have dp/dy=0 in the flow core and you would impose that conditions also on the wall boundary, in order to have the flow "free" to develop it's internal behavior in the y direction.

Hey Blanco, thanks for that thorough explanation!
When I simplify the viscid momentum equation I get:

0= -dp/dx + ny (d2u/dy2) (x-direction) and
0= -dp/dy + ny (d2v/dy2) (y-direction).

I see that for your example (with v=0) dp/dy is clearly zero. However, in most flows, we have more than one velocity component and still seem to use this boundary condition - even for low Reynolds numbers, where the viscid terms cannot be neglected. Does anybody have an explanation for that?

I am not using OpenFOAM. Did anybody of you ever check if the pressure boundary condition is actually applied at the wall? Since the pressure is calculated in the pressure-velocity-coupling, we do not solve for a differential pressure equation and so there is not really the need for a boundary condition, or am I wrong?

[francisco.angel]

Hi murx, I think in your case, the reasoning is pretty similar to the one presented by Blanco, this is, based on the continuity equation. At the wall (x=0, horizontal wall) we have:

du/dx = 0 (because u is zero at every point on the wall). So from the continuity equation, we deduce:

du/dx+dv/dy=dv/dy=0

From here, we deduce:

d2v/dy2 = d ( dv/dy ) / dy = d(0) / dy = 0

And from the y-momentum equation presented by you we obtain the zero-gradient at the wall boundary condition for pressure.

0=-dp/dy+nuy*d2vdy2=-dp/dy=0

For a vertical wall a similar reasoning can be applied. You can also check this document:
https://www.dropbox.com/s/k292ol82ya...ation.pdf?dl=0

The document wasn't written by me, I can't credit the author now, because I had it in my computer and honestly I don't remember how I obtained it.
Also, I remember to have seen developments similar to this one, regarding boundary conditions on walls in material dedicated to boundary layer, so I think you can find more material if you look that topic.
Regarding OpenFoam, in the cavity case, these are the boundary conditions for pressure field:

boundaryField
{
    movingWall
    {
        type            zeroGradient;
    }

    fixedWalls
    {
        type            zeroGradient;
    }

    frontAndBack
    {
        type            empty;
    }
}

So, zero gradient is used in the not empty walls.

Hope this helps.

[murx]

Originally Posted by francisco.angel

du/dx+dv/dy=dv/dy=0

From here, we deduce:

d2v/dy2 = d ( dv/dy ) / dy = d(0) / dy = 0

Hi Francesco,
I've seen that pdf document in another thread on cfd-online. Maybe that's where you got it But I think the deduction in the document and in your post is not correct. Because dv/dy = 0 does not mean that d2v/dy2 = 0. There is no need for the second derivative to be zero just because the first one is.

Originally Posted by francisco.angel

Regarding OpenFoam, in the cavity case, these are the boundary conditions for pressure field:

boundaryField
{
    movingWall
    {
        type            zeroGradient;
    }

    fixedWalls
    {
        type            zeroGradient;
    }

    frontAndBack
    {
        type            empty;
    }
}

So, zero gradient is used in the not empty walls.

Hope this helps.

Thanks. What I meant is to check the results and see if the pressure gradient is really zero. I think I remember a collegue of mine saying that when he uses a higher value for the orthogonal correctors, than this boundary condition is apperently not applied and the gradient is not zero.

[Blanco]

Hi murx,
My example was for a 2d flow (where v(y)=0 because of boundary conditions...) but nothing changes if we go to 3d.
In any case, the viscid term in my example is not present just because we have that v(y) = constant from the continuity equation, therefore the first AND the second derivative are 0. if you prefer another point of view is that continuity equation is valide everywhere in the fluid, therefore dv/dy = 0 hold everywhere, and the second derivative is null as a consequence (you're right, a null first derivative in one point wouldn't be sufficient to say that the second derivative is null too in that point, but it became a sufficient condition if the first derivative is null everywhere ).

Hope this help,
Best

[murx]

Originally Posted by Blanco

Hi murx,
My example was for a 2d flow (where v(y)=0 because of boundary conditions...) but nothing changes if we go to 3d.
In any case, the viscid term in my example is not present just because we have that v(y) = constant from the continuity equation, therefore the first AND the second derivative are 0. if you prefer another point of view is that continuity equation is valide everywhere in the fluid, therefore dv/dy = 0 hold everywhere, and the second derivative is null as a consequence (you're right, a null first derivative in one point wouldn't be sufficient to say that the second derivative is null too in that point, but it became a sufficient condition if the first derivative is null everywhere ).

Hope this help,
Best

Ok, I understand that in your case of a developed flow (du/dx = 0), v must be constant - and therefore always 0 since you have v(y=0)=0. And therefore the deduction of dp/dy = 0 is correct.

However, this is a very special case. Usually we use CFD to examine undeveloped flows, as I can analytically obtain the developed flow field for a lot of geometries. So for those undeveloped flows, I use a pressure boundary condition which can only be derived under the assumption of a developed flow, right?

I'm so interested in that because I am investigation the force that is exerted on an object by a certain flow and I am concerned that an unphysical pressure boundary condition might cause a wrong pressure field and therefore unphysical results.

[Blanco]

Ok, I see your concern, however the boundary conditions are not "derived" for a particular case, but they must be choosen to solve each particular problem.

What I would stress is that we are not imposing dp/dx or whatever because we have seen that it works for a particular flow...we impose just what we already know about a boundary. If we have constant pressure on that boundary, we impose it numerically, if we have wall, we simply impose dp/dn =0 because from thw physics we know that pressure is "similar" to a normal stress applied on a wall in normal direction. We know that the wall will exert a reaction force on the fluid equal to the pressure (per unit area), indeed that's why we have pressure in the momentum equations...because we have already used this condition (which is indees valid even if you don't have a wall but another fluid volume) to obtains the NS equations.

[murx]

Originally Posted by Blanco

Ok, I see your concern, however the boundary conditions are not "derived" for a particular case, but they must be choosen to solve each particular problem.

What I would stress is that we are not imposing dp/dx or whatever because we have seen that it works for a particular flow...we impose just what we already know about a boundary. If we have constant pressure on that boundary, we impose it numerically, if we have wall, we simply impose dp/dn =0 because from thw physics we know that pressure is "similar" to a normal stress applied on a wall in normal direction. We know that the wall will exert a reaction force on the fluid equal to the pressure (per unit area), indeed that's why we have pressure in the momentum equations...because we have already used this condition (which is indees valid even if you don't have a wall but another fluid volume) to obtains the NS equations.

Ok, I see. But lets use a simple example: A container filled with fluid, with the bottom in x-z-plane. No velocities. Here, your argumentation is of course right. The hydrostatic pressure times the area will be the force exerted back onto the fluid by the bottom of the container. But pressure gradient will be dp/dy = rho*g everywhere and dp/dy is physically not right. So I can't see that this actio = reactio argumentation proves that dp/dn must be zero.

[Blanco]

Ok, let's consider your example and consider what we get if we remove the bottom of the container and we substitute it with a virtual fuid layer...I call it virtual because in some way we can assure that the pressure inside this fluid layer is constant and that gravity is not working in this layer...this just to simply, anyway, what is the pressure value we need to provide at the interface of the two fluids, in order to balance the hydrostatic pressure coming from above? the answer is...the same value of the hydrostatic pressure coming from the above fluid, calculated at the interface. Therefore at the interface between the two fluids you have dp/dn = dp/dy = 0...the virtual fluid is indeed acting as a wall for the fluid above at rest. The boundary condition has to be considered on the boundary only, not in the whole domain as you did in your notes...it is indeed just a condition for the fluid touching the boundary. In this case we will have a discontinuity of the pressure derivative at the bottom boundary, because it is 0 at that boundary

Hope this clarify everything.
Best

[Blanco]

Btw, since we need dp/dn =0 just at the fluids interface, we can have dp/dn = ro*g*h in the above fluid core, the boundary condition needs to be valid only at the interface

[francisco.angel]

Hi murx:

Originally posted by murx
Because dv/dy = 0 does not mean that d2v/dy2 = 0. There is no need for the second derivative to be zero just because the first one is.

I understand, after reading other replies, that if dv/dy=0 at one point, this not imply that d2v/dy2=0 and I agree. I we move "along the surface" (in this case parallel to x), then u=0, du/dx=0 and dv/dy =0, but if we try to calculate d2v/dy2=0 we need to consider the values of dv/dy just outside the surface, and there dv/dy is not necessarily zero. Now I also have doubts about the physical reasoning behind this boundary condition
I also agree with you in that if we consider the pressure, then dp/dy=rho*g, sadly I can't see clearly the argument of Blanco on action-reaction pairs to prove dp/dn is zero at interfaces.

Originally posted by murx
Thanks. What I meant is to check the results and see if the pressure gradient is really zero. I think I remember a collegue of mine saying that when he uses a higher value for the orthogonal correctors, than this boundary condition is apparently not applied and the gradient is not zero.

I remember having read somewhere that in certain cases, boundary conditions were included by "completeness" but not used at the end, but my knowledge of OpenFoam is not that deep to say if this is true or not.

Best regards.

[Blanco]

Hi Francisco, in short:

For your first point - please consider that continuity is valid everywhere, therefore since du/dx =0 is valid everywhere, also dv/dy = = 0 is the results you get everywhere... this should solve your doubts.

For the second point concerning dp/dn = 0...my example was done considering an interface between a real fluid and a "virtual fluid" used to substitute a wall.. it does not matter indeed if you have a wall or another fluid instead, as long as it is at rest in the previous example. As you can see, imposing dp/dy = 0 is the only way to balance the forces acting on the fluids interfaces, therefore since it does not matter if we have fluid or a wall below the interface...the same conclusion is still valid also for the wall.

For the third point: boundary conditions are requested by the equations, in each specific problem we have to consider what are the real boundaries needed...it is up to the user to set up correctly the boundaries, avoiding unused entries. To see if a bounday is needed...just try to delete it!

Best

[Blanco]

Another hint : think about two fluids at rest inside a container, one above the other since we have ro1<ro2, where roN is the density of fluid N. We know from experience that fluid interface is planar, horizontal and parallel to the x axis. The y axis start from the free surface of the top fluids and points down the bottom of the container. We know That relative pressure is 0 at y = 0, try to figure out what is dp/dy at the interface between fluids 1 and 2...

[Blanco]

(Free surface is the x axis in my assumptions...)

[murx]

Originally Posted by Blanco

Another hint : think about two fluids at rest inside a container, one above the other since we have ro1<ro2, where roN is the density of fluid N. We know from experience that fluid interface is planar, horizontal and parallel to the x axis. The y axis start from the free surface of the top fluids and points down the bottom of the container. We know That relative pressure is 0 at y = 0, try to figure out what is dp/dy at the interface between fluids 1 and 2...

I would say that the pressure is discontinuous at the surface and the pressure gradient is either rho1*g or rho2*p, depending on which side you come from.

But lets take the argument for dp/dp=0 and apply it to the velocity of a wall in a pipe flow. The first velocity layer at the wall must have a velocity of zero, just as the wall. So following this argumentation, du/dn must be zero, which is of course wrong. What I am saying is that your argumentation would explain a Dirichlet boundary condition p_wall = pfluid, but not a Neumann condition. However, these are just my thoughts... Feel free to convince me that I am wrong

[Blanco]

Hi murx,

I'not here to convince you of course :-P I just would like to underline my point of view and my understanding... I think you are doing some confusion because you are mixing pressure and velocity boundary conditions, but they are not the same at all!
if we consider the pipe flow as you said, we have obviously v=0 at the wall (the first fluid layer will have 0 velocity if we are considering a continua, otherwise in a cfd code with finite volume you have that the first fluid cell centroid velocity is not 0...the boundary condition at the wall side of the cell tells that v is 0, but At the boundary of the cell!) and dv/dn will be very different from 0 as we know from experience...but What about pressure? What are boundaries condition on the pressure at the wall in your opinion?

I will continue with my previous example in my next post, as soon as I will have time
Best

[Blanco]

Hi again Murx,
I think I misunderstood your previous post and now I get what you mean...
The problem is that dirichlet and neumann conditions are quite similar if looked from "far away" but they have specific important details that make them distinguished. This is why I think you derive a dirichlet condition from my previous example, while I was pointing to a neumann condition at the end...let's try another path instead of continuing the example.

In my opinion:

-Dirichlet condition is used everywhere we already know a variable value at a boundary, therefore we want to impose that variable value EXTERNALLY to the fluid domain interior (from the boundary to the fluid domain). Examples are v at a wall (V=0), p on a constant pressure boundary, v on a constant velocity boundary, etc. This boundary condition force an external value of a variable on the fluid domain interior "touching" the boundary, leaving free the variable derivatives on the boundary. We are not worried about what are the derivatives,we just want to impose a variable value to the fluid surface.

-Neumann condition is used everywhere we really don't know what is the value of a certain variable on a boundary, but we want to get 0 flux of that variable across the boundary...we want the boundary to act as a mirror for the specific variable considered. We want the boundary to "react" to any information coming INTERNALLY from the fluid domain, assuming a boundary value for the variable equal to the internal fluid domain value of the same variable (information goes from fluid domain to the boundary). Examples are v on a constant pressure boundary, p and v on a symmetry plane, and p on a wall in my opinion.
To obtain what I described we impose a null first derivative of the variable across the boundary, and this condition (as all the boundary conditions) is valid only at the boundary. Indeed dp/dn =0 means that p is constant if we move from the boundary in its normal direction, but this derivative is again valid on the boundary only, it is no more valid at a dn far away from the boundary, with dn ->0. We are not forcing the boundary value of the variable, we are just saying that the specific variable will not vary between the boundary and the very first fluid layer touching it. If we integrate the relation dp/dn =0 we get p=const in normal direction at the boundary...what is the value of this constant? P has not a value at the boundary at the very first iteration, because it is not present in boundary initial conditions since we defined only its first derivative on the boundary...but p is defined in initial condition for the fluid domain interior! Therefore the p=const value is set from the fluid interior at the first iteration and we obtain pbound=pint, and this will happen at each iteration.

The difference in my opinion is that dirichlet condition is used to impose a known value of a variable to the fluid domain periphery, while neunamm condition is used to impose an unknown value of variable on the fluid periphery (imposing it to be equal to the fluid variable value touching that periphery, that will be solved from ns eq.).

Best

[Blanco]

Originally Posted by francisco.angel

Hi murx:
I understand, after reading other replies, that if dv/dy=0 at one point, this not imply that d2v/dy2=0 and I agree. I we move "along the surface" (in this case parallel to x), then u=0, du/dx=0 and dv/dy =0, but if we try to calculate d2v/dy2=0 we need to consider the values of dv/dy just outside the surface, and there dv/dy is not necessarily zero. Now I also have doubts about the physical reasoning behind this boundary condition

I think the document cited is correct, even if it doesn't underline why d2v/dy2 = 0 in the most clear way.

We have to consider that from the continuity equation in an incompressible flow we have that:

du/dx+dv/dy = 0

therefore

dv/dy = -du/dx

and this holds also on the wall boundary of course. This means that

d2v/dy2 = d/dy(dv/dy) = -d/dy (du/dx) = -d/dy (0) = 0

because du/dx holds everywhere on the wall surface (u(x,0)=v(x,0)=0 for every x).
Therefore d2v/dy2 can be computed without "going out from the wall surface" but linking the calculus to du/dx thanks to the 0 divergence condition on the velocity field.

This proof that the document is correct since other 0 derivatives are fully justified.

NOTE: Murx please note that the momentum equation you cited for the x-dir is not correct, the correct one is:

v*du/dy+dp/dx=ny*d2u/dy2

otherwise we wouldn't have to solve for anything but the pressure
It can be noticed that the equations we get from the doc are the same obtained in post n°2. Indeed example n° 2 can be used with a single wall, removing the constraints related to u=u(y) and v=v(y) but allowing a more general condition u=u(x,y) and v=v(x,y) and adding a boundary condition for the wall, u(x,0)=v(x,0)=0. The result is correctly shown in the cited doc.

Best regards

[murx]

I would like to come back to this topic, since I still didn't find a physical explanaition for the zero-gradient pressure condition on walls.

Originally Posted by Blanco

d2v/dy2 = d/dy(dv/dy) = -d/dy (du/dx) = -d/dy (0) = 0

because du/dx holds everywhere on the wall surface (u(x,0)=v(x,0)=0 for every x).
Therefore d2v/dy2 can be computed without "going out from the wall surface" but linking the calculus to du/dx thanks to the 0 divergence condition on the velocity field.

I agree that you can link the momentum equations to the continuity equations of course. However, you are implying that du/dx is not a function of y because du/dx(y=0)= 0. And that is a false conclusion. You cannot derive the dependence of du/dx on y by saying it is zero at y=0 and the derivative of 0 is 0. du/dx will definitely change with y unless you assume a developed flow in the given example, but that's obviously not the case we are discussing here.

Originally Posted by Blanco

NOTE: Murx please note that the momentum equation you cited for the x-dir is not correct, the correct one is:

v*du/dy+dp/dx=ny*d2u/dy2

otherwise we wouldn't have to solve for anything but the pressure
Best regards

The momentum equations I cited were simplified, that means I already left out the terms that can be omitted due to the no-slip condition.

You can drop all convective terms and all derivates along the wall because of the no-slip condition. If we stick to the example of a flat wall along the x-direction, the momentum equation in y-direction can be written:

dp/dy = my*d2v/dy2

or in dimensionless form:

dP/dY = 1/Re * d2V/dY

This shows, that for high Reynolds numbers, the zero-gradient pressure is physically meaningful but is not a good approximation for low Reynolds numbers. I think I have seen sombody else stating the latter, but I did not see a mathematical explanation so far. However, this seems to make sense and I doubt that the derivation of dp/dy = 0 in the document linked by francisco.angelis correct beause the make the same conclosion of d/dy (du/dx) = -d/dy (0) = 0 without giving understandable reason why du/dx is not a function of y.