Skip to content

Commit 194692c

Browse files
TobiKattmannTobiKattmann
committed
Const-ing vars in SetTau_Wall_WF func for NS and IncNS.
Plus some other minor cosmetic things.
1 parent dde0a9f commit 194692c

2 files changed

Lines changed: 57 additions & 58 deletions

File tree

‎SU2_CFD/src/solvers/CIncNSSolver.cpp‎

Lines changed: 31 additions & 31 deletions
Original file line numberDiff line numberDiff line change
@@ -643,21 +643,17 @@ void CIncNSSolver::BC_ConjugateHeat_Interface(CGeometry *geometry, CSolver **sol
643643
}
644644

645645
void CIncNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_container, const CConfig *config) {
646-
/*---
647-
The wall function implemented herein is based on Nichols and Nelson, AIAA J. v32 n6 2004.
648-
---*/
646+
/*--- The wall function implemented herein is based on Nichols and Nelson, AIAA J. v32 n6 2004. ---*/
649647

650-
unsigned long notConvergedCounter = 0; /*--- counts the number of wall cells that are not converged ---*/
651-
unsigned long smallYPlusCounter = 0; /*--- counts the number of wall cells where y+ < 5 ---*/
648+
unsigned long notConvergedCounter = 0; /*--- Counts the number of wall cells that are not converged ---*/
649+
unsigned long smallYPlusCounter = 0; /*--- Counts the number of wall cells where y+ < 5 ---*/
652650

653651
const su2double Gas_Constant = config->GetGas_ConstantND();
654652
const su2double Cp = (Gamma / Gamma_Minus_One) * Gas_Constant;
655-
const su2double tol = 1e-12; /*!< \brief convergence criterium for the Newton solver, note that 1e-10 is too large */
656-
const unsigned short max_iter = config->GetwallModel_MaxIter();
653+
const auto max_iter = config->GetwallModel_MaxIter();
657654
const su2double relax = config->GetwallModel_RelFac();
658655

659-
/*--- Compute the recovery factor
660-
* use Molecular (Laminar) Prandtl number (see Nichols & Nelson, nomenclature ) ---*/
656+
/*--- Compute the recovery factor use Molecular (Laminar) Prandtl number (see Nichols & Nelson, nomenclature ) ---*/
661657

662658
const su2double Recovery = pow(config->GetPrandtl_Lam(), (1.0/3.0));
663659

@@ -701,7 +697,7 @@ void CIncNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_containe
701697

702698
const auto Normal = geometry->vertex[iMarker][iVertex]->GetNormal();
703699

704-
su2double Area = GeometryToolbox::Norm(nDim, Normal);
700+
const su2double Area = GeometryToolbox::Norm(nDim, Normal);
705701

706702
su2double UnitNormal[MAXNDIM] = {0.0};
707703
for (auto iDim = 0u; iDim < nDim; iDim++)
@@ -716,13 +712,13 @@ void CIncNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_containe
716712

717713
/*--- Compute the wall-parallel velocity at first point off the wall ---*/
718714

719-
su2double VelNormal = GeometryToolbox::DotProduct(int(MAXNDIM), Vel, UnitNormal);
715+
const su2double VelNormal = GeometryToolbox::DotProduct(int(MAXNDIM), Vel, UnitNormal);
720716

721717
su2double VelTang[MAXNDIM] = {0.0};
722718
for (auto iDim = 0u; iDim < nDim; iDim++)
723719
VelTang[iDim] = Vel[iDim] - VelNormal*UnitNormal[iDim];
724720

725-
su2double VelTangMod = GeometryToolbox::Norm(int(MAXNDIM), VelTang);
721+
const su2double VelTangMod = GeometryToolbox::Norm(int(MAXNDIM), VelTang);
726722

727723
/*--- Compute normal distance of the interior point from the wall ---*/
728724

@@ -735,9 +731,9 @@ void CIncNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_containe
735731
* In incompressible flows, we can assume that there is no velocity-related
736732
* temperature change Prandtl: T+ = Pr*y+ ---*/
737733
su2double T_Wall = nodes->GetTemperature(iPoint);
738-
su2double Conductivity_Wall = nodes->GetThermalConductivity(iPoint);
734+
const su2double Conductivity_Wall = nodes->GetThermalConductivity(iPoint);
739735

740-
/*--- if a wall temperature was given, we compute the local heat flux using k*dT/dn ---*/
736+
/*--- If a wall temperature was given, we compute the local heat flux using k*dT/dn ---*/
741737

742738
su2double q_w = 0.0;
743739

@@ -746,32 +742,33 @@ void CIncNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_containe
746742
q_w = Conductivity_Wall * (T_Wall - T_n) / (WallDistMod);
747743
}
748744
else if (config->GetMarker_All_KindBC(iMarker) == HEAT_FLUX) {
749-
q_w = config->GetWall_HeatFlux(Marker_Tag)/config->GetHeat_Flux_Ref();
745+
q_w = config->GetWall_HeatFlux(Marker_Tag) / config->GetHeat_Flux_Ref();
750746
}
751747
else if (config->GetMarker_All_KindBC(iMarker) == HEAT_TRANSFER) {
752-
const su2double Transfer_Coefficient = config->GetWall_HeatTransfer_Coefficient(Marker_Tag) * config->GetTemperature_Ref()/config->GetHeat_Flux_Ref();
753-
const su2double Tinfinity = config->GetWall_HeatTransfer_Temperature(Marker_Tag)/config->GetTemperature_Ref();
748+
const su2double Transfer_Coefficient = config->GetWall_HeatTransfer_Coefficient(Marker_Tag) *
749+
config->GetTemperature_Ref()/config->GetHeat_Flux_Ref();
750+
const su2double Tinfinity = config->GetWall_HeatTransfer_Temperature(Marker_Tag) / config->GetTemperature_Ref();
754751
q_w = Transfer_Coefficient * (Tinfinity - T_Wall);
755752
}
756753

