SST omega is wrong in sharp corners

[bigfooted]

Describe the bug

In a 3D channel with square cross section and inflation layers, omega can be wrong in the corner. Below is a cross section of such a channel. Omega is computed using:
$$\omega_{wall} = 60 \nu / (\beta d^2)$$,
with d the "distance to the next point away from the wall" [Menter, 1994]. We take the point that is most normal to the wall, computed using FindNormal_Neighbor. However, in the case below, the returned points are the points on the wall. Since they are very far away, the corner omega is very small, order 1e3 instead of 1e6 in the points surrounding the corner point, which causes convergence issues.

Menter mentions that the interior point should actually have y+ < 3. So we should favor nodes that are closer. Maybe prefer nodes based on a weight factor using the distance and the angle?

[bigfooted]

In any case since the boundary condition is based on an y+ value, I think points on the wall should be excluded, if possible. There is an 'edge case' (pun intended) where you can have a triangle in the corner, so the corner node is connected to 2 nodes on the wall and no internal nodes.

[pcarruscag]

There are two versions of the w boundary condition, the "60" version at the face/wall and the "6" version on the first wall CV. We can switch to the latter, it's consistent with wall functions.

[rois1995]

Is the problem related to the normal of the corner point being computed wrong?

[pcarruscag]

The problem is evaluating what is the wall distance to use, it cannot be the wall distance of the wall point because that is 0, so you want the smallest non-zero wall distance of a neighbor of the wall point.
We are incorrectly using the distance to the "most normal neighbor" (a figure we invented in SU2) but like Nijso shows, for corner points that neighbor can be a wall point.

[pcarruscag]

Of course, all of this is brought about by the quirks of vertex centered discretizations :)

[bigfooted]

@pcarruscag do you mean like the approach of Menter and Esch (What Fluent uses):
https://www.researchgate.net/file.PostFileLoader.html?id=58500561b0366d121d756004&assetKey=AS:438841040478208@1481639265395

This is the paper where they introduce a blending function for the wall bc $\omega(y^+) = (\omega_{visc}(y^+)^2 + \omega_{log}(y^+)^2)^{0.5}$
and they use $\omega_{visc}(y) = \frac{6 \nu}{\beta y^2}$

[bigfooted]

..and they use the location of the cell center for everything, do you mean we should also take the cell center?

[pcarruscag]

They are both cell centers for us, but one coincides with the wall and the other is the first node of the wall, there are some testing benefits to having a consistent treatment for no-slip and wall-modelled walls. But I don't know if the "60" version at the wall is more robust.

[rois1995]

Maybe I am missing something, but in this particular case, if the wall normal on a boundary point is taken as the mean of the boundary edge's normal that are connected to that point, shouldn't the first normal point be the correct one?

[pcarruscag]

if the vertical and horizontal walls are different markers there are two vertices that share that node

[bigfooted]

See image below, for corner node P, we check the connected nodes {1,2,3} to see which one is 'most normal' to the corner node, which is point 1 if we compute it for the nodes on edge (P-3). For the nodes on edge (P-1), it ill be point 3. One value will be overwritten by the other, depending on the ordering of the markers. Both values are wrong.

[rois1995]

I see. So the problem only appears if adjacent wall surfaces have different markers. Otherwise the normal should be the averaged one and the first normal point should be correct.

[bigfooted]

Points 1 and 3 can be on the same marker for the problem to occur. We do not average. Also, 1 and 3 are on the wall and not in the interior and we need an interior point, which is 2. If we do not allow wall nodes then we will end up with point 2.

[bigfooted]

I have a quick fix in the above mentioned PR. What I also noticed is that the wall omega strongly depends on the local mesh. Here is omega on the 90 degree bend:

[rois1995]

In this case I see that the first cell height is not the same for every point, that is why you can see that difference. I don't think that there's much you can do in that case since the boundary condition on w has directly depends on the first cell height.