Possible bug in develop compared to master branch.

[EduardoMolina]

Describe the bug

Hi SU2 developers,

I will consider this as a possible bug until I hear the opinion of the SU2 community.

As a validation test case for turbulence models, I am running the NASA Hump test case (https://turbmodels.larc.nasa.gov/nasahump_val.html). However, the results obtained from the develop and master branch are quite different, specially the pressure distribution over the hump. The same difference also appears when other numerical scheme is used (SLAU2 instead of Roe). On the other hand, the skin friction seems not be affected. Please see attached figures.

Also, I attached the plot of the residual over iterations, showing that all simulations converged nicely.

Firstly, I thought that this was related to pull request #670 but I ran the same case with an older version of the dev branch and the results were identical. This suggest that the difference started from pull requests that were merged before #670.

I also attached the config and mesh files if someone wants to reproduce.

Thanks to Danah Park (PhD student) for pointed me out.

I appreciate your time and feedback.

Best,

Eduardo

To Reproduce
Please provide a config file and mesh for a minimal example that recreates the problem.
hump_ddes_full.su2.gz
turb_SA_hump.cfg.txt

Desktop (please complete the following information):

• OS: [Ubuntu 16.04]
• C++ compiler and version: [g++ (GCC) 5.4.0]
• MPI implementation and version: [OpenMPI 1.10.2]
• SU2 Version: [v6.2.0]

[jayantmukho]

Any chance you could run on both branches with SST to narrow down where the bug might be?

[economon]

Yep, it's not an elegant task, but prob best to bisect the git commit history until you can isolate the commit with the change (more or less).

[clarkpede]

I can provide some ideas as to where the change might have occurred. If you look at the git log for the regression tests, you can see every time that someone changed the expected results for a test. Specifically, we can see all the times that the SA regression tests have been changed, which is a necessary (but not sufficient) condition for a significant change in SA results. (It's possible a sneaky change slipped in with a small enough change in residuals so as to not be picked up by the regression tests.)

I looked at the following test cases:

• turb_SA_flatplate
• turb_SA_RAE2822
• turb_ONERAM6
• turb_NACA0012_sa

I only found one commit where the SA regression test values changed since Feb 18, 2019 (when v6.2.0 was released): fa7ef0a . That commit was part of PR #657, where a bug in the symmetry boundary condition was fixed. You do have symmetry boundary conditions in your cfg file.

Why don't you try re-running the test case with the versions before and after PR #657 was merged?

[EduardoMolina]

@jayantmukho : The difference is also observed when using SST.

@economon and @clarkpede : Thanks for the feedback. @clarkpede , I will re-run using the commit that you pointed and I will post the results here.

[EduardoMolina]

So, I compared results from ae64215 with e15ffa1

@clarkpede Thanks for the quick but huge help! I really appreciate it!

Now that we found were the possible bug is located. What should we do? Return to the previous MARKER_SYM as Euler BC and change the current MARKER_SYM to a new MARKER_SYM-FIX for example?

Best,

Eduardo

[economon]

Maybe this will help.. can you please try setting the viscosity and eddy viscosity inside V_reflected in the if (viscous) {} section of the BC_Sym routine with something like:

/*--- Set laminar and eddy viscosity at the reflected state ---*/

V_reflected[nDim+5] = node[iPoint]->GetLaminarViscosity();
V_reflected[nDim+6] = node[iPoint]->GetEddyViscosity();

[EduardoMolina]

Hi @economon

Setting both laminar and eddy viscosities in V_reflected didn't change the results.

Eduardo

[TobiKattmann]

Hi Eduardo,

As I wrote #657 here my bits (thanks @economon to notify me about this issue):
The top wall, which is the symmetry wall in your case, is not a straight line (i.e. line in 2D or plane in 3D) but has a bump in it. Therefore I would say the symmetry BC is inappropriate for that geometry/mesh.
Implementation-wise the BC_sym routine in develop assumes (but not checks) a constant unit normal over all boundary elements in a symmetry marker and I would say that this a necessary condition for a symmetry-marker.

To cite a comment from @economon from that PR #657

As we should always have a plane for this BC type, is this [tangent vector computation to the unit normal, Editors note] something you can check outside the loop, or complete as a preprocessing and store? Looks a little expensive to perform at every vertex during runtime

Because of that preprocessed computation/storage of the unit normal (and a tangent to this normal) the current implementation of BC_Sym cannot be "exploited" to serve as an "inviscid (slip) wall" for non-straight markers.

Using MARKER_EULER= ( sym ) is not correct either as the gradients of the velocity components (used in the viscous stress tensor of the momentum equations) are not computed correctly for the reflected state. Although using MARKER_EULER in the current develop will give you the same results as the master (or older develop before #657 ) with MARKER_SYM (... as it uses the same implementation).

If you pull the "Preprocessing" part from the develop Sym_BC into the loop over all vertices of the boundary marker, the unit normal (and tangential) gets recomputed for each vertex. That gives you results which should be a lot closer to what happens in the master. Plus handling the velocity gradient correctly. From my point of view this would be the correct way.
Apply this patch in the root of the latest develop code using patch -p1 < SymAsInviscidSlipWall.diff to get the code I mean.

SymAsInviscidSlipWall.diff.txt

The implementation of the Symmetry_BC could serve double as 'Inviscid (slip) wall' (which is prescribed in the NASA documents) for viscous flow if the 'Preprocessing' part gets recomputed for each boundary vertex. But how to integrate that best is s.th. to discuss. One could just pull it in at some higher computational cost, put a conditional together with a boolean (e.g. SYMMETRY_AS_SLIP_WALL= YES) or fixing and using MARKER_EULER for that case.

I hope that clarifies what happens with that boundary. Please let me know if that works for you and fixes the issue.

When I fixed the Regression Testcases for #657 i think that this one was not part of it, right? Maybe this validation case goes into the su2code/VandV repo :)

Cheers,Tobi

Edit: I think we can all agree, that it is Toms fault due to his requested changes in the PR 😉

[TobiKattmann]

And here a quick image of the Cp over the whole surface_flow.dat. The master and the descried fix are lying pretty much on top of each other. I reproduced the current dev results which you also found.

Additionally density and momentum_x for a master (sym, euler) + develop (sym, euler, sym with fix) Note that the 'dev sym fix' is close but not identical to the old implementation

[EduardoMolina]

Hi @TobiKattmann ,

Thanks for the detailed explanation!

The only thing that I don't agree is about the assumption that a symmetry BC is a plane or a line.

In my point of view, this is not always true and can lead to wrong interpretation of the results, specially if users and not developers are running the code. I suggest that the "symmetry as the slip wall" (as you mentioned) should be our default and the "symmetry plane" an option.

Please let me know your thoughts.
Thanks again.

Best,
Eduardo

[TobiKattmann]

Hi Eduardo,

both, 'symmetry plane' and 'symmetry as the slip wall' should work in the code (upon fixing this) and it should be made clear (e.g. in the config_template.cfg) what the BC is suitable for, independent of how it is named in the end.
The simplest solution is probably the patch I provided together with a short note in the config_template (maybe even a reference to this issue), as both 'options' will work as expected (using the very same code). It will basically 'restore' the behaviour of MARKER_SYM of before the fix in #657. The downside is some more computational cost (in case of a plane) which has to be done for each vertex on the marker. I cant say to which extend this takes more time but I would guess it is not too costly. Depends on the case as well. I can do a little check with your case here.

That's it for the diplomatic part :) (Please take the following with a grain)

For me a symmetry has to be a line or plane. To cite Jiri Blazek's book 'Computational Fluid Dynamics'(...) 3rd edition using a screenshot:

The chapter is called 'symmetry plane' which already is kinda biased. But there is no 'slip wall for viscous flow' mentioned to my knowledge. In another book (from Ferziger&Peric, which I only have in german) symmetry is only used in the context of symmetry planes as well.
To me, symmetry BC's are also linked to visualization using mirroring. And mirroring (as the word is commonly meant) is only possible along planes and its little brother, the straight line.

I can be wrong or it might be a question of research background. So how to proceed? Trial by combat might have some unpleasant aftermath for the winner so we might go the democratic way and do a poll :)

Cheers,
Tobi

[EduardoMolina]

Hi @TobiKattmann .

I am OK if we clearly explain the changes in the config template as my main concern is with the users side.
Thanks!

Eduardo

[economon]

I am also mostly concerned with making this less error-prone for the user, which means we should put in an error check or force a proper behavior.

Since the Euler bc does not have the viscous terms, maybe we can throw an error if folks try to use it with viscous problems for slip walls and require the Sym bc in that case (with the patch above). What do you both think?

[economon]

Oh, and this one is definitely my fault 😀