757-
/*--- incompressible formulation ---*/
754+
/*--- Incompressible formulation ---*/
758755

759756
su2double Density_Wall = nodes->GetDensity(iPoint);
760-
su2double Lam_Visc_Normal = nodes->GetLaminarViscosity(Point_Normal);
757+
const su2double Lam_Visc_Normal = nodes->GetLaminarViscosity(Point_Normal);
761758

762759
/*--- Compute the shear stress at the wall in the regular fashion
763760
* by using the stress tensor on the surface ---*/
764761

765762
su2double tau[MAXNDIM][MAXNDIM] = {{0.0}};
766-
su2double Lam_Visc_Wall = nodes->GetLaminarViscosity(iPoint);
763+
const su2double Lam_Visc_Wall = nodes->GetLaminarViscosity(iPoint);
767764
su2double Eddy_Visc_Wall = nodes->GetEddyViscosity(iPoint);
768765

769766
CNumerics::ComputeStressTensor(nDim, tau, nodes->GetGradient_Primitive(iPoint)+1, Lam_Visc_Wall);
770767

771768
su2double TauTangent[MAXNDIM] = {0.0};
772769
GeometryToolbox::TangentProjection(nDim, tau, UnitNormal, TauTangent);
773770

774-
su2double WallShearStress = GeometryToolbox::Norm(int(MAXNDIM), TauTangent);
771+
const su2double WallShearStress = GeometryToolbox::Norm(int(MAXNDIM), TauTangent);
775772

776773
/*--- Calculate the quantities from boundary layer theory and
777774
* iteratively solve for a new wall shear stress. Use the current wall
@@ -782,37 +779,40 @@ void CIncNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_containe
782779
su2double U_Tau = max(1.0e-6,sqrt(WallShearStress/Density_Wall));
783780
su2double Y_Plus = 0.99*config->GetwallModel_MinYPlus(); // use clipping value as minimum
784781

785-
su2double Y_Plus_Start = Density_Wall * U_Tau * WallDistMod / Lam_Visc_Wall;
782+
const su2double Y_Plus_Start = Density_Wall * U_Tau * WallDistMod / Lam_Visc_Wall;
786783

787784
/*--- Automatic switch off when y+ < "limit" according to Nichols & Nelson (2004) ---*/
788785

