From a0aa071d88a1ceea8da0a1ccd608e18130a65b28 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:09:09 +0200 Subject: [PATCH 01/14] Add the incompressible grid-motion Jacobian once, with the variable-density terms FDS added the grid-motion term inside a loop over the variables (nVar times). The term is now added once, as the exact derivative with respect to (p, u, h), including d(rho)/dh; the same terms are added to the incompressible JST and Lax-Friedrichs schemes. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- .../src/numerics/flow/convection/centered.cpp | 48 ++++++++++++------- SU2_CFD/src/numerics/flow/convection/fds.cpp | 24 ++++++---- 2 files changed, 46 insertions(+), 26 deletions(-) diff --git a/SU2_CFD/src/numerics/flow/convection/centered.cpp b/SU2_CFD/src/numerics/flow/convection/centered.cpp index b0130279f76f..54fdbdfdd8f9 100644 --- a/SU2_CFD/src/numerics/flow/convection/centered.cpp +++ b/SU2_CFD/src/numerics/flow/convection/centered.cpp @@ -161,17 +161,25 @@ CNumerics::ResidualType<> CCentLaxInc_Flow::ComputeResidual(const CConfig* confi ProjVelocity += 0.5*(GridVel_i[iDim]+GridVel_j[iDim])*Normal[iDim]; /*--- Residual contributions ---*/ - for (iVar = 0; iVar < nVar; iVar++) { + for (iVar = 0; iVar < nVar; iVar++) ProjFlux[iVar] -= ProjVelocity * 0.5*(U_i[iVar]+U_j[iVar]); - /*--- Jacobian contributions ---*/ - /*--- Implicit terms ---*/ - if (implicit) { - for (iDim = 0; iDim < nDim+1; iDim++){ - Jacobian_i[iDim+1][iDim+1] -= 0.5*ProjVelocity*DensityInc_i; - Jacobian_j[iDim+1][iDim+1] -= 0.5*ProjVelocity*DensityInc_j; + /*--- Jacobian of -ProjVelocity * 0.5*U, with U = (rho, rho*u, rho*h), w.r.t. the primitives (p, u, h), added + once. The density depends on h (variable density models) and not on the pressure. ---*/ + if (implicit) { + auto addGridJacobian = [&](su2double **Jac, su2double rho, const su2double *vel, su2double h, su2double T, + su2double cp) { + const su2double dRhodh = variable_density ? su2double(-rho / (T * cp)) : su2double(0.0); + const su2double factor = 0.5*ProjVelocity; + Jac[0][nDim+1] -= factor*dRhodh; + for (iDim = 0; iDim < nDim; iDim++) { + Jac[iDim+1][iDim+1] -= factor*rho; + Jac[iDim+1][nDim+1] -= factor*vel[iDim]*dRhodh; } - } + Jac[nDim+1][nDim+1] -= factor*(rho + h*dRhodh); + }; + addGridJacobian(Jacobian_i, DensityInc_i, Velocity_i, Enthalpy_i, Temperature_i, Cp_i); + addGridJacobian(Jacobian_j, DensityInc_j, Velocity_j, Enthalpy_j, Temperature_j, Cp_j); } } @@ -401,17 +409,25 @@ CNumerics::ResidualType<> CCentJSTInc_Flow::ComputeResidual(const CConfig* confi ProjVelocity += 0.5*(GridVel_i[iDim]+GridVel_j[iDim])*Normal[iDim]; /*--- Residual contributions ---*/ - for (iVar = 0; iVar < nVar; iVar++) { + for (iVar = 0; iVar < nVar; iVar++) ProjFlux[iVar] -= ProjVelocity * 0.5*(U_i[iVar]+U_j[iVar]); - /*--- Jacobian contributions ---*/ - /*--- Implicit terms ---*/ - if (implicit) { - for (iDim = 0; iDim < nDim+1; iDim++){ - Jacobian_i[iDim+1][iDim+1] -= 0.5*ProjVelocity*DensityInc_i; - Jacobian_j[iDim+1][iDim+1] -= 0.5*ProjVelocity*DensityInc_j; + /*--- Jacobian of -ProjVelocity * 0.5*U, with U = (rho, rho*u, rho*h), w.r.t. the primitives (p, u, h), added + once. The density depends on h (variable density models) and not on the pressure. ---*/ + if (implicit) { + auto addGridJacobian = [&](su2double **Jac, su2double rho, const su2double *vel, su2double h, su2double T, + su2double cp) { + const su2double dRhodh = variable_density ? su2double(-rho / (T * cp)) : su2double(0.0); + const su2double factor = 0.5*ProjVelocity; + Jac[0][nDim+1] -= factor*dRhodh; + for (iDim = 0; iDim < nDim; iDim++) { + Jac[iDim+1][iDim+1] -= factor*rho; + Jac[iDim+1][nDim+1] -= factor*vel[iDim]*dRhodh; } - } + Jac[nDim+1][nDim+1] -= factor*(rho + h*dRhodh); + }; + addGridJacobian(Jacobian_i, DensityInc_i, Velocity_i, Enthalpy_i, Temperature_i, Cp_i); + addGridJacobian(Jacobian_j, DensityInc_j, Velocity_j, Enthalpy_j, Temperature_j, Cp_j); } } diff --git a/SU2_CFD/src/numerics/flow/convection/fds.cpp b/SU2_CFD/src/numerics/flow/convection/fds.cpp index e806127da203..57c7ed9f6f5a 100644 --- a/SU2_CFD/src/numerics/flow/convection/fds.cpp +++ b/SU2_CFD/src/numerics/flow/convection/fds.cpp @@ -249,19 +249,23 @@ CNumerics::ResidualType<> CUpwFDSInc_Flow::ComputeResidual(const CConfig *config ProjVelocity += 0.5*(GridVel_i[iDim]+GridVel_j[iDim])*Normal[iDim]; /*--- Residual contributions ---*/ - for (iVar = 0; iVar < nVar; iVar++) { + for (iVar = 0; iVar < nVar; iVar++) Flux[iVar] -= ProjVelocity * 0.5*(U_i[iVar]+U_j[iVar]); - /*--- Jacobian contributions ---*/ - /*--- Implicit terms ---*/ - if (implicit) { - for (iDim = 0; iDim < nDim; iDim++){ - Jacobian_i[iDim+1][iDim+1] -= 0.5*ProjVelocity*DensityInc_i; - Jacobian_j[iDim+1][iDim+1] -= 0.5*ProjVelocity*DensityInc_j; + /*--- Jacobian of -ProjVelocity * 0.5*U, with U = (rho, rho*u, rho*h), w.r.t. the primitives (p, u, h), added + once. The density depends on h (variable density models) and not on the pressure. ---*/ + if (implicit) { + auto addGridJacobian = [&](su2double **Jac, su2double rho, const su2double *vel, su2double h, su2double dRhodh) { + const su2double factor = 0.5*ProjVelocity; + Jac[0][nDim+1] -= factor*dRhodh; + for (iDim = 0; iDim < nDim; iDim++) { + Jac[iDim+1][iDim+1] -= factor*rho; + Jac[iDim+1][nDim+1] -= factor*vel[iDim]*dRhodh; } - Jacobian_i[nDim+1][nDim+1] -= 0.5*ProjVelocity*DensityInc_i; - Jacobian_j[nDim+1][nDim+1] -= 0.5*ProjVelocity*DensityInc_j; - } + Jac[nDim+1][nDim+1] -= factor*(rho + h*dRhodh); + }; + addGridJacobian(Jacobian_i, DensityInc_i, Velocity_i, Enthalpy_i, dRhodh_i); + addGridJacobian(Jacobian_j, DensityInc_j, Velocity_j, Enthalpy_j, dRhodh_j); } } From d45339581c65b8ee01057d897f9ab67db43a1820 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:09:09 +0200 Subject: [PATCH 02/14] Reduce the minimum pressure with MPI_MIN SetRangePressure reduced the local minima with MPI_MAX, so with more than one rank MinP was the largest of the local minima (used by the mass-flow outlet). Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- SU2_CFD/src/solvers/CIncEulerSolver.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/SU2_CFD/src/solvers/CIncEulerSolver.cpp b/SU2_CFD/src/solvers/CIncEulerSolver.cpp index 512a99a97dc6..8c6ee1719a01 100644 --- a/SU2_CFD/src/solvers/CIncEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CIncEulerSolver.cpp @@ -2235,7 +2235,7 @@ void CIncEulerSolver::SetRangePressure(CGeometry *geometry, CSolver **solver_con BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS { minP = MinP; - SU2_MPI::Allreduce(&minP, &MinP, 1, MPI_DOUBLE, MPI_MAX, SU2_MPI::GetComm()); + SU2_MPI::Allreduce(&minP, &MinP, 1, MPI_DOUBLE, MPI_MIN, SU2_MPI::GetComm()); maxP = MaxP; SU2_MPI::Allreduce(&maxP, &MaxP, 1, MPI_DOUBLE, MPI_MAX, SU2_MPI::GetComm()); From 62d12060d152902b3461ebdf839d9c0332612a6c Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:09:09 +0200 Subject: [PATCH 03/14] Reject combinations of source terms that share one numerics object The flow solvers build one numerics object per source slot, so two active sources of the same slot overwrote each other (e.g. body force + gravity gave twice the body force). First slot: body force, rotating frame, axisymmetric, gravity, Boussinesq, streamwise periodic. Second slot: radiation, periodic outlet, vorticity confinement. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- Common/src/CConfig.cpp | 30 ++++++++++++++++++++++++++++++ 1 file changed, 30 insertions(+) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index fd4207b0c5fe..bfbd230b1f1d 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -4598,6 +4598,36 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i Rotating_Frame = (Kind_GridMovement == ROTATING_FRAME); + /*--- The finite volume flow solvers apply these source terms with a single numerics object, which is built + for only one of them: combinations would apply one source several times and ignore the others. Only these + solvers (and their discrete adjoints) are checked: e.g. elasticity uses BODY_FORCE and GRAVITY_FORCE together. ---*/ + + if (GetFluidProblem() && !GetNEMOProblem()) { + const bool inc = (Kind_Regime == ENUM_REGIME::INCOMPRESSIBLE); + const int nSources = int(Body_Force) + int(Rotating_Frame) + int(Axisymmetric) + int(GravityForce && !inc) + + int(inc && Kind_DensityModel == INC_DENSITYMODEL::BOUSSINESQ) + + int(inc && Kind_Streamwise_Periodic != ENUM_STREAMWISE_PERIODIC::NONE); + if (nSources > 1) { + SU2_MPI::Error("Only one of BODY_FORCE, GRID_MOVEMENT= ROTATING_FRAME, AXISYMMETRIC, GRAVITY_FORCE (compressible), " + "INC_DENSITY_MODEL= BOUSSINESQ and KIND_STREAMWISE_PERIODIC (incompressible) can be used at a time.", + CURRENT_FUNCTION); + } + } + + /*--- The second source object, for radiation, the energy correction of incompressible streamwise periodic flow + (without periodic temperature) or vorticity confinement, is also built for only one of them. ---*/ + + if (GetFluidProblem() && !GetNEMOProblem()) { + const bool inc = (Kind_Regime == ENUM_REGIME::INCOMPRESSIBLE); + const bool periodicOutlet = inc && Kind_Streamwise_Periodic != ENUM_STREAMWISE_PERIODIC::NONE && Energy_Equation && + !Streamwise_Periodic_Temperature; + if (int(Radiation) + int(periodicOutlet) + int(VorticityConfinement) > 1) { + SU2_MPI::Error("Only one of RADIATION, the energy correction of KIND_STREAMWISE_PERIODIC (incompressible, " + "STREAMWISE_PERIODIC_TEMPERATURE= NO) and VORTICITY_CONFINEMENT can be used at a time.", + CURRENT_FUNCTION); + } + } + /*--- In case the grid movement parameters have not been declared in the config file, set them equal to zero for safety. Also check to make sure that for each option, a value has been declared for each moving marker. ---*/ From f390fbda77e648dcd990e58de6dacd96be1a0759 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:09:09 +0200 Subject: [PATCH 04/14] Read the Riemann boundary kind in BC_TurboRiemann It read the Giles kind of the marker (copied from BC_Giles), which reads past the end of the Giles list and crashes when the case has no Giles marker. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- 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 32c265e55656..e107d1ed7072 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -5797,7 +5797,7 @@ void CEulerSolver::BC_TurboRiemann(CGeometry *geometry, CSolver **solver_contain ProjVelocity_i += Velocity_i[iDim]*UnitNormal[iDim]; su2double donorAverages[5] = {0.0}; - switch (config->GetKind_Data_Giles(Marker_Tag)){ + switch (config->GetKind_Data_Riemann(Marker_Tag)){ case MIXING_IN: case MIXING_IN_1D: case MIXING_OUT: case MIXING_OUT_1D: for (auto mixVar = 0u; mixVar < 5; mixVar++) donorAverages[mixVar] = GetMixingState(val_marker, iSpan, mixVar); break; From 7b32ea93c61b798cd8dd950958cbeb5bb6293038 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:09:09 +0200 Subject: [PATCH 05/14] Check the number of variables in the Green-Gauss gradient The boundary flux array had 20 entries and no bound check (NEMO with 10 or more species in 3D). Use 32, as the limiters, and stop with an error above it. The error of both checks is raised by the master thread only, since these functions run inside OpenMP parallel regions (MPI_THREAD_FUNNELED). Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- SU2_CFD/include/gradients/computeGradientsGreenGauss.hpp | 9 ++++++++- SU2_CFD/include/limiters/computeLimiters_impl.hpp | 6 +++++- 2 files changed, 13 insertions(+), 2 deletions(-) diff --git a/SU2_CFD/include/gradients/computeGradientsGreenGauss.hpp b/SU2_CFD/include/gradients/computeGradientsGreenGauss.hpp index f60c0422f0cc..89eb028c6aeb 100644 --- a/SU2_CFD/include/gradients/computeGradientsGreenGauss.hpp +++ b/SU2_CFD/include/gradients/computeGradientsGreenGauss.hpp @@ -69,7 +69,14 @@ void computeGradientsGreenGauss(CSolver* solver, MPI_QUANTITIES kindMpiComm, PER const auto chunkSize = computeStaticChunkSize(nPointDomain, omp_get_max_threads(), OMP_MAX_CHUNK); #endif - static constexpr size_t MAXNVAR = 20; + static constexpr size_t MAXNVAR = 32; + + /*--- Called by all threads of a parallel region: only the master may call MPI (MPI_THREAD_FUNNELED). ---*/ + if (varEnd > MAXNVAR) { + BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS + SU2_MPI::Error("Number of variables is too large, increase MAXNVAR.", CURRENT_FUNCTION); + END_SU2_OMP_SAFE_GLOBAL_ACCESS + } /*--- For each (non-halo) volume integrate over its faces (edges). ---*/ diff --git a/SU2_CFD/include/limiters/computeLimiters_impl.hpp b/SU2_CFD/include/limiters/computeLimiters_impl.hpp index 7d8beafb03e0..0cdc60207a27 100644 --- a/SU2_CFD/include/limiters/computeLimiters_impl.hpp +++ b/SU2_CFD/include/limiters/computeLimiters_impl.hpp @@ -80,8 +80,12 @@ void computeLimiters_impl(CSolver* solver, { constexpr size_t MAXNVAR = 32; - if (varEnd > MAXNVAR) + /*--- Called by all threads of a parallel region: only the master may call MPI (MPI_THREAD_FUNNELED). ---*/ + if (varEnd > MAXNVAR) { + BEGIN_SU2_OMP_SAFE_GLOBAL_ACCESS SU2_MPI::Error("Number of variables is too large, increase MAXNVAR.", CURRENT_FUNCTION); + END_SU2_OMP_SAFE_GLOBAL_ACCESS + } const size_t nPointDomain = geometry.GetnPointDomain(); const size_t nPoint = geometry.GetnPoint(); From 9eb34f11a4d8f2891f1f7bd4a99f3e5118c3c69b Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:09:10 +0200 Subject: [PATCH 06/14] Clamp MGLEVEL to 10 The multigrid arrays hold at most 10 coarse levels (MAX_MG_LEVELS); MGLEVEL= 12 with MG_MIN_MESHSIZE= 1 crashed with stack smashing. Clamp with a warning, MG_MIN_MESHSIZE can still reduce the levels. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- Common/src/CConfig.cpp | 7 +++++++ 1 file changed, 7 insertions(+) diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index bfbd230b1f1d..261002df1f4f 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -4987,6 +4987,13 @@ void CConfig::SetPostprocessing(SU2_COMPONENT val_software, unsigned short val_i if (Restart || ((Kind_MGCycle == MG_CYCLE::FULL) && ContinuousAdjoint)) Kind_MGCycle = MG_CYCLE::V; + /*--- The multigrid arrays hold at most 10 coarse levels (CMultiGridIntegration::MAX_MG_LEVELS). Clamp instead of + stopping, MG_MIN_MESHSIZE can still reduce the number of levels during the agglomeration. ---*/ + if (nMGLevels > 10) { + if (rank == MASTER_NODE) cout << "WARNING: MGLEVEL > 10 is not supported, using MGLEVEL= 10." << endl; + nMGLevels = 10; + } + FinestMesh = MESH_0; if (Kind_MGCycle == MG_CYCLE::FULL) FinestMesh = nMGLevels; From ab210631a4277f224f3fecc7783ff60f8caf4ab6 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:09:19 +0200 Subject: [PATCH 07/14] Clear the whole Jacobian in the non-physical guard of the general-gas Roe scheme Only the diagonal was zeroed (index [iVar][iVar] inside the jVar loop), so the off-diagonal entries kept the values of the previous edge. The guard of CUpwRoeBase_Flow already clears all entries. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- SU2_CFD/src/numerics/flow/convection/roe.cpp | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/SU2_CFD/src/numerics/flow/convection/roe.cpp b/SU2_CFD/src/numerics/flow/convection/roe.cpp index 73114576d409..371e23059d20 100644 --- a/SU2_CFD/src/numerics/flow/convection/roe.cpp +++ b/SU2_CFD/src/numerics/flow/convection/roe.cpp @@ -811,8 +811,8 @@ CNumerics::ResidualType<> CUpwGeneralRoe_Flow::ComputeResidual(const CConfig* co for (iVar = 0; iVar < nVar; iVar++) { Flux[iVar] = 0.0; for (jVar = 0; jVar < nVar; jVar++) { - Jacobian_i[iVar][iVar] = 0.0; - Jacobian_j[iVar][iVar] = 0.0; + Jacobian_i[iVar][jVar] = 0.0; + Jacobian_j[iVar][jVar] = 0.0; } } AD::SetPreaccOut(Flux, nVar); From 4d81f81ba7bd2afe9866b85c82faab6eaceaf0a4 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:09:20 +0200 Subject: [PATCH 08/14] Use 0.5 instead of the integer division 1/2 in the axisymmetric source Jacobian The term 1/2*(Gamma-1) of d(v (rho E + p))/d(rho v) was 0 (implicit Jacobian only). Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- SU2_CFD/src/numerics/flow/flow_sources.cpp | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/SU2_CFD/src/numerics/flow/flow_sources.cpp b/SU2_CFD/src/numerics/flow/flow_sources.cpp index 8cf9642b4251..0a3a2d0954b9 100644 --- a/SU2_CFD/src/numerics/flow/flow_sources.cpp +++ b/SU2_CFD/src/numerics/flow/flow_sources.cpp @@ -101,7 +101,7 @@ CNumerics::ResidualType<> CSourceAxisymmetric_Flow::ComputeResidual(const CConfi jacobian[3][0] = -Gamma*U_i[2]*U_i[3]/(U_i[0]*U_i[0]) + (Gamma-1)*U_i[2]*(U_i[1]*U_i[1]+U_i[2]*U_i[2])/(U_i[0]*U_i[0]*U_i[0]); jacobian[3][1] = -(Gamma-1)*U_i[2]*U_i[1]/(U_i[0]*U_i[0]); - jacobian[3][2] = Gamma*U_i[3]/U_i[0] - 1/2*(Gamma-1)*( (U_i[1]*U_i[1]+U_i[2]*U_i[2])/(U_i[0]*U_i[0]) + 2*U_i[2]*U_i[2]/(U_i[0]*U_i[0]) ); + jacobian[3][2] = Gamma*U_i[3]/U_i[0] - 0.5*(Gamma-1)*( (U_i[1]*U_i[1]+U_i[2]*U_i[2])/(U_i[0]*U_i[0]) + 2*U_i[2]*U_i[2]/(U_i[0]*U_i[0]) ); jacobian[3][3] = Gamma*U_i[2]/U_i[0]; for (iVar=0; iVar < nVar; iVar++) From cf8719ab47321c6d3958f46f7fde1739ec979151 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:09:55 +0200 Subject: [PATCH 09/14] Use wave strengths that match the 3D P matrix in L2Roe, LMRoe and the explicit general-gas Roe In 3D the first three columns of the P matrix combine the entropy wave, weighted by n_j, with the shear waves. Their strengths are n_j (Delta rho - Delta p / c^2) + (Delta u x n)_j; the entropy jump was used unweighted in the first one only, so a contact or shear layer on a face not aligned with x got a wrong flux. The standard Roe scheme already projects the jump with P^-1 since #299. Unit test roe_contact_tests.cpp: a contact and a shear layer give the upwind flux in 2D and 3D. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- SU2_CFD/src/numerics/flow/convection/roe.cpp | 27 ++-- .../SU2_CFD/numerics/roe_contact_tests.cpp | 125 ++++++++++++++++++ UnitTests/meson.build | 1 + 3 files changed, 144 insertions(+), 9 deletions(-) create mode 100644 UnitTests/SU2_CFD/numerics/roe_contact_tests.cpp diff --git a/SU2_CFD/src/numerics/flow/convection/roe.cpp b/SU2_CFD/src/numerics/flow/convection/roe.cpp index 371e23059d20..37c697515cfb 100644 --- a/SU2_CFD/src/numerics/flow/convection/roe.cpp +++ b/SU2_CFD/src/numerics/flow/convection/roe.cpp @@ -334,9 +334,12 @@ void CUpwL2Roe_Flow::FinalizeResidual(su2double *val_residual, su2double **val_J delta_wave[2] = proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); delta_wave[3] = -proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); } else { - delta_wave[0] = delta_rho - delta_p/RoeSoundSpeed2; - delta_wave[1] = (UnitNormal[0]*delta_vel[2]-UnitNormal[2]*delta_vel[0])*zeta; - delta_wave[2] = (UnitNormal[1]*delta_vel[0]-UnitNormal[0]*delta_vel[1])*zeta; + /*--- The first three columns of the 3D P matrix combine the entropy wave (weighted by the normal) with the + shear waves, so their strengths are n_j (Delta rho - Delta p / c^2) + (Delta u x n)_j. ---*/ + const su2double delta_entropy = delta_rho - delta_p/RoeSoundSpeed2; + delta_wave[0] = UnitNormal[0]*delta_entropy + (delta_vel[1]*UnitNormal[2]-delta_vel[2]*UnitNormal[1])*zeta; + delta_wave[1] = UnitNormal[1]*delta_entropy + (delta_vel[2]*UnitNormal[0]-delta_vel[0]*UnitNormal[2])*zeta; + delta_wave[2] = UnitNormal[2]*delta_entropy + (delta_vel[0]*UnitNormal[1]-delta_vel[1]*UnitNormal[0])*zeta; delta_wave[3] = proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); delta_wave[4] = -proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); } @@ -407,9 +410,12 @@ void CUpwLMRoe_Flow::FinalizeResidual(su2double *val_residual, su2double **val_J delta_wave[2] = proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); delta_wave[3] = -proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); } else { - delta_wave[0] = delta_rho - delta_p/RoeSoundSpeed2; - delta_wave[1] = (UnitNormal[0]*delta_vel[2]-UnitNormal[2]*delta_vel[0]); - delta_wave[2] = (UnitNormal[1]*delta_vel[0]-UnitNormal[0]*delta_vel[1]); + /*--- The first three columns of the 3D P matrix combine the entropy wave (weighted by the normal) with the + shear waves, so their strengths are n_j (Delta rho - Delta p / c^2) + (Delta u x n)_j. ---*/ + const su2double delta_entropy = delta_rho - delta_p/RoeSoundSpeed2; + delta_wave[0] = UnitNormal[0]*delta_entropy + (delta_vel[1]*UnitNormal[2]-delta_vel[2]*UnitNormal[1]); + delta_wave[1] = UnitNormal[1]*delta_entropy + (delta_vel[2]*UnitNormal[0]-delta_vel[0]*UnitNormal[2]); + delta_wave[2] = UnitNormal[2]*delta_entropy + (delta_vel[0]*UnitNormal[1]-delta_vel[1]*UnitNormal[0]); delta_wave[3] = proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); delta_wave[4] = -proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); } @@ -901,9 +907,12 @@ CNumerics::ResidualType<> CUpwGeneralRoe_Flow::ComputeResidual(const CConfig* co delta_wave[2] = proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); delta_wave[3] = -proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); } else { - delta_wave[0] = delta_rho - delta_p/(RoeSoundSpeed*RoeSoundSpeed); - delta_wave[1] = UnitNormal[0]*delta_vel[2]-UnitNormal[2]*delta_vel[0]; - delta_wave[2] = UnitNormal[1]*delta_vel[0]-UnitNormal[0]*delta_vel[1]; + /*--- The first three columns of the 3D P matrix combine the entropy wave (weighted by the normal) with the + shear waves, so their strengths are n_j (Delta rho - Delta p / c^2) + (Delta u x n)_j. ---*/ + const su2double delta_entropy = delta_rho - delta_p/(RoeSoundSpeed*RoeSoundSpeed); + delta_wave[0] = UnitNormal[0]*delta_entropy + delta_vel[1]*UnitNormal[2]-delta_vel[2]*UnitNormal[1]; + delta_wave[1] = UnitNormal[1]*delta_entropy + delta_vel[2]*UnitNormal[0]-delta_vel[0]*UnitNormal[2]; + delta_wave[2] = UnitNormal[2]*delta_entropy + delta_vel[0]*UnitNormal[1]-delta_vel[1]*UnitNormal[0]; delta_wave[3] = proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); delta_wave[4] = -proj_delta_vel + delta_p/(RoeDensity*RoeSoundSpeed); } diff --git a/UnitTests/SU2_CFD/numerics/roe_contact_tests.cpp b/UnitTests/SU2_CFD/numerics/roe_contact_tests.cpp new file mode 100644 index 000000000000..9c99a461f218 --- /dev/null +++ b/UnitTests/SU2_CFD/numerics/roe_contact_tests.cpp @@ -0,0 +1,125 @@ +/*! + * \file roe_contact_tests.cpp + * \brief Roe-type schemes must resolve a contact discontinuity and a shear layer exactly (upwind flux), in 2D + * and 3D, on faces not aligned with the axes. + * \version 8.5.0 "Harrier" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2026, SU2 Contributors (cf. AUTHORS.md) + * + * SU2 is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public + * License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * + * SU2 is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with SU2. If not, see . + */ + +#include "catch.hpp" +#include +#include +#include +#include "../../../SU2_CFD/include/numerics/flow/convection/roe.hpp" + +namespace { + +constexpr su2double gamma = 1.4; + +std::unique_ptr MakeConfig(bool implicit) { + std::stringstream options; + options << "SOLVER= EULER\nCONV_NUM_METHOD_FLOW= ROE\nENTROPY_FIX_COEFF= 0.0\nROE_KAPPA= 0.5\n"; + options << "TIME_DISCRE_FLOW= " << (implicit ? "EULER_IMPLICIT" : "EULER_EXPLICIT") << "\n"; + return std::make_unique(options, SU2_COMPONENT::SU2_CFD, false); +} + +/*--- Primitive (T, u, v, [w], p, rho, h, c) and secondary (dp/drho_e, dp/de_rho) variables of an ideal gas. ---*/ +void State(unsigned short nDim, su2double rho, const su2double* vel, su2double p, su2double* V, su2double* S) { + su2double q2 = 0.0; + for (unsigned short iDim = 0; iDim < nDim; iDim++) { + V[iDim + 1] = vel[iDim]; + q2 += vel[iDim] * vel[iDim]; + } + V[0] = p / (rho * 287.058); + V[nDim + 1] = p; + V[nDim + 2] = rho; + V[nDim + 3] = gamma / (gamma - 1) * p / rho + 0.5 * q2; + V[nDim + 4] = sqrt(gamma * p / rho); + S[0] = p / rho; + S[1] = (gamma - 1) * rho; +} + +template +std::unique_ptr MakeScheme(unsigned short nDim, const CConfig* config) { + return std::make_unique(nDim, nDim + 2, config); +} +template <> +std::unique_ptr MakeScheme(unsigned short nDim, const CConfig* config) { + return std::make_unique(nDim, nDim + 2, config, false); +} + +/*--- Two states with the same pressure and normal velocity, different density and (with shear) tangential + velocity: a contact discontinuity, whose exact (upwind) flux is the flux of the upstream state. ---*/ +template +void CheckContact(const char* name, unsigned short nDim, const su2double* vel, bool implicit, bool shear) { + auto config = MakeConfig(implicit); + auto numerics = MakeScheme(nDim, config.get()); + + const su2double p = 1.0, rho_i = 1.4, rho_j = 0.7; + const su2double n2[2] = {0.6, 0.8}, n3[3] = {2.0 / 7, 3.0 / 7, 6.0 / 7}; + const su2double* n = (nDim == 2) ? n2 : n3; + const su2double t2[2] = {0.16, -0.12}, t3[3] = {0.15, -0.1, 0.0}; // tangential velocity jumps + const su2double* t = (nDim == 2) ? t2 : t3; + su2double vel_j[3] = {0.0}; + for (unsigned short iDim = 0; iDim < nDim; iDim++) vel_j[iDim] = vel[iDim] + (shear ? t[iDim] : 0.0); + + su2double Vi[10] = {0.0}, Vj[10] = {0.0}, Si[2], Sj[2]; + State(nDim, rho_i, vel, p, Vi, Si); + State(nDim, rho_j, vel_j, p, Vj, Sj); + numerics->SetPrimitive(Vi, Vj); + numerics->SetSecondary(Si, Sj); + numerics->SetNormal(n); + const auto res = numerics->ComputeResidual(config.get()); + + su2double un = 0.0; + for (unsigned short iDim = 0; iDim < nDim; iDim++) un += vel[iDim] * n[iDim]; + const su2double rho_up = (un >= 0.0) ? rho_i : rho_j; + const su2double* V_up = (un >= 0.0) ? Vi : Vj; + + CAPTURE(name, nDim, un, implicit, shear); + CHECK(res.residual[0] == Approx(rho_up * un).margin(1e-12)); + for (unsigned short iDim = 0; iDim < nDim; iDim++) + CHECK(res.residual[iDim + 1] == Approx(rho_up * V_up[iDim + 1] * un + p * n[iDim]).margin(1e-12)); + CHECK(res.residual[nDim + 1] == Approx(rho_up * V_up[nDim + 3] * un).margin(1e-12)); +} + +template +void CheckContacts(const char* name, bool implicit = true, bool shear = true) { + const su2double zero[3] = {0.0, 0.0, 0.0}, moving[3] = {0.3, -0.1, 0.2}, reverse[3] = {-0.3, 0.1, -0.2}; + for (const unsigned short nDim : {2, 3}) { + for (const auto* vel : {zero, moving, reverse}) { + CheckContact(name, nDim, vel, implicit, false); + if (shear) CheckContact(name, nDim, vel, implicit, true); + } + } +} + +} // namespace + +TEST_CASE("Roe schemes resolve a contact and a shear layer exactly", "[Roe]") { + CheckContacts("Roe"); + CheckContacts("general Roe", true); + CheckContacts("general Roe", false); + CheckContacts("LMRoe"); + /*--- L2Roe scales the tangential velocity jump at low Mach number (not upwind by design). ---*/ + CheckContacts("L2Roe", true, false); +} diff --git a/UnitTests/meson.build b/UnitTests/meson.build index fc99bc101b3a..3dbbede4a513 100644 --- a/UnitTests/meson.build +++ b/UnitTests/meson.build @@ -15,6 +15,7 @@ su2_cfd_tests = files(['Common/CConfig_tests.cpp', 'Common/containers/CLookupTable_tests.cpp', 'Common/toolboxes/multilayer_perceptron/CLookUp_ANN_tests.cpp', 'SU2_CFD/numerics/CNumerics_tests.cpp', + 'SU2_CFD/numerics/roe_contact_tests.cpp', 'SU2_CFD/fluid/CFluidModel_tests.cpp', 'SU2_CFD/gradients.cpp', 'SU2_CFD/solvers/CNEMOEulerSolver_tests.cpp', From 26c15468d79ddf4b953e0080b5cfd8a034a50474 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:10:31 +0200 Subject: [PATCH 10/14] Fix the internal energy of the under-relaxation factor The factor compares the old and new internal energy, computed as rho*E - 0.5*|rho*u|^2 without the division by rho (the code comment already refers to the static energy). For dimensional cases this value can be negative, and those points got a zero update: in the axisymmetric air nozzle 28 % of the points were frozen, and the case could not re-converge from a restart. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- SU2_CFD/src/solvers/CEulerSolver.cpp | 7 +++++-- 1 file changed, 5 insertions(+), 2 deletions(-) diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index e107d1ed7072..15658e5a7459 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -2551,9 +2551,12 @@ void CEulerSolver::ComputeUnderRelaxationFactor(const CConfig* config) { su2double ratio = fabs(LinSysSol(iPoint, 0)) / max(nodes->GetSolution(iPoint, 0), EPS); su2double e_old = nodes->GetSolution(iPoint, nVar - 1); su2double e_new = e_old + LinSysSol(iPoint, nVar - 1); + /*--- rho*e = rho*E - |rho*u|^2 / (2*rho), with the old and the proposed new density. ---*/ + const su2double rho_old = max(nodes->GetSolution(iPoint, 0), EPS); + const su2double rho_new = max(rho_old + LinSysSol(iPoint, 0), EPS); for (unsigned short jVar = 1; jVar <= nDim; jVar++) { - e_old -= 0.5 * pow(nodes->GetSolution(iPoint, jVar), 2); - e_new -= 0.5 * pow(nodes->GetSolution(iPoint, jVar) + LinSysSol(iPoint, jVar), 2); + e_old -= 0.5 * pow(nodes->GetSolution(iPoint, jVar), 2) / rho_old; + e_new -= 0.5 * pow(nodes->GetSolution(iPoint, jVar) + LinSysSol(iPoint, jVar), 2) / rho_new; } ratio = fmax(ratio, fabs(e_new - e_old) / max(e_old, EPS)); From 611948a9079621d5d83d5315c8cb360e9b83c7e4 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:10:40 +0200 Subject: [PATCH 11/14] Check the MUSCL face states before evaluating the fluid model The consistent extrapolation and the low-Mach correction called the fluid model with the reconstructed pressure and density before checking that they are positive (possible NaN or fluid-model errors with real-gas models at start-up). Check them first, also rejecting zero and NaN, and check again after the low-Mach correction. Valid states give the same result as before. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- SU2_CFD/src/solvers/CEulerSolver.cpp | 42 ++++++++++++++++++---------- 1 file changed, 27 insertions(+), 15 deletions(-) diff --git a/SU2_CFD/src/solvers/CEulerSolver.cpp b/SU2_CFD/src/solvers/CEulerSolver.cpp index 15658e5a7459..6b8784f228aa 100644 --- a/SU2_CFD/src/solvers/CEulerSolver.cpp +++ b/SU2_CFD/src/solvers/CEulerSolver.cpp @@ -1947,27 +1947,39 @@ void CEulerSolver::Upwind_Residual(CGeometry *geometry, CSolver **solver_contain } - /*--- Recompute the reconstructed quantities in a thermodynamically consistent way. ---*/ + /*--- Check for non-physical solutions after reconstruction. If found, use the + cell-average value of the solution. This is a locally 1st order approximation, + which is typically only active during the start-up of a calculation. + Pressure and density are checked before the fluid model is evaluated with them. ---*/ - if (!ideal_gas || low_mach_corr) { - ComputeConsistentExtrapolation(GetFluidModel(), nDim, Primitive_i, Secondary_i); - ComputeConsistentExtrapolation(GetFluidModel(), nDim, Primitive_j, Secondary_j); - } + auto nonPositive = [](su2double value) { return !(value > 0.0) || !std::isfinite(SU2_TYPE::GetValue(value)); }; + bool neg_pres_or_rho_i = nonPositive(Primitive_i[prim_idx.Pressure()]) || nonPositive(Primitive_i[prim_idx.Density()]); + bool neg_pres_or_rho_j = nonPositive(Primitive_j[prim_idx.Pressure()]) || nonPositive(Primitive_j[prim_idx.Density()]); - /*--- Low-Mach number correction. ---*/ + if (!neg_pres_or_rho_i && !neg_pres_or_rho_j) { - if (low_mach_corr) { - LowMachPrimitiveCorrection(GetFluidModel(), nDim, Primitive_i, Primitive_j); - } + /*--- Recompute the reconstructed quantities in a thermodynamically consistent way. ---*/ - /*--- Check for non-physical solutions after reconstruction. If found, use the - cell-average value of the solution. This is a locally 1st order approximation, - which is typically only active during the start-up of a calculation. ---*/ + if (!ideal_gas || low_mach_corr) { + ComputeConsistentExtrapolation(GetFluidModel(), nDim, Primitive_i, Secondary_i); + ComputeConsistentExtrapolation(GetFluidModel(), nDim, Primitive_j, Secondary_j); + } + + /*--- Low-Mach number correction. ---*/ + + if (low_mach_corr) { + LowMachPrimitiveCorrection(GetFluidModel(), nDim, Primitive_i, Primitive_j); + } + + /*--- The correction changes the face states, check them again. ---*/ - bool neg_pres_or_rho_i = (Primitive_i[prim_idx.Pressure()] < 0.0) || (Primitive_i[prim_idx.Density()] < 0.0); - bool neg_pres_or_rho_j = (Primitive_j[prim_idx.Pressure()] < 0.0) || (Primitive_j[prim_idx.Density()] < 0.0); + neg_pres_or_rho_i = nonPositive(Primitive_i[prim_idx.Pressure()]) || nonPositive(Primitive_i[prim_idx.Density()]); + neg_pres_or_rho_j = nonPositive(Primitive_j[prim_idx.Pressure()]) || nonPositive(Primitive_j[prim_idx.Density()]); + } - su2double R = sqrt(fabs(Primitive_j[prim_idx.Density()]/Primitive_i[prim_idx.Density()])); + su2double R = 1.0; + if (!neg_pres_or_rho_i && !neg_pres_or_rho_j) + R = sqrt(Primitive_j[prim_idx.Density()]/Primitive_i[prim_idx.Density()]); su2double sq_vel = 0.0; for (auto iDim = 0u; iDim < nDim; iDim++) { su2double RoeVelocity = (R * Primitive_j[iDim + prim_idx.Velocity()] + From 62ea3ca08b3a97725b1be0890bea33f48e36bcb9 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 14:21:56 +0200 Subject: [PATCH 12/14] Compute HLLC on moving grids in the frame of the face, fix its Jacobians - Moving grids: the wave speeds, the pressure work of the moving face in the supersonic branches and the term V_g (p* - p)/(s - S_M) of the star-state energy were missing or inconsistent, so HLLC was not Galilean invariant (NACA0012 M 0.8, fixed grid vs still air and a translating grid: CL -12 %). The flux is now computed in the frame of the face (u' = u - w n, H' = H - w u_n + w^2/2) and transformed back (F = B F', J = B J' T). - Jacobians: sign of the velocity term of the general-gas dp/drho. The star-state Jacobians keep by default the density of the other state in dp*/dU (inexact, #2198), which was more robust with inexact linear solves in our tests; USE_ACCURATE_FLUX_JACOBIANS= YES now gives the exact derivative for fixed wave speeds. - Unit test hllc_jacobian_tests.cpp: finite differences in 2D and 3D, fixed and moving faces, both sides of the contact, all branches, default and exact. Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- Common/src/CConfig.cpp | 2 +- .../include/numerics/flow/convection/hllc.hpp | 30 +++ SU2_CFD/src/numerics/flow/convection/hllc.cpp | 180 +++++++++----- .../SU2_CFD/numerics/hllc_jacobian_tests.cpp | 220 ++++++++++++++++++ UnitTests/meson.build | 1 + config_template.cfg | 1 + 6 files changed, 381 insertions(+), 53 deletions(-) create mode 100644 UnitTests/SU2_CFD/numerics/hllc_jacobian_tests.cpp diff --git a/Common/src/CConfig.cpp b/Common/src/CConfig.cpp index 261002df1f4f..3a5bfee8f713 100644 --- a/Common/src/CConfig.cpp +++ b/Common/src/CConfig.cpp @@ -2171,7 +2171,7 @@ void CConfig::SetConfig_Options() { addDoubleArrayOption("JST_SENSOR_COEFF", 2, false, jst_coeff); /*!\brief LAX_SENSOR_COEFF \n DESCRIPTION: 1st order artificial dissipation coefficients for the Lax-Friedrichs method. \ingroup Config*/ addDoubleOption("LAX_SENSOR_COEFF", Kappa_1st_Flow, 0.15); - /*!\brief USE_ACCURATE_FLUX_JACOBIANS \n DESCRIPTION: Use numerically computed Jacobians for AUSM+up(2) and SLAU(2) \ingroup Config*/ + /*!\brief USE_ACCURATE_FLUX_JACOBIANS \n DESCRIPTION: Use numerically computed Jacobians for AUSM+up(2) and SLAU(2), exact Jacobians (for fixed wave speeds) for HLLC \ingroup Config*/ addBoolOption("USE_ACCURATE_FLUX_JACOBIANS", Use_Accurate_Jacobians, false); /*!\brief CENTRAL_JACOBIAN_FIX_FACTOR \n DESCRIPTION: Improve the numerical properties (diagonal dominance) of the global Jacobian matrix, 3 to 4 is "optimum" (central schemes) \ingroup Config*/ addDoubleOption("CENTRAL_JACOBIAN_FIX_FACTOR", Cent_Jac_Fix_Factor, 4.0); diff --git a/SU2_CFD/include/numerics/flow/convection/hllc.hpp b/SU2_CFD/include/numerics/flow/convection/hllc.hpp index b525df4bf9cf..29f7ce727a34 100644 --- a/SU2_CFD/include/numerics/flow/convection/hllc.hpp +++ b/SU2_CFD/include/numerics/flow/convection/hllc.hpp @@ -39,6 +39,7 @@ class CUpwHLLC_Flow final : public CNumerics { private: bool implicit, dynamic_grid; + bool useAccurateJacobian; /*!< \brief Jacobians computed exactly for fixed wave speeds (USE_ACCURATE_FLUX_JACOBIANS). */ unsigned short iDim, jDim, iVar, jVar; su2double *IntermediateState; @@ -55,6 +56,8 @@ class CUpwHLLC_Flow final : public CNumerics { su2double *dSm_dU, *dPI_dU, *drhoStar_dU, *dpStar_dU, *dEStar_dU; su2double* Flux; /*!< \brief The flux accross the face. */ + bool fixedWaveSpeeds = false; /*!< \brief Use fixed_sL, fixed_sR as wave speeds (verification of the Jacobians). */ + su2double fixed_sL = 0.0, fixed_sR = 0.0; su2double** Jacobian_i; /*!< \brief The Jacobian w.r.t. point i after computation. */ su2double** Jacobian_j; /*!< \brief The Jacobian w.r.t. point j after computation. */ @@ -79,6 +82,18 @@ class CUpwHLLC_Flow final : public CNumerics { */ ResidualType<> ComputeResidual(const CConfig* config) override; + /*! + * \brief Fix the left and right wave speeds (relative to the face), for the verification of the Jacobians, which + * treat them as constants. + * \param[in] active - Use the fixed values. + * \param[in] val_sL, val_sR - Left and right wave speeds. + */ + inline void SetFixedWaveSpeeds(bool active, su2double val_sL, su2double val_sR) { + fixedWaveSpeeds = active; + fixed_sL = val_sL; + fixed_sR = val_sR; + } + }; /*! @@ -91,6 +106,7 @@ class CUpwHLLC_Flow final : public CNumerics { class CUpwGeneralHLLC_Flow final : public CNumerics { private: bool implicit, dynamic_grid; + bool useAccurateJacobian; /*!< \brief Jacobians computed exactly for fixed wave speeds (USE_ACCURATE_FLUX_JACOBIANS). */ unsigned short iDim, jDim, iVar, jVar; su2double *IntermediateState; @@ -108,6 +124,8 @@ class CUpwGeneralHLLC_Flow final : public CNumerics { su2double *dSm_dU, *dPI_dU, *drhoStar_dU, *dpStar_dU, *dEStar_dU; su2double* Flux; /*!< \brief The flux accross the face. */ + bool fixedWaveSpeeds = false; /*!< \brief Use fixed_sL, fixed_sR as wave speeds (verification of the Jacobians). */ + su2double fixed_sL = 0.0, fixed_sR = 0.0; su2double** Jacobian_i; /*!< \brief The Jacobian w.r.t. point i after computation. */ su2double** Jacobian_j; /*!< \brief The Jacobian w.r.t. point j after computation. */ @@ -132,6 +150,18 @@ class CUpwGeneralHLLC_Flow final : public CNumerics { */ ResidualType<> ComputeResidual(const CConfig* config) override; + /*! + * \brief Fix the left and right wave speeds (relative to the face), for the verification of the Jacobians, which + * treat them as constants. + * \param[in] active - Use the fixed values. + * \param[in] val_sL, val_sR - Left and right wave speeds. + */ + inline void SetFixedWaveSpeeds(bool active, su2double val_sL, su2double val_sR) { + fixedWaveSpeeds = active; + fixed_sL = val_sL; + fixed_sR = val_sR; + } + /*! * \brief Compute the Average quantities for a general fluid flux between two nodes i and j. * Using the approach of Vinokur and Montagne' diff --git a/SU2_CFD/src/numerics/flow/convection/hllc.cpp b/SU2_CFD/src/numerics/flow/convection/hllc.cpp index ef6a2a5a41e9..18b210e72381 100644 --- a/SU2_CFD/src/numerics/flow/convection/hllc.cpp +++ b/SU2_CFD/src/numerics/flow/convection/hllc.cpp @@ -28,12 +28,59 @@ #include "../../../../include/numerics/flow/convection/hllc.hpp" #include "../../../../../Common/include/toolboxes/geometry_toolbox.hpp" +namespace { +/*--- Moving grids. The flux across a face that moves with normal velocity w is the flux of the states seen from the + face (velocity u' = u - c, c = w n, same thermodynamic state) transformed back: with U' = T U, + F_ALE(U) = B F'(U'), where F' is the flux for a fixed face. Mass: unchanged; momentum: + c F'_0; + energy: + c.F'_m + |c|^2/2 F'_0. The Jacobians are B J' T. ---*/ + +void ShiftToFaceFrame(unsigned short nDim, su2double w, const su2double* unitNormal, su2double* velocity, + su2double& enthalpy) { + const su2double projVel = GeometryToolbox::DotProduct(nDim, velocity, unitNormal); + for (unsigned short iDim = 0; iDim < nDim; iDim++) velocity[iDim] -= w * unitNormal[iDim]; + enthalpy += w * (0.5 * w - projVel); // H' = H - u.c + |c|^2/2 +} + +void FaceFrameToMovingFlux(unsigned short nDim, su2double w, const su2double* unitNormal, su2double* flux) { + const unsigned short nVar = nDim + 2; + su2double cDotMom = 0.0; + for (unsigned short iDim = 0; iDim < nDim; iDim++) cDotMom += w * unitNormal[iDim] * flux[iDim+1]; + flux[nVar-1] += cDotMom + 0.5 * w * w * flux[0]; + for (unsigned short iDim = 0; iDim < nDim; iDim++) flux[iDim+1] += w * unitNormal[iDim] * flux[0]; +} + +void FaceFrameToMovingJacobian(unsigned short nDim, su2double w, const su2double* unitNormal, su2double** jac) { + const unsigned short nVar = nDim + 2; + constexpr unsigned short maxDim = 3, maxVar = maxDim + 2; + su2double c[maxDim] = {0.0}, T[maxVar][maxVar] = {{0.0}}, B[maxVar][maxVar] = {{0.0}}, tmp[maxVar][maxVar]; + su2double c2 = 0.0; + for (unsigned short iDim = 0; iDim < nDim; iDim++) { c[iDim] = w * unitNormal[iDim]; c2 += c[iDim] * c[iDim]; } + for (unsigned short iVar = 0; iVar < nVar; iVar++) { T[iVar][iVar] = 1.0; B[iVar][iVar] = 1.0; } + for (unsigned short iDim = 0; iDim < nDim; iDim++) { + T[iDim+1][0] = -c[iDim]; T[nVar-1][iDim+1] = -c[iDim]; + B[iDim+1][0] = c[iDim]; B[nVar-1][iDim+1] = c[iDim]; + } + T[nVar-1][0] = 0.5 * c2; B[nVar-1][0] = 0.5 * c2; + for (unsigned short iVar = 0; iVar < nVar; iVar++) + for (unsigned short jVar = 0; jVar < nVar; jVar++) { + tmp[iVar][jVar] = 0.0; + for (unsigned short kVar = 0; kVar < nVar; kVar++) tmp[iVar][jVar] += jac[iVar][kVar] * T[kVar][jVar]; + } + for (unsigned short iVar = 0; iVar < nVar; iVar++) + for (unsigned short jVar = 0; jVar < nVar; jVar++) { + jac[iVar][jVar] = 0.0; + for (unsigned short kVar = 0; kVar < nVar; kVar++) jac[iVar][jVar] += B[iVar][kVar] * tmp[kVar][jVar]; + } +} +} // namespace + CUpwHLLC_Flow::CUpwHLLC_Flow(unsigned short val_nDim, unsigned short val_nVar, const CConfig* config) : CNumerics(val_nDim, val_nVar, config) { implicit = (config->GetKind_TimeIntScheme_Flow() == EULER_IMPLICIT); /* A grid is defined as dynamic if there's rigid grid movement or grid deformation AND the problem is time domain */ dynamic_grid = config->GetDynamic_Grid(); + useAccurateJacobian = config->GetUse_Accurate_Jacobians(); Gamma = config->GetGamma(); @@ -114,6 +161,16 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) Density_j = V_j[nDim+2]; Enthalpy_j = V_j[nDim+3]; + /*--- Moving grid: states seen from the face (see ShiftToFaceFrame). ---*/ + + ProjInterfaceVel = 0.0; + if (dynamic_grid) { + for (iDim = 0; iDim < nDim; iDim++) + ProjInterfaceVel += 0.5 * (GridVel_i[iDim] + GridVel_j[iDim]) * UnitNormal[iDim]; + ShiftToFaceFrame(nDim, ProjInterfaceVel, UnitNormal, Velocity_i, Enthalpy_i); + ShiftToFaceFrame(nDim, ProjInterfaceVel, UnitNormal, Velocity_j, Enthalpy_j); + } + sq_vel_i = GeometryToolbox::SquaredNorm(nDim, Velocity_i); sq_vel_j = GeometryToolbox::SquaredNorm(nDim, Velocity_j); @@ -128,22 +185,6 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) ProjVelocity_i = GeometryToolbox::DotProduct(nDim, Velocity_i, UnitNormal); ProjVelocity_j = GeometryToolbox::DotProduct(nDim, Velocity_j, UnitNormal); - /*--- Projected Grid Velocity ---*/ - - ProjInterfaceVel = 0; - - if (dynamic_grid) { - - for (iDim = 0; iDim < nDim; iDim++) - ProjInterfaceVel += 0.5 * (GridVel_i[iDim] + GridVel_j[iDim]) * UnitNormal[iDim]; - - SoundSpeed_i -= ProjInterfaceVel; - SoundSpeed_j += ProjInterfaceVel; - - ProjVelocity_i -= ProjInterfaceVel; - ProjVelocity_j -= ProjInterfaceVel; - } - /*--- Roe's averaging ---*/ Rrho = ( sqrt(Density_i) + sqrt(Density_j) ); @@ -152,7 +193,7 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) RoeVelocity[iDim] = ( Velocity_i[iDim] * sqrt(Density_i) + Velocity_j[iDim] * sqrt(Density_j) ) / Rrho; } sq_velRoe = GeometryToolbox::SquaredNorm(nDim, RoeVelocity); - RoeProjVelocity = GeometryToolbox::DotProduct(nDim, RoeVelocity, UnitNormal) - ProjInterfaceVel; + RoeProjVelocity = GeometryToolbox::DotProduct(nDim, RoeVelocity, UnitNormal); /*--- Mean Roe variables iPoint and jPoint ---*/ @@ -161,13 +202,15 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) /*--- Roe-averaged speed of sound ---*/ - RoeSoundSpeed = sqrt( Gamma_Minus_One * ( RoeEnthalpy - 0.5 * sq_velRoe ) ) - ProjInterfaceVel; + RoeSoundSpeed = sqrt( Gamma_Minus_One * ( RoeEnthalpy - 0.5 * sq_velRoe ) ); /*--- Speed of sound at L and R ---*/ sL = min( RoeProjVelocity - RoeSoundSpeed, ProjVelocity_i - SoundSpeed_i); sR = max( RoeProjVelocity + RoeSoundSpeed, ProjVelocity_j + SoundSpeed_j); + if (fixedWaveSpeeds) { sL = fixed_sL; sR = fixed_sR; } + /*--- speed of contact surface ---*/ RHO = Density_j * (sR - ProjVelocity_j) - Density_i * (sL - ProjVelocity_i); @@ -205,7 +248,7 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) Flux[0] = sM * IntermediateState[0]; for (iDim = 0; iDim < nDim; iDim++) Flux[iDim+1] = sM * IntermediateState[iDim+1] + pStar * UnitNormal[iDim]; - Flux[nVar-1] = sM * ( IntermediateState[nVar-1] + pStar ) + pStar * ProjInterfaceVel; + Flux[nVar-1] = sM * ( IntermediateState[nVar-1] + pStar ); } } else { @@ -234,16 +277,29 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) Flux[0] = sM * IntermediateState[0]; for (iDim = 0; iDim < nDim; iDim++) Flux[iDim+1] = sM * IntermediateState[iDim+1] + pStar * UnitNormal[iDim]; - Flux[nVar-1] = sM * (IntermediateState[nVar-1] + pStar ) + pStar * ProjInterfaceVel; + Flux[nVar-1] = sM * (IntermediateState[nVar-1] + pStar ); } } + if (dynamic_grid) FaceFrameToMovingFlux(nDim, ProjInterfaceVel, UnitNormal, Flux); + for (iVar = 0; iVar < nVar; iVar++) Flux[iVar] *= Area; /*--- Return early if the Jacobians do not need to be computed. ---*/ if (!implicit) return ResidualType<>(Flux, Jacobian_i, Jacobian_j); + /*--- Jacobians for the wave speeds held fixed (in the frame of the face). p* = rho_j (u_j - sR)(u_j - sM) + p_j + = rho_i (u_i - sL)(u_i - sM) + p_i, so its derivative contains rho_j (sR - u_j) dsM/dU (or rho_i (sL - u_i) dsM/dU). + By default the density of the other state is kept, as in the previous implementation (an inexact transcription + of the reference, see #2198): in our tests it was more robust with inexact linear solves than the exact + derivative, which USE_ACCURATE_FLUX_JACOBIANS= YES gives. ---*/ + + const su2double rhoPStar_R = useAccurateJacobian ? Density_j : Density_i; + const su2double rhoPStar_L = useAccurateJacobian ? Density_i : Density_j; + + + if (sM > 0.0) { if (sL > 0.0) { @@ -295,7 +351,7 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) /*--- Computing d/dU_L (pStar) ---*/ for (iVar = 0; iVar < nVar; iVar++) - dpStar_dU[iVar] = Density_i * (sR - ProjVelocity_j) * dSm_dU[iVar]; + dpStar_dU[iVar] = rhoPStar_R * (sR - ProjVelocity_j) * dSm_dU[iVar]; /*--- Computing d/dU_L (EStar) ---*/ @@ -352,7 +408,7 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) /*--- Computing d/dU_R (pStar) ---*/ for (iVar = 0; iVar < nVar; iVar++) - dpStar_dU[iVar] = Density_j * (sL - ProjVelocity_i) * dSm_dU[iVar]; + dpStar_dU[iVar] = rhoPStar_L * (sL - ProjVelocity_i) * dSm_dU[iVar]; /*--- Computing d/dU_R (EStar) ---*/ @@ -414,7 +470,7 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) /*--- Computing d/dU_L (pStar) ---*/ for (iVar = 0; iVar < nVar; iVar++) - dpStar_dU[iVar] = Density_i * (sR - ProjVelocity_j) * dSm_dU[iVar]; + dpStar_dU[iVar] = rhoPStar_R * (sR - ProjVelocity_j) * dSm_dU[iVar]; /*--- Computing d/dU_L (EStar) ---*/ @@ -466,7 +522,7 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) /*--- Computing d/dU_R (pStar) ---*/ for (iVar = 0; iVar < nVar; iVar++) - dpStar_dU[iVar] = Density_j * (sL - ProjVelocity_i) * dSm_dU[iVar]; + dpStar_dU[iVar] = rhoPStar_L * (sL - ProjVelocity_i) * dSm_dU[iVar]; /*--- Computing d/dU_R (rhoStar) ---*/ @@ -519,6 +575,11 @@ CNumerics::ResidualType<> CUpwHLLC_Flow::ComputeResidual(const CConfig* config) } } + if (dynamic_grid) { + FaceFrameToMovingJacobian(nDim, ProjInterfaceVel, UnitNormal, Jacobian_i); + FaceFrameToMovingJacobian(nDim, ProjInterfaceVel, UnitNormal, Jacobian_j); + } + /*--- Scale Jacobians by area (from Flux *= Area). ---*/ for (iVar = 0; iVar < nVar; iVar++) { for (jVar = 0; jVar < nVar; jVar++) { @@ -535,6 +596,7 @@ CUpwGeneralHLLC_Flow::CUpwGeneralHLLC_Flow(unsigned short val_nDim, unsigned sho /* A grid is defined as dynamic if there's rigid grid movement or grid deformation AND the problem is time domain */ dynamic_grid = config->GetDynamic_Grid(); + useAccurateJacobian = config->GetUse_Accurate_Jacobians(); Gamma = config->GetGamma(); @@ -613,6 +675,16 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c Density_j = V_j[nDim+2]; Enthalpy_j = V_j[nDim+3]; + /*--- Moving grid: states seen from the face (see ShiftToFaceFrame). ---*/ + + ProjInterfaceVel = 0.0; + if (dynamic_grid) { + for (iDim = 0; iDim < nDim; iDim++) + ProjInterfaceVel += 0.5 * (GridVel_i[iDim] + GridVel_j[iDim]) * UnitNormal[iDim]; + ShiftToFaceFrame(nDim, ProjInterfaceVel, UnitNormal, Velocity_i, Enthalpy_i); + ShiftToFaceFrame(nDim, ProjInterfaceVel, UnitNormal, Velocity_j, Enthalpy_j); + } + sq_vel_i = 0.0; sq_vel_j = 0.0; @@ -650,28 +722,12 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c } - /*--- Projected Grid Velocity ---*/ - - ProjInterfaceVel = 0; - - if (dynamic_grid) { - - for (iDim = 0; iDim < nDim; iDim++) - ProjInterfaceVel += 0.5 * ( GridVel_i[iDim] + GridVel_j[iDim] )*UnitNormal[iDim]; - - SoundSpeed_i -= ProjInterfaceVel; - SoundSpeed_j += ProjInterfaceVel; - - ProjVelocity_i -= ProjInterfaceVel; - ProjVelocity_j -= ProjInterfaceVel; - } - /*--- Roe's averaging ---*/ Rrho = ( sqrt(Density_i) + sqrt(Density_j) ); sq_velRoe = 0.0; - RoeProjVelocity = - ProjInterfaceVel; + RoeProjVelocity = 0.0; for (iDim = 0; iDim < nDim; iDim++) { RoeVelocity[iDim] = ( Velocity_i[iDim] * sqrt(Density_i) + Velocity_j[iDim] * sqrt(Density_j) ) / Rrho; @@ -691,13 +747,15 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c /*--- Roe-averaged speed of sound ---*/ //RoeSoundSpeed2 = RoeChi + RoeKappa * ( RoeEnthalpy - 0.5 * sq_velRoe ); - RoeSoundSpeed = sqrt( RoeChi + RoeKappa * ( RoeEnthalpy - 0.5 * sq_velRoe ) ) - ProjInterfaceVel; + RoeSoundSpeed = sqrt( RoeChi + RoeKappa * ( RoeEnthalpy - 0.5 * sq_velRoe ) ); /*--- Speed of sound at L and R ---*/ sL = min( RoeProjVelocity - RoeSoundSpeed, ProjVelocity_i - SoundSpeed_i ); sR = max( RoeProjVelocity + RoeSoundSpeed, ProjVelocity_j + SoundSpeed_j ); + if (fixedWaveSpeeds) { sL = fixed_sL; sR = fixed_sR; } + /*--- speed of contact surface ---*/ RHO = Density_j * (sR - ProjVelocity_j) - Density_i * (sL - ProjVelocity_i); @@ -734,7 +792,7 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c Flux[0] = sM * IntermediateState[0]; for (iDim = 0; iDim < nDim; iDim++) Flux[iDim+1] = sM * IntermediateState[iDim+1] + pStar * UnitNormal[iDim]; - Flux[nVar-1] = sM * ( IntermediateState[nVar-1] + pStar ) + pStar * ProjInterfaceVel; + Flux[nVar-1] = sM * ( IntermediateState[nVar-1] + pStar ); } } else { @@ -763,16 +821,29 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c Flux[0] = sM * IntermediateState[0]; for (iDim = 0; iDim < nDim; iDim++) Flux[iDim+1] = sM * IntermediateState[iDim+1] + pStar * UnitNormal[iDim]; - Flux[nVar-1] = sM * (IntermediateState[nVar-1] + pStar ) + pStar * ProjInterfaceVel; + Flux[nVar-1] = sM * (IntermediateState[nVar-1] + pStar ); } } + if (dynamic_grid) FaceFrameToMovingFlux(nDim, ProjInterfaceVel, UnitNormal, Flux); + for (iVar = 0; iVar < nVar; iVar++) Flux[iVar] *= Area; /*--- Return early if the Jacobians do not need to be computed. ---*/ if (!implicit) return ResidualType<>(Flux, Jacobian_i, Jacobian_j); + /*--- Jacobians for the wave speeds held fixed (in the frame of the face). p* = rho_j (u_j - sR)(u_j - sM) + p_j + = rho_i (u_i - sL)(u_i - sM) + p_i, so its derivative contains rho_j (sR - u_j) dsM/dU (or rho_i (sL - u_i) dsM/dU). + By default the density of the other state is kept, as in the previous implementation (an inexact transcription + of the reference, see #2198): in our tests it was more robust with inexact linear solves than the exact + derivative, which USE_ACCURATE_FLUX_JACOBIANS= YES gives. ---*/ + + const su2double rhoPStar_R = useAccurateJacobian ? Density_j : Density_i; + const su2double rhoPStar_L = useAccurateJacobian ? Density_i : Density_j; + + + if (sM > 0.0) { if (sL > 0.0) { @@ -800,7 +871,7 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c /*--- Computing pressure derivatives d/dU_L (PI) ---*/ - dPI_dU[0] = Chi_i - 0.5 * Kappa_i * sq_vel_i; + dPI_dU[0] = Chi_i + 0.5 * Kappa_i * sq_vel_i; for (iDim = 0; iDim < nDim; iDim++) dPI_dU[iDim+1] = - Kappa_i * Velocity_i[iDim]; dPI_dU[nVar-1] = Kappa_i; @@ -825,7 +896,7 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c /*--- Computing d/dU_L (pStar) ---*/ for (iVar = 0; iVar < nVar; iVar++) - dpStar_dU[iVar] = Density_i * (sR - ProjVelocity_j) * dSm_dU[iVar]; + dpStar_dU[iVar] = rhoPStar_R * (sR - ProjVelocity_j) * dSm_dU[iVar]; /*--- Computing d/dU_L (EStar) ---*/ @@ -874,7 +945,7 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c /*--- Computing pressure derivatives d/dU_R (PI) ---*/ - dPI_dU[0] = Chi_j - 0.5 * Kappa_j * sq_vel_j; + dPI_dU[0] = Chi_j + 0.5 * Kappa_j * sq_vel_j; for (iDim = 0; iDim < nDim; iDim++) dPI_dU[iDim+1] = - Kappa_j * Velocity_j[iDim]; dPI_dU[nVar-1] = Kappa_j; @@ -891,7 +962,7 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c /*--- Computing d/dU_R (pStar) ---*/ for (iVar = 0; iVar < nVar; iVar++) - dpStar_dU[iVar] = Density_j * (sL - ProjVelocity_i) * dSm_dU[iVar]; + dpStar_dU[iVar] = rhoPStar_L * (sL - ProjVelocity_i) * dSm_dU[iVar]; /*--- Computing d/dU_R (EStar) ---*/ @@ -944,7 +1015,7 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c /*--- Computing pressure derivatives d/dU_L (PI) ---*/ - dPI_dU[0] = Chi_i - 0.5 * Kappa_i * sq_vel_i; + dPI_dU[0] = Chi_i + 0.5 * Kappa_i * sq_vel_i; for (iDim = 0; iDim < nDim; iDim++) dPI_dU[iDim+1] = - Kappa_i * Velocity_i[iDim]; dPI_dU[nVar-1] = Kappa_i; @@ -961,7 +1032,7 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c /*--- Computing d/dU_L (pStar) ---*/ for (iVar = 0; iVar < nVar; iVar++) - dpStar_dU[iVar] = Density_i * (sR - ProjVelocity_j) * dSm_dU[iVar]; + dpStar_dU[iVar] = rhoPStar_R * (sR - ProjVelocity_j) * dSm_dU[iVar]; /*--- Computing d/dU_L (EStar) ---*/ @@ -995,7 +1066,7 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c /*--- Computing pressure derivatives d/dU_R (PI) ---*/ - dPI_dU[0] = Chi_j - 0.5 * Kappa_j * sq_vel_j; + dPI_dU[0] = Chi_j + 0.5 * Kappa_j * sq_vel_j; for (iDim = 0; iDim < nDim; iDim++) dPI_dU[iDim+1] = - Kappa_j * Velocity_j[iDim]; dPI_dU[nVar-1] = Kappa_j; @@ -1012,7 +1083,7 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c /*--- Computing d/dU_R (pStar) ---*/ for (iVar = 0; iVar < nVar; iVar++) - dpStar_dU[iVar] = Density_j * (sL - ProjVelocity_i) * dSm_dU[iVar]; + dpStar_dU[iVar] = rhoPStar_L * (sL - ProjVelocity_i) * dSm_dU[iVar]; /*--- Computing d/dU_R (rhoStar) ---*/ @@ -1064,6 +1135,11 @@ CNumerics::ResidualType<> CUpwGeneralHLLC_Flow::ComputeResidual(const CConfig* c } } + if (dynamic_grid) { + FaceFrameToMovingJacobian(nDim, ProjInterfaceVel, UnitNormal, Jacobian_i); + FaceFrameToMovingJacobian(nDim, ProjInterfaceVel, UnitNormal, Jacobian_j); + } + /*--- Scale Jacobians by area (from Flux *= Area). ---*/ for (iVar = 0; iVar < nVar; iVar++) { for (jVar = 0; jVar < nVar; jVar++) { diff --git a/UnitTests/SU2_CFD/numerics/hllc_jacobian_tests.cpp b/UnitTests/SU2_CFD/numerics/hllc_jacobian_tests.cpp new file mode 100644 index 000000000000..0a14682612ed --- /dev/null +++ b/UnitTests/SU2_CFD/numerics/hllc_jacobian_tests.cpp @@ -0,0 +1,220 @@ +/*! + * \file hllc_jacobian_tests.cpp + * \brief Finite-difference verification of the HLLC Jacobians (ideal and general gas, 2D and 3D), on fixed and moving + * faces, exact (USE_ACCURATE_FLUX_JACOBIANS) or default. + * \version 8.5.0 "Harrier" + * + * SU2 Project Website: https://su2code.github.io + * + * The SU2 Project is maintained by the SU2 Foundation + * (http://su2foundation.org) + * + * Copyright 2012-2026, SU2 Contributors (cf. AUTHORS.md) + * + * SU2 is free software; you can redistribute it and/or + * modify it under the terms of the GNU Lesser General Public + * License as published by the Free Software Foundation; either + * version 2.1 of the License, or (at your option) any later version. + * + * SU2 is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU + * Lesser General Public License for more details. + * + * You should have received a copy of the GNU Lesser General Public + * License along with SU2. If not, see . + */ + +#include "catch.hpp" +#include +#include +#include +#include +#include "../../../SU2_CFD/include/numerics/flow/convection/hllc.hpp" + +namespace { + +constexpr su2double gamma = 1.4; +constexpr unsigned short maxVar = 5; + +using State = std::array; // conservative variables (rho, rho u, rho v, [rho w], rho E) + +/*--- Primitive (T, u, ..., p, rho, h, c) and secondary (dp/drho_e, dp/de_rho) variables of an ideal gas whose total + energy contains the turbulent kinetic energy k. ---*/ +void Primitives(unsigned short nDim, const State& U, su2double k, su2double* V, su2double* S) { + const su2double rho = U[0]; + su2double q2 = 0.0; + for (unsigned short iDim = 0; iDim < nDim; iDim++) { + V[iDim + 1] = U[iDim + 1] / rho; + q2 += V[iDim + 1] * V[iDim + 1]; + } + const su2double e = U[nDim + 1] / rho - 0.5 * q2 - k; + const su2double p = (gamma - 1) * rho * e; + V[0] = p / (rho * 287.058); + V[nDim + 1] = p; + V[nDim + 2] = rho; + V[nDim + 3] = (U[nDim + 1] + p) / rho; + V[nDim + 4] = sqrt(gamma * p / rho); + S[0] = (gamma - 1) * e; + S[1] = (gamma - 1) * rho; +} + +/*--- Conservative state from density, velocity, pressure and k. ---*/ +State Conservative(unsigned short nDim, su2double rho, const su2double* vel, su2double p, su2double k) { + State U{}; + su2double q2 = 0.0; + U[0] = rho; + for (unsigned short iDim = 0; iDim < nDim; iDim++) { + U[iDim + 1] = rho * vel[iDim]; + q2 += vel[iDim] * vel[iDim]; + } + U[nDim + 1] = p / (gamma - 1) + rho * (0.5 * q2 + k); + return U; +} + +struct Result { + State flux; + su2double jac_i[maxVar][maxVar], jac_j[maxVar][maxVar]; +}; + +Result Evaluate(unsigned short nDim, CNumerics& numerics, const CConfig* config, const State& Ui, const State& Uj, + su2double k_i, su2double k_j) { + su2double Vi[20] = {0.0}, Vj[20] = {0.0}, Si[2], Sj[2]; + Primitives(nDim, Ui, k_i, Vi, Si); + Primitives(nDim, Uj, k_j, Vj, Sj); + numerics.SetPrimitive(Vi, Vj); + numerics.SetSecondary(Si, Sj); + numerics.SetTurbKineticEnergy(k_i, k_j); + const auto res = numerics.ComputeResidual(config); + Result r; + for (unsigned short iVar = 0; iVar < nDim + 2; iVar++) { + r.flux[iVar] = res.residual[iVar]; + for (unsigned short jVar = 0; jVar < nDim + 2; jVar++) { + r.jac_i[iVar][jVar] = res.jacobian_i[iVar][jVar]; + r.jac_j[iVar][jVar] = res.jacobian_j[iVar][jVar]; + } + } + return r; +} + +/*--- Largest difference between the Jacobians and central differences of the flux, relative to the largest entry. ---*/ +su2double JacobianError(unsigned short nDim, CNumerics& numerics, const CConfig* config, const State& Ui, + const State& Uj, su2double k_i, su2double k_j) { + const unsigned short nVar = nDim + 2; + const auto ref = Evaluate(nDim, numerics, config, Ui, Uj, k_i, k_j); + su2double maxErr = 0.0, maxJac = 0.0; + for (int side = 0; side < 2; side++) { + for (unsigned short jVar = 0; jVar < nVar; jVar++) { + State Up = side ? Uj : Ui, Um = Up; + const su2double h = 1e-6 * std::max(std::abs(Up[jVar]), su2double(1.0)); + Up[jVar] += h; + Um[jVar] -= h; + const auto fp = side ? Evaluate(nDim, numerics, config, Ui, Up, k_i, k_j) + : Evaluate(nDim, numerics, config, Up, Uj, k_i, k_j); + const auto fm = side ? Evaluate(nDim, numerics, config, Ui, Um, k_i, k_j) + : Evaluate(nDim, numerics, config, Um, Uj, k_i, k_j); + for (unsigned short iVar = 0; iVar < nVar; iVar++) { + const su2double fd = (fp.flux[iVar] - fm.flux[iVar]) / (2 * h); + const su2double jac = side ? ref.jac_j[iVar][jVar] : ref.jac_i[iVar][jVar]; + maxErr = std::max(maxErr, std::abs(jac - fd)); + maxJac = std::max(maxJac, std::abs(jac)); + } + } + } + return maxErr / maxJac; +} + +/*--- Contact speed sM of the scheme, in the frame of the face (w = grid velocity along the unit normal n). ---*/ +su2double ContactSpeed(unsigned short nDim, const State& Ui, const State& Uj, su2double k_i, su2double k_j, + const su2double* n, su2double w, su2double sL, su2double sR) { + su2double Vi[20] = {0.0}, Vj[20] = {0.0}, S[2]; + Primitives(nDim, Ui, k_i, Vi, S); + Primitives(nDim, Uj, k_j, Vj, S); + su2double u_i = -w, u_j = -w; + for (unsigned short iDim = 0; iDim < nDim; iDim++) { + u_i += Vi[iDim + 1] * n[iDim]; + u_j += Vj[iDim + 1] * n[iDim]; + } + const su2double p_i = Vi[nDim + 1], p_j = Vj[nDim + 1], rho_i = Vi[nDim + 2], rho_j = Vj[nDim + 2]; + return (p_i - p_j - rho_i * u_i * (sL - u_i) + rho_j * u_j * (sR - u_j)) / (rho_j * (sR - u_j) - rho_i * (sL - u_i)); +} + +std::unique_ptr MakeConfig(bool moving, bool accurate) { + std::stringstream options; + options << "SOLVER= EULER\nCONV_NUM_METHOD_FLOW= HLLC\nTIME_DISCRE_FLOW= EULER_IMPLICIT\n"; + options << "USE_ACCURATE_FLUX_JACOBIANS= " << (accurate ? "YES" : "NO") << "\n"; + if (moving) options << "GRID_MOVEMENT= STEADY_TRANSLATION\nTRANSLATION_RATE= 0.0 0.0 0.0\n"; + return std::make_unique(options, SU2_COMPONENT::SU2_CFD, false); +} + +/*--- With USE_ACCURATE_FLUX_JACOBIANS= YES the Jacobians are exact for fixed wave speeds in all branches. By default + the derivative of p* is approximated (star states), the supersonic branches remain exact. ---*/ +template +void CheckAllBranches(unsigned short nDim, bool moving, su2double k, bool accurate) { + auto config = MakeConfig(moving, accurate); + HLLC numerics(nDim, nDim + 2, config.get()); + + /*--- Oblique unit normal (face area 1.3) and, on moving faces, a grid velocity with a normal component. ---*/ + const su2double n2[3] = {0.8, 0.6, 0.0}, n3[3] = {2.0 / 7, 3.0 / 7, 6.0 / 7}; + const su2double* n = (nDim == 2) ? n2 : n3; + su2double normal[3] = {0.0}, gridVel[3] = {0.0}; + for (unsigned short iDim = 0; iDim < nDim; iDim++) normal[iDim] = 1.3 * n[iDim]; + if (moving) { + gridVel[0] = 0.3; + gridVel[1] = -0.1; + gridVel[2] = 0.2; + } + numerics.SetNormal(normal); + numerics.SetGridVel(gridVel, gridVel); + su2double w = 0.0; + for (unsigned short iDim = 0; iDim < nDim; iDim++) w += gridVel[iDim] * n[iDim]; + + /*--- Two state pairs with different density, pressure and velocity (also tangential): flow along +n (sM > 0, + left branches) and along -n (sM < 0, right branches). ---*/ + const su2double tangent[3] = {0.1, -0.2, 0.15}; + for (const su2double dir : {1.0, -1.0}) { + su2double vel_i[3], vel_j[3]; + for (unsigned short iDim = 0; iDim < nDim; iDim++) { + vel_i[iDim] = dir * 0.5 * n[iDim] + tangent[iDim]; + vel_j[iDim] = dir * 0.25 * n[iDim] - tangent[iDim]; + } + const State Ui = Conservative(nDim, 1.2, vel_i, 1.1, k); + const State Uj = Conservative(nDim, 0.9, vel_j, 0.8, 0.8 * k); + + /*--- Wave speeds relative to the face, held fixed (the Jacobians treat them as constants): the supersonic and + the star branch on the side of the flow. ---*/ + const su2double supersonicSpeeds[2] = {dir > 0 ? 0.02 : -2.0, dir > 0 ? 2.0 : -0.02}; + const su2double starSpeeds[2] = {-2.0, 2.0}; + for (const auto* s : {supersonicSpeeds, starSpeeds}) { + numerics.SetFixedWaveSpeeds(true, s[0], s[1]); + const bool supersonic = (s[0] > 0.0) || (s[1] < 0.0); + const su2double sM = ContactSpeed(nDim, Ui, Uj, k, 0.8 * k, n, w, s[0], s[1]); + CAPTURE(nDim, moving, k, accurate, dir, s[0], s[1], sM); + + /*--- Make sure the intended branch is the one selected. ---*/ + if (supersonic) + REQUIRE(((dir > 0) ? s[0] : -s[1]) > 0.0); + else + REQUIRE(sM * dir > 0.0); + + const su2double error = JacobianError(nDim, numerics, config.get(), Ui, Uj, k, 0.8 * k); + if (accurate || supersonic) + CHECK(error < 1e-6); + else + CHECK(error > 1e-3); // the approximation is active + } + } +} + +} // namespace + +TEST_CASE("HLLC Jacobians match finite differences", "[HLLC]") { + for (const unsigned short nDim : {2, 3}) { + for (const bool moving : {false, true}) { + for (const bool accurate : {false, true}) { + CheckAllBranches(nDim, moving, 0.0, accurate); + CheckAllBranches(nDim, moving, 0.0, accurate); + } + } + } +} diff --git a/UnitTests/meson.build b/UnitTests/meson.build index 3dbbede4a513..8c2e61614e74 100644 --- a/UnitTests/meson.build +++ b/UnitTests/meson.build @@ -16,6 +16,7 @@ su2_cfd_tests = files(['Common/CConfig_tests.cpp', 'Common/toolboxes/multilayer_perceptron/CLookUp_ANN_tests.cpp', 'SU2_CFD/numerics/CNumerics_tests.cpp', 'SU2_CFD/numerics/roe_contact_tests.cpp', + 'SU2_CFD/numerics/hllc_jacobian_tests.cpp', 'SU2_CFD/fluid/CFluidModel_tests.cpp', 'SU2_CFD/gradients.cpp', 'SU2_CFD/solvers/CNEMOEulerSolver_tests.cpp', diff --git a/config_template.cfg b/config_template.cfg index 88533d33c689..5c7132467cc4 100644 --- a/config_template.cfg +++ b/config_template.cfg @@ -1856,6 +1856,7 @@ LOW_MACH_PREC= NO % % Use numerically computed Jacobians for AUSM+up(2) and SLAU(2) % Slower per iteration but potentialy more stable and capable of higher CFL +% For HLLC: exact Jacobians for fixed wave speeds (the default approximation is more robust with inexact linear solves) USE_ACCURATE_FLUX_JACOBIANS= NO % % Use the vectorized version of the selected numerical method (available for JST family and Roe). From 8e860f049c8fb0e9d264c1be341e003b4c40cc3c Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 21:38:48 +0200 Subject: [PATCH 13/14] Update the reference values of the regression tests from CI The under-relaxation fix changes the convergence history of the compressible implicit cases; the axisymmetric source Jacobian and the minimum pressure change their cases. All 95 failures of the CI run on 62ea3ca08b were value differences (no crash); edge_PPR differs at iteration 20 only (both versions reach -12.63 by iteration 50). Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- TestCases/hybrid_regression.py | 32 ++++++------- TestCases/parallel_regression.py | 74 ++++++++++++++--------------- TestCases/parallel_regression_AD.py | 6 +-- TestCases/serial_regression.py | 52 ++++++++++---------- TestCases/serial_regression_AD.py | 2 +- TestCases/tutorials.py | 20 ++++---- TestCases/vandv.py | 4 +- 7 files changed, 95 insertions(+), 95 deletions(-) diff --git a/TestCases/hybrid_regression.py b/TestCases/hybrid_regression.py index 149f4a921607..a574e2e1c4a3 100644 --- a/TestCases/hybrid_regression.py +++ b/TestCases/hybrid_regression.py @@ -72,7 +72,7 @@ def main(): naca0012.cfg_dir = "euler/naca0012" naca0012.cfg_file = "inv_NACA0012_Roe.cfg" naca0012.test_iter = 20 - naca0012.test_vals = [-4.031227, -3.524077, 0.272155, 0.009374] + naca0012.test_vals = [-4.029205, -3.521852, 0.272095, 0.009425] test_list.append(naca0012) # Supersonic wedge @@ -80,7 +80,7 @@ def main(): wedge.cfg_dir = "euler/wedge" wedge.cfg_file = "inv_wedge_HLLC.cfg" wedge.test_iter = 20 - wedge.test_vals = [-0.480786, 5.329367, -0.244869, 0.043126] + wedge.test_vals = [-4.223565, 1.506875, -0.249526, 0.043953] test_list.append(wedge) # ONERA M6 Wing @@ -104,7 +104,7 @@ def main(): bluntbody.cfg_dir = "euler/bluntbody" bluntbody.cfg_file = "blunt.cfg" bluntbody.test_iter = 20 - bluntbody.test_vals = [0.666182, 7.055173, -0.000631, 3.917770] + bluntbody.test_vals = [0.565173, 6.945846, 0.000171, 1.970669] test_list.append(bluntbody) ########################## @@ -170,7 +170,7 @@ def main(): rae2822_sa.cfg_dir = "rans/rae2822" rae2822_sa.cfg_file = "turb_SA_RAE2822.cfg" rae2822_sa.test_iter = 20 - rae2822_sa.test_vals = [-2.639525, -4.027096, 0.792680, 0.020549, 1000.000000] + rae2822_sa.test_vals = [-2.858565, -4.906734, 0.782692, 0.018639, 1000.000000] test_list.append(rae2822_sa) # RAE2822 SST @@ -178,7 +178,7 @@ def main(): rae2822_sst.cfg_dir = "rans/rae2822" rae2822_sst.cfg_file = "turb_SST_RAE2822.cfg" rae2822_sst.test_iter = 20 - rae2822_sst.test_vals = [-1.501958, 5.889330, 0.635453, 0.021771, 100.000000] + rae2822_sst.test_vals = [-1.461774, 5.883885, 0.630267, 0.020404, 100.000000] test_list.append(rae2822_sst) # RAE2822 SST_SUST @@ -202,7 +202,7 @@ def main(): turb_oneram6.cfg_dir = "rans/oneram6" turb_oneram6.cfg_file = "turb_ONERAM6.cfg" turb_oneram6.test_iter = 10 - turb_oneram6.test_vals = [-2.418707, -6.631575, 0.238586, 0.159599, 0.000000] + turb_oneram6.test_vals = [-2.419971, -6.631992, 0.238576, 0.159694, 0.000000] test_list.append(turb_oneram6) # NACA0012 (SA, FUN3D finest grid results: CL=1.0983, CD=0.01242) @@ -265,7 +265,7 @@ def main(): axi_rans_air_nozzle_restart.cfg_dir = "axisymmetric_rans/air_nozzle" axi_rans_air_nozzle_restart.cfg_file = "air_nozzle_restart.cfg" axi_rans_air_nozzle_restart.test_iter = 10 - axi_rans_air_nozzle_restart.test_vals = [-11.083068, -5.374686, -8.880093, -4.073548, 0.000000] + axi_rans_air_nozzle_restart.test_vals = [-12.066180, -7.429820, -8.817907, -3.734612, 0.000000] axi_rans_air_nozzle_restart.test_vals_aarch64 = [-14.140441, -9.154674, -10.886121, -5.806594, 0.000000] test_list.append(axi_rans_air_nozzle_restart) @@ -348,7 +348,7 @@ def main(): hb_rans_preconditioning.cfg_dir = "harmonic_balance/hb_rans_preconditioning" hb_rans_preconditioning.cfg_file = "davis.cfg" hb_rans_preconditioning.test_iter = 25 - hb_rans_preconditioning.test_vals = [-1.905220, 0.481910, 0.598990, 3.605349, -5.945852] + hb_rans_preconditioning.test_vals = [-1.905523, 0.482128, 0.598472, 3.605108, -5.946081] test_list.append(hb_rans_preconditioning) ############################# @@ -558,7 +558,7 @@ def main(): unst_inc_turb_naca0015_sa.cfg_dir = "unsteady/pitching_naca0015_rans_inc" unst_inc_turb_naca0015_sa.cfg_file = "config_incomp_turb_sa.cfg" unst_inc_turb_naca0015_sa.test_iter = 1 - unst_inc_turb_naca0015_sa.test_vals = [-3.008640, -6.889022, 1.435181, 0.433567] + unst_inc_turb_naca0015_sa.test_vals = [-3.007258, -6.887901, 1.435511, 0.433955] unst_inc_turb_naca0015_sa.unsteady = True test_list.append(unst_inc_turb_naca0015_sa) @@ -581,7 +581,7 @@ def main(): edge_VW.cfg_dir = "nicf/edge" edge_VW.cfg_file = "edge_VW.cfg" edge_VW.test_iter = 25 - edge_VW.test_vals = [-11.670757, -5.483973, -0.000009, 0.000000] + edge_VW.test_vals = [-11.723086, -5.536517, -0.000009, 0.000000] test_list.append(edge_VW) # Rarefaction shock wave edge_PPR @@ -589,7 +589,7 @@ def main(): edge_PPR.cfg_dir = "nicf/edge" edge_PPR.cfg_file = "edge_PPR.cfg" edge_PPR.test_iter = 20 - edge_PPR.test_vals = [-12.338051, -6.157966, -0.000034, 0.000000] + edge_PPR.test_vals = [-10.993286, -4.827527, -0.000034, 0.000000] edge_PPR.test_vals_aarch64 = [ -7.139211, -0.980821, -0.000034, 0.000000] test_list.append(edge_PPR) @@ -659,7 +659,7 @@ def main(): channel_2D.cfg_dir = "sliding_interface/channel_2D" channel_2D.cfg_file = "channel_2D_WA.cfg" channel_2D.test_iter = 2 - channel_2D.test_vals = [2.000000, 0.000000, 0.466194, 0.350087, 0.398979] + channel_2D.test_vals = [2.000000, 0.000000, 0.466339, 0.346858, 0.397486] channel_2D.unsteady = True channel_2D.multizone = True test_list.append(channel_2D) @@ -669,7 +669,7 @@ def main(): channel_3D.cfg_dir = "sliding_interface/channel_3D" channel_3D.cfg_file = "channel_3D_WA.cfg" channel_3D.test_iter = 2 - channel_3D.test_vals = [2.000000, 0.000000, 0.632251, 0.534213, 0.431803] + channel_3D.test_vals = [2.000000, 0.000000, 0.632749, 0.522745, 0.423269] channel_3D.test_vals_aarch64 = [2.000000, 0.000000, 0.629112, 0.524948, 0.422396] channel_3D.unsteady = True channel_3D.multizone = True @@ -681,7 +681,7 @@ def main(): pipe.cfg_dir = "sliding_interface/pipe" pipe.cfg_file = "pipe_NN.cfg" pipe.test_iter = 2 - pipe.test_vals = [0.092415, 0.568971, 0.692864, 0.989451, 1.048245] + pipe.test_vals = [0.070225, 0.528891, 0.684177, 0.970350, 1.065268] pipe.unsteady = True pipe.multizone = True test_list.append(pipe) @@ -691,7 +691,7 @@ def main(): rotating_cylinders.cfg_dir = "sliding_interface/rotating_cylinders" rotating_cylinders.cfg_file = "rot_cylinders_WA.cfg" rotating_cylinders.test_iter = 3 - rotating_cylinders.test_vals = [3.000000, 0.000000, 0.664817, 1.125803, 1.117608] + rotating_cylinders.test_vals = [3.000000, 0.000000, 0.686432, 1.159463, 1.151135] rotating_cylinders.unsteady = True rotating_cylinders.multizone = True test_list.append(rotating_cylinders) @@ -701,7 +701,7 @@ def main(): supersonic_vortex_shedding.cfg_dir = "sliding_interface/supersonic_vortex_shedding" supersonic_vortex_shedding.cfg_file = "sup_vor_shed_WA.cfg" supersonic_vortex_shedding.test_iter = 5 - supersonic_vortex_shedding.test_vals = [5.000000, 0.000000, 0.899645, 1.076181] + supersonic_vortex_shedding.test_vals = [5.000000, 0.000000, 1.126036, 1.081765] supersonic_vortex_shedding.unsteady = True supersonic_vortex_shedding.multizone = True test_list.append(supersonic_vortex_shedding) diff --git a/TestCases/parallel_regression.py b/TestCases/parallel_regression.py index 36766da22528..d60bfa8de4cb 100755 --- a/TestCases/parallel_regression.py +++ b/TestCases/parallel_regression.py @@ -266,7 +266,7 @@ def main(): naca0012.cfg_dir = "euler/naca0012" naca0012.cfg_file = "inv_NACA0012_Roe_su2bin.cfg" naca0012.test_iter = 20 - naca0012.test_vals = [-4.026356, -3.521920, 0.274929, 0.009852] + naca0012.test_vals = [-4.024996, -3.520626, 0.275033, 0.009889] test_list.append(naca0012) # Supersonic wedge @@ -274,7 +274,7 @@ def main(): wedge.cfg_dir = "euler/wedge" wedge.cfg_file = "inv_wedge_HLLC.cfg" wedge.test_iter = 20 - wedge.test_vals = [-4.142624, 1.585638, -0.249527, 0.043953] + wedge.test_vals = [-4.018855, 1.708167, -0.249527, 0.043953] test_list.append(wedge) # ONERA M6 Wing @@ -300,7 +300,7 @@ def main(): polar_naca0012.cfg_file = "inv_NACA0012.cfg" polar_naca0012.polar = True polar_naca0012.test_iter = 10 - polar_naca0012.test_vals = [-1.315764, 4.158361, 0.002646, 0.114559] + polar_naca0012.test_vals = [-1.315882, 4.158313, 0.002698, 0.114644] polar_naca0012.test_vals_aarch64 = [-1.083394, 4.386134, 0.001588, 0.033513] polar_naca0012.command = TestCase.Command(exec = "compute_polar.py", param = "-i 11") # flaky test on arm64 @@ -312,7 +312,7 @@ def main(): bluntbody.cfg_dir = "euler/bluntbody" bluntbody.cfg_file = "blunt.cfg" bluntbody.test_iter = 20 - bluntbody.test_vals = [0.666133, 7.054751, -0.000051, 3.912894] + bluntbody.test_vals = [0.565207, 6.945909, -0.000001, 1.970304] test_list.append(bluntbody) # Equivalent area NACA64-206 @@ -320,7 +320,7 @@ def main(): ea_naca64206.cfg_dir = "optimization_euler/equivalentarea_naca64206" ea_naca64206.cfg_file = "NACA64206.cfg" ea_naca64206.test_iter = 10 - ea_naca64206.test_vals = [-1.125734, -0.473873, -0.002416, 67775.000000] + ea_naca64206.test_vals = [-1.126282, -0.474499, -0.002405, 67775.000000] test_list.append(ea_naca64206) # SUPERSONIC FLOW PAST A RAMP IN A CHANNEL @@ -345,7 +345,7 @@ def main(): MFR_coupling.cfg_dir = "euler/turbofan_MFR_coupling" MFR_coupling.cfg_file = "MFR_coupling.cfg" MFR_coupling.test_iter = 100 - MFR_coupling.test_vals = [-211.240000, 150.030000, 20.151000] + MFR_coupling.test_vals = [-211.250000, 150.030000, 20.151000] test_list.append(MFR_coupling) ########################## @@ -390,7 +390,7 @@ def main(): cylinder.cfg_dir = "navierstokes/cylinder" cylinder.cfg_file = "lam_cylinder.cfg" cylinder.test_iter = 25 - cylinder.test_vals = [-8.557685, -3.072812, -0.018904, 1.610320, 0.000000] + cylinder.test_vals = [-8.557718, -3.072847, -0.018912, 1.610317, 0.000000] test_list.append(cylinder) # Laminar cylinder (low Mach correction) @@ -430,7 +430,7 @@ def main(): rae2822_sa.cfg_dir = "rans/rae2822" rae2822_sa.cfg_file = "turb_SA_RAE2822.cfg" rae2822_sa.test_iter = 20 - rae2822_sa.test_vals = [-2.888183, -5.136146, 0.788144, 0.019032, 1000.000000] + rae2822_sa.test_vals = [-2.998634, -4.043296, 0.794248, 0.019260, 1000.000000] test_list.append(rae2822_sa) # RAE2822 SST @@ -438,7 +438,7 @@ def main(): rae2822_sst.cfg_dir = "rans/rae2822" rae2822_sst.cfg_file = "turb_SST_RAE2822.cfg" rae2822_sst.test_iter = 20 - rae2822_sst.test_vals = [-1.495621, 5.875629, 0.638353, 0.020410, 100.000000] + rae2822_sst.test_vals = [-1.485957, 5.891491, 0.675200, 0.023427, 100.000000] test_list.append(rae2822_sst) # RAE2822 SST_SUST @@ -470,7 +470,7 @@ def main(): turb_flatplate_CC_Wilcox.cfg_dir = "rans/flatplate" turb_flatplate_CC_Wilcox.cfg_file = "turb_SST_flatplate_compressibility_Wilcox.cfg" turb_flatplate_CC_Wilcox.test_iter = 20 - turb_flatplate_CC_Wilcox.test_vals = [-1.169750, 2.167973, 1.488567, 5.244807, -3.795767, 9.699598] + turb_flatplate_CC_Wilcox.test_vals = [-1.204997, 2.149390, 1.396996, 5.216302, -3.821981, 9.707535] test_list.append(turb_flatplate_CC_Wilcox) # Flat plate SST compressibility correction Sarkar @@ -478,7 +478,7 @@ def main(): turb_flatplate_CC_Sarkar.cfg_dir = "rans/flatplate" turb_flatplate_CC_Sarkar.cfg_file = "turb_SST_flatplate_compressibility_Sarkar.cfg" turb_flatplate_CC_Sarkar.test_iter = 20 - turb_flatplate_CC_Sarkar.test_vals = [-1.169750, 2.167973, 1.488567, 5.244807, -3.795765, 9.699598] + turb_flatplate_CC_Sarkar.test_vals = [-1.204997, 2.149390, 1.396996, 5.216302, -3.821978, 9.707535] test_list.append(turb_flatplate_CC_Sarkar) # FLAT PLATE, ROUGHNESS BC KNOPP SST @@ -502,7 +502,7 @@ def main(): turb_oneram6.cfg_dir = "rans/oneram6" turb_oneram6.cfg_file = "turb_ONERAM6.cfg" turb_oneram6.test_iter = 10 - turb_oneram6.test_vals = [-2.418711, -6.631576, 0.238587, 0.159599, 0.000000] + turb_oneram6.test_vals = [-2.419976, -6.631992, 0.238577, 0.159694, 0.000000] turb_oneram6.timeout = 3200 test_list.append(turb_oneram6) @@ -511,7 +511,7 @@ def main(): turb_oneram6_vc.cfg_dir = "rans/oneram6" turb_oneram6_vc.cfg_file = "turb_ONERAM6_vc.cfg" turb_oneram6_vc.test_iter = 15 - turb_oneram6_vc.test_vals = [-2.294311, -6.572307, 0.234416, 0.144830, 0.000000] + turb_oneram6_vc.test_vals = [-2.295077, -6.572734, 0.234428, 0.145063, 0.000000] turb_oneram6_vc.timeout = 3200 test_list.append(turb_oneram6_vc) @@ -520,7 +520,7 @@ def main(): turb_oneram6_nk.cfg_dir = "rans/oneram6" turb_oneram6_nk.cfg_file = "turb_ONERAM6_nk.cfg" turb_oneram6_nk.test_iter = 20 - turb_oneram6_nk.test_vals = [-4.848809, -4.451213, -11.426556, 0.221554, 0.049197, 2.000000, -0.880072, 10.000000] + turb_oneram6_nk.test_vals = [-4.848121, -4.450658, -11.426212, 0.221522, 0.049358, 2.000000, -0.879182, 10.000000] turb_oneram6_nk.timeout = 600 turb_oneram6_nk.tol = 0.0001 test_list.append(turb_oneram6_nk) @@ -619,7 +619,7 @@ def main(): axi_rans_air_nozzle_restart.cfg_dir = "axisymmetric_rans/air_nozzle" axi_rans_air_nozzle_restart.cfg_file = "air_nozzle_restart.cfg" axi_rans_air_nozzle_restart.test_iter = 10 - axi_rans_air_nozzle_restart.test_vals = [-11.056236, -5.334116, -8.842310, -4.067917, 0.000000] + axi_rans_air_nozzle_restart.test_vals = [-12.069274, -7.508838, -8.813357, -3.732828, 0.000000] axi_rans_air_nozzle_restart.test_vals_aarch64 = [-14.143310, -9.163287, -10.858232, -5.787715, 0.000000] axi_rans_air_nozzle_restart.tol = 0.0001 test_list.append(axi_rans_air_nozzle_restart) @@ -664,7 +664,7 @@ def main(): inc_nozzle.cfg_dir = "incomp_euler/nozzle" inc_nozzle.cfg_file = "inv_nozzle.cfg" inc_nozzle.test_iter = 20 - inc_nozzle.test_vals = [-6.051005, -5.361765, 0.006999, 0.127295] + inc_nozzle.test_vals = [-6.357399, -5.618580, 0.004548, 0.126243] test_list.append(inc_nozzle) # Laminar wall mounted cylinder, Euler walls, cylinder wall diagonally split @@ -883,7 +883,7 @@ def main(): turbmod_sa_bsl_rae2822.cfg_dir = "turbulence_models/sa/rae2822" turbmod_sa_bsl_rae2822.cfg_file = "turb_SA_BSL_RAE2822.cfg" turbmod_sa_bsl_rae2822.test_iter = 20 - turbmod_sa_bsl_rae2822.test_vals = [-3.182935, -0.447793, -0.875861, -3.286652, 0.769925, 0.017709] + turbmod_sa_bsl_rae2822.test_vals = [-3.191783, -0.473872, -0.892710, -5.169793, 0.777160, 0.017964] test_list.append(turbmod_sa_bsl_rae2822) # SA Negative @@ -891,7 +891,7 @@ def main(): turbmod_sa_neg_rae2822.cfg_dir = "turbulence_models/sa/rae2822" turbmod_sa_neg_rae2822.cfg_file = "turb_SA_NEG_RAE2822.cfg" turbmod_sa_neg_rae2822.test_iter = 10 - turbmod_sa_neg_rae2822.test_vals = [1.570952, 1.355242, -0.869320, 1.244343, 0.452923, 0.000000] + turbmod_sa_neg_rae2822.test_vals = [1.429980, 1.127884, 4.865536, 1.277981, 0.472937, 0.000000] turbmod_sa_neg_rae2822.test_vals_aarch64 = [-1.345593, 1.448310, 1.208721, -0.846597, 1.248410, 0.489117, 0.000000] test_list.append(turbmod_sa_neg_rae2822) @@ -900,7 +900,7 @@ def main(): turbmod_sa_comp_rae2822.cfg_dir = "turbulence_models/sa/rae2822" turbmod_sa_comp_rae2822.cfg_file = "turb_SA_COMP_RAE2822.cfg" turbmod_sa_comp_rae2822.test_iter = 20 - turbmod_sa_comp_rae2822.test_vals = [-3.192776, -0.471890, -0.893254, -5.190537, 0.777261, 0.017904] + turbmod_sa_comp_rae2822.test_vals = [-3.190862, -0.473847, -0.893269, -5.195424, 0.777296, 0.017920] test_list.append(turbmod_sa_comp_rae2822) # SA Edwards @@ -908,7 +908,7 @@ def main(): turbmod_sa_edw_rae2822.cfg_dir = "turbulence_models/sa/rae2822" turbmod_sa_edw_rae2822.cfg_file = "turb_SA_EDW_RAE2822.cfg" turbmod_sa_edw_rae2822.test_iter = 20 - turbmod_sa_edw_rae2822.test_vals = [-3.198033, -0.439446, -0.883563, -3.899428, 0.786167, 0.017883] + turbmod_sa_edw_rae2822.test_vals = [-3.189014, -0.465457, -0.890004, -5.206434, 0.781697, 0.017765] test_list.append(turbmod_sa_edw_rae2822) # SA Compressibility and Edwards @@ -916,7 +916,7 @@ def main(): turbmod_sa_comp_edw_rae2822.cfg_dir = "turbulence_models/sa/rae2822" turbmod_sa_comp_edw_rae2822.cfg_file = "turb_SA_COMP_EDW_RAE2822.cfg" turbmod_sa_comp_edw_rae2822.test_iter = 20 - turbmod_sa_comp_edw_rae2822.test_vals = [-2.943295, -0.168387, -0.632480, -5.127730, 0.802507, 0.020196] + turbmod_sa_comp_edw_rae2822.test_vals = [-3.056328, -0.293497, -0.766206, -5.230932, 0.798848, 0.019838] test_list.append(turbmod_sa_comp_edw_rae2822) # SA QCR @@ -924,7 +924,7 @@ def main(): turbmod_sa_qcr_rae2822.cfg_dir = "turbulence_models/sa/rae2822" turbmod_sa_qcr_rae2822.cfg_file = "turb_SA_QCR_RAE2822.cfg" turbmod_sa_qcr_rae2822.test_iter = 20 - turbmod_sa_qcr_rae2822.test_vals = [-3.056565, -0.419621, -0.805532, -3.857890, 0.760796, 0.016908] + turbmod_sa_qcr_rae2822.test_vals = [-3.149635, -0.387144, -0.827745, -2.226502, 0.780050, 0.015548] test_list.append(turbmod_sa_qcr_rae2822) ############################ @@ -1050,7 +1050,7 @@ def main(): turb_naca0012_2c.cfg_dir = "rans_uq/naca0012" turb_naca0012_2c.cfg_file = "turb_NACA0012_uq_2c.cfg" turb_naca0012_2c.test_iter = 10 - turb_naca0012_2c.test_vals = [-5.482691, 1.262918, 0.529450, -0.023537] + turb_naca0012_2c.test_vals = [-5.482691, 1.262918, 0.529442, -0.023584] turb_naca0012_2c.test_vals_aarch64 = [-5.484365, 1.264701, 0.501741, -0.033109] test_list.append(turb_naca0012_2c) @@ -1097,7 +1097,7 @@ def main(): hb_rans_preconditioning.cfg_file = "davis.cfg" hb_rans_preconditioning.test_iter = 25 hb_rans_preconditioning.tol = 0.00001 - hb_rans_preconditioning.test_vals = [-1.905206, 0.481900, 0.599004, 3.605361, -5.945837] + hb_rans_preconditioning.test_vals = [-1.905509, 0.482119, 0.598486, 3.605120, -5.946067] test_list.append(hb_rans_preconditioning) ###################################### @@ -1109,7 +1109,7 @@ def main(): rot_naca0012.cfg_dir = "rotating/naca0012" rot_naca0012.cfg_file = "rot_NACA0012.cfg" rot_naca0012.test_iter = 25 - rot_naca0012.test_vals = [-0.023432, 0.083323, 2.000000, -1.181593] + rot_naca0012.test_vals = [-0.023360, 0.082942, 2.000000, -1.182142] test_list.append(rot_naca0012) # Lid-driven cavity @@ -1125,7 +1125,7 @@ def main(): spinning_cylinder.cfg_dir = "moving_wall/spinning_cylinder" spinning_cylinder.cfg_file = "spinning_cylinder.cfg" spinning_cylinder.test_iter = 25 - spinning_cylinder.test_vals = [-7.471334, -2.007097, 2.183550, 1.847845] + spinning_cylinder.test_vals = [-7.471322, -2.007082, 2.183565, 1.847848] test_list.append(spinning_cylinder) ###################################### @@ -1175,7 +1175,7 @@ def main(): unst_inc_turb_naca0015_sa.cfg_dir = "unsteady/pitching_naca0015_rans_inc" unst_inc_turb_naca0015_sa.cfg_file = "config_incomp_turb_sa.cfg" unst_inc_turb_naca0015_sa.test_iter = 1 - unst_inc_turb_naca0015_sa.test_vals = [-3.004011, -6.876257, 1.487888, 0.421869] + unst_inc_turb_naca0015_sa.test_vals = [-3.002783, -6.875096, 1.487783, 0.421690] unst_inc_turb_naca0015_sa.unsteady = True test_list.append(unst_inc_turb_naca0015_sa) @@ -1197,7 +1197,7 @@ def main(): edge_VW.cfg_dir = "nicf/edge" edge_VW.cfg_file = "edge_VW.cfg" edge_VW.test_iter = 20 - edge_VW.test_vals = [-4.444040, 1.757269, -0.000009, 0.000000] + edge_VW.test_vals = [-4.478306, 1.723006, -0.000009, 0.000000] test_list.append(edge_VW) # Rarefaction shock wave edge_PPR @@ -1205,7 +1205,7 @@ def main(): edge_PPR.cfg_dir = "nicf/edge" edge_PPR.cfg_file = "edge_PPR.cfg" edge_PPR.test_iter = 20 - edge_PPR.test_vals = [-9.533869, -3.385210, -0.000034, 0.000000] + edge_PPR.test_vals = [-8.395110, -2.236894, -0.000034, 0.000000] test_list.append(edge_PPR) # Rarefaction Q1D nozzle, include CoolProp fluid model @@ -1231,7 +1231,7 @@ def main(): datadriven_fluidModel.cfg_dir = "nicf/datadriven" datadriven_fluidModel.cfg_file = "datadriven_nozzle.cfg" datadriven_fluidModel.test_iter = 50 - datadriven_fluidModel.test_vals = [-4.879421, -2.353319, -2.344004, 0.390613, -1.559800, 1.322761] + datadriven_fluidModel.test_vals = [-6.283594, -3.253961, -4.172935, -0.949415, -2.046039, 1.546933] test_list.append(datadriven_fluidModel) ###################################### @@ -1260,7 +1260,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.167491, 1.598507, -2.928578, 2.573644, -2.527392, 3.016170, 106370.000000, 106370.000000, 5.726800, 64.383000] + axial_stage2D.test_vals = [1.167500, 1.598546, -2.928573, 2.573647, -2.527331, 3.016254, 106370.000000, 106370.000000, 5.726800, 64.383000] test_list.append(axial_stage2D) # 2D transonic stator restart @@ -1308,7 +1308,7 @@ def main(): channel_2D.cfg_dir = "sliding_interface/channel_2D" channel_2D.cfg_file = "channel_2D_WA.cfg" channel_2D.test_iter = 2 - channel_2D.test_vals = [2.000000, 0.000000, 0.466215, 0.350091, 0.398996] + channel_2D.test_vals = [2.000000, 0.000000, 0.466408, 0.346806, 0.397328] channel_2D.timeout = 100 channel_2D.unsteady = True channel_2D.multizone = True @@ -1319,7 +1319,7 @@ def main(): channel_3D.cfg_dir = "sliding_interface/channel_3D" channel_3D.cfg_file = "channel_3D_WA.cfg" channel_3D.test_iter = 2 - channel_3D.test_vals = [2.000000, 0.000000, 0.632256, 0.534230, 0.431926] + channel_3D.test_vals = [2.000000, 0.000000, 0.632752, 0.522763, 0.423355] channel_3D.test_vals_aarch64 = [2.000000, 0.000000, 0.629119, 0.524959, 0.422390] channel_3D.unsteady = True channel_3D.multizone = True @@ -1330,7 +1330,7 @@ def main(): pipe.cfg_dir = "sliding_interface/pipe" pipe.cfg_file = "pipe_NN.cfg" pipe.test_iter = 2 - pipe.test_vals = [0.092415, 0.568973, 0.692859, 0.989459, 1.048237] + pipe.test_vals = [0.070225, 0.528893, 0.684173, 0.970355, 1.065262] pipe.unsteady = True pipe.multizone = True test_list.append(pipe) @@ -1340,7 +1340,7 @@ def main(): rotating_cylinders.cfg_dir = "sliding_interface/rotating_cylinders" rotating_cylinders.cfg_file = "rot_cylinders_WA.cfg" rotating_cylinders.test_iter = 3 - rotating_cylinders.test_vals = [3.000000, 0.000000, 0.664821, 1.125807, 1.117610] + rotating_cylinders.test_vals = [3.000000, 0.000000, 0.686431, 1.159459, 1.151139] rotating_cylinders.unsteady = True rotating_cylinders.multizone = True test_list.append(rotating_cylinders) @@ -1350,7 +1350,7 @@ def main(): supersonic_vortex_shedding.cfg_dir = "sliding_interface/supersonic_vortex_shedding" supersonic_vortex_shedding.cfg_file = "sup_vor_shed_WA.cfg" supersonic_vortex_shedding.test_iter = 5 - supersonic_vortex_shedding.test_vals = [5.000000, 0.000000, 0.899643, 1.076189] + supersonic_vortex_shedding.test_vals = [5.000000, 0.000000, 1.126036, 1.081775] supersonic_vortex_shedding.unsteady = True supersonic_vortex_shedding.multizone = True test_list.append(supersonic_vortex_shedding) @@ -1551,7 +1551,7 @@ def main(): pywrapper_naca0012.cfg_dir = "euler/naca0012" pywrapper_naca0012.cfg_file = "inv_NACA0012_Roe.cfg" pywrapper_naca0012.test_iter = 80 - pywrapper_naca0012.test_vals = [-5.699819, -5.186257, 0.334856, 0.023343] + pywrapper_naca0012.test_vals = [-5.701280, -5.186480, 0.334856, 0.023342] pywrapper_naca0012.command = TestCase.Command("mpirun -np 2", "SU2_CFD.py", "--parallel -f") test_list.append(pywrapper_naca0012) diff --git a/TestCases/parallel_regression_AD.py b/TestCases/parallel_regression_AD.py index 0b5d8724f62b..f5126ace6feb 100644 --- a/TestCases/parallel_regression_AD.py +++ b/TestCases/parallel_regression_AD.py @@ -150,7 +150,7 @@ def main(): discadj_axisymmetric_rans_nozzle.cfg_dir = "axisymmetric_rans/air_nozzle" discadj_axisymmetric_rans_nozzle.cfg_file = "air_nozzle_restart.cfg" discadj_axisymmetric_rans_nozzle.test_iter = 10 - discadj_axisymmetric_rans_nozzle.test_vals = [9.909657, 5.078045, 7.129068, 2.490955] + discadj_axisymmetric_rans_nozzle.test_vals = [9.733989, 5.141750, 7.108464, 2.490402] discadj_axisymmetric_rans_nozzle.no_restart = True test_list.append(discadj_axisymmetric_rans_nozzle) @@ -239,7 +239,7 @@ def main(): discadj_sliding_interface.cfg_dir = "disc_adj_ffi/sliding_interface" discadj_sliding_interface.cfg_file = "circles.cfg" discadj_sliding_interface.test_iter = 10 - discadj_sliding_interface.test_vals = [10.000000, -4.395153, -4.419435] + discadj_sliding_interface.test_vals = [10.000000, -4.403529, -4.419435] test_list.append(discadj_sliding_interface) ####################################################### @@ -315,7 +315,7 @@ def main(): discadj_fsi2.cfg_dir = "disc_adj_fsi/Airfoil_2d" discadj_fsi2.cfg_file = "config.cfg" discadj_fsi2.test_iter = 8 - discadj_fsi2.test_vals = [-3.824633, 1.979516, -3.863368, 0.295450, 3.839800] + discadj_fsi2.test_vals = [-3.824591, 1.979595, -3.863368, 0.295450, 3.839800] discadj_fsi2.test_vals_aarch64 = [-3.824870, 1.979160, -3.863368, 0.295450, 3.839800] discadj_fsi2.tol = 0.00001 test_list.append(discadj_fsi2) diff --git a/TestCases/serial_regression.py b/TestCases/serial_regression.py index ab1de7fe1431..967a4eaf2d97 100755 --- a/TestCases/serial_regression.py +++ b/TestCases/serial_regression.py @@ -103,7 +103,7 @@ def main(): naca0012.cfg_dir = "euler/naca0012" naca0012.cfg_file = "inv_NACA0012_Roe.cfg" naca0012.test_iter = 20 - naca0012.test_vals = [-4.032065, -3.525196, 0.273051, 0.009680] + naca0012.test_vals = [-4.031446, -3.524530, 0.273056, 0.009703] test_list.append(naca0012) # NACA0012 - FMG test @@ -111,7 +111,7 @@ def main(): naca0012_FMG.cfg_dir = "euler/naca0012" naca0012_FMG.cfg_file = "inv_NACA0012.cfg" naca0012_FMG.test_iter = 20 - naca0012_FMG.test_vals = [-4.160277, -3.553897, 0.296156, 0.016889] + naca0012_FMG.test_vals = [-4.160287, -3.553906, 0.296154, 0.016890] test_list.append(naca0012_FMG) # Supersonic wedge @@ -119,7 +119,7 @@ def main(): wedge.cfg_dir = "euler/wedge" wedge.cfg_file = "inv_wedge_HLLC.cfg" wedge.test_iter = 20 - wedge.test_vals = [-4.399745, 1.331048, -0.249527, 0.043953] + wedge.test_vals = [-4.359580, 1.370530, -0.249526, 0.043953] test_list.append(wedge) # ONERA M6 Wing @@ -127,7 +127,7 @@ def main(): oneram6.cfg_dir = "euler/oneram6" oneram6.cfg_file = "inv_ONERAM6.cfg" oneram6.test_iter = 10 - oneram6.test_vals = [-4.825465, -4.205219, 0.277033, 0.006614] + oneram6.test_vals = [-4.825509, -4.205271, 0.277086, 0.006619] oneram6.timeout = 9600 test_list.append(oneram6) @@ -145,7 +145,7 @@ def main(): polar_naca0012.cfg_file = "inv_NACA0012.cfg" polar_naca0012.polar = True polar_naca0012.test_iter = 10 - polar_naca0012.test_vals = [-1.311660, 4.165589, 0.000067, 0.113157] + polar_naca0012.test_vals = [-1.311796, 4.165508, 0.000092, 0.113306] polar_naca0012.test_vals_aarch64 = [-1.063447, 4.401847, 0.000291, 0.031696] polar_naca0012.command = TestCase.Command(exec = "compute_polar.py", param = "-n 1 -i 11") # flaky test on arm64 @@ -157,7 +157,7 @@ def main(): bluntbody.cfg_dir = "euler/bluntbody" bluntbody.cfg_file = "blunt.cfg" bluntbody.test_iter = 20 - bluntbody.test_vals = [0.666162, 7.055023, -0.000093, 3.916825] + bluntbody.test_vals = [0.565185, 6.945875, 0.000007, 1.970309] test_list.append(bluntbody) ########################## @@ -184,7 +184,7 @@ def main(): cylinder.cfg_dir = "navierstokes/cylinder" cylinder.cfg_file = "lam_cylinder.cfg" cylinder.test_iter = 25 - cylinder.test_vals = [-8.725123, -3.234819, -0.028022, 1.605527, 0.000000] + cylinder.test_vals = [-8.725171, -3.234867, -0.028017, 1.605526, 0.000000] test_list.append(cylinder) # Laminar cylinder (low Mach correction) @@ -228,7 +228,7 @@ def main(): rae2822_sa.cfg_dir = "rans/rae2822" rae2822_sa.cfg_file = "turb_SA_RAE2822.cfg" rae2822_sa.test_iter = 20 - rae2822_sa.test_vals = [-3.131221, -5.161756, 0.778275, 0.019566, 1000.000000] + rae2822_sa.test_vals = [-3.056253, -5.155475, 0.783190, 0.020552, 1000.000000] test_list.append(rae2822_sa) # RAE2822 SST @@ -236,7 +236,7 @@ def main(): rae2822_sst.cfg_dir = "rans/rae2822" rae2822_sst.cfg_file = "turb_SST_RAE2822.cfg" rae2822_sst.test_iter = 20 - rae2822_sst.test_vals = [-1.745363, -1.484660, 5.883796, 0.579910, 0.017299, 100.000000] + rae2822_sst.test_vals = [-1.916371, -1.461911, 5.894621, 0.611793, 0.016873, 100.000000] test_list.append(rae2822_sst) # RAE2822 SST_SUST @@ -244,7 +244,7 @@ def main(): rae2822_sst_sust.cfg_dir = "rans/rae2822" rae2822_sst_sust.cfg_file = "turb_SST_SUST_RAE2822.cfg" rae2822_sst_sust.test_iter = 20 - rae2822_sst_sust.test_vals = [-2.464319, 5.852571, 0.490309, 0.042521] + rae2822_sst_sust.test_vals = [-2.464309, 5.852571, 0.490309, 0.042519] test_list.append(rae2822_sst_sust) # Flat plate @@ -284,7 +284,7 @@ def main(): turb_oneram6.cfg_dir = "rans/oneram6" turb_oneram6.cfg_file = "turb_ONERAM6.cfg" turb_oneram6.test_iter = 10 - turb_oneram6.test_vals = [-2.418713, -6.631577, 0.238587, 0.159599, 0.000000] + turb_oneram6.test_vals = [-2.419977, -6.631993, 0.238577, 0.159694, 0.000000] turb_oneram6.timeout = 3200 test_list.append(turb_oneram6) @@ -367,7 +367,7 @@ def main(): axi_rans_air_nozzle_restart.cfg_dir = "axisymmetric_rans/air_nozzle" axi_rans_air_nozzle_restart.cfg_file = "air_nozzle_restart.cfg" axi_rans_air_nozzle_restart.test_iter = 10 - axi_rans_air_nozzle_restart.test_vals = [-11.054281, -5.328905, -8.835591, -4.056830, 0.000000] + axi_rans_air_nozzle_restart.test_vals = [-12.067607, -7.525507, -8.816792, -3.735401, 0.000000] axi_rans_air_nozzle_restart.test_vals_aarch64 = [-14.143715, -9.170705, -10.848554, -5.776746, 0.000000] axi_rans_air_nozzle_restart.tol = 0.0001 test_list.append(axi_rans_air_nozzle_restart) @@ -377,7 +377,7 @@ def main(): axi_rans_air_nozzle_species.cfg_dir = "axisymmetric_rans/air_nozzle" axi_rans_air_nozzle_species.cfg_file = "air_nozzle_species.cfg" axi_rans_air_nozzle_species.test_iter = 10 - axi_rans_air_nozzle_species.test_vals = [-1.690665, 3.882506, -2.928702, 5.760933, -3.144560, 0.000000] + axi_rans_air_nozzle_species.test_vals = [-2.078338, 3.412490, -1.554413, 5.882895, -2.377710, 0.000000] axi_rans_air_nozzle_species.tol = 0.0001 test_list.append(axi_rans_air_nozzle_species) @@ -744,7 +744,7 @@ def main(): turb_naca0012_1c.cfg_dir = "rans_uq/naca0012" turb_naca0012_1c.cfg_file = "turb_NACA0012_uq_1c.cfg" turb_naca0012_1c.test_iter = 10 - turb_naca0012_1c.test_vals = [-4.981181, 1.343583, 0.597121, 0.017223] + turb_naca0012_1c.test_vals = [-4.981181, 1.343583, 0.597111, 0.017244] turb_naca0012_1c.test_vals_aarch64 = [-4.992791, 1.342873, 0.557941, 0.003269] test_list.append(turb_naca0012_1c) @@ -798,7 +798,7 @@ def main(): hb_rans_preconditioning.cfg_dir = "harmonic_balance/hb_rans_preconditioning" hb_rans_preconditioning.cfg_file = "davis.cfg" hb_rans_preconditioning.test_iter = 25 - hb_rans_preconditioning.test_vals = [-1.905206, 0.481900, 0.599004, 3.605363, -5.945837] + hb_rans_preconditioning.test_vals = [-1.905509, 0.482119, 0.598486, 3.605121, -5.946067] test_list.append(hb_rans_preconditioning) ###################################### @@ -810,7 +810,7 @@ def main(): rot_naca0012.cfg_dir = "rotating/naca0012" rot_naca0012.cfg_file = "rot_NACA0012.cfg" rot_naca0012.test_iter = 25 - rot_naca0012.test_vals = [-1.242808, 4.266616, -0.022207, 0.084206, 2.000000, -1.252601] + rot_naca0012.test_vals = [-1.240352, 4.268727, -0.022164, 0.083891, 2.000000, -1.253214] test_list.append(rot_naca0012) # Lid-driven cavity @@ -874,7 +874,7 @@ def main(): unst_inc_turb_naca0015_sa.cfg_dir = "unsteady/pitching_naca0015_rans_inc" unst_inc_turb_naca0015_sa.cfg_file = "config_incomp_turb_sa.cfg" unst_inc_turb_naca0015_sa.test_iter = 1 - unst_inc_turb_naca0015_sa.test_vals = [-3.007635, -6.879816, 1.445300, 0.419281] + unst_inc_turb_naca0015_sa.test_vals = [-3.005435, -6.878657, 1.445291, 0.420526] unst_inc_turb_naca0015_sa.unsteady = True test_list.append(unst_inc_turb_naca0015_sa) # unsteady pitching NACA64A010, RANS, SA @@ -928,7 +928,7 @@ def main(): edge_VW.cfg_dir = "nicf/edge" edge_VW.cfg_file = "edge_VW.cfg" edge_VW.test_iter = 20 - edge_VW.test_vals = [-8.395792, -2.198370, -0.000009, 0.000000] + edge_VW.test_vals = [-8.901513, -2.707051, -0.000009, 0.000000] test_list.append(edge_VW) # Rarefaction shock wave edge_PPR @@ -936,7 +936,7 @@ def main(): edge_PPR.cfg_dir = "nicf/edge" edge_PPR.cfg_file = "edge_PPR.cfg" edge_PPR.test_iter = 20 - edge_PPR.test_vals = [-12.584167, -6.372084, -0.000034, 0.000000] + edge_PPR.test_vals = [-12.051551, -5.877885, -0.000034, 0.000000] test_list.append(edge_PPR) @@ -966,7 +966,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.167512, 1.598496, -2.928579, 2.573642, -2.527392, 3.016170, 106370.000000, 106370.000000, 5.726800, 64.383000] + axial_stage2D.test_vals = [1.167521, 1.598536, -2.928575, 2.573646, -2.527330, 3.016253, 106370.000000, 106370.000000, 5.726800, 64.383000] test_list.append(axial_stage2D) # 2D transonic stator restart @@ -1025,7 +1025,7 @@ def main(): channel_2D.cfg_dir = "sliding_interface/channel_2D" channel_2D.cfg_file = "channel_2D_WA.cfg" channel_2D.test_iter = 2 - channel_2D.test_vals = [2.000000, 0.000000, 0.466202, 0.350093, 0.399014] + channel_2D.test_vals = [2.000000, 0.000000, 0.466389, 0.346857, 0.397472] channel_2D.timeout = 100 channel_2D.unsteady = True channel_2D.multizone = True @@ -1036,7 +1036,7 @@ def main(): channel_3D.cfg_dir = "sliding_interface/channel_3D" channel_3D.cfg_file = "channel_3D_WA.cfg" channel_3D.test_iter = 1 - channel_3D.test_vals = [1.000000, 0.000000, 0.617237, 0.798808, 0.692608] + channel_3D.test_vals = [1.000000, 0.000000, 0.620573, 0.810294, 0.718984] channel_3D.test_vals_aarch64 = [1.000000, 0.000000, 0.611996, 0.798988, 0.702357] channel_3D.unsteady = True channel_3D.multizone = True @@ -1048,7 +1048,7 @@ def main(): pipe.cfg_dir = "sliding_interface/pipe" pipe.cfg_file = "pipe_NN.cfg" pipe.test_iter = 2 - pipe.test_vals = [0.568972, 0.692858, 0.989456, 1.048237] + pipe.test_vals = [0.528897, 0.684175, 0.970353, 1.065262] pipe.unsteady = True pipe.multizone = True test_list.append(pipe) @@ -1058,7 +1058,7 @@ def main(): rotating_cylinders.cfg_dir = "sliding_interface/rotating_cylinders" rotating_cylinders.cfg_file = "rot_cylinders_WA.cfg" rotating_cylinders.test_iter = 3 - rotating_cylinders.test_vals = [3.000000, 0.000000, 0.664822, 1.125812, 1.117605] + rotating_cylinders.test_vals = [3.000000, 0.000000, 0.686432, 1.159461, 1.151134] rotating_cylinders.unsteady = True rotating_cylinders.multizone = True test_list.append(rotating_cylinders) @@ -1068,7 +1068,7 @@ def main(): supersonic_vortex_shedding.cfg_dir = "sliding_interface/supersonic_vortex_shedding" supersonic_vortex_shedding.cfg_file = "sup_vor_shed_WA.cfg" supersonic_vortex_shedding.test_iter = 5 - supersonic_vortex_shedding.test_vals = [5.000000, 0.000000, 0.899641, 1.076189] + supersonic_vortex_shedding.test_vals = [5.000000, 0.000000, 1.126036, 1.081775] supersonic_vortex_shedding.unsteady = True supersonic_vortex_shedding.multizone = True test_list.append(supersonic_vortex_shedding) @@ -1632,7 +1632,7 @@ def main(): pywrapper_naca0012.cfg_dir = "euler/naca0012" pywrapper_naca0012.cfg_file = "inv_NACA0012_Roe.cfg" pywrapper_naca0012.test_iter = 20 - pywrapper_naca0012.test_vals = [-4.032065, -3.525196, 0.273051, 0.009680] + pywrapper_naca0012.test_vals = [-4.031446, -3.524530, 0.273056, 0.009703] pywrapper_naca0012.command = TestCase.Command(exec = "SU2_CFD.py", param = "-f") pywrapper_naca0012.timeout = 1600 pywrapper_naca0012.tol = 0.00001 diff --git a/TestCases/serial_regression_AD.py b/TestCases/serial_regression_AD.py index 8a0277294406..288f8a47efd3 100644 --- a/TestCases/serial_regression_AD.py +++ b/TestCases/serial_regression_AD.py @@ -188,7 +188,7 @@ def main(): discadj_sliding_interface.cfg_dir = "disc_adj_ffi/sliding_interface" discadj_sliding_interface.cfg_file = "circles.cfg" discadj_sliding_interface.test_iter = 10 - discadj_sliding_interface.test_vals = [10.000000, -4.395143, -4.419431] + discadj_sliding_interface.test_vals = [10.000000, -4.403519, -4.419431] test_list.append(discadj_sliding_interface) ####################################################### diff --git a/TestCases/tutorials.py b/TestCases/tutorials.py index d59aad17bf96..f2e16e6a3dbb 100644 --- a/TestCases/tutorials.py +++ b/TestCases/tutorials.py @@ -187,7 +187,7 @@ def main(): tutorial_inv_bump.cfg_dir = "../Tutorials/compressible_flow/Inviscid_Bump" tutorial_inv_bump.cfg_file = "inv_channel.cfg" tutorial_inv_bump.test_iter = 0 - tutorial_inv_bump.test_vals = [-1.437425, 4.075857, 0.086038, 0.127546] + tutorial_inv_bump.test_vals = [-1.437425, 4.075857, 0.085674, 0.123798] test_list.append(tutorial_inv_bump) # Inviscid Wedge @@ -195,7 +195,7 @@ def main(): tutorial_inv_wedge.cfg_dir = "../Tutorials/compressible_flow/Inviscid_Wedge" tutorial_inv_wedge.cfg_file = "inv_wedge_HLLC.cfg" tutorial_inv_wedge.test_iter = 0 - tutorial_inv_wedge.test_vals = [-0.481460, 5.253008, -0.215274, 0.037861] + tutorial_inv_wedge.test_vals = [-0.481460, 5.253008, -0.191841, 0.033699] tutorial_inv_wedge.no_restart = True test_list.append(tutorial_inv_wedge) @@ -204,7 +204,7 @@ def main(): tutorial_inv_onera.cfg_dir = "../Tutorials/compressible_flow/Inviscid_ONERAM6" tutorial_inv_onera.cfg_file = "inv_ONERAM6.cfg" tutorial_inv_onera.test_iter = 0 - tutorial_inv_onera.test_vals = [-5.204928, -4.597762, 0.256743, 0.113097] + tutorial_inv_onera.test_vals = [-5.204928, -4.597762, 0.243530, 0.111160] tutorial_inv_onera.no_restart = True test_list.append(tutorial_inv_onera) @@ -213,7 +213,7 @@ def main(): tutorial_lam_cylinder.cfg_dir = "../Tutorials/compressible_flow/Laminar_Cylinder" tutorial_lam_cylinder.cfg_file = "lam_cylinder.cfg" tutorial_lam_cylinder.test_iter = 0 - tutorial_lam_cylinder.test_vals = [-6.162141, -0.699617, -0.080518, 43.028531] + tutorial_lam_cylinder.test_vals = [-6.162141, -0.699617, -0.080602, 42.990503] tutorial_lam_cylinder.no_restart = True test_list.append(tutorial_lam_cylinder) @@ -287,7 +287,7 @@ def main(): tutorial_turb_oneram6.cfg_dir = "../Tutorials/compressible_flow/Turbulent_ONERAM6" tutorial_turb_oneram6.cfg_file = "turb_ONERAM6.cfg" tutorial_turb_oneram6.test_iter = 0 - tutorial_turb_oneram6.test_vals = [-4.564441, -11.537596, 0.293674, 0.087628] + tutorial_turb_oneram6.test_vals = [-4.564441, -11.534937, 0.250153, 0.080951] test_list.append(tutorial_turb_oneram6) # NICFD Nozzle @@ -295,7 +295,7 @@ def main(): tutorial_nicfd_nozzle.cfg_dir = "../Tutorials/compressible_flow/NICFD_nozzle" tutorial_nicfd_nozzle.cfg_file = "NICFD_nozzle.cfg" tutorial_nicfd_nozzle.test_iter = 20 - tutorial_nicfd_nozzle.test_vals = [-1.800072, -6.002691, 4.635204, 0.000000, 0.000000] + tutorial_nicfd_nozzle.test_vals = [-1.799851, -6.002738, 4.635397, 0.000000, 0.000000] tutorial_nicfd_nozzle.no_restart = True test_list.append(tutorial_nicfd_nozzle) @@ -314,7 +314,7 @@ def main(): tutorial_unst_naca0012.cfg_dir = "../Tutorials/compressible_flow/Unsteady_NACA0012" tutorial_unst_naca0012.cfg_file = "unsteady_naca0012.cfg" tutorial_unst_naca0012.test_iter = 520 - tutorial_unst_naca0012.test_vals = [520.000000, 0.000000, -5.292139, 0.000000, 0.306930, 0.787006, 0.001794, 0.011018] + tutorial_unst_naca0012.test_vals = [520.000000, 0.000000, -5.292186, 0.000000, 0.304255, 0.792893, 0.001907, 0.010888] tutorial_unst_naca0012.test_vals_aarch64 = [520, 0, -5.292359, 0, 0.284720, 0.766329, 0.000954, 0.007565] tutorial_unst_naca0012.unsteady = True test_list.append(tutorial_unst_naca0012) @@ -324,7 +324,7 @@ def main(): propeller_var_load.cfg_dir = "../Tutorials/compressible_flow/ActuatorDisk_VariableLoad" propeller_var_load.cfg_file = "propeller_variable_load.cfg" propeller_var_load.test_iter = 20 - propeller_var_load.test_vals = [-1.831126, -4.534989, -0.000323, 0.171579] + propeller_var_load.test_vals = [-1.831998, -4.534991, -0.000323, 0.171494] propeller_var_load.timeout = 3200 test_list.append(propeller_var_load) @@ -335,7 +335,7 @@ def main(): tutorial_design_inv_naca0012.cfg_dir = "../Tutorials/design/Inviscid_2D_Unconstrained_NACA0012" tutorial_design_inv_naca0012.cfg_file = "inv_NACA0012_basic.cfg" tutorial_design_inv_naca0012.test_iter = 0 - tutorial_design_inv_naca0012.test_vals = [-3.585391, -2.989014, 0.146155, 0.179374] + tutorial_design_inv_naca0012.test_vals = [-3.585391, -2.989014, 0.141838, 0.172252] tutorial_design_inv_naca0012.no_restart = True test_list.append(tutorial_design_inv_naca0012) @@ -344,7 +344,7 @@ def main(): tutorial_design_turb_rae2822.cfg_dir = "../Tutorials/design/Turbulent_2D_Constrained_RAE2822" tutorial_design_turb_rae2822.cfg_file = "turb_SA_RAE2822.cfg" tutorial_design_turb_rae2822.test_iter = 0 - tutorial_design_turb_rae2822.test_vals = [-1.700114, -4.941829, 0.205309, 0.184990] + tutorial_design_turb_rae2822.test_vals = [-1.700114, -4.941830, 0.161066, 0.170510] tutorial_design_turb_rae2822.no_restart = True test_list.append(tutorial_design_turb_rae2822) diff --git a/TestCases/vandv.py b/TestCases/vandv.py index 516f51cc5ebd..aa9fb4398d52 100644 --- a/TestCases/vandv.py +++ b/TestCases/vandv.py @@ -81,7 +81,7 @@ def main(): swbli_sa.cfg_dir = "vandv/rans/swbli" swbli_sa.cfg_file = "config_sa.cfg" swbli_sa.test_iter = 5 - swbli_sa.test_vals = [-11.502718, -10.939184, -12.034284, -10.581169, -16.088844, 0.002242, -1.664946, 1.258100] + swbli_sa.test_vals = [-11.504424, -10.941741, -12.049925, -10.586263, -16.090385, 0.002242, -1.614365, 1.340100] swbli_sa.test_vals_aarch64 = [-11.504424, -10.941741, -12.049925, -10.586263, -16.090385, 0.002242, -1.614365, 1.340100] test_list.append(swbli_sa) @@ -91,7 +91,7 @@ def main(): swbli_sst.cfg_dir = "vandv/rans/swbli" swbli_sst.cfg_file = "config_sst.cfg" swbli_sst.test_iter = 5 - swbli_sst.test_vals = [-11.319578, -10.641523, -11.224600, -10.150215, -11.407538, -2.637660, 0.001816, -1.839484, -3.514593, 11.136000] + swbli_sst.test_vals = [-11.569218, -10.909086, -11.607984, -10.431163, -11.407588, -2.637660, 0.001816, -1.305818, -3.514590, 13.399000] test_list.append(swbli_sst) # DSMA661 - SA From 09415a28dccc616c0e1c662b0c1518a96169d643 Mon Sep 17 00:00:00 2001 From: rois1995 Date: Tue, 29 Sep 2026 22:42:31 +0200 Subject: [PATCH 14/14] Update the reference files of the python-wrapper and gradient regression tests from CI Co-Authored-By: Claude Opus 5.5 (1M context) Claude-Session: https://claude.ai/code/session_01LEL91DW5WPbPwgFtCvHga6 --- .../naca0012/of_grad_directdiff.dat.ref | 6 +- .../multiple_ffd/naca0012/of_grad_cd.dat.ref | 4 +- .../naca0012/of_grad_directdiff.dat.ref | 4 +- .../translating_NACA0012/forces_0.csv.ref | 400 +++++++++--------- .../forces_0.csv.ref | 400 +++++++++--------- 5 files changed, 407 insertions(+), 407 deletions(-) diff --git a/TestCases/cont_adj_euler/naca0012/of_grad_directdiff.dat.ref b/TestCases/cont_adj_euler/naca0012/of_grad_directdiff.dat.ref index c87c1c9a2a49..2e6cd4903a0d 100644 --- a/TestCases/cont_adj_euler/naca0012/of_grad_directdiff.dat.ref +++ b/TestCases/cont_adj_euler/naca0012/of_grad_directdiff.dat.ref @@ -1,4 +1,4 @@ VARIABLES="VARIABLE" , "DRAG" , "EFFICIENCY" , "FORCE_X" , "FORCE_Y" , "FORCE_Z" , "LIFT" , "MOMENT_X" , "MOMENT_Y" , "MOMENT_Z" , "SIDEFORCE" - 0 , -0.4704786435 , -33.40446141 , -0.3888715474 , -3.745132593 , 0.0 , -3.735758164 , 0.0 , 0.0 , 0.4834798526 , 0.0 - 1 , 0.412031843 , -39.78218403 , 0.4880368175 , -3.47876392 , 0.0 , -3.488582535 , 0.0 , 0.0 , 0.2021305183 , 0.0 - 2 , 0.8908446367 , -37.475187 , 0.9514902759 , -2.769632304 , 0.0 , -2.789729858 , 0.0 , 0.0 , -0.1182204408 , 0.0 + 0 , -0.4834688689 , -33.38043366 , -0.4014654728 , -3.763436475 , 0.0 , -3.753782956 , 0.0 , 0.0 , 0.4865376034 , 0.0 + 1 , 0.4087842111 , -39.69099805 , 0.4848740806 , -3.482690031 , 0.0 , -3.492438718 , 0.0 , 0.0 , 0.2039858924 , 0.0 + 2 , 0.8921290437 , -37.32658951 , 0.9526325615 , -2.763104959 , 0.0 , -2.783228985 , 0.0 , 0.0 , -0.117214663 , 0.0 diff --git a/TestCases/multiple_ffd/naca0012/of_grad_cd.dat.ref b/TestCases/multiple_ffd/naca0012/of_grad_cd.dat.ref index 2f2e345e4221..e98387350993 100644 --- a/TestCases/multiple_ffd/naca0012/of_grad_cd.dat.ref +++ b/TestCases/multiple_ffd/naca0012/of_grad_cd.dat.ref @@ -1,3 +1,3 @@ VARIABLES="VARIABLE" , "GRADIENT" , "FINDIFF_STEP" - 0 , 0.0440381 , 0.001 - 1 , -0.0848429 , 0.001 + 0 , 0.0439624 , 0.001 + 1 , -0.0848386 , 0.001 diff --git a/TestCases/multiple_ffd/naca0012/of_grad_directdiff.dat.ref b/TestCases/multiple_ffd/naca0012/of_grad_directdiff.dat.ref index ef3b2d9f463a..dfad81314570 100644 --- a/TestCases/multiple_ffd/naca0012/of_grad_directdiff.dat.ref +++ b/TestCases/multiple_ffd/naca0012/of_grad_directdiff.dat.ref @@ -1,3 +1,3 @@ VARIABLES="VARIABLE" , "DRAG" , "EFFICIENCY" , "FORCE_X" , "FORCE_Y" , "FORCE_Z" , "LIFT" , "MOMENT_X" , "MOMENT_Y" , "MOMENT_Z" , "SIDEFORCE" - 0 , 0.04790586765 , -0.4797536912 , 0.04825151062 , -0.01531800012 , 0.0 , -0.01636695601 , 0.0 , 0.0 , 0.03036927931 , 0.0 - 1 , -0.09320474204 , 0.3762440542 , -0.0926027588 , -0.02860524743 , 0.0 , -0.02657832162 , 0.0 , 0.0 , 0.06626460628 , 0.0 + 0 , 0.04764411115 , -0.463384101 , 0.04796070834 , -0.01398970686 , 0.0 , -0.01503263502 , 0.0 , 0.0 , 0.03019600866 , 0.0 + 1 , -0.09347608871 , 0.3646170468 , -0.0928475452 , -0.02982544755 , 0.0 , -0.02779289138 , 0.0 , 0.0 , 0.06643204278 , 0.0 diff --git a/TestCases/py_wrapper/translating_NACA0012/forces_0.csv.ref b/TestCases/py_wrapper/translating_NACA0012/forces_0.csv.ref index 3aa219b7120e..0c9091d233ae 100644 --- a/TestCases/py_wrapper/translating_NACA0012/forces_0.csv.ref +++ b/TestCases/py_wrapper/translating_NACA0012/forces_0.csv.ref @@ -1,200 +1,200 @@ -199, -0.99, -0.00, 0.00 -0, -2.91, 20.02, 0.00 -1, -4.27, 29.42, 0.00 -2, -5.43, 37.48, 0.00 -3, -6.49, 44.88, 0.00 -4, -7.43, 51.52, 0.00 -5, -8.23, 57.28, 0.00 -6, -8.93, 62.42, 0.00 -7, -9.50, 66.76, 0.00 -8, -9.96, 70.36, 0.00 -9, -10.31, 73.26, 0.00 -10, -10.55, 75.47, 0.00 -11, -10.71, 77.10, 0.00 -12, -10.75, 78.01, 0.00 -13, -10.71, 78.35, 0.00 -14, -10.57, 78.03, 0.00 -15, -10.36, 77.16, 0.00 -16, -10.07, 75.79, 0.00 -17, -9.73, 73.92, 0.00 -18, -9.31, 71.49, 0.00 -19, -8.85, 68.73, 0.00 -20, -8.32, 65.37, 0.00 -21, -7.77, 61.78, 0.00 -22, -7.15, 57.64, 0.00 -23, -6.53, 53.27, 0.00 -24, -5.85, 48.43, 0.00 -25, -5.18, 43.45, 0.00 -26, -4.46, 37.94, 0.00 -27, -3.74, 32.28, 0.00 -28, -2.99, 26.23, 0.00 -29, -2.26, 20.13, 0.00 -30, -1.50, 13.62, 0.00 -31, -0.77, 7.14, 0.00 -32, -0.03, 0.27, 0.00 -33, 0.64, -6.16, 0.00 -34, 1.33, -12.94, 0.00 -35, 1.88, -18.69, 0.00 -36, 2.50, -25.43, 0.00 -37, 2.94, -30.60, 0.00 -38, 3.61, -38.39, 0.00 -39, 4.12, -44.94, 0.00 -40, 5.20, -58.26, 0.00 -41, 6.45, -74.32, 0.00 -42, 8.48, -100.75, 0.00 -43, 10.40, -127.62, 0.00 -44, 12.32, -156.48, 0.00 -45, 14.32, -188.86, 0.00 -46, 15.99, -219.70, 0.00 -47, 17.77, -255.19, 0.00 -48, 19.05, -287.25, 0.00 -49, 20.43, -325.17, 0.00 -50, 21.21, -358.33, 0.00 -51, 21.60, -390.04, 0.00 -52, 20.80, -404.89, 0.00 -53, 19.46, -412.52, 0.00 -54, 17.55, -410.28, 0.00 -55, 15.38, -403.24, 0.00 -56, 12.80, -384.55, 0.00 -57, 10.01, -354.85, 0.00 -58, 7.33, -320.56, 0.00 -59, 5.62, -325.28, 0.00 -60, 4.54, -398.71, 0.00 -61, 2.33, -444.80, 0.00 -62, -0.55, -444.61, 0.00 -63, -3.45, -430.49, 0.00 -64, -6.33, -418.34, 0.00 -65, -9.16, -405.22, 0.00 -66, -11.93, -391.89, 0.00 -67, -14.56, -376.47, 0.00 -68, -17.08, -360.79, 0.00 -69, -19.39, -343.61, 0.00 -70, -21.54, -326.20, 0.00 -71, -23.38, -307.04, 0.00 -72, -25.04, -288.42, 0.00 -73, -26.35, -268.52, 0.00 -74, -27.44, -249.13, 0.00 -75, -28.26, -229.94, 0.00 -76, -28.70, -210.22, 0.00 -77, -28.76, -190.36, 0.00 -78, -28.57, -171.29, 0.00 -79, -27.84, -151.55, 0.00 -80, -26.53, -131.24, 0.00 -81, -24.47, -110.11, 0.00 -82, -22.10, -90.41, 0.00 -83, -17.63, -65.56, 0.00 -84, -13.82, -46.59, 0.00 -85, -9.00, -27.47, 0.00 -86, -3.97, -10.94, 0.00 -87, 2.04, 5.04, 0.00 -88, 8.66, 19.11, 0.00 -89, 17.22, 33.67, 0.00 -90, 26.22, 44.98, 0.00 -91, 38.13, 56.68, 0.00 -92, 50.18, 63.55, 0.00 -93, 64.64, 68.16, 0.00 -94, 81.14, 68.93, 0.00 -95, 92.61, 60.32, 0.00 -96, 98.61, 45.20, 0.00 -97, 91.85, 24.86, 0.00 -98, 60.95, 8.18, 0.00 -99, 54.21, -0.00, 0.00 -100, 87.82, -11.78, 0.00 -101, 110.00, -29.77, 0.00 -102, 92.95, -42.60, 0.00 -103, 73.45, -47.84, 0.00 -104, 54.95, -46.68, 0.00 -105, 39.22, -41.36, 0.00 -106, 29.77, -37.70, 0.00 -107, 13.13, -19.52, 0.00 -108, 1.46, -2.51, 0.00 -109, -12.03, 23.53, 0.00 -110, -21.14, 46.66, 0.00 -111, -30.42, 75.21, 0.00 -112, -36.34, 100.05, 0.00 -113, -42.50, 129.72, 0.00 -114, -45.24, 152.58, 0.00 -115, -47.68, 177.27, 0.00 -116, -48.11, 196.85, 0.00 -117, -49.30, 221.84, 0.00 -118, -49.86, 246.67, 0.00 -119, -50.26, 273.55, 0.00 -120, -49.99, 299.73, 0.00 -121, -49.24, 325.87, 0.00 -122, -47.80, 350.18, 0.00 -123, -46.06, 374.79, 0.00 -124, -43.92, 398.75, 0.00 -125, -41.45, 422.34, 0.00 -126, -38.71, 445.76, 0.00 -127, -35.65, 468.20, 0.00 -128, -32.39, 490.54, 0.00 -129, -28.90, 511.99, 0.00 -130, -25.24, 533.17, 0.00 -131, -21.41, 553.63, 0.00 -132, -17.46, 573.42, 0.00 -133, -13.40, 592.80, 0.00 -134, -9.25, 611.37, 0.00 -135, -5.04, 629.20, 0.00 -136, -0.80, 646.64, 0.00 -137, 3.47, 662.94, 0.00 -138, 7.73, 678.10, 0.00 -139, 11.96, 692.64, 0.00 -140, 16.15, 706.26, 0.00 -141, 20.27, 718.38, 0.00 -142, 24.31, 730.12, 0.00 -143, 28.27, 741.00, 0.00 -144, 32.10, 750.49, 0.00 -145, 35.79, 758.77, 0.00 -146, 39.36, 766.12, 0.00 -147, 42.77, 772.47, 0.00 -148, 46.02, 777.49, 0.00 -149, 49.12, 781.76, 0.00 -150, 52.03, 784.73, 0.00 -151, 54.76, 786.49, 0.00 -152, 57.27, 786.70, 0.00 -153, 59.65, 786.61, 0.00 -154, 62.04, 787.97, 0.00 -155, 63.93, 784.13, 0.00 -156, 65.88, 782.40, 0.00 -157, 67.06, 772.81, 0.00 -158, 69.82, 782.23, 0.00 -159, 35.12, 383.21, 0.00 -160, -13.10, -139.35, 0.00 -161, -9.96, -103.51, 0.00 -162, -7.72, -78.42, 0.00 -163, -7.29, -72.50, 0.00 -164, -7.05, -68.73, 0.00 -165, -6.60, -63.13, 0.00 -166, -6.57, -61.71, 0.00 -167, -6.59, -60.74, 0.00 -168, -6.81, -61.72, 0.00 -169, -7.08, -63.08, 0.00 -170, -7.42, -65.08, 0.00 -171, -7.76, -67.01, 0.00 -172, -8.15, -69.40, 0.00 -173, -8.55, -71.76, 0.00 -174, -8.94, -73.94, 0.00 -175, -9.30, -75.91, 0.00 -176, -9.67, -77.93, 0.00 -177, -10.01, -79.65, 0.00 -178, -10.32, -81.15, 0.00 -179, -10.61, -82.44, 0.00 -180, -10.85, -83.37, 0.00 -181, -11.06, -84.05, 0.00 -182, -11.21, -84.33, 0.00 -183, -11.31, -84.26, 0.00 -184, -11.35, -83.81, 0.00 -185, -11.34, -82.98, 0.00 -186, -11.25, -81.65, 0.00 -187, -11.11, -80.01, 0.00 -188, -10.89, -77.90, 0.00 -189, -10.54, -74.91, 0.00 -190, -10.09, -71.22, 0.00 -191, -9.54, -66.99, 0.00 -192, -8.89, -62.12, 0.00 -193, -8.13, -56.59, 0.00 -194, -7.28, -50.53, 0.00 -195, -6.31, -43.67, 0.00 -196, -5.23, -36.08, 0.00 -197, -4.05, -27.89, 0.00 -198, -2.67, -18.38, 0.00 +199, -0.94, -0.00, 0.00 +0, -2.79, 19.21, 0.00 +1, -4.07, 28.05, 0.00 +2, -5.12, 35.35, 0.00 +3, -6.07, 41.99, 0.00 +4, -6.89, 47.80, 0.00 +5, -7.57, 52.73, 0.00 +6, -8.16, 57.03, 0.00 +7, -8.63, 60.64, 0.00 +8, -8.99, 63.49, 0.00 +9, -9.25, 65.69, 0.00 +10, -9.40, 67.22, 0.00 +11, -9.47, 68.21, 0.00 +12, -9.44, 68.52, 0.00 +13, -9.32, 68.22, 0.00 +14, -9.13, 67.36, 0.00 +15, -8.84, 65.83, 0.00 +16, -8.50, 63.93, 0.00 +17, -8.06, 61.28, 0.00 +18, -7.60, 58.37, 0.00 +19, -7.04, 54.67, 0.00 +20, -6.48, 50.91, 0.00 +21, -5.83, 46.37, 0.00 +22, -5.21, 41.95, 0.00 +23, -4.50, 36.73, 0.00 +24, -3.85, 31.89, 0.00 +25, -3.13, 26.30, 0.00 +26, -2.48, 21.12, 0.00 +27, -1.81, 15.63, 0.00 +28, -1.16, 10.19, 0.00 +29, -0.56, 5.01, 0.00 +30, -0.01, 0.07, 0.00 +31, 0.30, -2.76, 0.00 +32, 0.83, -7.82, 0.00 +33, 0.85, -8.14, 0.00 +34, 1.47, -14.33, 0.00 +35, 1.50, -14.91, 0.00 +36, 2.34, -23.79, 0.00 +37, 2.58, -26.78, 0.00 +38, 3.60, -38.33, 0.00 +39, 4.15, -45.32, 0.00 +40, 5.35, -59.93, 0.00 +41, 6.34, -73.02, 0.00 +42, 7.53, -89.47, 0.00 +43, 8.62, -105.75, 0.00 +44, 9.52, -120.93, 0.00 +45, 10.46, -137.99, 0.00 +46, 10.89, -149.59, 0.00 +47, 11.06, -158.86, 0.00 +48, 10.91, -164.51, 0.00 +49, 10.20, -162.29, 0.00 +50, 10.41, -175.95, 0.00 +51, 10.69, -193.01, 0.00 +52, 11.35, -220.88, 0.00 +53, 12.09, -256.22, 0.00 +54, 11.70, -273.60, 0.00 +55, 11.18, -293.17, 0.00 +56, 9.73, -292.14, 0.00 +57, 8.95, -317.19, 0.00 +58, 8.19, -357.93, 0.00 +59, 7.34, -424.98, 0.00 +60, 5.17, -454.05, 0.00 +61, 2.26, -432.87, 0.00 +62, -0.52, -418.79, 0.00 +63, -3.27, -407.56, 0.00 +64, -5.99, -395.49, 0.00 +65, -8.62, -381.54, 0.00 +66, -11.16, -366.54, 0.00 +67, -13.55, -350.35, 0.00 +68, -15.82, -334.23, 0.00 +69, -17.87, -316.57, 0.00 +70, -19.66, -297.72, 0.00 +71, -21.14, -277.56, 0.00 +72, -22.40, -258.00, 0.00 +73, -23.24, -236.77, 0.00 +74, -23.93, -217.26, 0.00 +75, -23.95, -194.85, 0.00 +76, -23.92, -175.24, 0.00 +77, -23.06, -152.62, 0.00 +78, -22.31, -133.79, 0.00 +79, -20.38, -110.93, 0.00 +80, -18.51, -91.56, 0.00 +81, -14.90, -67.05, 0.00 +82, -11.41, -46.70, 0.00 +83, -4.97, -18.48, 0.00 +84, -0.30, -1.01, 0.00 +85, 9.93, 30.30, 0.00 +86, 16.77, 46.17, 0.00 +87, 32.33, 79.93, 0.00 +88, 43.53, 96.08, 0.00 +89, 60.42, 118.16, 0.00 +90, 71.17, 122.11, 0.00 +91, 98.75, 146.81, 0.00 +92, 112.94, 143.03, 0.00 +93, 162.93, 171.80, 0.00 +94, 168.48, 143.12, 0.00 +95, 162.93, 106.12, 0.00 +96, 161.49, 74.02, 0.00 +97, 226.14, 61.20, 0.00 +98, 173.55, 23.29, 0.00 +99, 81.85, -0.00, 0.00 +100, 101.93, -13.68, 0.00 +101, 111.85, -30.27, 0.00 +102, 83.84, -38.43, 0.00 +103, 59.18, -38.55, 0.00 +104, 39.46, -33.52, 0.00 +105, 20.52, -21.64, 0.00 +106, 9.26, -11.73, 0.00 +107, -5.26, 7.82, 0.00 +108, -14.72, 25.25, 0.00 +109, -24.47, 47.86, 0.00 +110, -30.46, 67.22, 0.00 +111, -36.98, 91.44, 0.00 +112, -40.61, 111.80, 0.00 +113, -44.24, 135.04, 0.00 +114, -47.23, 159.29, 0.00 +115, -50.25, 186.80, 0.00 +116, -52.01, 212.78, 0.00 +117, -49.79, 224.02, 0.00 +118, -49.93, 247.03, 0.00 +119, -50.10, 272.69, 0.00 +120, -49.50, 296.79, 0.00 +121, -48.43, 320.54, 0.00 +122, -46.68, 341.98, 0.00 +123, -45.16, 367.44, 0.00 +124, -43.12, 391.55, 0.00 +125, -41.04, 418.22, 0.00 +126, -38.38, 441.96, 0.00 +127, -35.65, 468.10, 0.00 +128, -32.39, 490.58, 0.00 +129, -29.12, 515.86, 0.00 +130, -25.43, 537.19, 0.00 +131, -21.72, 561.44, 0.00 +132, -17.69, 581.14, 0.00 +133, -13.62, 602.77, 0.00 +134, -9.39, 620.72, 0.00 +135, -5.13, 639.76, 0.00 +136, -0.81, 656.52, 0.00 +137, 3.52, 673.60, 0.00 +138, 7.84, 688.01, 0.00 +139, 12.15, 703.61, 0.00 +140, 16.38, 716.38, 0.00 +141, 20.56, 728.93, 0.00 +142, 24.62, 739.26, 0.00 +143, 28.61, 749.93, 0.00 +144, 32.43, 758.20, 0.00 +145, 36.13, 765.81, 0.00 +146, 39.66, 772.00, 0.00 +147, 43.08, 777.93, 0.00 +148, 46.29, 782.12, 0.00 +149, 49.44, 786.83, 0.00 +150, 52.29, 788.58, 0.00 +151, 54.86, 787.96, 0.00 +152, 57.31, 787.22, 0.00 +153, 59.58, 785.72, 0.00 +154, 62.23, 790.29, 0.00 +155, 62.41, 765.50, 0.00 +156, 64.45, 765.33, 0.00 +157, 54.69, 630.25, 0.00 +158, 48.23, 540.32, 0.00 +159, -12.34, -134.66, 0.00 +160, -29.22, -310.97, 0.00 +161, -7.18, -74.61, 0.00 +162, -3.54, -35.94, 0.00 +163, -3.39, -33.73, 0.00 +164, -3.11, -30.30, 0.00 +165, -2.96, -28.33, 0.00 +166, -2.85, -26.76, 0.00 +167, -2.79, -25.71, 0.00 +168, -2.88, -26.13, 0.00 +169, -3.27, -29.13, 0.00 +170, -3.57, -31.34, 0.00 +171, -3.86, -33.32, 0.00 +172, -4.21, -35.80, 0.00 +173, -4.70, -39.42, 0.00 +174, -5.09, -42.12, 0.00 +175, -5.57, -45.47, 0.00 +176, -6.00, -48.33, 0.00 +177, -6.47, -51.51, 0.00 +178, -6.92, -54.40, 0.00 +179, -7.37, -57.26, 0.00 +180, -7.77, -59.72, 0.00 +181, -8.16, -62.01, 0.00 +182, -8.49, -63.90, 0.00 +183, -8.78, -65.39, 0.00 +184, -9.01, -66.52, 0.00 +185, -9.18, -67.19, 0.00 +186, -9.30, -67.46, 0.00 +187, -9.30, -67.00, 0.00 +188, -9.21, -65.84, 0.00 +189, -9.06, -64.37, 0.00 +190, -8.81, -62.23, 0.00 +191, -8.46, -59.42, 0.00 +192, -8.00, -55.90, 0.00 +193, -7.42, -51.65, 0.00 +194, -6.74, -46.78, 0.00 +195, -5.94, -41.07, 0.00 +196, -5.00, -34.54, 0.00 +197, -3.98, -27.42, 0.00 +198, -2.72, -18.72, 0.00 diff --git a/TestCases/py_wrapper/updated_moving_frame_NACA12/forces_0.csv.ref b/TestCases/py_wrapper/updated_moving_frame_NACA12/forces_0.csv.ref index 1c87fcffcb67..5168d86660a2 100644 --- a/TestCases/py_wrapper/updated_moving_frame_NACA12/forces_0.csv.ref +++ b/TestCases/py_wrapper/updated_moving_frame_NACA12/forces_0.csv.ref @@ -1,200 +1,200 @@ -199, -0.98, -0.00, 0.00 -0, -2.89, 19.89, 0.00 -1, -4.12, 28.36, 0.00 -2, -5.17, 35.67, 0.00 -3, -6.10, 42.20, 0.00 -4, -6.93, 48.05, 0.00 -5, -7.63, 53.11, 0.00 -6, -8.21, 57.37, 0.00 -7, -8.68, 60.97, 0.00 -8, -8.99, 63.51, 0.00 -9, -9.20, 65.36, 0.00 -10, -9.28, 66.37, 0.00 -11, -9.27, 66.73, 0.00 -12, -9.17, 66.55, 0.00 -13, -8.99, 65.76, 0.00 -14, -8.76, 64.65, 0.00 -15, -8.46, 63.02, 0.00 -16, -8.12, 61.09, 0.00 -17, -7.73, 58.74, 0.00 -18, -7.26, 55.80, 0.00 -19, -6.76, 52.48, 0.00 -20, -6.17, 48.47, 0.00 -21, -5.53, 44.02, 0.00 -22, -4.86, 39.16, 0.00 -23, -4.13, 33.74, 0.00 -24, -3.43, 28.38, 0.00 -25, -2.68, 22.45, 0.00 -26, -1.99, 16.92, 0.00 -27, -1.26, 10.86, 0.00 -28, -0.62, 5.40, 0.00 -29, 0.05, -0.48, 0.00 -30, 0.57, -5.20, 0.00 -31, 1.09, -10.01, 0.00 -32, 1.48, -13.85, 0.00 -33, 1.88, -17.94, 0.00 -34, 2.18, -21.27, 0.00 -35, 2.56, -25.45, 0.00 -36, 2.72, -27.60, 0.00 -37, 2.84, -29.52, 0.00 -38, 2.54, -27.06, 0.00 -39, 2.39, -26.05, 0.00 -40, 2.12, -23.72, 0.00 -41, 2.85, -32.88, 0.00 -42, 0.07, -0.80, 0.00 -43, -7.64, 93.76, 0.00 -44, -15.56, 197.66, 0.00 -45, -18.11, 238.88, 0.00 -46, 4.57, -62.81, 0.00 -47, 28.46, -408.83, 0.00 -48, 30.15, -454.80, 0.00 -49, 29.55, -470.30, 0.00 -50, 27.36, -462.29, 0.00 -51, 25.32, -457.28, 0.00 -52, 23.08, -449.36, 0.00 -53, 21.91, -464.51, 0.00 -54, 21.44, -501.19, 0.00 -55, 21.98, -576.18, 0.00 -56, 20.44, -613.98, 0.00 -57, 16.27, -576.74, 0.00 -58, 12.57, -549.58, 0.00 -59, 9.29, -537.86, 0.00 -60, 6.03, -529.16, 0.00 -61, 2.72, -519.13, 0.00 -62, -0.62, -505.94, 0.00 -63, -3.93, -490.73, 0.00 -64, -7.20, -475.71, 0.00 -65, -10.38, -459.45, 0.00 -66, -13.47, -442.43, 0.00 -67, -16.41, -424.36, 0.00 -68, -19.21, -405.79, 0.00 -69, -21.81, -386.39, 0.00 -70, -24.21, -366.71, 0.00 -71, -26.37, -346.30, 0.00 -72, -28.28, -325.64, 0.00 -73, -29.91, -304.75, 0.00 -74, -31.21, -283.42, 0.00 -75, -32.28, -262.70, 0.00 -76, -32.88, -240.86, 0.00 -77, -33.29, -220.29, 0.00 -78, -33.20, -199.09, 0.00 -79, -33.11, -180.21, 0.00 -80, -32.42, -160.41, 0.00 -81, -32.04, -144.14, 0.00 -82, -30.93, -126.55, 0.00 -83, -30.26, -112.49, 0.00 -84, -28.67, -96.70, 0.00 -85, -26.72, -81.57, 0.00 -86, -23.77, -65.45, 0.00 -87, -19.10, -47.23, 0.00 -88, -13.13, -28.97, 0.00 -89, -4.87, -9.52, 0.00 -90, 5.12, 8.78, 0.00 -91, 16.89, 25.11, 0.00 -92, 29.83, 37.77, 0.00 -93, 43.55, 45.92, 0.00 -94, 58.79, 49.94, 0.00 -95, 71.13, 46.33, 0.00 -96, 80.33, 36.82, 0.00 -97, 91.25, 24.70, 0.00 -98, 89.60, 12.02, 0.00 -99, 52.39, -0.00, 0.00 -100, 121.03, -16.24, 0.00 -101, 165.10, -44.68, 0.00 -102, 144.64, -66.30, 0.00 -103, 117.07, -76.25, 0.00 -104, 94.88, -80.60, 0.00 -105, 66.02, -69.62, 0.00 -106, 47.83, -60.58, 0.00 -107, 28.23, -41.96, 0.00 -108, 15.85, -27.20, 0.00 -109, 3.13, -6.12, 0.00 -110, -5.71, 12.59, 0.00 -111, -14.69, 36.31, 0.00 -112, -21.31, 58.67, 0.00 -113, -28.22, 86.14, 0.00 -114, -33.21, 112.01, 0.00 -115, -38.21, 142.05, 0.00 -116, -41.32, 169.06, 0.00 -117, -44.19, 198.82, 0.00 -118, -45.47, 224.95, 0.00 -119, -46.71, 254.22, 0.00 -120, -46.59, 279.39, 0.00 -121, -46.26, 306.19, 0.00 -122, -45.07, 330.17, 0.00 -123, -43.64, 355.12, 0.00 -124, -41.58, 377.56, 0.00 -125, -39.33, 400.75, 0.00 -126, -36.72, 422.86, 0.00 -127, -33.87, 444.82, 0.00 -128, -30.76, 465.82, 0.00 -129, -27.43, 486.02, 0.00 -130, -23.94, 505.83, 0.00 -131, -20.30, 524.79, 0.00 -132, -16.54, 543.40, 0.00 -133, -12.68, 561.07, 0.00 -134, -8.75, 578.23, 0.00 -135, -4.76, 594.19, 0.00 -136, -0.75, 610.20, 0.00 -137, 3.27, 624.70, 0.00 -138, 7.28, 638.58, 0.00 -139, 11.25, 651.52, 0.00 -140, 15.19, 664.02, 0.00 -141, 19.03, 674.72, 0.00 -142, 22.82, 685.17, 0.00 -143, 26.50, 694.58, 0.00 -144, 30.07, 703.02, 0.00 -145, 33.49, 710.00, 0.00 -146, 36.80, 716.33, 0.00 -147, 39.96, 721.72, 0.00 -148, 42.96, 725.88, 0.00 -149, 45.80, 729.01, 0.00 -150, 48.40, 730.04, 0.00 -151, 50.69, 728.06, 0.00 -152, 53.19, 730.63, 0.00 -153, 53.67, 707.80, 0.00 -154, 40.89, 519.28, 0.00 -155, -10.41, -127.67, 0.00 -156, -9.92, -117.86, 0.00 -157, -0.17, -2.01, 0.00 -158, 0.20, 2.24, 0.00 -159, -0.45, -4.95, 0.00 -160, 0.03, 0.33, 0.00 -161, 0.43, 4.50, 0.00 -162, 0.61, 6.20, 0.00 -163, 0.63, 6.25, 0.00 -164, 0.56, 5.43, 0.00 -165, 0.40, 3.85, 0.00 -166, 0.18, 1.72, 0.00 -167, -0.12, -1.09, 0.00 -168, -0.49, -4.46, 0.00 -169, -0.97, -8.64, 0.00 -170, -1.49, -13.10, 0.00 -171, -2.10, -18.13, 0.00 -172, -2.71, -23.05, 0.00 -173, -3.37, -28.27, 0.00 -174, -4.02, -33.22, 0.00 -175, -4.70, -38.39, 0.00 -176, -5.37, -43.29, 0.00 -177, -6.05, -48.13, 0.00 -178, -6.69, -52.54, 0.00 -179, -7.31, -56.81, 0.00 -180, -7.88, -60.56, 0.00 -181, -8.40, -63.82, 0.00 -182, -8.85, -66.57, 0.00 -183, -9.18, -68.42, 0.00 -184, -9.47, -69.89, 0.00 -185, -9.62, -70.38, 0.00 -186, -9.71, -70.43, 0.00 -187, -9.68, -69.72, 0.00 -188, -9.58, -68.53, 0.00 -189, -9.42, -66.89, 0.00 -190, -9.15, -64.59, 0.00 -191, -8.81, -61.89, 0.00 -192, -8.34, -58.33, 0.00 -193, -7.78, -54.16, 0.00 -194, -7.09, -49.18, 0.00 -195, -6.25, -43.20, 0.00 -196, -5.26, -36.32, 0.00 -197, -4.13, -28.49, 0.00 -198, -2.86, -19.72, 0.00 +199, -0.99, -0.00, 0.00 +0, -2.91, 20.03, 0.00 +1, -4.14, 28.50, 0.00 +2, -5.20, 35.86, 0.00 +3, -6.14, 42.47, 0.00 +4, -6.97, 48.37, 0.00 +5, -7.67, 53.42, 0.00 +6, -8.25, 57.65, 0.00 +7, -8.70, 61.14, 0.00 +8, -9.00, 63.59, 0.00 +9, -9.19, 65.26, 0.00 +10, -9.25, 66.17, 0.00 +11, -9.21, 66.36, 0.00 +12, -9.11, 66.11, 0.00 +13, -8.91, 65.21, 0.00 +14, -8.68, 64.05, 0.00 +15, -8.38, 62.39, 0.00 +16, -8.03, 60.38, 0.00 +17, -7.63, 58.02, 0.00 +18, -7.15, 54.94, 0.00 +19, -6.64, 51.57, 0.00 +20, -6.03, 47.40, 0.00 +21, -5.39, 42.89, 0.00 +22, -4.70, 37.86, 0.00 +23, -3.97, 32.41, 0.00 +24, -3.25, 26.89, 0.00 +25, -2.49, 20.93, 0.00 +26, -1.80, 15.28, 0.00 +27, -1.07, 9.26, 0.00 +28, -0.43, 3.74, 0.00 +29, 0.22, -2.00, 0.00 +30, 0.74, -6.72, 0.00 +31, 1.23, -11.32, 0.00 +32, 1.61, -15.15, 0.00 +33, 1.99, -19.02, 0.00 +34, 2.29, -22.32, 0.00 +35, 2.64, -26.26, 0.00 +36, 2.79, -28.35, 0.00 +37, 2.88, -29.98, 0.00 +38, 2.58, -27.41, 0.00 +39, 2.40, -26.23, 0.00 +40, 2.14, -23.97, 0.00 +41, 2.88, -33.18, 0.00 +42, 0.07, -0.88, 0.00 +43, -7.69, 94.38, 0.00 +44, -15.57, 197.73, 0.00 +45, -17.94, 236.64, 0.00 +46, 4.70, -64.62, 0.00 +47, 28.52, -409.61, 0.00 +48, 30.19, -455.33, 0.00 +49, 29.54, -470.12, 0.00 +50, 27.35, -462.05, 0.00 +51, 25.31, -457.07, 0.00 +52, 23.08, -449.24, 0.00 +53, 21.94, -465.13, 0.00 +54, 21.48, -502.18, 0.00 +55, 22.02, -577.36, 0.00 +56, 20.47, -614.69, 0.00 +57, 16.29, -577.49, 0.00 +58, 12.59, -550.58, 0.00 +59, 9.32, -539.47, 0.00 +60, 6.05, -531.03, 0.00 +61, 2.73, -521.29, 0.00 +62, -0.63, -508.34, 0.00 +63, -3.96, -493.40, 0.00 +64, -7.24, -478.67, 0.00 +65, -10.46, -462.56, 0.00 +66, -13.57, -445.82, 0.00 +67, -16.54, -427.72, 0.00 +68, -19.38, -409.38, 0.00 +69, -22.00, -389.70, 0.00 +70, -24.45, -370.21, 0.00 +71, -26.59, -349.18, 0.00 +72, -28.53, -328.51, 0.00 +73, -30.06, -306.33, 0.00 +74, -31.33, -284.50, 0.00 +75, -32.15, -261.58, 0.00 +76, -32.58, -238.65, 0.00 +77, -32.48, -214.97, 0.00 +78, -32.01, -191.95, 0.00 +79, -31.12, -169.41, 0.00 +80, -29.76, -147.23, 0.00 +81, -28.42, -127.89, 0.00 +82, -26.46, -108.25, 0.00 +83, -25.23, -93.79, 0.00 +84, -23.18, -78.17, 0.00 +85, -22.06, -67.33, 0.00 +86, -20.11, -55.36, 0.00 +87, -17.52, -43.33, 0.00 +88, -14.90, -32.88, 0.00 +89, -7.49, -14.65, 0.00 +90, -1.78, -3.05, 0.00 +91, 11.94, 17.75, 0.00 +92, 21.51, 27.24, 0.00 +93, 37.02, 39.03, 0.00 +94, 47.69, 40.51, 0.00 +95, 60.96, 39.71, 0.00 +96, 70.98, 32.53, 0.00 +97, 78.52, 21.25, 0.00 +98, 83.66, 11.23, 0.00 +99, 49.17, -0.00, 0.00 +100, 124.12, -16.66, 0.00 +101, 184.07, -49.82, 0.00 +102, 164.48, -75.39, 0.00 +103, 133.69, -87.08, 0.00 +104, 109.43, -92.96, 0.00 +105, 76.86, -81.04, 0.00 +106, 57.12, -72.34, 0.00 +107, 33.93, -50.44, 0.00 +108, 19.70, -33.80, 0.00 +109, 3.84, -7.50, 0.00 +110, -6.08, 13.41, 0.00 +111, -17.27, 42.71, 0.00 +112, -24.34, 67.01, 0.00 +113, -32.29, 98.57, 0.00 +114, -37.23, 125.57, 0.00 +115, -42.27, 157.13, 0.00 +116, -45.00, 184.10, 0.00 +117, -47.22, 212.46, 0.00 +118, -48.02, 237.54, 0.00 +119, -48.54, 264.19, 0.00 +120, -48.03, 288.01, 0.00 +121, -47.18, 312.25, 0.00 +122, -45.73, 334.99, 0.00 +123, -43.96, 357.74, 0.00 +124, -41.77, 379.29, 0.00 +125, -39.34, 400.82, 0.00 +126, -36.68, 422.40, 0.00 +127, -33.75, 443.26, 0.00 +128, -30.64, 464.03, 0.00 +129, -27.30, 483.59, 0.00 +130, -23.83, 503.39, 0.00 +131, -20.19, 522.07, 0.00 +132, -16.46, 540.78, 0.00 +133, -12.62, 558.35, 0.00 +134, -8.71, 575.65, 0.00 +135, -4.74, 591.58, 0.00 +136, -0.75, 607.73, 0.00 +137, 3.26, 622.30, 0.00 +138, 7.25, 636.34, 0.00 +139, 11.22, 649.44, 0.00 +140, 15.14, 662.10, 0.00 +141, 18.99, 673.05, 0.00 +142, 22.76, 683.64, 0.00 +143, 26.45, 693.25, 0.00 +144, 30.01, 701.79, 0.00 +145, 33.44, 708.94, 0.00 +146, 36.75, 715.33, 0.00 +147, 39.92, 720.85, 0.00 +148, 42.91, 725.05, 0.00 +149, 45.76, 728.29, 0.00 +150, 48.36, 729.34, 0.00 +151, 50.65, 727.41, 0.00 +152, 53.13, 729.86, 0.00 +153, 53.65, 707.51, 0.00 +154, 41.35, 525.16, 0.00 +155, -10.50, -128.78, 0.00 +156, -10.40, -123.56, 0.00 +157, -0.38, -4.42, 0.00 +158, 0.07, 0.74, 0.00 +159, -0.56, -6.11, 0.00 +160, -0.06, -0.63, 0.00 +161, 0.35, 3.66, 0.00 +162, 0.54, 5.47, 0.00 +163, 0.57, 5.63, 0.00 +164, 0.50, 4.88, 0.00 +165, 0.34, 3.29, 0.00 +166, 0.13, 1.18, 0.00 +167, -0.18, -1.65, 0.00 +168, -0.56, -5.04, 0.00 +169, -1.04, -9.24, 0.00 +170, -1.56, -13.68, 0.00 +171, -2.17, -18.73, 0.00 +172, -2.78, -23.65, 0.00 +173, -3.45, -28.94, 0.00 +174, -4.09, -33.84, 0.00 +175, -4.78, -39.02, 0.00 +176, -5.44, -43.83, 0.00 +177, -6.12, -48.67, 0.00 +178, -6.74, -52.99, 0.00 +179, -7.38, -57.30, 0.00 +180, -7.94, -61.00, 0.00 +181, -8.47, -64.36, 0.00 +182, -8.92, -67.13, 0.00 +183, -9.27, -69.09, 0.00 +184, -9.56, -70.59, 0.00 +185, -9.71, -71.08, 0.00 +186, -9.80, -71.10, 0.00 +187, -9.76, -70.26, 0.00 +188, -9.64, -68.96, 0.00 +189, -9.45, -67.16, 0.00 +190, -9.17, -64.76, 0.00 +191, -8.82, -61.96, 0.00 +192, -8.34, -58.34, 0.00 +193, -7.78, -54.14, 0.00 +194, -7.08, -49.13, 0.00 +195, -6.24, -43.17, 0.00 +196, -5.26, -36.30, 0.00 +197, -4.14, -28.53, 0.00 +198, -2.87, -19.79, 0.00