Skip to content
Open
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
6 changes: 6 additions & 0 deletions Common/include/CConfig.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -815,6 +815,7 @@ class CConfig {
unsigned short ActDisk_Jump; /*!< \brief Format of the output files. */
unsigned long StartWindowIteration; /*!< \brief Starting Iteration for long time Windowing apporach . */
unsigned short nCFL_AdaptParam; /*!< \brief Number of CFL parameters provided in config. */
unsigned long outlierMitigationParam[4]; /*!< \brief Parameters of outlier mitigation strategy. */
bool CFL_Adapt; /*!< \brief Use adaptive CFL number. */
bool HB_Precondition; /*!< \brief Flag to turn on harmonic balance source term preconditioning */
su2double RefArea, /*!< \brief Reference area for coefficient computation. */
Expand Down Expand Up @@ -1712,6 +1713,11 @@ class CConfig {
*/
bool GetCFL_Adapt(void) const { return CFL_Adapt; }

/*!
* \brief Get the outlier mitigation parameters.
*/
const unsigned long* GetOutlierMitigationParam() const { return outlierMitigationParam; }

/*!
* \brief Get the value of the limits for the sections.
* \return Value of the limits for the sections.
Expand Down
7 changes: 6 additions & 1 deletion Common/src/CConfig.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1819,6 +1819,11 @@ void CConfig::SetConfig_Options() {
addDoubleOption("CFL_NUMBER", CFLFineGrid, 1.25);
/* DESCRIPTION: Max time step in local time stepping simulations */
addDoubleOption("MAX_DELTA_TIME", Max_DeltaTime, 1000000);
/* !\brief OUTLIER_MITIGATION_PARAM
* DESCRIPTION: Parameters of the outlier mitigation strategy: start iteration, update frequency, print frequency,
* and number of standard deviations (N * sigma) to identify outliers statistically. \ingroup Config*/
outlierMitigationParam[0] = 999999; outlierMitigationParam[1] = 5; outlierMitigationParam[2] = 2; outlierMitigationParam[3] = 5;
addULongArrayOption("OUTLIER_MITIGATION_PARAM", 4, true, outlierMitigationParam);
/* DESCRIPTION: Activate The adaptive CFL number. */
addBoolOption("CFL_ADAPT", CFL_Adapt, false);
/* !\brief CFL_ADAPT_PARAM
Expand Down Expand Up @@ -2029,7 +2034,7 @@ void CConfig::SetConfig_Options() {
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,
/*! brief RAMP_MUSCL_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);
Expand Down
17 changes: 12 additions & 5 deletions SU2_CFD/include/numerics_simd/flow/convection/common.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -188,6 +188,7 @@
* \param[in] V1st - Pair of compressible flow primitives for nodes i,j.
* \param[in] vector_ij - Distance vector from i to j.
* \param[in] solution - Entire solution container (a derived CVariable).
* \param[out] nonPhysical - Signals that the edge is treated as non-physical.
* \return Pair of primitive variables.
*/
template<class ReconVarType, class PrimVarType, size_t nDim, class VariableType>
Expand All @@ -201,7 +202,8 @@
const LIMITER limiterType,
const CPair<PrimVarType>& V1st,
const VectorDbl<nDim>& vector_ij,
const VariableType& solution) {
const VariableType& solution,
Double& nonPhysical) {
static_assert(ReconVarType::nVar <= PrimVarType::nVar);

const auto& gradients = solution.GetGradient_Reconstruction();
Expand Down Expand Up @@ -262,15 +264,20 @@
const Double neg_sound_speed = enthalpy * (R+1) < 0.5 * v_squared;

/*--- Revert to first order if the state is non-physical. ---*/
Double bad_recon = fmax(neg_p_or_rho, neg_sound_speed);
Double nonPhysical = fmax(neg_p_or_rho, neg_sound_speed);

Check notice

Code scanning / CodeQL

Declaration hides parameter Note

Local variable 'nonPhysical' hides a
parameter of the same name
.
/*--- Handle SIMD dimensions 1 by 1. ---*/
for (size_t k = 0; k < Double::Size; ++k) {
bad_recon[k] = solution.UpdateNonPhysicalEdgeCounter(iEdge[k], bad_recon[k]);
nonPhysical[k] = solution.UpdateNonPhysicalEdgeCounter(iEdge[k], nonPhysical[k]);
nonPhysical[k] = fmax(nonPhysical[k],
fmax(solution.OutlierMitigation(iPoint[k]),
solution.OutlierMitigation(jPoint[k])) / VariableType::MAX_OUTLIER_MITIGATION);
}
for (size_t iVar = 0; iVar < ReconVarType::nVar; ++iVar) {
V.i.all(iVar) = bad_recon * V1st.i.all(iVar) + (1-bad_recon) * V.i.all(iVar);
V.j.all(iVar) = bad_recon * V1st.j.all(iVar) + (1-bad_recon) * V.j.all(iVar);
V.i.all(iVar) = nonPhysical * V1st.i.all(iVar) + (1-nonPhysical) * V.i.all(iVar);
V.j.all(iVar) = nonPhysical * V1st.j.all(iVar) + (1-nonPhysical) * V.j.all(iVar);
}
} else {
nonPhysical = 0;
}
return V;
}
Expand Down
17 changes: 10 additions & 7 deletions SU2_CFD/include/numerics_simd/flow/convection/upwind.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -118,8 +118,10 @@ class CUpwindBase : public Base {
V1st.j.all = gatherVariables<nPrimVar>(jPoint, solution.GetPrimitive());

/*--- Recompute density and enthalpy instead of reconstructing. ---*/
Double nonPhysical;
auto V = reconstructPrimitives<CCompressiblePrimitives<nDim,nPrimVarGrad> >(
iEdge, iPoint, jPoint, gamma, gasConst, muscl, umusclKappa, umusclRamp, typeLimiter, V1st, vector_ij, solution);
iEdge, iPoint, jPoint, gamma, gasConst, muscl, umusclKappa, umusclRamp,
typeLimiter, V1st, vector_ij, solution, nonPhysical);

