Illegal access to normal neighbor of halo vertex

[maxaehle]

Describe the bug
In CFVMFlowSolverBase<>::Friction_Forces, there is the following piece of code (permalink):

    /*--- Loop over the vertices to compute the forces ---*/
    
    for (iVertex = 0; iVertex < geometry->nVertex[iMarker]; iVertex++) {
      iPoint = geometry->vertex[iMarker][iVertex]->GetNode();
      iPointNormal = geometry->vertex[iMarker][iVertex]->GetNormal_Neighbor();
      ...
    }

However, a normal neighbor exists only for vertices on the domain (see CMultiGridGeometry::FindNormal_Neighbor and CPhysicalGeometry::FindNormal_Neighbor), not for halo points.
Edit: The link for CPhysicalGeometry was wrong. CPhysicalGeometry actually computes the normal neighbors for halo points, CMultigridGeometry doesn't.

iPointNormal is needed to retrieve Coord_Normal, which is needed to compute WallDist/WallDistMod, which are needed to compute YPlus. To fix the above problem, I would therefore suggest to compute YPlus only for domain points.

Looking for locations where YPlus could be accessed (to see whether YPlus is needed on non-domain points), I found e.g. CTurbSSTSolver::SetTurbVars_WF. However there again the normal neighbor is accessed without checking that we're on a domain point, so likely this should also be changed.

[pcarruscag]

Have you checked if you ever get halo point in normal markers?

[maxaehle]

Yes, in the testcase described here: #1362 (comment) I saw in the debugger that halo vertices occur in the loop.

[pcarruscag]

Right, in the viscous forces case, it should be ok to move the "if Domain" up and avoiding computations on halo points.

[pcarruscag]

I guess this is a problem for FULL MG only, because CPhysicalGeometry finds normal neighbors also for halo points.

[pcarruscag]

@bigfooted , I see we are skiping halo nodes when we compute TauWall via wall functions, probably also because of this PointNormal problem?
This may be a problem, because some edges might be missing one tau wall when computing the viscous fluxes.
If that is the case the easier solution is to communicate TauWall.

[maxaehle]

By the way, is there a reason why CMultiGridGeometry::FindNormal_Neighbor does not compute the normal neighbor for halo nodes (if anyone happens to know)?

[maxaehle]

in the viscous forces case, it should be ok to move the "if Domain" up and avoiding computations on halo points.

Corrected in commit f6f986e in #1362 similar to that, moving the code depending on the nearest neighbor access (i.e. YPlus calculation) down into an "if Domain" block. Note however that after this change YPlus is not computed on halo cells any more, which might be problematic e.g. for CTurbSSTSolver::SetTurbVars_WF.

[pcarruscag]

Maybe it has something to do with the weird shapes of agglomerated volumes, which could make it difficult to find a consistent normal neighbor across ranks.

When Wall Functions are used, y+ is computed elsewhere and not in Viscous forces, and that elsewhere function is already skiping halos, hence my earlier comment on potentially having to communicate tau wall, etc.

[stale]

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. If this is still a relevant issue please comment on it to restart the discussion. Thank you for your contributions.