From 83e080c60d2d3f628be884a320623b6319d2cc6c Mon Sep 17 00:00:00 2001 From: bigfooted Date: Tue, 16 Apr 2024 22:33:09 +0200 Subject: [PATCH 01/11] fix corner node --- SU2_CFD/src/solvers/CEulerSolver.cpp | 17 +++++++++++++++++ 1 file changed, 17 insertions(+) diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index 16cdfe59a4b..9385c9c874c 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -7111,6 +7111,23 @@ void CEulerSolver::BC_Inlet(CGeometry *geometry, CSolver **solver_container, break; } + /*--- Check if the inlet node is shared with a viscous wall. ---*/ + + if (geometry->nodes->GetViscousBoundary(iPoint)) { + + V_inlet[0] = nodes->GetTemperature(iPoint); + + /*--- Impose the wall velocity from the interior. ---*/ + for (iDim = 0; iDim < nDim; iDim++){ + V_inlet[iDim+1] = nodes->GetVelocity(iPoint,iDim); + } + + /*--- Match the pressure and energy at the wall. ---*/ + V_inlet[nDim+1] = nodes->GetPressure(iPoint); + V_inlet[nDim+3] = nodes->GetPressure(iPoint)/(Density*Gamma_Minus_One) + nodes->GetPressure(iPoint)/Density; // + Pressure/Density; + } + + /*--- Set various quantities in the solver class ---*/ conv_numerics->SetPrimitive(V_domain, V_inlet); From 8624e0df58e87694951b036375444a1e53eb095d Mon Sep 17 00:00:00 2001 From: bigfooted Date: Tue, 16 Apr 2024 23:38:05 +0200 Subject: [PATCH 02/11] fix corner node temp and dens --- SU2_CFD/src/solvers/CEulerSolver.cpp | 16 ++++++++++++---- 1 file changed, 12 insertions(+), 4 deletions(-) diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index 9385c9c874c..feaf719b308 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -7115,16 +7115,24 @@ void CEulerSolver::BC_Inlet(CGeometry *geometry, CSolver **solver_container, if (geometry->nodes->GetViscousBoundary(iPoint)) { - V_inlet[0] = nodes->GetTemperature(iPoint); - /*--- Impose the wall velocity from the interior. ---*/ + + Velocity2 = 0.0; for (iDim = 0; iDim < nDim; iDim++){ V_inlet[iDim+1] = nodes->GetVelocity(iPoint,iDim); + Velocity2 += V_inlet[iDim+1] * V_inlet[iDim+1]; } + Pressure = nodes->GetPressure(iPoint); + /*--- Match the pressure and energy at the wall. ---*/ - V_inlet[nDim+1] = nodes->GetPressure(iPoint); - V_inlet[nDim+3] = nodes->GetPressure(iPoint)/(Density*Gamma_Minus_One) + nodes->GetPressure(iPoint)/Density; // + Pressure/Density; + Density = Pressure/(Gas_Constant*Temperature); + Energy = Pressure/(Density*Gamma_Minus_One) + 0.5*Velocity2; + if (tkeNeeded) Energy += GetTke_Inf(); + + V_inlet[nDim+1] = Pressure; + V_inlet[nDim+2] = Density; + V_inlet[nDim+3] = Energy + Pressure/Density; } From f4f8a2bf74b2b69627f3ee2692e446a1acb0f440 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 17 Apr 2024 00:09:54 +0200 Subject: [PATCH 03/11] fix warning for temperature maybe used uninitialized --- SU2_CFD/src/solvers/CEulerSolver.cpp | 7 ++++--- 1 file changed, 4 insertions(+), 3 deletions(-) diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index feaf719b308..7bac11d2e4d 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -6880,10 +6880,11 @@ void CEulerSolver::BC_Inlet(CGeometry *geometry, CSolver **solver_container, CConfig *config, unsigned short val_marker) { unsigned short iDim; unsigned long iVertex, iPoint; - su2double P_Total, T_Total, Velocity[MAXNDIM], Velocity2, H_Total, Temperature, Riemann, + su2double P_Total, T_Total, Velocity[MAXNDIM], Velocity2, H_Total, Riemann, Pressure, Density, Energy, Flow_Dir[MAXNDIM], Mach2, SoundSpeed2, SoundSpeed_Total2, Vel_Mag, alpha, aa, bb, cc, dd, Area, UnitNormal[MAXNDIM], Normal[MAXNDIM]; su2double *V_inlet, *V_domain; + su2double Temperature = 288.15; const bool implicit = (config->GetKind_TimeIntScheme() == EULER_IMPLICIT); const su2double Two_Gamma_M1 = 2.0 / Gamma_Minus_One; @@ -7096,8 +7097,8 @@ void CEulerSolver::BC_Inlet(CGeometry *geometry, CSolver **solver_container, if (tkeNeeded) Energy += GetTke_Inf(); /*--- Primitive variables, using the derived quantities ---*/ - - V_inlet[0] = Pressure / ( Gas_Constant * Density); + Temperature = Pressure / ( Gas_Constant * Density); + V_inlet[0] = Temperature; for (iDim = 0; iDim < nDim; iDim++) V_inlet[iDim+1] = Vel_Mag*Flow_Dir[iDim]; V_inlet[nDim+1] = Pressure; From 4844ba99783a68a6f6789017f0023c3d22e540a6 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Wed, 17 Apr 2024 00:11:51 +0200 Subject: [PATCH 04/11] move comment --- SU2_CFD/src/solvers/CEulerSolver.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index 7bac11d2e4d..bf465b5cd14 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -7124,9 +7124,9 @@ void CEulerSolver::BC_Inlet(CGeometry *geometry, CSolver **solver_container, Velocity2 += V_inlet[iDim+1] * V_inlet[iDim+1]; } - Pressure = nodes->GetPressure(iPoint); + /*--- Match the pressure, density and energy at the wall. ---*/ - /*--- Match the pressure and energy at the wall. ---*/ + Pressure = nodes->GetPressure(iPoint); Density = Pressure/(Gas_Constant*Temperature); Energy = Pressure/(Density*Gamma_Minus_One) + 0.5*Velocity2; if (tkeNeeded) Energy += GetTke_Inf(); From 0cb2f38c5dc923836e193b9c7dcf1a525a51434c Mon Sep 17 00:00:00 2001 From: bigfooted Date: Sat, 20 Apr 2024 11:28:14 +0200 Subject: [PATCH 05/11] fix bounded scalar for compressible --- Common/src/CConfig.cpp | 6 +- SU2_CFD/include/solvers/CSpeciesSolver.hpp | 20 ++++ SU2_CFD/src/solvers/CEulerSolver.cpp | 114 ++++++++++++++------- SU2_CFD/src/solvers/CSpeciesSolver.cpp | 16 ++- TestCases/parallel_regression.py | 2 +- 5 files changed, 115 insertions(+), 43 deletions(-) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 18da76556c4..e4a9f562c7a 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -5538,16 +5538,12 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i /*--- Define some variables for flamelet model. ---*/ if (Kind_Species_Model == SPECIES_MODEL::FLAMELET) { /*--- The controlling variables are progress variable, total enthalpy, and optionally mixture fraction ---*/ - //n_control_vars = nSpecies - n_user_scalars; if (n_control_vars != (nSpecies - n_user_scalars)) - SU2_MPI::Error("Number of initial species incompatbile with number of controlling variables and user scalars.", CURRENT_FUNCTION); + SU2_MPI::Error("Number of initial species incompatible with number of controlling variables and user scalars.", CURRENT_FUNCTION); /*--- We can have additional user defined transported scalars ---*/ n_scalars = n_control_vars + n_user_scalars; } - if (Kind_Regime == ENUM_REGIME::COMPRESSIBLE && GetBounded_Scalar()) { - SU2_MPI::Error("BOUNDED_SCALAR discretization can only be used for incompressible problems.", CURRENT_FUNCTION); - } } void CConfig::SetMarkers(SU2_COMPONENT val_software) { diff --git a/SU2_CFD/include/solvers/CSpeciesSolver.hpp b/SU2_CFD/include/solvers/CSpeciesSolver.hpp index 1ee78d90286..717f855c59b 100644 --- a/SU2_CFD/include/solvers/CSpeciesSolver.hpp +++ b/SU2_CFD/include/solvers/CSpeciesSolver.hpp @@ -153,6 +153,26 @@ class CSpeciesSolver : public CScalarSolver { * (like BC_Isothermal_Wall) are implemented the respective CHeatSolver implementations can act as a starting * point ---*/ + /*! + * \brief Generic implementation of the isothermal wall also covering CHT cases, + * for which the wall temperature is given by GetConjugateHeatVariable. + */ + void BC_Isothermal_Wall_Generic(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, + CNumerics* visc_numerics, CConfig* config, unsigned short val_marker, + bool cht_mode = false); + + /*! + * \brief Impose the scalar wall boundary condition. + * \param[in] geometry - Geometrical definition of the problem. + * \param[in] solver_container - Container vector with all the solutions. + * \param[in] conv_numerics - Description of the numerical method. + * \param[in] visc_numerics - Description of the numerical method. + * \param[in] config - Definition of the particular problem. + * \param[in] val_marker - Surface marker where the boundary condition is applied. + */ + void BC_Isothermal_Wall(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, + CNumerics* visc_numerics, CConfig* config, unsigned short val_marker) override; + /*! * \brief Source term computation for axisymmetric flow. * \param[in] geometry - Geometrical definition of the problem. diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index bf465b5cd14..f49b19ef5a4 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -341,6 +341,10 @@ CEulerSolver::CEulerSolver(CGeometry *geometry, CConfig *config, CommunicateInitialState(geometry, config); + /*--- Sizing edge mass flux array ---*/ + if (config->GetBounded_Scalar()) + EdgeMassFluxes.resize(geometry->GetnEdge()) = su2double(0.0); + /*--- Add the solver name.. ---*/ SolverName = "C.FLOW"; @@ -1749,6 +1753,7 @@ void CEulerSolver::Centered_Residual(CGeometry *geometry, CSolver **solver_conta CConfig *config, unsigned short iMesh, unsigned short iRKStep) { EdgeFluxResidual(geometry, solver_container, config); + const bool bounded_scalar = config->GetBounded_Scalar(); } void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_container, @@ -1772,6 +1777,7 @@ void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_contain const bool muscl = (config->GetMUSCL_Flow() && (iMesh == MESH_0)); const bool limiter = (config->GetKind_SlopeLimit_Flow() != LIMITER::NONE); const bool van_albada = (config->GetKind_SlopeLimit_Flow() == LIMITER::VAN_ALBADA_EDGE); + const bool bounded_scalar = config->GetBounded_Scalar(); /*--- Non-physical counter. ---*/ unsigned long counter_local = 0; @@ -1937,6 +1943,8 @@ void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_contain auto residual = numerics->ComputeResidual(config); + if (bounded_scalar) EdgeMassFluxes[iEdge] = residual[0]; + /*--- Set the final value of the Roe dissipation coefficient ---*/ if ((kind_dissipation != NO_ROELOWDISS) && (MGLevel != MESH_0)) { @@ -5067,8 +5075,8 @@ void CEulerSolver::BC_Riemann(CGeometry *geometry, CSolver **solver_container, su2double Velocity_i[MAXNDIM]={0}; for (auto iDim=0u; iDim < nDim; iDim++) Velocity_i[iDim] = nodes->GetVelocity(iPoint,iDim); - - const auto Velocity2_i = GeometryToolbox::SquaredNorm(nDim, Velocity_i); + + const auto Velocity2_i = GeometryToolbox::SquaredNorm(nDim, Velocity_i); const auto Density_i = nodes->GetDensity(iPoint); @@ -5090,7 +5098,7 @@ void CEulerSolver::BC_Riemann(CGeometry *geometry, CSolver **solver_container, T_Total{0}, P_Total{0}, Density_e{0}, StaticEnthalpy_e{0}, StaticEnergy_e{0}; su2double Velocity2_e{0}, NormalVelocity{0}, TangVelocity{0}, VelMag_e{0}; - su2double Velocity_e[MAXNDIM] = {0}; + su2double Velocity_e[MAXNDIM] = {0}; const su2double * Flow_Dir, * Mach; switch(config->GetKind_Data_Riemann(Marker_Tag)) { @@ -5146,10 +5154,10 @@ void CEulerSolver::BC_Riemann(CGeometry *geometry, CSolver **solver_container, GetFluidModel()->SetTDState_PT(P_static, T_static); /* --- Compute the boundary state u_e --- */ - for (auto iDim = 0u; iDim < nDim; iDim++) + for (auto iDim = 0u; iDim < nDim; iDim++) Velocity_e[iDim] = Mach[iDim]*GetFluidModel()->GetSoundSpeed(); - Velocity2_e = GeometryToolbox::SquaredNorm(nDim, Velocity_e); + Velocity2_e = GeometryToolbox::SquaredNorm(nDim, Velocity_e); Density_e = GetFluidModel()->GetDensity(); StaticEnergy_e = GetFluidModel()->GetStaticEnergy(); Energy_e = StaticEnergy_e + 0.5 * Velocity2_e; @@ -5176,7 +5184,7 @@ void CEulerSolver::BC_Riemann(CGeometry *geometry, CSolver **solver_container, for (auto iDim = 0u; iDim < nDim; iDim++) Velocity_e[iDim] = Mach[iDim]*GetFluidModel()->GetSoundSpeed(); - Velocity2_e = GeometryToolbox::SquaredNorm(nDim, Velocity_e); + Velocity2_e = GeometryToolbox::SquaredNorm(nDim, Velocity_e); Density_e = GetFluidModel()->GetDensity(); StaticEnergy_e = GetFluidModel()->GetStaticEnergy(); Energy_e = StaticEnergy_e + 0.5 * Velocity2_e; @@ -5208,10 +5216,10 @@ void CEulerSolver::BC_Riemann(CGeometry *geometry, CSolver **solver_container, /* --- Compute the boundary state u_e --- */ GetFluidModel()->SetTDState_Prho(Pressure_e, Density_e); - for (auto iDim = 0u; iDim < nDim; iDim++) + for (auto iDim = 0u; iDim < nDim; iDim++) Velocity_e[iDim] = Velocity_i[iDim]; - Velocity2_e = GeometryToolbox::SquaredNorm(nDim, Velocity_e); + Velocity2_e = GeometryToolbox::SquaredNorm(nDim, Velocity_e); Energy_e = GetFluidModel()->GetStaticEnergy() + 0.5*Velocity2_e; break; @@ -5235,7 +5243,7 @@ void CEulerSolver::BC_Riemann(CGeometry *geometry, CSolver **solver_container, } /*--- Flow eigenvalues, boundary state u_e and u_i ---*/ - su2double Lambda_i[MAXNVAR] = {0}, u_e[MAXNVAR] = {0}, u_i[MAXNVAR]={0}, u_b[MAXNVAR]={0}, + su2double Lambda_i[MAXNVAR] = {0}, u_e[MAXNVAR] = {0}, u_i[MAXNVAR]={0}, u_b[MAXNVAR]={0}, dw[MAXNVAR]={0}; u_e[0] = Density_e; u_i[0] = Density_i; @@ -5302,7 +5310,7 @@ void CEulerSolver::BC_Riemann(CGeometry *geometry, CSolver **solver_container, LinSysRes.AddBlock(iPoint, Residual); if (implicit) { - su2double **Jacobian_b = new su2double*[nVar], + su2double **Jacobian_b = new su2double*[nVar], **Jacobian_i = new su2double*[nVar]; su2double **DubDu = new su2double*[nVar]; for (auto iVar = 0u; iVar < nVar; iVar++){ @@ -5333,7 +5341,7 @@ void CEulerSolver::BC_Riemann(CGeometry *geometry, CSolver **solver_container, if (dynamic_grid){ const auto gridVel = geometry->nodes->GetGridVel(iPoint); const auto projVelocity = GeometryToolbox::DotProduct(nDim, gridVel, Normal); - for (auto iVar = 0u; iVar < nVar; iVar++) + for (auto iVar = 0u; iVar < nVar; iVar++) Jacobian_b[iVar][iVar] -= projVelocity; } @@ -5394,7 +5402,7 @@ void CEulerSolver::BC_Riemann(CGeometry *geometry, CSolver **solver_container, visc_numerics->SetPrimVarGradient(nodes->GetGradient_Primitive(iPoint), nodes->GetGradient_Primitive(iPoint)); /*--- Secondary variables ---*/ - + auto S_domain = nodes->GetSecondary(iPoint); /*--- Compute secondary thermodynamic properties (partial derivatives...) ---*/ @@ -6922,7 +6930,7 @@ void CEulerSolver::BC_Inlet(CGeometry *geometry, CSolver **solver_container, V_domain = nodes->GetPrimitive(iPoint); - /*--- Build the fictitious intlet state based on characteristics ---*/ + /*--- Build the fictitious inlet state based on characteristics ---*/ /*--- Subsonic inflow: there is one outgoing characteristic (u-c), @@ -7116,27 +7124,61 @@ void CEulerSolver::BC_Inlet(CGeometry *geometry, CSolver **solver_container, if (geometry->nodes->GetViscousBoundary(iPoint)) { - /*--- Impose the wall velocity from the interior. ---*/ + switch (Kind_Inlet) { - Velocity2 = 0.0; - for (iDim = 0; iDim < nDim; iDim++){ - V_inlet[iDim+1] = nodes->GetVelocity(iPoint,iDim); - Velocity2 += V_inlet[iDim+1] * V_inlet[iDim+1]; - } + /*--- Total properties have been specified at the inlet. ---*/ - /*--- Match the pressure, density and energy at the wall. ---*/ + case INLET_TYPE::TOTAL_CONDITIONS: { - Pressure = nodes->GetPressure(iPoint); - Density = Pressure/(Gas_Constant*Temperature); - Energy = Pressure/(Density*Gamma_Minus_One) + 0.5*Velocity2; - if (tkeNeeded) Energy += GetTke_Inf(); + /*--- Impose the wall velocity from the interior wall node. ---*/ - V_inlet[nDim+1] = Pressure; - V_inlet[nDim+2] = Density; - V_inlet[nDim+3] = Energy + Pressure/Density; - } + Velocity2 = 0.0; + for (iDim = 0; iDim < nDim; iDim++) { + V_inlet[iDim+1] = nodes->GetVelocity(iPoint,iDim); + Velocity2 += V_inlet[iDim+1] * V_inlet[iDim+1]; + } + + /*--- Match the pressure, density and energy at the wall. ---*/ + + Pressure = nodes->GetPressure(iPoint); + Density = Pressure / (Gas_Constant * Temperature); + Energy = Pressure / (Density*Gamma_Minus_One) + 0.5 * Velocity2; + if (tkeNeeded) Energy += GetTke_Inf(); + + V_inlet[nDim+1] = Pressure; + V_inlet[nDim+2] = Density; + V_inlet[nDim+3] = Energy + Pressure/Density; + break; + } + case INLET_TYPE::MASS_FLOW: { + /*--- Impose the wall velocity from the interior wall node. ---*/ + + Velocity2 = 0.0; + for (iDim = 0; iDim < nDim; iDim++) { + V_inlet[iDim+1] = nodes->GetVelocity(iPoint,iDim); + Velocity2 += V_inlet[iDim+1] * V_inlet[iDim+1]; + } + + Pressure = nodes->GetPressure(iPoint); + Density = nodes->GetDensity(iPoint); + + Energy = Pressure / (Density * Gamma_Minus_One) + 0.5 * Velocity2; + + if (tkeNeeded) Energy += GetTke_Inf(); + + V_inlet[nDim+3] = Energy + Pressure/Density; + + break; + } + + default: + SU2_MPI::Error("Unsupported INLET_TYPE.", CURRENT_FUNCTION); + break; + + } + } /*--- Set various quantities in the solver class ---*/ conv_numerics->SetPrimitive(V_domain, V_inlet); @@ -8821,9 +8863,9 @@ void CEulerSolver::PreprocessAverage(CSolver **solver, CGeometry *geometry, CCon const auto nSpanWiseSections = config->GetnSpanWiseSections(); const auto iZone = config->GetiZone(); - + for (auto iSpan= 0u; iSpan < nSpanWiseSections; iSpan++){ - su2double TotalAreaVelocity[MAXNDIM]={0.0}, + su2double TotalAreaVelocity[MAXNDIM]={0.0}, TotalAreaPressure{0}, TotalAreaDensity{0}; for (auto iMarker = 0u; iMarker < config->GetnMarker_All(); iMarker++){ @@ -9013,7 +9055,7 @@ void CEulerSolver::TurboAverageProcess(CSolver **solver, CGeometry *geometry, CC TotalAreaVelocity[iDim] += Area*Velocity[iDim]; TotalMassVelocity[iDim] += Area*(Density*TurboVelocity[0] )*Velocity[iDim]; } - + TotalFluxes[0] += Area*(Density*TurboVelocity[0]); TotalFluxes[1] += Area*(Density*TurboVelocity[0]*TurboVelocity[0] + Pressure); for (auto iDim = 2; iDim < nDim+1; iDim++) @@ -9062,7 +9104,7 @@ void CEulerSolver::TurboAverageProcess(CSolver **solver, CGeometry *geometry, CC UpdateTotalQuantities(iMarker, jSpan, iVertex); } } - } + } } // marker_flag match } // iMarkerTP match @@ -9277,7 +9319,7 @@ void CEulerSolver::TurboAverageProcess(CSolver **solver, CGeometry *geometry, CC OmegaOut[iMarkerTP - 1][iSpan] = AverageOmega[iMarker][iSpan]; NuOut[iMarkerTP - 1][iSpan] = AverageNu[iMarker][iSpan]; } - + auto TurboVel = (marker_flag == INFLOW) ? TurboVelocityIn[iMarkerTP - 1][iSpan] : TurboVelocityOut[iMarkerTP - 1][iSpan]; if (performance_average_process == MIXEDOUT) { @@ -9307,7 +9349,7 @@ void CEulerSolver::TurboAverageProcess(CSolver **solver, CGeometry *geometry, CC if(config->GetKind_Data_Giles(Marker_Tag) == RADIAL_EQUILIBRIUM){ RadialEquilibriumPressure[iMarker][nSpanWiseSections/2] = config->GetGiles_Var1(Marker_Tag)/config->GetPressure_Ref(); } - } + } for (auto iSpan= nSpanWiseSections/2; iSpan < nSpanWiseSections-1; iSpan++){ const auto Radius2 = geometry->GetTurboRadius(iMarker,iSpan+1); const auto Radius1 = geometry->GetTurboRadius(iMarker,iSpan); @@ -9351,9 +9393,9 @@ void CEulerSolver::MixedOut_Average(CConfig *config, su2double val_init_pressure su2double vel[MAXNDIM] = {0}; vel[0] = (val_Averaged_Flux[1] - pressure_mix) / val_Averaged_Flux[0]; - for (auto iDim = 1u; iDim < nDim; iDim++) + for (auto iDim = 1u; iDim < nDim; iDim++) vel[iDim] = val_Averaged_Flux[iDim+1] / val_Averaged_Flux[0]; - + const su2double velsq = GeometryToolbox::DotProduct(nDim, vel, vel); diff --git a/SU2_CFD/src/solvers/CSpeciesSolver.cpp b/SU2_CFD/src/solvers/CSpeciesSolver.cpp index f5f5ed16603..144c2afc4a4 100644 --- a/SU2_CFD/src/solvers/CSpeciesSolver.cpp +++ b/SU2_CFD/src/solvers/CSpeciesSolver.cpp @@ -345,7 +345,7 @@ void CSpeciesSolver::BC_Inlet(CGeometry* geometry, CSolver** solver_container, C /*--- Identify the boundary by string name ---*/ string Marker_Tag = config->GetMarker_All_TagBound(val_marker); - + if (config->GetMarker_StrongBC(Marker_Tag)==true) { nodes->SetSolution_Old(iPoint, Inlet_SpeciesVars[val_marker][iVertex]); @@ -406,6 +406,20 @@ void CSpeciesSolver::BC_Inlet(CGeometry* geometry, CSolver** solver_container, C END_SU2_OMP_FOR } +void CSpeciesSolver::BC_Isothermal_Wall_Generic(CGeometry* geometry, CSolver** solver_container, CNumerics* conv_numerics, + CNumerics* visc_numerics, CConfig* config, unsigned short val_marker, + bool cht_mode) { + +} + +void CSpeciesSolver::BC_Isothermal_Wall(CGeometry *geometry, CSolver **solver_container, CNumerics *conv_numerics, + CNumerics *visc_numerics, CConfig *config, unsigned short val_marker) { + + BC_Isothermal_Wall_Generic(geometry, solver_container, conv_numerics, visc_numerics, config, val_marker); + +} + + void CSpeciesSolver::SetInletAtVertex(const su2double *val_inlet, unsigned short iMarker, unsigned long iVertex) { diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index 51b18e48741..567e9adc016 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -302,7 +302,7 @@ def main(): flatplate_udobj.cfg_dir = "user_defined_functions" flatplate_udobj.cfg_file = "lam_flatplate.cfg" flatplate_udobj.test_iter = 20 - flatplate_udobj.test_vals = [-6.653802, -1.181430, -0.794887, 0.000611, -0.000369, 0.000736, -0.001104, 596.690000, 299.800000, 296.890000, 21.492000, 0.563990, 37.148000, 2.278700] + flatplate_udobj.test_vals = [-6.653278, -1.180894, -0.794611, 0.000613, -0.000370, 0.000737, -0.001107, 596.680000, 299.790000, 296.890000, 21.488000, 0.564630, 37.153000, 2.284900] test_list.append(flatplate_udobj) # Laminar cylinder (steady) From 819b46ca14d2a3114f9957575bb6abcdb9294a36 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Sat, 20 Apr 2024 16:19:26 +0200 Subject: [PATCH 06/11] fix bounded scalar for compressible --- SU2_CFD/src/solvers/CEulerSolver.cpp | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index f49b19ef5a4..5e24a1ad039 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -1753,7 +1753,6 @@ void CEulerSolver::Centered_Residual(CGeometry *geometry, CSolver **solver_conta CConfig *config, unsigned short iMesh, unsigned short iRKStep) { EdgeFluxResidual(geometry, solver_container, config); - const bool bounded_scalar = config->GetBounded_Scalar(); } void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_container, @@ -1777,8 +1776,7 @@ void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_contain const bool muscl = (config->GetMUSCL_Flow() && (iMesh == MESH_0)); const bool limiter = (config->GetKind_SlopeLimit_Flow() != LIMITER::NONE); const bool van_albada = (config->GetKind_SlopeLimit_Flow() == LIMITER::VAN_ALBADA_EDGE); - const bool bounded_scalar = config->GetBounded_Scalar(); - + //const bool bounded_scalar = config->GetBounded_Scalar(); /*--- Non-physical counter. ---*/ unsigned long counter_local = 0; SU2_OMP_MASTER @@ -1943,7 +1941,9 @@ void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_contain auto residual = numerics->ComputeResidual(config); - if (bounded_scalar) EdgeMassFluxes[iEdge] = residual[0]; + // this does not seem to be necessary + cout << "bounded" << endl; + //if (bounded_scalar) EdgeMassFluxes[iEdge] = residual[0]; /*--- Set the final value of the Roe dissipation coefficient ---*/ From 7a5104d5627c001f2975b6e2e750cd8d796a58ea Mon Sep 17 00:00:00 2001 From: bigfooted Date: Sat, 10 Aug 2024 16:08:41 +0200 Subject: [PATCH 07/11] update regression tests --- TestCases/hybrid_regression.py | 2 +- TestCases/parallel_regression.py | 4 ++-- TestCases/parallel_regression_AD.py | 2 +- TestCases/serial_regression.py | 4 ++-- 4 files changed, 6 insertions(+), 6 deletions(-) diff --git a/TestCases/hybrid_regression.py b/TestCases/hybrid_regression.py index bd7635e1413..3b93ec68f68 100644 --- a/TestCases/hybrid_regression.py +++ b/TestCases/hybrid_regression.py @@ -139,7 +139,7 @@ def main(): poiseuille_profile.cfg_dir = "navierstokes/poiseuille" poiseuille_profile.cfg_file = "profile_poiseuille.cfg" poiseuille_profile.test_iter = 10 - poiseuille_profile.test_vals = [-12.494728, -7.712546, -0.000000, 2.085796] + poiseuille_profile.test_vals = [-12.494736, -7.712091, -0.000000, 2.085796] poiseuille_profile.test_vals_aarch64 = [-12.494717, -7.711274, -0.000000, 2.085796] test_list.append(poiseuille_profile) diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index 6ca9e503c0d..e2dd41e0c18 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -344,7 +344,7 @@ def main(): poiseuille_profile.cfg_dir = "navierstokes/poiseuille" poiseuille_profile.cfg_file = "profile_poiseuille.cfg" poiseuille_profile.test_iter = 10 - poiseuille_profile.test_vals = [-12.492939, -7.672950, -0.000000, 2.085796] + poiseuille_profile.test_vals = [-12.492865, -7.672526, -0.000000, 2.085796] poiseuille_profile.test_vals_aarch64 = [-12.492864, -7.671632, -0.000000, 2.085796] test_list.append(poiseuille_profile) @@ -1015,7 +1015,7 @@ def main(): flatplate_unsteady.cfg_dir = "navierstokes/flatplate" flatplate_unsteady.cfg_file = "lam_flatplate_unst.cfg" flatplate_unsteady.test_iter = 3 - flatplate_unsteady.test_vals = [7.9509e-06, -8.868859, -8.231652, -6.283262, -5.466675, -3.391163, 0.002078, -0.343642] + flatplate_unsteady.test_vals = [0.000008, -8.869318, -8.232471, -6.283031, -5.467151, -3.391524, 0.002078, -0.343471] flatplate_unsteady.unsteady = True test_list.append(flatplate_unsteady) diff --git a/TestCases/parallel_regression_AD.py b/TestCases/parallel_regression_AD.py index 59576fdc410..9c69ebf4792 100644 --- a/TestCases/parallel_regression_AD.py +++ b/TestCases/parallel_regression_AD.py @@ -524,7 +524,7 @@ def main(): pywrapper_wavy_wall_steady.cfg_dir = "py_wrapper/wavy_wall" pywrapper_wavy_wall_steady.cfg_file = "run_steady.py" pywrapper_wavy_wall_steady.test_iter = 100 - pywrapper_wavy_wall_steady.test_vals = [-1.360044, 2.580709, -2.892473] + pywrapper_wavy_wall_steady.test_vals = [-1.359409, 2.580816, -2.892697] pywrapper_wavy_wall_steady.command = TestCase.Command("mpirun -n 2", "python", "run_steady.py") pywrapper_wavy_wall_steady.timeout = 1600 pywrapper_wavy_wall_steady.tol = 0.00001 diff --git a/TestCases/serial_regression.py b/TestCases/serial_regression.py index fcdacb032a7..70e937959a9 100644 --- a/TestCases/serial_regression.py +++ b/TestCases/serial_regression.py @@ -198,7 +198,7 @@ def main(): poiseuille_profile.cfg_dir = "navierstokes/poiseuille" poiseuille_profile.cfg_file = "profile_poiseuille.cfg" poiseuille_profile.test_iter = 10 - poiseuille_profile.test_vals = [-12.494681, -7.711642, -0.000000, 2.085796] + poiseuille_profile.test_vals = [-12.494668, -7.710678, -0.000000, 2.085796] poiseuille_profile.test_vals_aarch64 = [-12.494684, -7.711379, -0.000000, 2.085796] #last 4 columns test_list.append(poiseuille_profile) @@ -1605,7 +1605,7 @@ def main(): pywrapper_custom_inlet.cfg_dir = "py_wrapper/custom_inlet" pywrapper_custom_inlet.cfg_file = "lam_flatplate.cfg" pywrapper_custom_inlet.test_iter = 20 - pywrapper_custom_inlet.test_vals = [-4.124164, -1.544359, -3.808866, 1.338411, -0.752679, 0.161436, -1.2391e-02, 5.1662e-01, -5.2901e-01] + pywrapper_custom_inlet.test_vals = [-4.124228, -1.544435, -3.806720, 1.338337, -0.752645, 0.162022, -0.012134, 0.516870, -0.529000] pywrapper_custom_inlet.command = TestCase.Command(exec = "python", param = "run.py") pywrapper_custom_inlet.timeout = 1600 pywrapper_custom_inlet.tol = 0.0001 From d179a4e116bd99ddcc6e22aaa9facc51a22933fd Mon Sep 17 00:00:00 2001 From: bigfooted Date: Sat, 10 Aug 2024 16:14:20 +0200 Subject: [PATCH 08/11] remove comments on compressible bounded scalar --- SU2_CFD/src/solvers/CEulerSolver.cpp | 10 +--------- 1 file changed, 1 insertion(+), 9 deletions(-) diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index 73c2ead44fa..789c068d894 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -341,10 +341,6 @@ CEulerSolver::CEulerSolver(CGeometry *geometry, CConfig *config, CommunicateInitialState(geometry, config); - /*--- Sizing edge mass flux array ---*/ - if (config->GetBounded_Scalar()) - EdgeMassFluxes.resize(geometry->GetnEdge()) = su2double(0.0); - /*--- Add the solver name.. ---*/ SolverName = "C.FLOW"; @@ -1776,7 +1772,7 @@ void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_contain const bool muscl = (config->GetMUSCL_Flow() && (iMesh == MESH_0)); const bool limiter = (config->GetKind_SlopeLimit_Flow() != LIMITER::NONE); const bool van_albada = (config->GetKind_SlopeLimit_Flow() == LIMITER::VAN_ALBADA_EDGE); - //const bool bounded_scalar = config->GetBounded_Scalar(); + /*--- Non-physical counter. ---*/ unsigned long counter_local = 0; SU2_OMP_MASTER @@ -1941,10 +1937,6 @@ void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_contain auto residual = numerics->ComputeResidual(config); - // this does not seem to be necessary - cout << "bounded" << endl; - //if (bounded_scalar) EdgeMassFluxes[iEdge] = residual[0]; - /*--- Set the final value of the Roe dissipation coefficient ---*/ if ((kind_dissipation != NO_ROELOWDISS) && (MGLevel != MESH_0)) { From 1757060262add279f321123d48d8c3f298288f44 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Sat, 10 Aug 2024 18:13:34 +0200 Subject: [PATCH 09/11] update regression --- TestCases/parallel_regression.py | 4 ++-- TestCases/serial_regression.py | 2 +- 2 files changed, 3 insertions(+), 3 deletions(-) diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index 02815976ef2..9dc96e30afd 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -311,7 +311,7 @@ def main(): flatplate_udobj.cfg_dir = "user_defined_functions" flatplate_udobj.cfg_file = "lam_flatplate.cfg" flatplate_udobj.test_iter = 20 - flatplate_udobj.test_vals = [-6.664134, -1.190073, -0.954366, 0.000641, -0.000633, 0.000548, -0.001181, 596.940000, 300.020000, 296.920000, 22.201000, 0.525750, 37.278000, 2.347900] + flatplate_udobj.test_vals = [-6.663514, -1.189436, -0.954938, 0.000644, -0.000635, 0.000549, -0.001183, 596.930000, 300.010000, 296.920000, 22.187000, 0.527090, 37.284000, 2.354200] test_list.append(flatplate_udobj) # Laminar cylinder (steady) @@ -1023,7 +1023,7 @@ def main(): flatplate_unsteady.cfg_dir = "navierstokes/flatplate" flatplate_unsteady.cfg_file = "lam_flatplate_unst.cfg" flatplate_unsteady.test_iter = 3 - flatplate_unsteady.test_vals = [0.000008, -8.869318, -8.232471, -6.283031, -5.467151, -3.391524, 0.002078, -0.343471] + flatplate_unsteady.test_vals = [-8.876870, -8.250064, -6.293918, -5.469416, -3.398952, 0.002075, -0.324126] flatplate_unsteady.unsteady = True test_list.append(flatplate_unsteady) diff --git a/TestCases/serial_regression.py b/TestCases/serial_regression.py index e1b0b85f0e8..6804176eac4 100644 --- a/TestCases/serial_regression.py +++ b/TestCases/serial_regression.py @@ -1605,7 +1605,7 @@ def main(): pywrapper_custom_inlet.cfg_dir = "py_wrapper/custom_inlet" pywrapper_custom_inlet.cfg_file = "lam_flatplate.cfg" pywrapper_custom_inlet.test_iter = 20 - pywrapper_custom_inlet.test_vals = [-4.124228, -1.544435, -3.806720, 1.338337, -0.752645, 0.162022, -0.012134, 0.516870, -0.529000] + pywrapper_custom_inlet.test_vals = [-4.120437, -1.540129, -3.563086, 1.342567, -0.748783, 0.161976, -0.012959, 0.516320, -0.529280] pywrapper_custom_inlet.command = TestCase.Command(exec = "python", param = "run.py") pywrapper_custom_inlet.timeout = 1600 pywrapper_custom_inlet.tol = 0.0001 From 0faa67e9c9f5dd7a9210ce00abd0892c6a8faad5 Mon Sep 17 00:00:00 2001 From: bigfooted Date: Sat, 10 Aug 2024 23:15:41 +0200 Subject: [PATCH 10/11] update regression profile --- TestCases/hybrid_regression.py | 2 +- TestCases/parallel_regression.py | 2 +- TestCases/serial_regression.py | 2 +- 3 files changed, 3 insertions(+), 3 deletions(-) diff --git a/TestCases/hybrid_regression.py b/TestCases/hybrid_regression.py index a1ae9811bcb..ad6eb0b13ab 100644 --- a/TestCases/hybrid_regression.py +++ b/TestCases/hybrid_regression.py @@ -139,7 +139,7 @@ def main(): poiseuille_profile.cfg_dir = "navierstokes/poiseuille" poiseuille_profile.cfg_file = "profile_poiseuille.cfg" poiseuille_profile.test_iter = 10 - poiseuille_profile.test_vals = [-12.494736, -7.712091, -0.000000, 2.085796] + poiseuille_profile.test_vals = [-12.485991, -7.613038, -0.000000, 2.085796] poiseuille_profile.test_vals_aarch64 = [-12.494717, -7.711274, -0.000000, 2.085796] test_list.append(poiseuille_profile) diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index e5b06814dd1..c21f7282fb0 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -344,7 +344,7 @@ def main(): poiseuille_profile.cfg_dir = "navierstokes/poiseuille" poiseuille_profile.cfg_file = "profile_poiseuille.cfg" poiseuille_profile.test_iter = 10 - poiseuille_profile.test_vals = [-12.492865, -7.672526, -0.000000, 2.085796] + poiseuille_profile.test_vals = [-12.483988, -7.577838, -0.000000, 2.085796] poiseuille_profile.test_vals_aarch64 = [-12.492864, -7.671632, -0.000000, 2.085796] test_list.append(poiseuille_profile) diff --git a/TestCases/serial_regression.py b/TestCases/serial_regression.py index 6804176eac4..17350f42bae 100644 --- a/TestCases/serial_regression.py +++ b/TestCases/serial_regression.py @@ -198,7 +198,7 @@ def main(): poiseuille_profile.cfg_dir = "navierstokes/poiseuille" poiseuille_profile.cfg_file = "profile_poiseuille.cfg" poiseuille_profile.test_iter = 10 - poiseuille_profile.test_vals = [-12.494668, -7.710678, -0.000000, 2.085796] + poiseuille_profile.test_vals = [-12.485982, -7.613005, -0.000000, 2.085796] poiseuille_profile.test_vals_aarch64 = [-12.494684, -7.711379, -0.000000, 2.085796] #last 4 columns test_list.append(poiseuille_profile) From 48077f0c07e558542cf969f2e6416bfed8aa5749 Mon Sep 17 00:00:00 2001 From: Nijso Date: Tue, 27 Aug 2024 18:59:49 +0200 Subject: [PATCH 11/11] Update SU2_CFD/src/solvers/CEulerSolver.cpp Co-authored-by: Pedro Gomes <38071223+pcarruscag@users.noreply.github.com> --- SU2_CFD/src/solvers/CEulerSolver.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index 29212a68775..b0253ed592a 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -7080,7 +7080,7 @@ void CEulerSolver::BC_Inlet(CGeometry *geometry, CSolver **solver_container, if (tkeNeeded) Energy += GetTke_Inf(); /*--- Primitive variables, using the derived quantities ---*/ - Temperature = Pressure / ( Gas_Constant * Density); + const su2double Temperature = Pressure / ( Gas_Constant * Density); V_inlet[0] = Temperature; for (iDim = 0; iDim < nDim; iDim++) V_inlet[iDim+1] = Vel_Mag*Flow_Dir[iDim];