789786
if (Y_Plus_Start < config->GetwallModel_MinYPlus()) {
790787
smallYPlusCounter++;
791788
continue;
792789
}
790+
791+
/*--- Convergence criterium for the Newton solver, note that 1e-10 is too large ---*/
792+
const su2double tol = 1e-12;
793793
while (fabs(diff) > tol) {
794794

795795
/*--- Friction velocity and u+ ---*/
796796

797-
su2double U_Plus = VelTangMod/U_Tau;
797+
const su2double U_Plus = VelTangMod/U_Tau;
798798

799799
/*--- Gamma, Beta, Q, and Phi, defined by Nichols & Nelson (2004) page 1110 ---*/
800800

801-
su2double Gam = Recovery*U_Tau*U_Tau/(2.0*Cp*T_Wall);
802-
su2double Beta = q_w*Lam_Visc_Wall/(Density_Wall*T_Wall*Conductivity_Wall*U_Tau);
803-
su2double Q = sqrt(Beta*Beta + 4.0*Gam);
804-
su2double Phi = asin(-1.0*Beta/Q);
801+
const su2double Gam = Recovery*U_Tau*U_Tau/(2.0*Cp*T_Wall);
802+
const su2double Beta = q_w*Lam_Visc_Wall/(Density_Wall*T_Wall*Conductivity_Wall*U_Tau);
803+
const su2double Q = sqrt(Beta*Beta + 4.0*Gam);
804+
const su2double Phi = asin(-1.0*Beta/Q);
805805

806806
/*--- Y+ defined by White & Christoph (compressibility and heat transfer) negative value for (2.0*Gam*U_Plus - Beta)/Q ---*/
807807

808-
su2double Y_Plus_White = exp((kappa/sqrt(Gam))*(asin((2.0*Gam*U_Plus - Beta)/Q) - Phi))*exp(-1.0*kappa*B);
808+
const su2double Y_Plus_White = exp((kappa/sqrt(Gam))*(asin((2.0*Gam*U_Plus - Beta)/Q) - Phi))*exp(-1.0*kappa*B);
809809

810810
/*--- Spalding's universal form for the BL velocity with the
811811
* outer velocity form of White & Christoph above. ---*/
812812
const su2double kUp = kappa*U_Plus;
813813
Y_Plus = U_Plus + Y_Plus_White - (exp(-1.0*kappa*B)* (1.0 + kUp + 0.5*kUp*kUp + kUp*kUp*kUp/6.0));
814814

815-
su2double dypw_dyp = 2.0*Y_Plus_White*(kappa*sqrt(Gam)/Q)*sqrt(1.0 - pow(2.0*Gam*U_Plus - Beta,2.0)/(Q*Q));
815+
const su2double dypw_dyp = 2.0*Y_Plus_White*(kappa*sqrt(Gam)/Q)*sqrt(1.0 - pow(2.0*Gam*U_Plus - Beta,2.0)/(Q*Q));
816816

817817
Eddy_Visc_Wall = Lam_Visc_Wall*(1.0 + dypw_dyp - kappa*exp(-1.0*kappa*B)*
818818
(1.0 + kappa*U_Plus + kappa*kappa*U_Plus*U_Plus/2.0)
@@ -825,7 +825,7 @@ void CIncNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_containe
825825

826826
/* --- Gradient of function defined above --- */
827827

828-
su2double grad_diff = Density_Wall * WallDistMod / Lam_Visc_Wall + VelTangMod / (U_Tau * U_Tau) +
828+
const su2double grad_diff = Density_Wall * WallDistMod / Lam_Visc_Wall + VelTangMod / (U_Tau * U_Tau) +
829829
kappa /(U_Tau * sqrt(Gam)) * asin(U_Plus * sqrt(Gam)) * Y_Plus_White -
830830
exp(-1.0 * B * kappa) * (0.5 * pow(VelTangMod * kappa / U_Tau, 3) +
831831
pow(VelTangMod * kappa / U_Tau, 2) + VelTangMod * kappa / U_Tau) / U_Tau;
@@ -852,7 +852,7 @@ void CIncNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_containe
852852
EddyViscWall[iMarker][iVertex] = Eddy_Visc_Wall;
853853
UTau[iMarker][iVertex] = U_Tau;
854854

855-
su2double Tau_Wall = (1.0/Density_Wall)*pow(Y_Plus*Lam_Visc_Wall/WallDistMod,2.0);
855+
const su2double Tau_Wall = (1.0/Density_Wall)*pow(Y_Plus*Lam_Visc_Wall/WallDistMod,2.0);
856856

857857
/*--- Store this value for the wall shear stress at the node. ---*/
858858
nodes->SetTau_Wall(iPoint, Tau_Wall);

‎SU2_CFD/src/solvers/CNSSolver.cpp‎

Lines changed: 26 additions & 27 deletions
Original file line numberDiff line numberDiff line change
@@ -843,7 +843,7 @@ void CNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_container,
843843

844844
const auto Normal = geometry->vertex[iMarker][iVertex]->GetNormal();
845845

846-
su2double Area = GeometryToolbox::Norm(nDim, Normal);
846+
const su2double Area = GeometryToolbox::Norm(nDim, Normal);
847847

848848
su2double UnitNormal[MAXNDIM] = {0.0};
849849
for (auto iDim = 0u; iDim < nDim; iDim++)
@@ -856,60 +856,59 @@ void CNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_container,
856856
for (auto iDim = 0u; iDim < nDim; iDim++)
857857
Vel[iDim] = nodes->GetVelocity(Point_Normal,iDim);
858858

859-
su2double P_Normal = nodes->GetPressure(Point_Normal);
860-
su2double T_Normal = nodes->GetTemperature(Point_Normal);
859+
const su2double P_Normal = nodes->GetPressure(Point_Normal);
860+
const su2double T_Normal = nodes->GetTemperature(Point_Normal);
861861

862862
/*--- Compute the wall-parallel velocity at first point off the wall ---*/
863863

864-
su2double VelNormal = GeometryToolbox::DotProduct(int(MAXNDIM), Vel, UnitNormal);
864+
const su2double VelNormal = GeometryToolbox::DotProduct(int(MAXNDIM), Vel, UnitNormal);
865865

866866
su2double VelTang[MAXNDIM] = {0.0};
867867
for (auto iDim = 0u; iDim < nDim; iDim++)
868868
VelTang[iDim] = Vel[iDim] - VelNormal*UnitNormal[iDim];
869869

870-
su2double VelTangMod = GeometryToolbox::Norm(int(MAXNDIM), VelTang);
870+
const su2double VelTangMod = GeometryToolbox::Norm(int(MAXNDIM), VelTang);
871871

872872
/*--- Compute normal distance of the interior point from the wall ---*/
873873

874874
su2double WallDist[MAXNDIM] = {0.0};
875875
GeometryToolbox::Distance(nDim, Coord, Coord_Normal, WallDist);
876876

877-
su2double WallDistMod = GeometryToolbox::Norm(int(MAXNDIM), WallDist);
877+
const su2double WallDistMod = GeometryToolbox::Norm(int(MAXNDIM), WallDist);
878878

879-
/*--- initial value for wall temperature ---*/
879+
/*--- Initial value for wall temperature ---*/
880880

881-
882881
su2double q_w = 0.0;
883882

884883
if (config->GetMarker_All_KindBC(iMarker) == HEAT_FLUX) {
885-
q_w = config->GetWall_HeatFlux(Marker_Tag)/config->GetHeat_Flux_Ref();
884+
q_w = config->GetWall_HeatFlux(Marker_Tag) / config->GetHeat_Flux_Ref();
886885
}
887886

888887
/*--- Extrapolate the pressure from the interior & compute the
889888
wall density using the equation of state ---*/
890889

891-
/*--- compressible formulation ---*/
890+
/*--- Compressible formulation ---*/
892891

893892
su2double T_Wall = nodes->GetTemperature(iPoint);
894-
su2double P_Wall = P_Normal;
895-
su2double Density_Wall = P_Wall/(Gas_Constant*T_Wall);
896-
su2double Lam_Visc_Normal = nodes->GetLaminarViscosity(Point_Normal);
897-
su2double Conductivity_Wall = nodes->GetThermalConductivity(iPoint);
893+
const su2double P_Wall = P_Normal;
894+
su2double Density_Wall = P_Wall / (Gas_Constant * T_Wall);
895+
const su2double Lam_Visc_Normal = nodes->GetLaminarViscosity(Point_Normal);
896+
const su2double Conductivity_Wall = nodes->GetThermalConductivity(iPoint);
898897

899898

900899
/*--- Compute the shear stress at the wall in the regular fashion
901900
* by using the stress tensor on the surface ---*/
902901

903902
su2double tau[MAXNDIM][MAXNDIM] = {{0.0}};
904-
su2double Lam_Visc_Wall = nodes->GetLaminarViscosity(iPoint);
903+
const su2double Lam_Visc_Wall = nodes->GetLaminarViscosity(iPoint);
905904
su2double Eddy_Visc_Wall = nodes->GetEddyViscosity(iPoint);
906905

907906
CNumerics::ComputeStressTensor(nDim, tau, nodes->GetGradient_Primitive(iPoint)+1, Lam_Visc_Wall);
908907

909908
su2double TauTangent[MAXNDIM] = {0.0};
910909
GeometryToolbox::TangentProjection(nDim, tau, UnitNormal, TauTangent);
911910

912-
su2double WallShearStress = GeometryToolbox::Norm(int(MAXNDIM), TauTangent);
911+
const su2double WallShearStress = GeometryToolbox::Norm(int(MAXNDIM), TauTangent);
913912

914913
/*--- Calculate the quantities from boundary layer theory and
915914
* iteratively solve for a new wall shear stress. Use the current wall
@@ -920,7 +919,7 @@ void CNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_container,
920919
su2double U_Tau = max(1.0e-6,sqrt(WallShearStress/Density_Wall));
921920
su2double Y_Plus = 0.99*config->GetwallModel_MinYPlus(); // use clipping value as minimum
922921

923-
su2double Y_Plus_Start = Density_Wall * U_Tau * WallDistMod / Lam_Visc_Wall;
922+
const su2double Y_Plus_Start = Density_Wall * U_Tau * WallDistMod / Lam_Visc_Wall;
924923

925924
/*--- Automatic switch off when y+ < "limit" according to Nichols & Nelson (2004) ---*/
926925

@@ -932,20 +931,20 @@ void CNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_container,
932931

933932
/*--- Friction velocity and u+ ---*/
934933

935-
su2double U_Plus = VelTangMod/U_Tau;
934+
const su2double U_Plus = VelTangMod/U_Tau;
936935

937936
/*--- Gamma, Beta, Q, and Phi, defined by Nichols & Nelson (2004) page 1110 ---*/
938937

939-
su2double Gam = Recovery*U_Tau*U_Tau/(2.0*Cp*T_Wall);
940-
su2double Beta = q_w*Lam_Visc_Wall/(Density_Wall*T_Wall*Conductivity_Wall*U_Tau);
941-
su2double Q = sqrt(Beta*Beta + 4.0*Gam);
942-
su2double Phi = asin(-1.0*Beta/Q);
938+
const su2double Gam = Recovery*U_Tau*U_Tau/(2.0*Cp*T_Wall);
939+
const su2double Beta = q_w*Lam_Visc_Wall/(Density_Wall*T_Wall*Conductivity_Wall*U_Tau);
940+
const su2double Q = sqrt(Beta*Beta + 4.0*Gam);
941+
const su2double Phi = asin(-1.0*Beta/Q);
943942

944943
/*--- Crocco-Busemann equation for wall temperature (eq. 11 of Nichols and Nelson) ---*/
945944
/*--- update T_Wall due to aerodynamic heating, unless the wall is isothermal ---*/
946945

947946
if (config->GetMarker_All_KindBC(iMarker) != ISOTHERMAL) {
948-
su2double denum = (1.0 + Beta*U_Plus - Gam*U_Plus*U_Plus);
947+
const su2double denum = (1.0 + Beta*U_Plus - Gam*U_Plus*U_Plus);
949948
if (denum > EPS){
950949
T_Wall = T_Normal / denum;
951950
nodes->SetTemperature(iPoint,T_Wall);
@@ -960,14 +959,14 @@ void CNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_container,
960959

961960
/*--- Y+ defined by White & Christoph (compressibility and heat transfer) negative value for (2.0*Gam*U_Plus - Beta)/Q ---*/
962961

963-
su2double Y_Plus_White = exp((kappa/sqrt(Gam))*(asin((2.0*Gam*U_Plus - Beta)/Q) - Phi))*exp(-1.0*kappa*B);
962+
const su2double Y_Plus_White = exp((kappa/sqrt(Gam))*(asin((2.0*Gam*U_Plus - Beta)/Q) - Phi))*exp(-1.0*kappa*B);
964963

965964
/*--- Spalding's universal form for the BL velocity with the
966965
* outer velocity form of White & Christoph above. ---*/
967966
const su2double kUp = kappa*U_Plus;
968967
Y_Plus = U_Plus + Y_Plus_White - (exp(-1.0*kappa*B)* (1.0 + kUp + 0.5*kUp*kUp + kUp*kUp*kUp/6.0));
969968

970-
su2double dypw_dyp = 2.0*Y_Plus_White*(kappa*sqrt(Gam)/Q)*sqrt(1.0 - pow(2.0*Gam*U_Plus - Beta,2.0)/(Q*Q));
969+
const su2double dypw_dyp = 2.0*Y_Plus_White*(kappa*sqrt(Gam)/Q)*sqrt(1.0 - pow(2.0*Gam*U_Plus - Beta,2.0)/(Q*Q));
971970

972971
Eddy_Visc_Wall = Lam_Visc_Wall*(1.0 + dypw_dyp - kappa*exp(-1.0*kappa*B)*
973972
(1.0 + kappa*U_Plus + kappa*kappa*U_Plus*U_Plus/2.0)
@@ -980,7 +979,7 @@ void CNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_container,
980979

981980
/* --- Gradient of function defined above --- */
982981

983-
su2double grad_diff = Density_Wall * WallDistMod / Lam_Visc_Wall + VelTangMod / (U_Tau * U_Tau) +
982+
const su2double grad_diff = Density_Wall * WallDistMod / Lam_Visc_Wall + VelTangMod / (U_Tau * U_Tau) +
984983
kappa /(U_Tau * sqrt(Gam)) * asin(U_Plus * sqrt(Gam)) * Y_Plus_White -
985984
exp(-1.0 * B * kappa) * (0.5 * pow(VelTangMod * kappa / U_Tau, 3) +
986985
pow(VelTangMod * kappa / U_Tau, 2) + VelTangMod * kappa / U_Tau) / U_Tau;
@@ -1008,7 +1007,7 @@ void CNSSolver::SetTau_Wall_WF(CGeometry *geometry, CSolver **solver_container,
10081007
EddyViscWall[iMarker][iVertex] = Eddy_Visc_Wall;
10091008
UTau[iMarker][iVertex] = U_Tau;
10101009

1011-
su2double Tau_Wall = (1.0/Density_Wall)*pow(Y_Plus*Lam_Visc_Wall/WallDistMod,2.0);
1010+
const su2double Tau_Wall = (1.0/Density_Wall)*pow(Y_Plus*Lam_Visc_Wall/WallDistMod,2.0);
10121011

10131012
/*--- Store this value for the wall shear stress at the node. ---*/
10141013

0 commit comments

Comments
 (0)