Pressure solution is not unique in splitting scheme for Navier Stokes equations

General
1

Hello,
In this Fenics example for the Navier Stokes equations, in the second step of the splitting scheme the boundary conditions (BCs) for the pressure correction \phi are

  • a Dirichlet BC in the boundary \partial \Omega_N. This results from the fact that the pressure is constrained by a Dirichlet BC on \partial \Omega_N.
  • a Neumann BC in the boundary \partial \Omega_D (where Dirichlet BCs are imposed on the velocity).

Now consider a slightly different problem: Navier Stokes equations in a disk, in which I impose, on the disk boundary \partial \Omega

  • a traction BC \sigma_{ij}n_j = ... where \sigma is the stress tensor and \hat n the unit normal
  • a BC on the normal velocity v_i n_i = ...

Proceeding along the lines as in the Fenics example, the second step is
\nabla^2 \phi = \frac{\rho}{\delta t} \nabla \cdot u^* \qquad\text{in } \Omega
with BC
\hat{n} \cdot \vec{\nabla} \phi = 0 on \partial \Omega.

The solution is not unique, and defined modulo an additive constant. How can I get the correct pressure field and get around this non-uniqueness?

Thank you

2

So this is your Neumann boundary. Then there would be a Dirichlet condition there, right?

You just write out \sigma_{ij} n_j = 2\mu[\nabla^s u]_{ij} \cdot n_j - p n_i.

(This is what you should always do, really. Often one can assume the first part is small, which is why the “free outflow” condition is often misunderstood as a pressure condition. )

3

Yes, like I said in my original post, ‘a BC on the normal velocity …’’

What do you mean by ‘you just out’? Also, \nabla^s is not defined.

The variational problem that I posted determines \phi up to an additive constant, and I would like to know how to fix that constant.

Thank you!

4

Typo, fixed now. \nabla^s is the symmetric gradient (the strain rate operator). \sigma = 2\mu \nabla^s u - p I is how the stress tensor is defined for an incompressible newtonian fluid.

This is not related to the normal component boundary condition though… As you indicated in your initial post, the Dirichlet condition for the pressure problem occurs at the Neumann (or really, the natural) boundaries for the Navier-Stokes problem. For the Navier-Stokes equations, the natural boundary condition is:

\boldsymbol{\sigma} \cdot \boldsymbol{n} = \boldsymbol{t}_N

i.e.,

2\mu \nabla^s \boldsymbol{u} \cdot \boldsymbol{n} - p \boldsymbol{n} = \boldsymbol{t}_N

Or, by dotting with \boldsymbol{n}

p = (2\mu \nabla^s \boldsymbol{u} \cdot \boldsymbol{n}) \cdot \boldsymbol{n} - \boldsymbol{t} \cdot \boldsymbol{n}

Which gives you your Dirichlet expression for p. That Dirichlet expression fixes the constant you’re referring to.

5

Yes. So what if I

  • solve for \phi by adding a temporary Dirichlet condition on a vertex on the mesh, setting for example \phi = 0 on that vertex. This temporary boundary condition selects one out of the many degenerate solutions, and prevents the solve from going to ‘large’ solutions.
  • The resulting pressure is p_{n-1/2} = p_{n-3/2} + \phi + C, and I determine C from the condition 𝑝 =(2𝜇∇𝑠𝒖 ⋅𝒏) ⋅𝒏 −𝒕_N ⋅𝒏?