diff --git a/Common/include/CConfig.hpp b/Common/include/CConfig.hpp index 794622bb1a1e..db178a5d83c5 100644 --- a/Common/include/CConfig.hpp +++ b/Common/include/CConfig.hpp @@ -1097,6 +1097,7 @@ class CConfig { bool SpatialFourier; /*!< \brief option for computing the fourier transforms for subsonic non-reflecting BC. */ bool RampMotionFrame; /*!< \brief option for ramping up or down the motion Frame values */ bool RampOutlet; /*!< \brief option for ramping up or down the outlet values */ + bool RampMUSCL; bool RampRotatingFrame; /*!< \brief option for ramping up or down the motion Frame values */ bool RampTranslationFrame; /*!< \brief option for ramping up or down the outlet values */ bool RampOutletMassFlow; /*!< \brief option for ramping up or down the motion Frame values */ @@ -1113,6 +1114,13 @@ class CConfig { array kt_polycoeffs{{0.0}}; /*!< \brief Array for thermal conductivity polynomial coefficients. */ bool Body_Force; /*!< \brief Flag to know if a body force is included in the formulation. */ + struct CMUSCLRampParam { + su2double RampMUSCLPower; /*!< \brief Exponent by which to raise the MUSCL ramp to the power of */ + MUSCL_RAMP_TYPE Kind_MUSCLRamp; /*!< \brief The kind of MUSCL ramp */ + unsigned long rampMUSCLCoeff[3]; /*!< \brief ramp MUSCL value coefficients for the COption class. */ + } RampMUSCLParam; + su2double rampMUSCLValue; /*!< \brief Current value of the MUSCL ramp */ + ENUM_STREAMWISE_PERIODIC Kind_Streamwise_Periodic; /*!< \brief Kind of Streamwise periodic flow (pressure drop or massflow) */ bool Streamwise_Periodic_Temperature; /*!< \brief Use real periodicity for Energy equation or otherwise outlet source term. */ su2double Streamwise_Periodic_PressureDrop; /*!< \brief Value of prescribed pressure drop [Pa] which results in an artificial body force vector. */ @@ -1343,6 +1351,8 @@ class CConfig { void addUShortArrayOption(const string& name, int size, bool allow_fewer, unsigned short* option_field); + void addULongArrayOption(const string& name, int size, bool allow_fewer, unsigned long* option_field); + void addDoubleListOption(const string& name, unsigned short & size, su2double * & option_field); void addShortListOption(const string& name, unsigned short & size, short * & option_field); @@ -5179,6 +5189,12 @@ class CConfig { */ bool GetRampOutflow(void) const { return RampOutlet; } + /*! + * \brief Get MUSCL ramp option. + * \return Ramp MUSCL option + */ + bool GetMUSCLRamp(void) const { return RampMUSCL; } + /*! * \brief General interface for accessing ramp coefficient information * \return coeff for ramps @@ -5189,6 +5205,23 @@ class CConfig { else return 0; }; + /*! + * \brief Set MUSCL ramp value. + */ + void SetMUSCLRampValue(su2double ramp_value) { rampMUSCLValue = ramp_value; } + + /*! + * \brief Get MUSCL ramp value. + * \return Ramp MUSCL value + */ + su2double GetMUSCLRampValue(void) const { return rampMUSCLValue; } + + /*! + * \brief Get MUSCL ramp paramaters. + * \return Ramp MUSCL kind + */ + const CMUSCLRampParam& GetMUSCLRampParam(void) const { return RampMUSCLParam; } + /*! * \brief Generic interface for setting monitor outlet values for the ramp. */ diff --git a/Common/include/option_structure.hpp b/Common/include/option_structure.hpp index bf99257b1387..4a8fe76aa9c1 100644 --- a/Common/include/option_structure.hpp +++ b/Common/include/option_structure.hpp @@ -1931,7 +1931,18 @@ enum TURBO_MARKER_TYPE{ enum class RAMP_TYPE{ GRID, /*!< \brief flag for rotational/translational ramps */ - BOUNDARY /*!< \brief flag for pressure/mass flow ramps*/ + BOUNDARY, /*!< \brief flag for pressure/mass flow ramps*/ + MUSCL /*!< \brief flag for MUSCL ramps */ +}; + +enum class MUSCL_RAMP_TYPE{ + ITERATION, /*!< \brief flag for linear iteration-based ramp */ + SMOOTH_FUNCTION /*!< \brief flag for smooth cosine ramp */ +}; + +static const MapType MUSCLRamp_Map = { + MakePair("ITERATION", MUSCL_RAMP_TYPE::ITERATION) + MakePair("SMOOTH_FUNCTION", MUSCL_RAMP_TYPE::SMOOTH_FUNCTION) }; /*! diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index c37a471be58e..013609ef7d9c 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -409,6 +409,13 @@ void CConfig::addULongListOption(const string& name, unsigned short & size, unsi option_map.insert(pair(name, val)); } +void CConfig::addULongArrayOption(const string& name, const int size, const bool allow_fewer, unsigned long* option_field) { + assert(option_map.find(name) == option_map.end()); + all_options.insert(pair(name, true)); + COptionBase* val = new COptionArray(name, size, allow_fewer, option_field); + option_map.insert(pair(name, val)); +} + void CConfig::addStringListOption(const string& name, unsigned short & num_marker, string* & option_field) { assert(option_map.find(name) == option_map.end()); all_options.insert(pair(name, true)); @@ -1986,6 +1993,16 @@ void CConfig::SetConfig_Options() { addBoolOption("MUSCL_FLOW", MUSCL_Flow, true); /*!\brief MUSCL_KAPPA_FLOW \n DESCRIPTION: Blending coefficient for the U-MUSCL scheme \ingroup Config*/ addDoubleOption("MUSCL_KAPPA_FLOW", MUSCL_Kappa_Flow, 0.0); + /*!\brief RAMP_MUSCL \n DESCRIPTION: Enable ramping of the MUSCL scheme from 1st to 2nd order using specified method*/ + addBoolOption("RAMP_MUSCL", RampMUSCL, false); + /*! brief RAMP_OUTLET_COEFF \n DESCRIPTION: the 1st coeff is the ramp start iteration, + * the 2nd coeff is the iteration update frequenct, 3rd coeff is the total number of iterations */ + RampMUSCLParam.rampMUSCLCoeff[0] = 0.0; RampMUSCLParam.rampMUSCLCoeff[1] = 1.0; RampMUSCLParam.rampMUSCLCoeff[2] = 500.0; + addULongArrayOption("RAMP_MUSCL_COEFF", 3, false, RampMUSCLParam.rampMUSCLCoeff); + /*!\brief RAMP_MUSCL_POWER \n DESRCIPTION: Exponent of the MUSCL ramp formulation */ + addDoubleOption("RAMP_MUSCL_POWER", RampMUSCLParam.RampMUSCLPower, 1.0); + /*!\brief KIND_MUSCL_RAMP \n DESCRIPTION: The kind of MUSCL Ramp to be applied */ + addEnumOption("KIND_MUSCL_RAMP", RampMUSCLParam.Kind_MUSCLRamp, MUSCLRamp_Map, MUSCL_RAMP_TYPE::ITERATION); /*!\brief SLOPE_LIMITER_FLOW * DESCRIPTION: Slope limiter for the direct solution. \n OPTIONS: See \link Limiter_Map \endlink \n DEFAULT VENKATAKRISHNAN \ingroup Config*/ addEnumOption("SLOPE_LIMITER_FLOW", Kind_SlopeLimit_Flow, Limiter_Map, LIMITER::VENKATAKRISHNAN); @@ -4479,6 +4496,13 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i } } + if(RampMUSCL && !DiscreteAdjoint){ + if (RampMUSCLParam.RampMUSCLPower <= 0.0) SU2_MPI::Error("RAMP_MUSCL_POWER cannot be less than or equal to zero!", CURRENT_FUNCTION); + rampMUSCLValue = 0.0; + } else { + rampMUSCLValue = 1.0; + } + /*--- Check on extra Relaxation factor for Giles---*/ if(extrarelfac[1] > 0.5){ extrarelfac[1] = 0.5; @@ -6968,6 +6992,21 @@ void CConfig::SetOutput(SU2_COMPONENT val_software, unsigned short val_izone) { auto PrintLimiterInfo = [&](const LIMITER kind_limiter, const su2double kappa) { cout << "Second order integration in space, with slope limiter.\n"; + if (RampMUSCL) { + cout << "Ramping MUSCL sheme from first to second order starting at iter " << RampMUSCLParam.rampMUSCLCoeff[RAMP_COEFF::INITIAL_VALUE] + << ", ending at iter " << RampMUSCLParam.rampMUSCLCoeff[RAMP_COEFF::FINAL_ITER] + RampMUSCLParam.rampMUSCLCoeff[RAMP_COEFF::INITIAL_VALUE] + << ", updating every " << RampMUSCLParam.rampMUSCLCoeff[RAMP_COEFF::UPDATE_FREQ] << " iterations." << endl; + string MUSCLRampType; + switch (RampMUSCLParam.Kind_MUSCLRamp) { + case MUSCL_RAMP_TYPE::ITERATION: + MUSCLRampType = "linear"; + break; + case MUSCL_RAMP_TYPE::SMOOTH_FUNCTION: + MUSCLRampType = "cosine"; + break; + } + cout << "Ramp applied according to a " << MUSCLRampType << " function, raised to the power " << RampMUSCLParam.RampMUSCLPower << "." << endl; + } if (kappa != 0.0) cout << "U-MUSCL reconstruction, with coefficient: " << kappa << ".\n"; switch (kind_limiter) { case LIMITER::NONE: diff --git a/SU2_CFD/include/numerics_simd/flow/convection/common.hpp b/SU2_CFD/include/numerics_simd/flow/convection/common.hpp index 8ef67caebec8..08c1d67511c0 100644 --- a/SU2_CFD/include/numerics_simd/flow/convection/common.hpp +++ b/SU2_CFD/include/numerics_simd/flow/convection/common.hpp @@ -71,9 +71,10 @@ FORCEINLINE Double musclReconstruction(const GradType& grad, const Double delta, size_t iVar, Double kappa, - Double relax) { + Double relax, + Double ramp_val) { const Double proj = dot(grad[iVar], vector_ij); - return relax * umusclProjection(proj, delta, kappa); + return ramp_val * relax * umusclProjection(proj, delta, kappa); } /*! @@ -86,7 +87,8 @@ FORCEINLINE void musclUnlimited(Int iPoint, const Gradient_t& gradient, CPair& V, Double kappa, - Double relax) { + Double relax, + Double ramp_val) { constexpr auto nVarGrad = nVarGrad_ > 0 ? nVarGrad_ : VarType::nVar; auto grad_i = gatherVariables(iPoint, gradient); @@ -97,8 +99,8 @@ FORCEINLINE void musclUnlimited(Int iPoint, const Double delta_ij = V.j.all(iVar) - V.i.all(iVar); /*--- U-MUSCL reconstructed variables ---*/ - const Double proj_i = musclReconstruction(grad_i, vector_ij, delta_ij, iVar, kappa, relax); - const Double proj_j = musclReconstruction(grad_j, vector_ij, delta_ij, iVar, kappa, relax); + const Double proj_i = musclReconstruction(grad_i, vector_ij, delta_ij, iVar, kappa, relax, ramp_val); + const Double proj_j = musclReconstruction(grad_j, vector_ij, delta_ij, iVar, kappa, relax, ramp_val); /*--- Apply reconstruction: V_L = V_i + 0.5 * dV_ij^kap ---*/ V.i.all(iVar) += 0.5 * proj_i; @@ -117,7 +119,8 @@ FORCEINLINE void musclPointLimited(Int iPoint, const Gradient_t& gradient, CPair& V, Double kappa, - Double relax) { + Double relax, + Double ramp_val) { constexpr auto nVarGrad = nVarGrad_ > 0 ? nVarGrad_ : VarType::nVar; auto lim_i = gatherVariables(iPoint, limiter); @@ -131,8 +134,8 @@ FORCEINLINE void musclPointLimited(Int iPoint, const Double delta_ij = V.j.all(iVar) - V.i.all(iVar); /*--- U-MUSCL reconstructed variables ---*/ - const Double proj_i = musclReconstruction(grad_i, vector_ij, delta_ij, iVar, kappa, relax); - const Double proj_j = musclReconstruction(grad_j, vector_ij, delta_ij, iVar, kappa, relax); + const Double proj_i = musclReconstruction(grad_i, vector_ij, delta_ij, iVar, kappa, relax, ramp_val); + const Double proj_j = musclReconstruction(grad_j, vector_ij, delta_ij, iVar, kappa, relax, ramp_val); /*--- Apply reconstruction: V_L = V_i + 0.5 * lim * dV_ij^kap ---*/ V.i.all(iVar) += 0.5 * lim_i(iVar) * proj_i; @@ -150,7 +153,8 @@ FORCEINLINE void musclEdgeLimited(Int iPoint, const Gradient_t& gradient, CPair& V, Double kappa, - Double relax) { + Double relax, + Double ramp_val) { constexpr auto nVarGrad = nVarGrad_ > 0 ? nVarGrad_ : VarType::nVar; auto grad_i = gatherVariables(iPoint, gradient); @@ -162,8 +166,8 @@ FORCEINLINE void musclEdgeLimited(Int iPoint, const Double delta_ij_2 = pow(delta_ij, 2) + 1e-6; /*--- U-MUSCL reconstructed variables ---*/ - const Double proj_i = musclReconstruction(grad_i, vector_ij, delta_ij, iVar, kappa, relax); - const Double proj_j = musclReconstruction(grad_j, vector_ij, delta_ij, iVar, kappa, relax); + const Double proj_i = musclReconstruction(grad_i, vector_ij, delta_ij, iVar, kappa, relax, ramp_val); + const Double proj_j = musclReconstruction(grad_j, vector_ij, delta_ij, iVar, kappa, relax, ramp_val); /// TODO: Customize the limiter function. const Double lim_i = (delta_ij_2 + proj_i*delta_ij) / (pow(proj_i,2) + delta_ij_2); @@ -200,7 +204,8 @@ FORCEINLINE CPair reconstructPrimitives(Int iEdge, Int iPoint, Int LIMITER limiterType, const CPair& V1st, const VectorDbl& vector_ij, - const VariableType& solution) { + const VariableType& solution, + const su2double& ramp_val) { static_assert(ReconVarType::nVar <= PrimVarType::nVar); const auto& gradients = solution.GetGradient_Reconstruction(); @@ -218,13 +223,13 @@ FORCEINLINE CPair reconstructPrimitives(Int iEdge, Int iPoint, Int constexpr auto nVarGrad = ReconVarType::nVar - 2; switch (limiterType) { case LIMITER::NONE: - musclUnlimited(iPoint, jPoint, vector_ij, gradients, V, kappa, relax); + musclUnlimited(iPoint, jPoint, vector_ij, gradients, V, kappa, relax, ramp_val); break; case LIMITER::VAN_ALBADA_EDGE: - musclEdgeLimited(iPoint, jPoint, vector_ij, gradients, V, kappa, relax); + musclEdgeLimited(iPoint, jPoint, vector_ij, gradients, V, kappa, relax, ramp_val); break; default: - musclPointLimited(iPoint, jPoint, vector_ij, limiters, gradients, V, kappa, relax); + musclPointLimited(iPoint, jPoint, vector_ij, limiters, gradients, V, kappa, relax, ramp_val); break; } /*--- Recompute density using the reconstructed pressure and temperature. ---*/ diff --git a/SU2_CFD/include/numerics_simd/flow/convection/roe.hpp b/SU2_CFD/include/numerics_simd/flow/convection/roe.hpp index 6b1d09e8b9e1..597ddc7efedb 100644 --- a/SU2_CFD/include/numerics_simd/flow/convection/roe.hpp +++ b/SU2_CFD/include/numerics_simd/flow/convection/roe.hpp @@ -123,7 +123,7 @@ class CRoeBase : public Base { /*--- Recompute density and enthalpy instead of reconstructing. ---*/ auto V = reconstructPrimitives >( - iEdge, iPoint, jPoint, gamma, gasConst, muscl, umusclKappa, nkRelax, typeLimiter, V1st, vector_ij, solution); + iEdge, iPoint, jPoint, gamma, gasConst, muscl, umusclKappa, nkRelax, typeLimiter, V1st, vector_ij, solution, config.GetMUSCLRampValue()); /*--- Compute conservative variables. ---*/ diff --git a/SU2_CFD/include/solvers/CScalarSolver.inl b/SU2_CFD/include/solvers/CScalarSolver.inl index bcf9fae33602..62029bd26edc 100644 --- a/SU2_CFD/include/solvers/CScalarSolver.inl +++ b/SU2_CFD/include/solvers/CScalarSolver.inl @@ -222,8 +222,8 @@ void CScalarSolver::Upwind_Residual(CGeometry* geometry, CSolver** for (auto iVar = 0u; iVar < solver_container[FLOW_SOL]->GetnPrimVarGrad(); iVar++) { const su2double V_ij = V_j[iVar] - V_i[iVar]; - su2double Project_Grad_i = MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, V_ij, kappaFlow); - su2double Project_Grad_j = MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, V_ij, kappaFlow); + su2double Project_Grad_i = MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, V_ij, kappaFlow, config->GetMUSCLRampValue()); + su2double Project_Grad_j = MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, V_ij, kappaFlow, config->GetMUSCLRampValue()); if (limiterFlow) { Project_Grad_i *= Limiter_i[iVar]; @@ -251,8 +251,8 @@ void CScalarSolver::Upwind_Residual(CGeometry* geometry, CSolver** for (auto iVar = 0u; iVar < nVar; iVar++) { const su2double U_ij = Scalar_j[iVar] - Scalar_i[iVar]; - su2double Project_Grad_i = MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, U_ij, kappa); - su2double Project_Grad_j = MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, U_ij, kappa); + su2double Project_Grad_i = MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, U_ij, kappa, config->GetMUSCLRampValue()); + su2double Project_Grad_j = MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, U_ij, kappa, config->GetMUSCLRampValue()); if (limiter) { Project_Grad_i *= Limiter_i[iVar]; diff --git a/SU2_CFD/include/solvers/CSolver.hpp b/SU2_CFD/include/solvers/CSolver.hpp index f0b3a4c493d1..884671abf369 100644 --- a/SU2_CFD/include/solvers/CSolver.hpp +++ b/SU2_CFD/include/solvers/CSolver.hpp @@ -587,11 +587,12 @@ class CSolver { * \param[in] vector_ij - Distance vector. * \param[in] delta_ij - Centered difference. * \param[in] kappa - Blending coefficient for U-MUSCL reconstruction. + * \param[in] ramp_val - Value of the ramp * \return - Projected variable. */ - inline su2double MUSCL_Reconstruction(const su2double* grad, const su2double* vector_ij, su2double delta_ij, su2double kappa) { + inline su2double MUSCL_Reconstruction(const su2double* grad, const su2double* vector_ij, su2double delta_ij, su2double kappa, su2double ramp_val) { su2double project_grad = GeometryToolbox::DotProduct(nDim, grad, vector_ij); - return LimiterHelpers<>::umusclProjection(project_grad, delta_ij, kappa); + return ramp_val * LimiterHelpers<>::umusclProjection(project_grad, delta_ij, kappa); } /*! diff --git a/SU2_CFD/src/iteration/CFluidIteration.cpp b/SU2_CFD/src/iteration/CFluidIteration.cpp index 529b3e12453c..fea25ba3b45d 100644 --- a/SU2_CFD/src/iteration/CFluidIteration.cpp +++ b/SU2_CFD/src/iteration/CFluidIteration.cpp @@ -241,6 +241,9 @@ bool CFluidIteration::Monitor(COutput* output, CIntegration**** integration, CGe if (config[val_iZone]->GetRampOutflow()) UpdateRamp(geometry, config, config[val_iZone]->GetInnerIter(), val_iZone, RAMP_TYPE::BOUNDARY); + if (config[val_iZone]->GetMUSCLRamp()) + UpdateRamp(geometry, config, config[val_iZone]->GetInnerIter(), val_iZone, RAMP_TYPE::MUSCL); + output->SetHistoryOutput(geometry[val_iZone][val_iInst][MESH_0], solver[val_iZone][val_iInst][MESH_0], config[val_iZone], config[val_iZone]->GetTimeIter(), config[val_iZone]->GetOuterIter(), config[val_iZone]->GetInnerIter()); @@ -330,6 +333,29 @@ void CFluidIteration::UpdateRamp(CGeometry**** geometry_container, CConfig** con } } } + + if (ramp_flag == RAMP_TYPE::MUSCL) { + const auto RampMUSCLParam = config->GetMUSCLRampParam(); + const long unsigned startIter = RampMUSCLParam.rampMUSCLCoeff[RAMP_COEFF::INITIAL_VALUE]; + const long unsigned updateFreq = RampMUSCLParam.rampMUSCLCoeff[RAMP_COEFF::UPDATE_FREQ]; + const long unsigned rampLength = RampMUSCLParam.rampMUSCLCoeff[RAMP_COEFF::FINAL_ITER]; + auto iterFrac = (static_cast(iter - startIter)/static_cast((rampLength + startIter) - startIter)); + if (iter < startIter) return; + if ((iter == startIter) && (rank == MASTER_NODE)) cout << "Beginning to ramp MUSCL scheme..." << endl; + if ((iter % updateFreq == 0 && iter < (rampLength + startIter)) || (iter == (rampLength + startIter))) { + switch (RampMUSCLParam.Kind_MUSCLRamp) { + case MUSCL_RAMP_TYPE::ITERATION: + config->SetMUSCLRampValue(std::pow(std::min(1.0, iterFrac), RampMUSCLParam.RampMUSCLPower)); + break; + case MUSCL_RAMP_TYPE::SMOOTH_FUNCTION: + config->SetMUSCLRampValue(std::pow((0.5 * (1 - cos(M_PI * std::min(1.0, iterFrac)))), RampMUSCLParam.RampMUSCLPower)); + break; + default: + break; + } + if (rank == MASTER_NODE) cout << "MUSCL Ramp value updated. New Value: " << config->GetMUSCLRampValue() << endl; + } + } } void CFluidIteration::ComputeTurboPerformance(CSolver***** solver, CGeometry**** geometry_container, CConfig** config_container, unsigned long ExtIter) { diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index dcca14edd19d..06ade8873507 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -1881,8 +1881,8 @@ void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_contain for (auto iVar = 0u; iVar < nPrimVarGrad; iVar++) { const su2double V_ij = V_j[iVar] - V_i[iVar]; - const su2double Project_Grad_i = nkRelax * MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, V_ij, kappa); - const su2double Project_Grad_j = nkRelax * MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, V_ij, kappa); + const su2double Project_Grad_i = nkRelax * MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, V_ij, kappa, config->GetMUSCLRampValue()); + const su2double Project_Grad_j = nkRelax * MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, V_ij, kappa, config->GetMUSCLRampValue()); su2double lim_i = 1.0; su2double lim_j = 1.0; diff --git a/SU2_CFD/src/solvers/CIncEulerSolver.cpp b/SU2_CFD/src/solvers/CIncEulerSolver.cpp index 3df1ffbe3dd8..501286ecf206 100644 --- a/SU2_CFD/src/solvers/CIncEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CIncEulerSolver.cpp @@ -1279,8 +1279,8 @@ void CIncEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_cont for (auto iVar = 0u; iVar < nPrimVarGrad; iVar++) { const su2double V_ij = V_j[iVar] - V_i[iVar]; - const su2double Project_Grad_i = nkRelax * MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, V_ij, kappa); - const su2double Project_Grad_j = nkRelax * MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, V_ij, kappa); + const su2double Project_Grad_i = nkRelax * MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, V_ij, kappa, config->GetMUSCLRampValue()); + const su2double Project_Grad_j = nkRelax * MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, V_ij, kappa, config->GetMUSCLRampValue()); su2double lim_i = 1.0; su2double lim_j = 1.0; diff --git a/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp b/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp index 65f3ea81f123..1e3ab602aa12 100644 --- a/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CNEMOEulerSolver.cpp @@ -542,8 +542,8 @@ void CNEMOEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_con for (auto iVar = 0ul; iVar < nPrimVarGrad; iVar++) { const su2double V_ij = V_j[iVar] - V_i[iVar]; - Project_Grad_i[iVar] = nkRelax * MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, V_ij, kappa); - Project_Grad_j[iVar] = nkRelax * MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, V_ij, kappa); + Project_Grad_i[iVar] = nkRelax * MUSCL_Reconstruction(Gradient_i[iVar], Vector_ij, V_ij, kappa, config->GetMUSCLRampValue()); + Project_Grad_j[iVar] = nkRelax * MUSCL_Reconstruction(Gradient_j[iVar], Vector_ij, V_ij, kappa, config->GetMUSCLRampValue()); if (limiter) { if (van_albada) { diff --git a/TestCases/hybrid_regression.py b/TestCases/hybrid_regression.py index c7ad4c989645..312b2b66cfb3 100644 --- a/TestCases/hybrid_regression.py +++ b/TestCases/hybrid_regression.py @@ -564,7 +564,7 @@ def main(): axial_stage2D.cfg_dir = "turbomachinery/axial_stage_2D" axial_stage2D.cfg_file = "Axial_stage2D.cfg" axial_stage2D.test_iter = 20 - axial_stage2D.test_vals = [1.084448, 1.526930, -2.895083, 2.607569, -2.479664, 3.063779, 106380.000000, 106380.000000, 5.733600, 64.737000] + axial_stage2D.test_vals = [1.065797, 1.519589, -2.928280, 2.573904, -2.526637, 3.017140, 106370.000000, 106370.000000, 5.726800, 64.383000] test_list.append(axial_stage2D) # 2D transonic stator restart diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index c96bc56ac810..fb89f7ba1655 100644 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -1108,7 +1108,7 @@ def main(): axial_stage2D.cfg_dir = "turbomachinery/axial_stage_2D" axial_stage2D.cfg_file = "Axial_stage2D.cfg" axial_stage2D.test_iter = 20 - axial_stage2D.test_vals = [1.084454, 1.526942, -2.895082, 2.607570, -2.479664, 3.063779, 106380.000000, 106380.000000, 5.733600, 64.737000] + axial_stage2D.test_vals = [1.065803, 1.519598, -2.928278, 2.573906, -2.526640, 3.017138, 106370.000000, 106370.000000, 5.726800, 64.383000] test_list.append(axial_stage2D) # 2D transonic stator restart diff --git a/TestCases/serial_regression.py b/TestCases/serial_regression.py index 94bc257ddb8c..f4005b3cb1b7 100644 --- a/TestCases/serial_regression.py +++ b/TestCases/serial_regression.py @@ -875,7 +875,7 @@ def main(): axial_stage2D.cfg_dir = "turbomachinery/axial_stage_2D" axial_stage2D.cfg_file = "Axial_stage2D.cfg" axial_stage2D.test_iter = 20 - axial_stage2D.test_vals = [1.084452, 1.526941, -2.895084, 2.607568, -2.479664, 3.063779, 106380.000000, 106380.000000, 5.733600, 64.737000] + axial_stage2D.test_vals = [1.065801, 1.519596, -2.928281, 2.573903, -2.526639, 3.017139, 106370.000000, 106370.000000, 5.726800, 64.383000] test_list.append(axial_stage2D) # 2D transonic stator restart diff --git a/TestCases/turbomachinery/axial_stage_2D/Axial_stage2D.cfg b/TestCases/turbomachinery/axial_stage_2D/Axial_stage2D.cfg index 042a04c01cbb..957d2fdb6e1f 100755 --- a/TestCases/turbomachinery/axial_stage_2D/Axial_stage2D.cfg +++ b/TestCases/turbomachinery/axial_stage_2D/Axial_stage2D.cfg @@ -118,6 +118,10 @@ LIMITER_ITER= 999999 % CONV_NUM_METHOD_FLOW= ROE MUSCL_FLOW= YES +RAMP_MUSCL= YES +% NOTE: These coefficients are not practical for real calculations, this is for the purpose of regression testing. +RAMP_MUSCL_COEFF= (0.0, 1.0, 20.0) +RAMP_MUSCL_POWER = 1 SLOPE_LIMITER_FLOW= VAN_ALBADA_EDGE ENTROPY_FIX_COEFF= 0.01 JST_SENSOR_COEFF= ( 0.5, 0.02 ) diff --git a/config_template.cfg b/config_template.cfg index aa7a1a97fc97..94e8541686fd 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -1463,6 +1463,22 @@ CUSTOM_OBJFUNC= 'DRAG + 10 * pow(fmax(0.4-LIFT, 0), 2)' % Required for 2nd order upwind schemes (NO, YES) MUSCL_FLOW= YES % +% Option to ramp the MUSCL reconstruction scheme gradually from first to second order +RAMP_MUSCL= NO +% +% Coefficients for the MUSCL ramp, 1st is ramp start iteration, +% 2nd is the ramp update frequency, 3rd is the total number of iterations of the ramp +RAMP_MUSCL_COEFF= (1000.0, 10.0, 500.0) +% +% Value for modifying how aggresively the ramp is applied (cannot be less than or equal to 0) +% < 1 = Aggresive ramp up early, slower later, > 1 slow initial ramping, more aggresive towards end +RAMP_MUSCL_POWER = 1 +% +% Type of MUSCL ramp function used +% ITERATION = [min(1.0, iter/endIter)]^RAMP_MUSCL_POWER +% SMOOTH_FUNCTION = [1/2 * ( 1 - cos ( pi * min(1.0, iter/endIter)))]^RAMP_MUSCL_POWER +KIND_MUSCL_RAMP= ITERATION +% % Coefficient for blending upwind and central differences in MUSCL scheme. % Values range from -1 (fully one-sided) to 1 (central difference). MUSCL_KAPPA_FLOW= 0.0