Skip to content

Commit

Permalink
Merge pull request #2271 from sun5k/feature_Compressibility_Correction
Browse files Browse the repository at this point in the history 
Compressibility correction for SST model
  • Loading branch information
@pcarruscag
pcarruscag committed Jul 21, 2024
2 parents 3cfd801 + c85ee1a commit 1044ccc
Show file tree
Hide file tree
Showing 7 changed files with 377 additions and 6 deletions.
19 changes: 19 additions & 0 deletions Common/include/option_structure.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -990,6 +990,8 @@ enum class SST_OPTIONS {
V, /*!< \brief Menter k-w SST model with vorticity production terms. */
KL, /*!< \brief Menter k-w SST model with Kato-Launder production terms. */
UQ, /*!< \brief Menter k-w SST model with uncertainty quantification modifications. */
COMP_Wilcox, /*!< \brief Menter k-w SST model with Compressibility correction of Wilcox. */
COMP_Sarkar, /*!< \brief Menter k-w SST model with Compressibility correction of Sarkar. */
};
static const MapType<std::string, SST_OPTIONS> SST_Options_Map = {
MakePair("NONE", SST_OPTIONS::NONE)
Expand All @@ -1002,6 +1004,8 @@ static const MapType<std::string, SST_OPTIONS> SST_Options_Map = {
MakePair("VORTICITY", SST_OPTIONS::V)
MakePair("KATO-LAUNDER", SST_OPTIONS::KL)
MakePair("UQ", SST_OPTIONS::UQ)
MakePair("COMPRESSIBILITY-WILCOX", SST_OPTIONS::COMP_Wilcox)
MakePair("COMPRESSIBILITY-SARKAR", SST_OPTIONS::COMP_Sarkar)
};

/*!
Expand All @@ -1013,6 +1017,8 @@ struct SST_ParsedOptions {
bool sust = false; /*!< \brief Bool for SST model with sustaining terms. */
bool uq = false; /*!< \brief Bool for using uncertainty quantification. */
bool modified = false; /*!< \brief Bool for modified (m) SST model. */
bool compWilcox = false; /*!< \brief Bool for compressibility correction of Wilcox. */
bool compSarkar = false; /*!< \brief Bool for compressibility correction of Sarkar. */
};

