-
Notifications
You must be signed in to change notification settings - Fork 915
Reduce flow solver triplication (comp, incomp, nemo) #1177
New issue
Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.
By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.
Already on GitHub? Sign in to your account
Conversation
| template<class SoundSpeedFunc, class LambdaViscFunc> | ||
| FORCEINLINE void SetTime_Step_impl(const SoundSpeedFunc& soundSpeed, | ||
| const LambdaViscFunc& lambdaVisc, | ||
| CGeometry *geometry, |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The idea is to generalize this type of function, the details of the implementation are provided by the derived solver when it instantiates (calls) the generic implementation.
There is no explicit "if incompressible" or "if nemo" or "if compressible" logic.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks for going through this...I made some attempts to update time integration schemes in another branch but didnt push yet.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
oops sorry about that, I was looking at this kind of routine and had a "what if" moment... and now it's done
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
No apologies! Im happy its done.
| /*--- Define an object to compute the speed of sound. ---*/ | ||
| struct SoundSpeed { | ||
| FORCEINLINE su2double operator() (const CEulerVariable& nodes, unsigned long iPoint, unsigned long jPoint) const { | ||
| return 0.5 * (nodes.GetSoundSpeed(iPoint) + nodes.GetSoundSpeed(jPoint)); | ||
| } |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The function objects are defined at function scope (SetTime_Step in this example) to make very clear what are the details for each solver.
Using virtual functions or implementing these details in the corresponding variables classes would work, but it would make finding things a lot harder.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think this looks good. Can you explain the FORCEINLINE?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Hmmm where to start... inline is mostly a mechanism to bypass ODR (one definition rule) if you define a function in an hpp (instead of just declaring), and include that hpp multiple times you are defining the same function multiple times which is not legal, declaring it inline makes it legal (for class methods the rules are less strict).
Because a function is defined inline, the compiler MAY decide (based on a bunch of heuristics I imagine) to copy paste its body into the places where it is used, instead of creating a function call (less indirection and more opportunity to optimize).
Most compilers will have a way for us to ask them to "just do it", i.e. forget about the heuristics and just copy the thing, in SU2 I called that FORCEINLINE.
Summary: FORCEINLINE is triple dog daring the compiler to inline the function.
| * \brief Update the solution using the classical RK4 scheme. | ||
| * \param[in] geometry - Geometrical definition of the problem. | ||
| * \param[in] solver_container - Container vector with all the solutions. | ||
| * \param[in] config - Definition of the particular problem. | ||
| * \param[in] iRKStep - Runge-Kutta step. | ||
| */ | ||
| void ClassicalRK4_Iteration(CGeometry *geometry, CSolver **solver_container, | ||
| CConfig *config, unsigned short iRKStep) final; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
This was missing from NEMO and the incompressible solver, any special reason?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thats interesting, there shouldn't be a reason to it is missing.
WallyMaier
left a comment
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
@pcarruscag this work looks great. Ill be on the lookout to reduce these dup/triplications in future PRs.
| /*--- Define an object to compute the speed of sound. ---*/ | ||
| struct SoundSpeed { | ||
| FORCEINLINE su2double operator() (const CEulerVariable& nodes, unsigned long iPoint, unsigned long jPoint) const { | ||
| return 0.5 * (nodes.GetSoundSpeed(iPoint) + nodes.GetSoundSpeed(jPoint)); | ||
| } |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I think this looks good. Can you explain the FORCEINLINE?
| conv_numerics->SetNormal(Normal); | ||
|
|
||
| for (iDim = 0; iDim < nDim; iDim++) | ||
| V_infty[iDim+1] = GetVelocity_Inf(iDim); |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
These calls V_infty[iDim+1], are one of the main differences between the NEMO solvers and everything. Would it be wasteful to use the the NEMO indices through the code? I haven't looked though enough of the code to know the full impact, but it seems that it could go a long way to reduce duplication.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It would be a little yes, because the velocity offset for nemo is not known at compile-time.
But I have something that might interest you in a subsequent PR.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Oh boy, I cant wait
TobiKattmann
left a comment
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Thanks for these efforts of Cut, Paste, Abstract, Survive 💐 . And yet again: sorry for being too late to the party
(I hope you don't follow all of Mr. Grylls advices, especially in the realm of refreshing beverages 😄 ) (Mr. Grylls's advice? Mr. Grylls' advice? i dont know...)
| */ | ||
| void SetResidual_DualTime(CGeometry *geometry, CSolver **solver_container, CConfig *config, unsigned short iRKStep, | ||
| unsigned short iMesh, unsigned short RunTime_EqSystem) override; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I wasn't aware that one can withhold the template parameters for the prototype and therefore can overload the function just by using different template parameters https://gcc.godbolt.org/z/9vsEfo . Makes sense though. I learned sth new today :)
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
That is not exactly what happens there, the method is not a template, the class is.
Each flow solver, incomp, nemo, comp, instantiates a different flavor of CFVMFlowSolverBase.
Each template instantiation creates a different SetResidual_DualTime that overrides the CSolver virtual function.
I don't think you can have method templates of virtual functions, as there would be no way to instantiate them, templates and overloading happen at compile time, and virtual functions at runtime.
That is why it is important to use the virtual specifiers (override, final) to guarantee that a method really overrides something (exact same signature, decide which method to call at runtime based on type) rather than overloading it (different signature, including template params, decide what to call at compile time based on arguments).
What you can have is class templates with virtual functions, each class template instantiation is like a totally different class.
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yep, thanks for the explanation. https://gcc.godbolt.org/z/T6555e this should better reflect the SU2 reality
| */ | ||
| inline virtual void SetBeta_Parameter(CGeometry *geometry, | ||
| CSolver **solver_container, |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
👍 for putting the inc specific functions where they belong
| /*--- Check if the node belongs to the domain (i.e., not a halo node) ---*/ | ||
|
|
||
| if (geometry->nodes->GetDomain(iPoint)) { | ||
| if (!geometry->nodes->GetDomain(iPoint)) continue; |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I guess because of this change (and following mass indentation ) the diffs sometimes confusing💫 although from what I see you didn't change a lot more
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It was just that, geometry toolbox, and avoiding those old arrays allocated at Solver level ("Residual", "Res_Conv", "Vector", "Primitive" and so on) in preparation for hybrid parallel stuff.
|
|
||
| void CIncEulerSolver::SetResidual_DualTime(CGeometry *geometry, CSolver **solver_container, CConfig *config, | ||
| unsigned short iRKStep, unsigned short iMesh, unsigned short RunTime_EqSystem) { |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Ok I guess here it was not straight forward to have that covered in the FVMBase class implementation?
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
It is a bit more complicated because of the "if energy" and the conversion from primitive to conservative variables
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
👍 ok I saw the V2U stuff and wondered whether that could be generalized in some way, but if the energy stuff comes on top it might not make too much sense
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Yeah it changes whether the loops run to nVar or to nVar-1... It could be done but it felt convoluted.
| void CIncNSSolver::BC_HeatFlux_Wall(CGeometry *geometry, CSolver**, CNumerics*, | ||
| CNumerics*, CConfig *config, unsigned short val_marker) { | ||
|
|
||
| void CIncNSSolver::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **solver_container, CNumerics *conv_numerics, | ||
| CConfig *config, unsigned short val_marker) { | ||
| BC_Wall_Generic(geometry, config, val_marker, HEAT_FLUX); | ||
| } | ||
|
|
||
| unsigned short iVar, jVar, iDim, Wall_Function; | ||
| unsigned long iVertex, iPoint, total_index, Point_Normal; | ||
| void CIncNSSolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver**, CNumerics*, | ||
| CNumerics*, CConfig *config, unsigned short val_marker) { | ||
|
|
||
| su2double *Coord_i, *Coord_j, dist_ij; | ||
| su2double *GridVel, There, Tconjugate, Twall= 0.0, Temperature_Ref, thermal_conductivity, HF_FactorHere, HF_FactorConjugate; | ||
| BC_Wall_Generic(geometry, config, val_marker, ISOTHERMAL); | ||
| } | ||
|
|
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
Nice simplification, especially useful if '!energy=true' ;)
| channel_2D.cfg_dir = "sliding_interface/channel_2D" | ||
| channel_2D.cfg_file = "channel_2D_WA.cfg" | ||
| channel_2D.test_iter = 2 | ||
| channel_2D.test_vals = [2.000000, 0.000000, 0.398053, 0.352779, 0.405462] |
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
I might have missed it in the diff but why did this Testcase change? It is only a minor change but I am curious where it was introduced
There was a problem hiding this comment.
Choose a reason for hiding this comment
The reason will be displayed to describe this comment to others. Learn more.
The sliding mesh cases are sensitive to very minor things
|
Never late to the party @TobiKattmann nothing is set in stone with git. I will have more incomp changes at some point and was a bit over eager to merge to keep myself from pilling them in just one PR. |
Proposed Changes
Move more routines to the FVM flow solver base class.
At the risk of sounding like a broken record, one function means:
One test, one optimization, one debug, one fix.
Implementing new features by copy-paste-modify is like eating your legs to get to the end of a marathon, it ain't gonna work in the long run xD
Cut, Paste, Abstract, Survive
Related Work
#1044
PR Checklist