/*--- Compute conservative variables. ---*/

Expand All @@ -132,8 +134,8 @@ class CUpwindBase : public Base {
const auto derived = static_cast<const Derived*>(this);
VectorDbl<nVar> flux;
MatrixDbl<nVar> jac_i, jac_j;
derived->finalizeFlux(flux, jac_i, jac_j, implicit, area, unitNormal,
normal, V, U, iPoint, jPoint, solution, geometry);
derived->finalizeFlux(flux, jac_i, jac_j, implicit, area, unitNormal, normal,
V, U, iPoint, jPoint, nonPhysical, solution, geometry);

/*--- Add the contributions from the base class (static decorator). ---*/

Expand Down Expand Up @@ -199,6 +201,7 @@ class CRoeScheme : public CUpwindBase<CRoeScheme<Decorator>, Decorator> {
const CPair<ConsVarType>& U,
const Int& iPoint,
const Int& jPoint,
const Double& nonPhysical,
const CEulerVariable& solution,
const CGeometry& geometry,
Ts&...) const {
Expand Down Expand Up @@ -227,10 +230,9 @@ class CRoeScheme : public CUpwindBase<CRoeScheme<Decorator>, Decorator> {

/*--- Apply Mavriplis' entropy correction to eigenvalues. ---*/

Double maxLambda = abs(projVel) + roeAvg.speedSound;

Double lambdaMin = fmax(entropyFix, nonPhysical) * (abs(projVel) + roeAvg.speedSound);
for (size_t iVar = 0; iVar < nVar; ++iVar) {
lambda(iVar) = fmax(abs(lambda(iVar)), entropyFix*maxLambda);
lambda(iVar) = fmax(abs(lambda(iVar)), lambdaMin);
}

/*--- Inviscid fluxes and Jacobians. ---*/
Expand Down Expand Up @@ -348,6 +350,7 @@ class CMSWScheme : public CUpwindBase<CMSWScheme<Decorator>, Decorator> {
const CPair<ConsVarType>& U,
const Int& iPoint,
const Int& jPoint,
const Double& nonPhysical,
const CEulerVariable& solution,
const CGeometry& geometry,
Ts&...) const {
Expand All @@ -358,7 +361,7 @@ class CMSWScheme : public CUpwindBase<CMSWScheme<Decorator>, Decorator> {
const auto sj = gatherVariables(jPoint, solution.GetSensor());

const Double dp = fmax(si, sj) - alpha * 0.06;
const Double w = 0.25 * (1 - sign(dp)) * (1 - exp(-100 * abs(dp)));
const Double w = 0.25 * (1 - sign(dp) * (1 - exp(-100 * abs(dp)))) * (1 - nonPhysical);
const Double onemw = 1 - w;

CPair<CCompressiblePrimitives<nDim, nPrimVarGrad>> Vweighted;
Expand Down
53 changes: 31 additions & 22 deletions SU2_CFD/include/solvers/CEulerSolver.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -116,6 +116,9 @@ class CEulerSolver : public CFVMFlowSolverBase<CEulerVariable, ENUM_REGIME::COMP

vector<CFluidModel*> FluidModel; /*!< \brief fluid model used in the solver. */

/*!< \brief Variables for outlier detection. */
su2double MeanTemperature, StdDevTemperature;

/*--- Turbomachinery Solver Variables ---*/

vector<su2activematrix> AverageFlux;
Expand Down Expand Up @@ -782,11 +785,17 @@ class CEulerSolver : public CFVMFlowSolverBase<CEulerVariable, ENUM_REGIME::COMP
void PrepareImplicitIteration(CGeometry *geometry, CSolver**, CConfig *config) final;

/*!
* \brief Complete an implicit iteration.
* \param[in] geometry - Geometrical definition of the problem.
* \brief Compute a suitable under-relaxation parameter to limit the change in the solution variables over
* a nonlinear iteration for stability.
* \param[in] config - Definition of the particular problem.
*/
void CompleteImplicitIteration(CGeometry *geometry, CSolver**, CConfig *config) final;
void ComputeUnderRelaxationFactor(const CConfig *config) final;

/*!
* \brief Identify points where the static temperature is an outlier to treat them
* as non-physical, e.g. force the reconstruction to be 1st order.
*/
void IdentifySolutionOutliers(const CConfig *config, unsigned long iter) final;

/*!
* \brief Provide the mass flow rate.
Expand Down Expand Up @@ -1127,20 +1136,21 @@ class CEulerSolver : public CFVMFlowSolverBase<CEulerVariable, ENUM_REGIME::COMP
unsigned short marker_flag,
TURBOMACHINERY_TYPE kind_turb){

if ((kind_turb == TURBOMACHINERY_TYPE::AXIAL && nDim == 3) || (kind_turb == TURBOMACHINERY_TYPE::CENTRIPETAL_AXIAL && marker_flag == OUTFLOW) || (kind_turb == TURBOMACHINERY_TYPE::AXIAL_CENTRIFUGAL && marker_flag == INFLOW) ){
turboVelocity[2] = turboNormal[0]*cartesianVelocity[0] + cartesianVelocity[1]*turboNormal[1];
turboVelocity[1] = turboNormal[0]*cartesianVelocity[1] - turboNormal[1]*cartesianVelocity[0];
if ((kind_turb == TURBOMACHINERY_TYPE::AXIAL && nDim == 3) ||
(kind_turb == TURBOMACHINERY_TYPE::CENTRIPETAL_AXIAL && marker_flag == OUTFLOW) ||
(kind_turb == TURBOMACHINERY_TYPE::AXIAL_CENTRIFUGAL && marker_flag == INFLOW)) {
turboVelocity[2] = turboNormal[0]*cartesianVelocity[0] + cartesianVelocity[1]*turboNormal[1];
turboVelocity[1] = turboNormal[0]*cartesianVelocity[1] - turboNormal[1]*cartesianVelocity[0];
turboVelocity[0] = cartesianVelocity[2];
}
else{
turboVelocity[0] = turboNormal[0]*cartesianVelocity[0] + cartesianVelocity[1]*turboNormal[1];
turboVelocity[1] = turboNormal[0]*cartesianVelocity[1] - turboNormal[1]*cartesianVelocity[0];
if (marker_flag == INFLOW){
} else {
turboVelocity[0] = turboNormal[0]*cartesianVelocity[0] + cartesianVelocity[1]*turboNormal[1];
turboVelocity[1] = turboNormal[0]*cartesianVelocity[1] - turboNormal[1]*cartesianVelocity[0];

if (marker_flag == INFLOW) {
turboVelocity[0] *= -1.0;
turboVelocity[1] *= -1.0;
}
if(nDim == 3)
turboVelocity[2] = cartesianVelocity[2];
if (nDim == 3) turboVelocity[2] = cartesianVelocity[2];
}
}

Expand All @@ -1156,22 +1166,21 @@ class CEulerSolver : public CFVMFlowSolverBase<CEulerVariable, ENUM_REGIME::COMP
unsigned short marker_flag,
TURBOMACHINERY_TYPE kind_turb){

if ((kind_turb == TURBOMACHINERY_TYPE::AXIAL && nDim == 3) || (kind_turb == TURBOMACHINERY_TYPE::CENTRIPETAL_AXIAL && marker_flag == OUTFLOW) || (kind_turb == TURBOMACHINERY_TYPE::AXIAL_CENTRIFUGAL && marker_flag == INFLOW)){
if ((kind_turb == TURBOMACHINERY_TYPE::AXIAL && nDim == 3) ||
(kind_turb == TURBOMACHINERY_TYPE::CENTRIPETAL_AXIAL && marker_flag == OUTFLOW) ||
(kind_turb == TURBOMACHINERY_TYPE::AXIAL_CENTRIFUGAL && marker_flag == INFLOW)) {
cartesianVelocity[0] = turboVelocity[2]*turboNormal[0] - turboVelocity[1]*turboNormal[1];
cartesianVelocity[1] = turboVelocity[2]*turboNormal[1] + turboVelocity[1]*turboNormal[0];
cartesianVelocity[2] = turboVelocity[0];
}
else{
cartesianVelocity[0] = turboVelocity[0]*turboNormal[0] - turboVelocity[1]*turboNormal[1];
cartesianVelocity[1] = turboVelocity[0]*turboNormal[1] + turboVelocity[1]*turboNormal[0];
} else {
cartesianVelocity[0] = turboVelocity[0]*turboNormal[0] - turboVelocity[1]*turboNormal[1];
cartesianVelocity[1] = turboVelocity[0]*turboNormal[1] + turboVelocity[1]*turboNormal[0];

if (marker_flag == INFLOW){
if (marker_flag == INFLOW) {
cartesianVelocity[0] *= -1.0;
cartesianVelocity[1] *= -1.0;
}

if(nDim == 3)
cartesianVelocity[2] = turboVelocity[2];
if (nDim == 3) cartesianVelocity[2] = turboVelocity[2];
}
}

Expand Down
44 changes: 10 additions & 34 deletions SU2_CFD/include/solvers/CFVMFlowSolverBase.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -345,7 +345,9 @@ class CFVMFlowSolverBase : public CSolver {
* \brief Compute a suitable under-relaxation parameter to limit the change in the solution variables over a nonlinear
* iteration for stability.
*/
virtual void ComputeUnderRelaxationFactor(const CConfig* config);
virtual void ComputeUnderRelaxationFactor(const CConfig* config) {
SU2_MPI::Error("Not implemented for this solver.", CURRENT_FUNCTION);
}

/*!
* \brief General implementation to load a flow solution from a restart file.
Expand Down Expand Up @@ -976,39 +978,6 @@ class CFVMFlowSolverBase : public CSolver {

}

/*!
* \brief Generic implementation to complete an implicit iteration, i.e. update the solution.
* \tparam compute_ur - Whether to use automatic under-relaxation for the update.
*/
template<bool compute_ur>
void CompleteImplicitIteration_impl(CGeometry *geometry, CConfig *config) {

if (compute_ur) ComputeUnderRelaxationFactor(config);

/*--- Update solution with under-relaxation and communicate it. ---*/

if (!config->GetContinuous_Adjoint()) {
SU2_OMP_FOR_STAT(omp_chunk_size)
for (unsigned long iPoint = 0; iPoint < nPointDomain; iPoint++) {
for (unsigned short iVar = 0; iVar < nVar; iVar++) {
nodes->AddSolution(iPoint, iVar, nodes->GetUnderRelaxation(iPoint)*LinSysSol[iPoint*nVar+iVar]);
}
}
END_SU2_OMP_FOR
}

for (unsigned short iPeriodic = 1; iPeriodic <= config->GetnMarker_Periodic()/2; iPeriodic++) {
InitiatePeriodicComms(geometry, config, iPeriodic, PERIODIC_IMPLICIT);
CompletePeriodicComms(geometry, config, iPeriodic, PERIODIC_IMPLICIT);
}

InitiateComms(geometry, config, MPI_QUANTITIES::SOLUTION);
CompleteComms(geometry, config, MPI_QUANTITIES::SOLUTION);

/*--- For verification cases, compute the global error metrics. ---*/
ComputeVerificationError(geometry, config);
}

/*!
* \brief Evaluate the vorticity and strain rate magnitude.
*/
Expand Down Expand Up @@ -1075,6 +1044,13 @@ class CFVMFlowSolverBase : public CSolver {
*/
void ImplicitEuler_Iteration(CGeometry *geometry, CSolver **solver_container, CConfig *config) final;

/*!
* \brief Complete an implicit iteration.
* \param[in] geometry - Geometrical definition of the problem.
* \param[in] config - Definition of the particular problem.
*/
void CompleteImplicitIteration(CGeometry *geometry, CSolver**, CConfig *config) final;

/*!
* \brief Set the total residual adding the term that comes from the Dual Time Strategy.
* \param[in] geometry - Geometrical definition of the problem.
Expand Down
45 changes: 18 additions & 27 deletions SU2_CFD/include/solvers/CFVMFlowSolverBase.inl
Original file line number Diff line number Diff line change
Expand Up @@ -596,42 +596,33 @@ void CFVMFlowSolverBase<V, R>::ComputeVerificationError(CGeometry* geometry, CCo
}

template <class V, ENUM_REGIME R>
void CFVMFlowSolverBase<V, R>::ComputeUnderRelaxationFactor(const CConfig* config) {
void CFVMFlowSolverBase<V, R>::CompleteImplicitIteration(CGeometry *geometry, CSolver**, CConfig *config) {
SU2_ZONE_SCOPED

/* Loop over the solution update given by relaxing the linear
system for this nonlinear iteration. */
if constexpr (R == ENUM_REGIME::COMPRESSIBLE) ComputeUnderRelaxationFactor(config);

const su2double allowableRatio = config->GetMaxUpdateFractionFlow();
/*--- Update solution with under-relaxation and communicate it. ---*/

SU2_OMP_FOR_STAT(omp_chunk_size)
for (unsigned long iPoint = 0; iPoint < nPointDomain; iPoint++) {
su2double localUnderRelaxation = 1.0;

for (unsigned short iVar = 0; iVar < nVar; iVar++) {
/* We impose a limit on the maximum percentage that the
density and energy can change over a nonlinear iteration. */

if ((iVar == 0) || (iVar == nVar - 1)) {
const unsigned long index = iPoint * nVar + iVar;
su2double ratio = fabs(LinSysSol[index]) / (fabs(nodes->GetSolution(iPoint, iVar)) + EPS);
if (ratio > allowableRatio) {
localUnderRelaxation = min(allowableRatio / ratio, localUnderRelaxation);
}
if (!config->GetContinuous_Adjoint()) {
SU2_OMP_FOR_STAT(omp_chunk_size)
for (unsigned long iPoint = 0; iPoint < nPointDomain; iPoint++) {
for (unsigned short iVar = 0; iVar < nVar; iVar++) {
nodes->AddSolution(iPoint, iVar, nodes->GetUnderRelaxation(iPoint) * LinSysSol(iPoint, iVar));
}
}
END_SU2_OMP_FOR
}

/* Threshold the relaxation factor in the event that there is
a very small value. This helps avoid catastrophic crashes due
to non-realizable states by canceling the update. */

if (localUnderRelaxation < 1e-10) localUnderRelaxation = 0.0;
for (unsigned short iPeriodic = 1; iPeriodic <= config->GetnMarker_Periodic()/2; iPeriodic++) {
InitiatePeriodicComms(geometry, config, iPeriodic, PERIODIC_IMPLICIT);
CompletePeriodicComms(geometry, config, iPeriodic, PERIODIC_IMPLICIT);
}

/* Store the under-relaxation factor for this point. */
InitiateComms(geometry, config, MPI_QUANTITIES::SOLUTION);
CompleteComms(geometry, config, MPI_QUANTITIES::SOLUTION);

nodes->SetUnderRelaxation(iPoint, localUnderRelaxation);
}
END_SU2_OMP_FOR
/*--- For verification cases, compute the global error metrics. ---*/
ComputeVerificationError(geometry, config);
}

template <class V, ENUM_REGIME R>
Expand Down
7 changes: 0 additions & 7 deletions SU2_CFD/include/solvers/CIncEulerSolver.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -365,13 +365,6 @@ class CIncEulerSolver : public CFVMFlowSolverBase<CIncEulerVariable, ENUM_REGIME
*/
void PrepareImplicitIteration(CGeometry *geometry, CSolver**, CConfig *config) final;

/*!
* \brief Complete an implicit iteration.
* \param[in] geometry - Geometrical definition of the problem.
* \param[in] config - Definition of the particular problem.
*/
void CompleteImplicitIteration(CGeometry *geometry, CSolver**, CConfig *config) final;

/*!
* \brief Set the total residual adding the term that comes from the Dual Time Strategy.
* \param[in] geometry - Geometrical definition of the problem.
Expand Down
7 changes: 0 additions & 7 deletions SU2_CFD/include/solvers/CNEMOEulerSolver.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -329,13 +329,6 @@ class CNEMOEulerSolver : public CFVMFlowSolverBase<CNEMOEulerVariable, ENUM_REGI
*/
void PrepareImplicitIteration(CGeometry *geometry, CSolver**, CConfig *config) final;

/*!
* \brief Complete an implicit iteration.
* \param[in] geometry - Geometrical definition of the problem.
* \param[in] config - Definition of the particular problem.
*/
void CompleteImplicitIteration(CGeometry *geometry, CSolver**, CConfig *config) final;

/*!
* \brief Print verification error to screen.
* \param[in] config - Definition of the particular problem.
Expand Down
Loading
Loading