/*!
Expand Down Expand Up @@ -1048,6 +1054,8 @@ inline SST_ParsedOptions ParseSSTOptions(const SST_OPTIONS *SST_Options, unsigne
const bool sst_v = IsPresent(SST_OPTIONS::V);
const bool sst_kl = IsPresent(SST_OPTIONS::KL);
const bool sst_uq = IsPresent(SST_OPTIONS::UQ);
const bool sst_compWilcox = IsPresent(SST_OPTIONS::COMP_Wilcox);
const bool sst_compSarkar = IsPresent(SST_OPTIONS::COMP_Sarkar);

if (sst_1994 && sst_2003) {
SU2_MPI::Error("Two versions (1994 and 2003) selected for SST_OPTIONS. Please choose only one.", CURRENT_FUNCTION);
Expand All @@ -1068,9 +1076,20 @@ inline SST_ParsedOptions ParseSSTOptions(const SST_OPTIONS *SST_Options, unsigne
SSTParsedOptions.production = SST_OPTIONS::UQ;
}

// Parse compressibility options
if (sst_compWilcox && sst_compSarkar) {
SU2_MPI::Error("Please select only one compressibility correction (COMPRESSIBILITY-WILCOX or COMPRESSIBILITY-SARKAR).", CURRENT_FUNCTION);
} else if (sst_compWilcox) {
SSTParsedOptions.production = SST_OPTIONS::COMP_Wilcox;
} else if (sst_compSarkar) {
SSTParsedOptions.production = SST_OPTIONS::COMP_Sarkar;
}

SSTParsedOptions.sust = sst_sust;
SSTParsedOptions.modified = sst_m;
SSTParsedOptions.uq = sst_uq;
SSTParsedOptions.compWilcox = sst_compWilcox;
SSTParsedOptions.compSarkar = sst_compSarkar;
return SSTParsedOptions;
}

Expand Down
14 changes: 14 additions & 0 deletions Common/src/CConfig.cpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -3475,6 +3475,14 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i
saParsedOptions = ParseSAOptions(SA_Options, nSA_Options, rank);
}

if (Kind_Solver == MAIN_SOLVER::INC_RANS && sstParsedOptions.compSarkar){
SU2_MPI::Error("COMPRESSIBILITY-SARKAR only supported for SOLVER= RANS", CURRENT_FUNCTION);
}

if (Kind_Solver == MAIN_SOLVER::INC_RANS && sstParsedOptions.compWilcox){
SU2_MPI::Error("COMPRESSIBILITY-WILCOX only supported for SOLVER= RANS", CURRENT_FUNCTION);
}

/*--- Check if turbulence model can be used for AXISYMMETRIC case---*/
if (Axisymmetric && Kind_Turb_Model != TURB_MODEL::NONE && Kind_Turb_Model != TURB_MODEL::SST){
SU2_MPI::Error("Axisymmetry is currently only supported for KIND_TURB_MODEL chosen as SST", CURRENT_FUNCTION);
Expand Down Expand Up @@ -6170,6 +6178,12 @@ void CConfig::SetOutput(SU2_COMPONENT val_software, unsigned short val_izone) {
cout << "\nperturbing the Reynold's Stress Matrix towards " << eig_val_comp << " component turbulence";
if (uq_permute) cout << " (permuting eigenvectors)";
break;
case SST_OPTIONS::COMP_Wilcox:
cout << " with compressibility correction of Wilcox";
break;
case SST_OPTIONS::COMP_Sarkar:
cout << " with compressibility correction of Sarkar";
break;
default:
cout << " with no production modification";
break;
Expand Down
32 changes: 27 additions & 5 deletions SU2_CFD/include/numerics/turbulent/turb_sources.hpp
View file
Original file line number Diff line number Diff line change
Expand Up @@ -740,6 +740,7 @@ class CSourcePieceWise_TurbSST final : public CNumerics {
AD::SetPreaccIn(Vorticity_i, 3);
AD::SetPreaccIn(V_i[idx.Density()], V_i[idx.LaminarViscosity()], V_i[idx.EddyViscosity()]);
AD::SetPreaccIn(V_i[idx.Velocity() + 1]);
AD::SetPreaccIn(V_i[idx.SoundSpeed()]);

Density_i = V_i[idx.Density()];
Laminar_Viscosity_i = V_i[idx.LaminarViscosity()];
Expand Down Expand Up @@ -778,6 +779,8 @@ class CSourcePieceWise_TurbSST final : public CNumerics {

const su2double VorticityMag = GeometryToolbox::Norm(3, Vorticity_i);
su2double P_Base = 0;
su2double zetaFMt = 0.0;
const su2double Mt = sqrt(2.0 * ScalarVar_i[0]) / V_i[idx.SoundSpeed()];

/*--- Apply production term modifications ---*/
switch (sstParsedOptions.production) {
Expand All @@ -795,6 +798,20 @@ class CSourcePieceWise_TurbSST final : public CNumerics {
P_Base = sqrt(StrainMag_i*VorticityMag);
break;

case SST_OPTIONS::COMP_Wilcox:
P_Base = StrainMag_i;
if (Mt >= 0.25) {
zetaFMt = 2.0 * (Mt * Mt - 0.25 * 0.25);
}
break;

case SST_OPTIONS::COMP_Sarkar:
P_Base = StrainMag_i;
if (Mt >= 0.25) {
zetaFMt = 0.5 * (Mt * Mt);
}
break;

default:
/*--- Base production term for SST-1994 and SST-2003 ---*/
P_Base = StrainMag_i;
Expand Down Expand Up @@ -831,10 +848,15 @@ class CSourcePieceWise_TurbSST final : public CNumerics {
pw = max(pw, sust_w);
}

if (sstParsedOptions.production == SST_OPTIONS::COMP_Sarkar) {
const su2double Dilatation_Sarkar = -0.15 * pk * Mt + 0.2 * beta_star * (1.0 +zetaFMt) * Density_i * ScalarVar_i[1] * ScalarVar_i[0] * Mt * Mt;
pk += Dilatation_Sarkar;
}

/*--- Dissipation ---*/

su2double dk = beta_star * Density_i * ScalarVar_i[1] * ScalarVar_i[0];
su2double dw = beta_blended * Density_i * ScalarVar_i[1] * ScalarVar_i[1];
su2double dk = beta_star * Density_i * ScalarVar_i[1] * ScalarVar_i[0] * (1.0 + zetaFMt);
su2double dw = beta_blended * Density_i * ScalarVar_i[1] * ScalarVar_i[1] * (1.0 - 0.09/beta_blended * zetaFMt);

/*--- LM model coupling with production and dissipation term for k transport equation---*/
if (config->GetKind_Trans_Model() == TURB_TRANS_MODEL::LM) {
Expand Down Expand Up @@ -862,10 +884,10 @@ class CSourcePieceWise_TurbSST final : public CNumerics {

/*--- Implicit part ---*/

Jacobian_i[0][0] = -beta_star * ScalarVar_i[1] * Volume;
Jacobian_i[0][1] = -beta_star * ScalarVar_i[0] * Volume;
Jacobian_i[0][0] = -beta_star * ScalarVar_i[1] * Volume * (1.0 + zetaFMt);
Jacobian_i[0][1] = -beta_star * ScalarVar_i[0] * Volume * (1.0 + zetaFMt);
Jacobian_i[1][0] = 0.0;
Jacobian_i[1][1] = -2.0 * beta_blended * ScalarVar_i[1] * Volume;
Jacobian_i[1][1] = -2.0 * beta_blended * ScalarVar_i[1] * Volume * (1.0 - 0.09/beta_blended * zetaFMt);
}

AD::SetPreaccOut(Residual, nVar);
Expand Down
16 changes: 16 additions & 0 deletions TestCases/parallel_regression.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -392,6 +392,22 @@ def main():
turb_flatplate_species.test_vals = [-4.147727, -0.634899, -1.770894, 1.334987, -3.250340, 9.000000, -6.700853, 5.000000, -6.991055, 10.000000, -6.033829, 0.996033, 0.996033]
test_list.append(turb_flatplate_species)

# Flat plate SST compressibility correction Wilcox
turb_flatplate_CC_Wilcox = TestCase('turb_flatplate_CC_Wilcox')
turb_flatplate_CC_Wilcox.cfg_dir = "rans/flatplate"
turb_flatplate_CC_Wilcox.cfg_file = "turb_SST_flatplate_compressibility_Wilcox.cfg"
turb_flatplate_CC_Wilcox.test_iter = 20
turb_flatplate_CC_Wilcox.test_vals = [-1.280847, 1.974242, 1.440510, 5.038429, -4.051126, 8.520857]
test_list.append(turb_flatplate_CC_Wilcox)

# Flat plate SST compressibility correction Sarkar
turb_flatplate_CC_Sarkar = TestCase('turb_flatplate_CC_Sarkar')
turb_flatplate_CC_Sarkar.cfg_dir = "rans/flatplate"
turb_flatplate_CC_Sarkar.cfg_file = "turb_SST_flatplate_compressibility_Sarkar.cfg"
turb_flatplate_CC_Sarkar.test_iter = 20
turb_flatplate_CC_Sarkar.test_vals = [-1.280847, 1.974242, 1.440510, 5.038429, -4.051129, 8.520857]
test_list.append(turb_flatplate_CC_Sarkar)

# ONERA M6 Wing
turb_oneram6 = TestCase('turb_oneram6')
turb_oneram6.cfg_dir = "rans/oneram6"
Expand Down
150 changes: 150 additions & 0 deletions TestCases/rans/flatplate/turb_SST_flatplate_compressibility_Sarkar.cfg
View file
Original file line number Diff line number Diff line change
@@ -0,0 +1,150 @@
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% %
% SU2 configuration file %
% Case description: Turbulent flow over flat plate(2DZPH) %
% Author: Sunoh. Kang %
% Institution: Pusan National University %
% Date: 2024.04.30 %
% File Version 8.0.1 "Harrier" %
% %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% ------------- DIRECT, ADJOINT, AND LINEARIZED PROBLEM DEFINITION ------------%
%
SOLVER= RANS
KIND_TURB_MODEL= SST
SST_OPTIONS= V2003m, COMPRESSIBILITY-SARKAR
MATH_PROBLEM= DIRECT
RESTART_SOL= NO

% ------------------------ COMPRESSIBLE FREE-STREAM DEFINITION ----------------%
%
MACH_NUMBER= 5.0
AOA= 0.000
INIT_OPTION= TD_CONDITIONS
FREESTREAM_OPTION= TEMPERATURE_FS
FREESTREAM_TEMPERATURE= 300.3333
FREESTREAM_PRESSURE= 13753.90558
FREESTREAM_TURBULENCEINTENSITY = 0.00002
FREESTREAM_TURB2LAMVISCRATIO = 0.009
REF_DIMENSIONALIZATION= DIMENSIONAL
REYNOLDS_NUMBER= 15E6

% ---- IDEAL GAS, POLYTPOPIC, VAN DER WALLS AND PENG ROBINSOPN CNSTANTS ------ %
%
FLUID_MODEL= IDEAL_GAS
GAMMA_VALUE= 1.4
GAS_CONSTANT= 287.058
ACENTRIC_FACTOR= 0.035
SPECIFIC_HEAT_CP= 1004.703
THERMAL_EXPANSION_COEFF= 0.00347
MOLECULAR_WEIGHT= 28.96

% --------------------------- VISCOSITY MODEL ---------------------------------%
%

VISCOSITY_MODEL= SUTHERLAND
MU_REF= 1.716E-5
MU_T_REF= 273.15
SUTHERLAND_CONSTANT= 110.4

% ------------------------ THERMAL CONDUCTIVITY MODEL -------------------------%
%

CONDUCTIVITY_MODEL= CONSTANT_PRANDTL
PRANDTL_LAM= 0.72
TURBULENT_CONDUCTIVITY_MODEL= CONSTANT_PRANDTL_TURB
PRANDTL_TURB= 0.9

% ---------------------- REFERENCE VALUE DEFINITION ---------------------------%
%

REF_ORIGIN_MOMENT_X = 0.00
REF_ORIGIN_MOMENT_Y = 0.00
REF_ORIGIN_MOMENT_Z = 0.00
REF_LENGTH= 1.0
REF_AREA= 1.0

% -------------------- BOUNDARY CONDITION DEFINITION --------------------------%
%

MARKER_PLOTTING = ( wall )
MARKER_SUPERSONIC_INLET= ( inlet, 300.3333 , 13753.90558 , 1737.082226 , 0.0 , 0.0 )
MARKER_FAR= ( farfield )
MARKER_SYM= ( symmetry )
MARKER_SUPERSONIC_OUTLET= (outlet, 34384.76396)
MARKER_ISOTHERMAL= ( wall, 327.36297)

% ------------- COMMON PARAMETERS DEFINING THE NUMERICAL METHOD ---------------%
%

NUM_METHOD_GRAD= WEIGHTED_LEAST_SQUARES
CFL_NUMBER= 0.5
CFL_ADAPT= YES
CFL_ADAPT_PARAM= ( 0.1, 1.5, 0.5, 50.0 )

% --------- SLOPE LIMITER DEFINITION, DISSIPATION SENSOR DEFINITION -----------------%
%

MUSCL_FLOW= YES
SLOPE_LIMITER_FLOW= VENKATAKRISHNAN_WANG
SLOPE_LIMITER_TURB= NONE
VENKAT_LIMITER_COEFF= 0.05
JST_SENSOR_COEFF= (0.5, 0.02)

% ------------------------ LINEAR SOLVER DEFINITION ---------------------------%
%

LINEAR_SOLVER= FGMRES
LINEAR_SOLVER_PREC= LU_SGS
LINEAR_SOLVER_ERROR= 1E-6
LINEAR_SOLVER_ITER= 5

% -------------------- FLOW NUMERICAL METHOD DEFINITION -----------------------%
%

CONV_NUM_METHOD_FLOW= AUSM
TIME_DISCRE_FLOW= EULER_IMPLICIT

% -------------------- TURBULENT NUMERICAL METHOD DEFINITION ------------------%
%

CONV_NUM_METHOD_TURB= SCALAR_UPWIND
TIME_DISCRE_TURB= EULER_IMPLICIT

% ------------------------- SCREEN/HISTORY VOLUME OUTPUT --------------------------%
%

SCREEN_OUTPUT= (INNER_ITER, RMS_DENSITY, RMS_MOMENTUM-X, RMS_MOMENTUM-Y, RMS_ENERGY, RMS_TKE, RMS_DISSIPATION)
SCREEN_WRT_FREQ_INNER= 1
SCREEN_WRT_FREQ_OUTER= 1
SCREEN_WRT_FREQ_TIME= 1
HISTORY_WRT_FREQ_INNER= 1
HISTORY_WRT_FREQ_OUTER= 1
HISTORY_WRT_FREQ_TIME= 1
OUTPUT_WRT_FREQ= 1000

% ------------------------- INPUT/OUTPUT FILE INFORMATION --------------------------%
%

MESH_FILENAME= mesh_flatplate_turb_137x97.su2
MESH_FORMAT= SU2
MESH_OUT_FILENAME= mesh_out.su2
SOLUTION_FILENAME= restart_flow.dat
SOLUTION_ADJ_FILENAME= solution_adj.dat
OUTPUT_FILES= (RESTART, PARAVIEW, SURFACE_PARAVIEW, TECPLOT, SURFACE_TECPLOT, CSV, SURFACE_CSV )
CONV_FILENAME= history
BREAKDOWN_FILENAME= forces_breakdown.dat
RESTART_FILENAME= restart_flow.dat
RESTART_ADJ_FILENAME= restart_adj.dat

% ------------------------------- SOLVER CONTROL ------------------------------%
%

ITER= 100
OUTER_ITER= 1
TIME_ITER= 100
CONV_RESIDUAL_MINVAL= -14
CONV_STARTITER= 10
CONV_CAUCHY_ELEMS= 100
CONV_CAUCHY_EPS= 1E-10
Loading

0 comments on commit 1044ccc

Please sign in to comment.