Consistent mass fluxes between bulk flow and turbulence solvers #726

pcarruscag · 2019-07-01T13:32:57Z

Proposed Changes

This is mainly to address #721 but it may have good side effects in other areas (see comments below).
Every turbulent testcase will fail and so far I only tested SA and SST on a simple airfoil case.

Related Work

Addresses #721
Should help a bit with #594
May help with #711

PR Checklist

I am submitting my contribution to the develop branch.
My contribution generates no new compiler warnings (try with the '-Wall -Wextra -Wno-unused-parameter -Wno-empty-body' compiler flags).
My contribution is commented and consistent with SU2 style.
I have added a test case that demonstrates my contribution, if necessary.

pcarruscag · 2019-07-01T13:36:38Z

SU2_CFD/src/numerics_direct_turbulent.cpp

-  a0 = 0.5*(q_ij+fabs(q_ij));
-  a1 = 0.5*(q_ij-fabs(q_ij));
+  a0 = max(0.0, MassFlux);
+  a1 = min(0.0, MassFlux);


Here it is possible to set a0 = -a1 on the assumption that at convergence the summation of fluxes is 0, with that diagonal dominance is guaranteed (at least for the convective part) which improves the condition number of the assembled Jacobian matrix.

pcarruscag · 2019-07-01T13:41:10Z

SU2_CFD/src/numerics_direct_turbulent.cpp

  AD::SetPreaccIn(TurbVar_i, nVar);  AD::SetPreaccIn(TurbVar_j, nVar);
-  if (grid_movement) {
-    AD::SetPreaccIn(GridVel_i, nDim); AD::SetPreaccIn(GridVel_j, nDim);
-  }


By not computing the mass flux again we greatly reduce the number of input variables for preaccumulation.

pcarruscag · 2019-07-01T13:48:20Z

SU2_CFD/src/numerics_direct_turbulent.cpp


-  val_residual[0] = a0*TurbVar_i[0]+a1*TurbVar_j[0];
+  su2double OneOnAveRho = 2.0/(V_i[nDim+2]+V_j[nDim+2]);
+


SA is funny and so density is still needed here...

jayantmukho · 2019-07-08T17:20:22Z

SU2_CFD/src/numerics_direct_turbulent.cpp

   val_Jacobian_i[1][0] = 0.0;
-   val_Jacobian_i[1][1] = -2.0*beta_blended*TurbVar_i[1]*Volume;
+   val_Jacobian_i[1][1] = -2.0*beta_blended*Density_i*TurbVar_i[1]*Volume;
  }


Damn, I've been staring at these lines of code for a while and I never realized that Density was missing from the Jacobian computation. Good catch!

It was not missing, before the Jacobians were w.r.t. rhok and rhow, that is why in the diffusion part they had to be divided by density, I changed it to k and w (only for the SST model) since the mass flux already contains rho. I am not 100% sure this is correct though... In fact CTurbSolver updates the solution via node[iPoint]->AddConservativeSolution. So a better catch by you actually!

jayantmukho · 2019-07-08T17:22:20Z

I can run this on some OneraM6 meshes to confirm that it works in 3D as well.

pcarruscag · 2019-07-08T18:15:19Z

Thanks @jayantmukho that would be great. If the meshes you have show grid convergence the results should converge to the same values, both ways of computing the fluxes are second order, but hopefully by using consistent fluxes the error constant is smaller.

clarkpede · 2019-07-08T18:15:03Z

