-
Notifications
You must be signed in to change notification settings - Fork 847
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
[WIP] Implementation of Amplification Factor Transport(AFT) 2019b tranistion model #2422
base: develop
Are you sure you want to change the base?
Conversation
Primary back up
Test for SAAFT2019b
Renew_Min_Max
@@ -1407,6 +1451,25 @@ | |||
AddVolumeOutput("TURB_INDEX", "Turb_index", "PRIMITIVE", "Turbulence index"); | |||
break; | |||
|
|||
case TURB_TRANS_MODEL::AFT: | |||
//nodes -> SetAFT_Wonder_Func(iPoint, HL, H12, dNdRet, Ret0, D_H12, l_H12, m_H12, kv, Rev, Rev0, F_growth, F_crit, PAF, Pg, Dg); |
Check notice
Code scanning / CodeQL
Commented-out code Note
AFT2017b
Test_Strain
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 👍 Looks good but we can re-use a lot of code from the LM solver and see the question about the aux var gradient which looks suspicious.
NONE, /*!< \brief No option / default. */ | ||
AFT2017b, /*!< \brief using AFT2017b model. */ | ||
AFT2019b /*!< \brief using AFT2019b model. */ |
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.
Please align the start of the comments and add some links to reference papers.
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.
Sure!
The reference paper is "https://doi.org/10.2514/6.2019-0039."
I'll add the above link.
* \param[in] rank - MPI rank. | ||
* \return Struct with AFT options. | ||
*/ | ||
inline AFT_ParsedOptions ParseAFTOptions(const AFT_OPTIONS *AFT_Options, unsigned short nAFT_Options, int rank) { |
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.
Rank is not used, can you remove it from this function and maybe LMOptions if it is also not used there.
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.
Okay. I remove the rank both AFT and LM.
switch (Kind_Turb_Model) { | ||
case TURB_MODEL::NONE: SU2_MPI::Error("No turbulence model has been selected but AFT transition model is active.", CURRENT_FUNCTION); break; | ||
case TURB_MODEL::SST: SU2_MPI::Error("k-w SST turbulence model has been selected but AFT transition model is active.", CURRENT_FUNCTION); break; | ||
} |
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 block of code appears twice, please refactor this a little to avoid duplicating it.
case AFT_CORRELATION::AFT2017b: { | ||
dNdRet = 0.028*(H12 - 1.0) - 0.0345 * exp(-pow(3.87/(H12 -1.0) - 2.52, 2.0)); | ||
break; | ||
} | ||
|
||
case AFT_CORRELATION::AFT2019b: { | ||
dNdRet = 0.028*(H12 - 1.0) - 0.0345 * exp(-pow(3.87/(H12 -1.0) - 2.52, 2.0)); | ||
break; | ||
} |
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.
See if there are other opportunities to avoid repeating formulas.
case AFT_CORRELATION::AFT2017b: { | |
dNdRet = 0.028*(H12 - 1.0) - 0.0345 * exp(-pow(3.87/(H12 -1.0) - 2.52, 2.0)); | |
break; | |
} | |
case AFT_CORRELATION::AFT2019b: { | |
dNdRet = 0.028*(H12 - 1.0) - 0.0345 * exp(-pow(3.87/(H12 -1.0) - 2.52, 2.0)); | |
break; | |
} | |
case AFT_CORRELATION::AFT2017b: | |
case AFT_CORRELATION::AFT2019b: { | |
dNdRet = 0.028*(H12 - 1.0) - 0.0345 * exp(-pow(3.87/(H12 -1.0) - 2.52, 2.0)); | |
break; | |
} |
void FinishResidualCalc(const CConfig* config) override { | ||
const bool implicit = config->GetKind_TimeIntScheme() == EULER_IMPLICIT; | ||
|
||
/*--- Compute mean effective dynamic viscosity ---*/ | ||
const su2double diff_i_AF = Laminar_Viscosity_i + Eddy_Viscosity_i; | ||
const su2double diff_j_AF = Laminar_Viscosity_j + Eddy_Viscosity_j; | ||
const su2double diff_i_Gamma = Laminar_Viscosity_i + Eddy_Viscosity_i; | ||
const su2double diff_j_Gamma = Laminar_Viscosity_j + Eddy_Viscosity_j; | ||
|
||
const su2double diff_AF = 0.5*(diff_i_AF + diff_j_AF); | ||
const su2double diff_Gamma = 0.5*(diff_i_Gamma + diff_j_Gamma); | ||
|
||
Flux[0] = diff_AF*Proj_Mean_GradScalarVar[0]; | ||
Flux[1] = diff_Gamma*Proj_Mean_GradScalarVar[1]; | ||
|
||
/*--- For Jacobians -> Use of TSL (Thin Shear Layer) approx. to compute derivatives of the gradients ---*/ | ||
if (implicit) { | ||
const su2double proj_on_rho_i = proj_vector_ij/Density_i; | ||
Jacobian_i[0][0] = -diff_AF*proj_on_rho_i; Jacobian_i[0][1] = 0.0; | ||
Jacobian_i[1][0] = 0.0; Jacobian_i[1][1] = -diff_Gamma*proj_on_rho_i; | ||
|
||
const su2double proj_on_rho_j = proj_vector_ij/Density_j; | ||
Jacobian_j[0][0] = diff_AF*proj_on_rho_j; Jacobian_j[0][1] = 0.0; | ||
Jacobian_j[1][0] = 0.0; Jacobian_j[1][1] = diff_Gamma*proj_on_rho_j; | ||
} | ||
} |
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 only differs from LM by some constant factors, refactor it to use the same 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.
You need to find a way to share the common parts between the LM solver and this one, this way there is too much duplication.
switch (config->GetKind_Gradient_Method()) { | ||
case GREEN_GAUSS: | ||
SetAuxVar_Gradient_GG(geometry, config); | ||
break; | ||
case WEIGHTED_LEAST_SQUARES: | ||
SetAuxVar_Gradient_LS(geometry, config); | ||
default: | ||
break; | ||
} |
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.
Why are you computing the AuxVarGradient after you use it above instead of before?
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 wanted to discuss this with someone.
On the SA-AFT2019b model, HL is defined as follows :
where, d is the distance from the wall.
The AuxVar[0] is for the d and AuxVarGradient[0] is for the . So, the d and
are the constant value when the geometry and mesh are not deformed.
And, the Auxvar[1] is for the and AuxVarGradient[1] is for the
.
The velocity values change with each iteration. So, the Auxvar[1] and it's gradient must be redefined and calculated each iteration. Therefore, I decided to do the gradient calculation in the Postprocessing.
On second thought, the Auxvar assign in "trans_sources.hpp", or would you have any other ideas?
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.
In the paper it is mentioned that grad d is the wall-normal unit vector, don't we have that stored already? I remember seeing it in another PR, maybe it was not merged 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.
In the paper it is mentioned that grad d is the wall-normal unit vector, don't we have that stored already? I remember seeing it in another PR, maybe it was not merged yet.
OK. I'll try to find it. :D
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 it is in the simplified 2015 PR
HLGradTerm = nodes->GetAuxVarGradient(iPoint, 1, 0) * nodes->GetAuxVarGradient(iPoint, 0, 0) + nodes->GetAuxVarGradient(iPoint, 1, 1) * nodes->GetAuxVarGradient(iPoint, 0, 1); | ||
|
||
if(nDim == 3) { | ||
HLGradTerm += nodes->GetAuxVarGradient(iPoint, 1, 2) * nodes->GetAuxVarGradient(iPoint, 0, 2); | ||
} |
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.
GeometryToolbox SquaredNorm
su2double Temp3 = flowNodes->GetVelocity(iPoint, 0) * nodes->GetAuxVarGradient(iPoint, 0, 0); | ||
su2double Temp4 = flowNodes->GetVelocity(iPoint, 1) * nodes->GetAuxVarGradient(iPoint, 0, 1); |
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.
Add some documentation please, if the initial aux var grad is 0 then you set the aux var to 0, then its gradient is 0 and this remains 0 forever? Looks like a bug.
|
||
/*-- destruction term of Intermeittency(Gamma) --*/ | ||
const su2double Dg = c_2 * Density_i * VorticityMag * F_turb * (c_3 * exp(lnIntermittency) - 1.0); | ||
nodes -> SetAFT_Wonder_Func(iPoint, HL, H12, dNdRet, Ret0, D_H12, l_H12, m_H12, kv, Rev, Rev0, F_growth, F_crit, PAF, Pg, Dg); |
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 is a lot of memory and these variables are only used for post processing (output) right?
If so these can be computed in the output class when setting the output values.
That saves doing this work every iteration and storing all these 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.
Thank you for your advise.
Current class, this function is for debugging. Later, I'll change to compute in the output class.
Thank you for comment @pcarruscag. I'll try to address each of your comments one by one. For now, the SA-AFT2017b model is for testing purposes, so I'll delete the relevant parts. |
case TURB_MODEL::NONE: SU2_MPI::Error("No turbulence model has been selected but AFT transition model is active.", CURRENT_FUNCTION); break; | ||
case TURB_MODEL::SST: SU2_MPI::Error("k-w SST turbulence model has been selected but AFT transition model is active.", CURRENT_FUNCTION); break; | ||
} | ||
cout << "-2017b" << endl; break; |
Check warning
Code scanning / CodeQL
Dead code due to goto or break statement Warning
if (...) ...
case TURB_MODEL::NONE: SU2_MPI::Error("No turbulence model has been selected but AFT transition model is active.", CURRENT_FUNCTION); break; | ||
case TURB_MODEL::SST: SU2_MPI::Error("k-w SST turbulence model has been selected but AFT transition model is active.", CURRENT_FUNCTION); break; | ||
} | ||
cout << "-2017b" << endl; break; |
Check warning
Code scanning / CodeQL
Dead code due to goto or break statement Warning
Proposed Changes
Hi, all.
This PR is for implementing well known subsonic transition model (SA-AFT2019b).
Coder's SA-AFT2019b model is composed of the amplification factor and the logarithmic intermittency transport equation as shown below.
![image](https://private-user-images.githubusercontent.com/50565314/403425697-9f3d41a7-dbd0-43d6-9327-4bd653ca4f73.png?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3Mzk0Mzk5MjgsIm5iZiI6MTczOTQzOTYyOCwicGF0aCI6Ii81MDU2NTMxNC80MDM0MjU2OTctOWYzZDQxYTctZGJkMC00M2Q2LTkzMjctNGJkNjUzY2E0ZjczLnBuZz9YLUFtei1BbGdvcml0aG09QVdTNC1ITUFDLVNIQTI1NiZYLUFtei1DcmVkZW50aWFsPUFLSUFWQ09EWUxTQTUzUFFLNFpBJTJGMjAyNTAyMTMlMkZ1cy1lYXN0LTElMkZzMyUyRmF3czRfcmVxdWVzdCZYLUFtei1EYXRlPTIwMjUwMjEzVDA5NDAyOFomWC1BbXotRXhwaXJlcz0zMDAmWC1BbXotU2lnbmF0dXJlPWU4ZmQ3YTMzNjk5M2RkZjg1OGI5MmRkNjFjMzQ3YTRlMmFmMDkxNTY4NmFlNzgzNTgwODI5ZmQ2NDk5YWJmYzAmWC1BbXotU2lnbmVkSGVhZGVycz1ob3N0In0.XnP89UCpFs_ZgTV_IiS13D9tcdTd-DkY8iKs2zabtmM)
Current State
The primary implementation of the model has been completed and is now in the validation. The validation problems are refer : https://doi.org/10.2514/6.2023-3530 and https://doi.org/10.2514/6.2019-0039.
The transition onset location and the amplification factor field are slightly different. I'm still debugging and checking various things.
Validation Problem of Current State
Firstly, I'm trying the validation problem for 2-D problems (flat plate and NLF-0416 airfoil).
The current results are shown below.
S&K Flat plate from https://doi.org/10.2514/6.2023-3530
NLF-0416 airfoil from https://doi.org/10.2514/6.2023-3530
![image](https://private-user-images.githubusercontent.com/50565314/403423239-baabe820-615c-42f3-b6f9-7485207965f3.png?jwt=eyJhbGciOiJIUzI1NiIsInR5cCI6IkpXVCJ9.eyJpc3MiOiJnaXRodWIuY29tIiwiYXVkIjoicmF3LmdpdGh1YnVzZXJjb250ZW50LmNvbSIsImtleSI6ImtleTUiLCJleHAiOjE3Mzk0Mzk5MjgsIm5iZiI6MTczOTQzOTYyOCwicGF0aCI6Ii81MDU2NTMxNC80MDM0MjMyMzktYmFhYmU4MjAtNjE1Yy00MmYzLWI2ZjktNzQ4NTIwNzk2NWYzLnBuZz9YLUFtei1BbGdvcml0aG09QVdTNC1ITUFDLVNIQTI1NiZYLUFtei1DcmVkZW50aWFsPUFLSUFWQ09EWUxTQTUzUFFLNFpBJTJGMjAyNTAyMTMlMkZ1cy1lYXN0LTElMkZzMyUyRmF3czRfcmVxdWVzdCZYLUFtei1EYXRlPTIwMjUwMjEzVDA5NDAyOFomWC1BbXotRXhwaXJlcz0zMDAmWC1BbXotU2lnbmF0dXJlPTRjZDg0YTg1NmFiODI4OTQwZjc2YjFiY2Y5N2IwNjBkN2JlZjEzOGRhZmUyNjYyMzEyYmUyZGUwYmY0ZTYxOTUmWC1BbXotU2lnbmVkSGVhZGVycz1ob3N0In0.ZWhsemveK5h6bus-Z2TEyLjgG6Ht6kvU6n30AfQvcMk)
PR Checklist
Put an X by all that apply. You can fill this out after submitting the PR. If you have any questions, don't hesitate to ask! We want to help. These are a guide for you to know what the reviewers will be looking for in your contribution.
pre-commit run --all
to format old commits.