SU2_CFD/src/numerics_direct_turbulent.cpp

                                               CConfig *config) {

-  val_residual[0] = a0*Density_i*TurbVar_i[0]+a1*Density_j*TurbVar_j[0];
-  val_residual[1] = a0*Density_i*TurbVar_i[1]+a1*Density_j*TurbVar_j[1];


What's your rationale behind making these density changes?

I understand the problems with reconstruction and why the value of density is questionable. But it seems to me like removing density is incorrect. There's two reasons I think that:

The SST variables are updated by the CVariable::AddConservativeSolution call. In that functions, the delta in the variables is assumed to be weighted by density, e.g. rho*k. Based on that you can conclude that the equations are equations for the density-weighted variables, and the residuals should also be density weighted.

In the source term for the SST variables, (CSourcePieceWise_TurbSST) there are terms like (2/3*rho*k*diverg*Volume). These residuals are density weighted. To be consistent, all the residuals should be density weighted.

For a converged, steady state solution, I think these changes will only change the conditioning. Similarly, I expect no impact for a low-Mach test case. Have you tested these changes on a case like the RAE 2822?

Edit: I just saw the comments you made while I was typing out my review. The mass flux used here (in the convection terms) explains the first question. I still question the density changes in the source terms though.

(I've been looking at the code a bit more to get an answer).
The rationale was: I did not want to divide the mass flux by an average density which would then break that "summation = 0" property we can use to improve conditioning a bit. Following that logic I linearized the other terms w.r.t. k and w instead of rho k and rho w, I did not change the source residual only its Jacobian. But I forgot to see how the solution was being updated so yeah in the end it is wrong.
Seeing how SA models need a volumetric flux anyway I will probably revert those changes to SST.

clarkpede · 2019-07-08T18:23:06Z

SU2_CFD/include/numerics_structure.hpp

   * \brief Adds any extra variables to AD
   */
-  void ExtraADPreaccIn();
+  inline void ExtraADPreaccIn() {}


The "inline" here is redundant. Functions defined inside a class definition are automatically "inline."

In common practice yes but apparently that is not part of the standard (which I have not read) and so some compilers might complain about multiple definition.

Ah ok. So no harm done in leaving it in.

pcarruscag · 2019-07-09T23:28:26Z

@clarkpede you are right this seems to make no difference to drag and lift regardless of Mach number and coarseness of the grid. But, it is absolutely terrible for convergence in the asymptotic region and therefore to the discrete adjoint, and so it's a garbage idea even if it saves a bit of memory/computation.
I'll post my results in #721 but in the meantime let no one else waste time with this, I owe you and @jayantmukho a beer after this...

clarkpede · 2019-07-10T01:34:26Z

@clarkpede you are right this seems to make no difference to drag and lift regardless of Mach number and coarseness of the grid. But, it is absolutely terrible for convergence in the asymptotic region and therefore to the discrete adjoint, and so it's a garbage idea even if it saves a bit of memory/computation.
I'll post my results in #721 but in the meantime let no one else waste time with this, I owe you and @jayantmukho a beer after this...

I wouldn't say it's a garbage idea. You correctly identified a shortcoming in the code (inconsistent mass fluxes). That is a problem that should be fixed, and you came closer than anyone else has to fixing it. But it turns out it's a trickier problem than you or I expected.

pcarruscag added 2 commits July 1, 2019 13:41

Euler/NS solver stores fluxes that CSolverTurb can access

17e5df0

forgot to comment a function

9c33db0

pcarruscag commented Jul 1, 2019

View reviewed changes

need to mark mass flux for preaccumulation

b434c61

pcarruscag commented Jul 1, 2019

View reviewed changes

jayantmukho reviewed Jul 8, 2019

View reviewed changes

clarkpede reviewed Jul 8, 2019

View reviewed changes

pcarruscag closed this Jul 9, 2019

pcarruscag mentioned this pull request Mar 25, 2020

Hybrid parallel coloring fallback strategies (better strong scaling and user friendliness) #908

Merged

5 tasks

bigfooted mentioned this pull request Jul 4, 2022

Consistent mass fluxes between bulk flow and scalar solvers and boundedness correction #1697

Merged

5 tasks


		val_residual[0] = a0TurbVar_i[0]+a1TurbVar_j[0];
		su2double OneOnAveRho = 2.0/(V_i[nDim+2]+V_j[nDim+2]);

Consistent mass fluxes between bulk flow and turbulence solvers #726

Consistent mass fluxes between bulk flow and turbulence solvers #726

Uh oh!

Conversation

pcarruscag commented Jul 1, 2019

Proposed Changes

Related Work

PR Checklist

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

jayantmukho commented Jul 8, 2019

Uh oh!

pcarruscag commented Jul 8, 2019

Uh oh!

clarkpede Jul 8, 2019 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

Choose a reason for hiding this comment

Uh oh!

pcarruscag commented Jul 9, 2019

Uh oh!

clarkpede commented Jul 10, 2019 • edited Loading Uh oh! There was an error while loading. Please reload this page.

Uh oh!

Uh oh!

Reviewers

Assignees

Labels

Projects

Milestone

Development

Uh oh!

3 participants

clarkpede Jul 8, 2019 •

edited

Loading

clarkpede commented Jul 10, 2019 •

edited

Loading